为什么在某些情况下变换一个匀称的多边形不起作用?

2024-09-28 05:20:14 发布

您现在位置:Python中文网/ 问答频道 /正文

我正在尝试使用shapely计算地理坐标多边形的大小,这似乎需要转换成合适的投影,以产生平方米的结果。我在网上找到了几个例子,但我无法让它在我的例子多边形中工作

因此,我尝试使用与我找到的代码片段相同的多边形示例,我注意到它适用于某些整体,而不适用于其他整体。为了重现结果,下面是最简单的示例代码:

import json
import pyproj

from shapely.ops import transform
from shapely.geometry import Polygon, mapping
from functools import partial

coords1 = [(-97.59238135821987, 43.47456565304017),
           (-97.59244690469288, 43.47962399877412),
           (-97.59191951546768, 43.47962728271748),
           (-97.59185396090983, 43.47456565304017),
           (-97.59238135821987, 43.47456565304017)]

coords1 = reversed(coords1) # Not sure if important, but https://geojsonlint.com says it's wrong handedness
                            # Doesn't seem to affect the error message though

coords2 = [(13.65374516425911, 52.38533382814119),
           (13.65239769133293, 52.38675829106993),
           (13.64970274383571, 52.38675829106993),
           (13.64835527090953, 52.38533382814119),
           (13.64970274383571, 52.38390931824483),
           (13.65239769133293, 52.38390931824483),
           (13.65374516425911, 52.38533382814119)]

coords = coords1 # DOES NOT WORK
#coords = coords2 # WORKS

polygon = Polygon(coords)

# Print GeoJON to check on https://geojsonlint.com
print(json.dumps(mapping(polygon)))

projection = partial(pyproj.transform, 
                     pyproj.Proj('epsg:4326'), 
                     pyproj.Proj('esri:54009'))

transform(projection, polygon)

coords1coords2都是从假定有效的代码片段复制而来的。然而,只有coords2对我有效。我使用了https://geojsonlint.com来查看这两个多边形之间是否存在差异,并且似乎多边形的利手/方向无效。我不知道shapely是否在乎,但是颠倒顺序-https://geojsonlint.com说它是有效的GeoJSON,然后它在地图上显示多边形-,并不会改变错误

因此,它与coords2一起工作,但是当我使用coords1时,我得到以下错误:

~/env/anaconda3/envs/py36/lib/python3.6/site-packages/shapely/geometry/base.py in _repr_svg_(self)
    398             if xmin == xmax and ymin == ymax:
    399                 # This is a point; buffer using an arbitrary size
--> 400                 xmin, ymin, xmax, ymax = self.buffer(1).bounds
    401             else:
    402                 # Expand bounds by a fraction of the data ranges

ValueError: not enough values to unpack (expected 4, got 0)

我假设coords1(以及我自己数据中的示例多边形)有一些不同的地方导致了这个问题,但我无法判断与coords2相比有什么不同

简言之,coords1coords2之间有什么区别,一个工作,另一个不工作

更新:我通过在投影的定义中添加always_xy=True来实现它。与shapely提供的较新语法一起,避免了partial,工作代码段如下所示:

    project = pyproj.Transformer.from_proj(
        pyproj.Proj('epsg:4326'), # source coordinate system
        pyproj.Proj('epsg:3857'),
        always_xy=True
    ) # destination coordinate system
 
    transform(project.transform, polygon)

老实说,即使在阅读了文档之后,我也不知道always_xy在做什么。因此,我不想提供答案


Tags: 代码fromhttpsimportcom示例transform多边形
1条回答
网友
1楼 · 发布于 2024-09-28 05:20:14

我认为您做得很好,只是相反的方法不会创建新的数据集

尝试使用此函数创建反向顺序列表:

def rev_slice(mylist):
    '''
    return a revered list
    mylist: is a list
    '''
    a = mylist[::-1]
    return a

按如下方式执行函数:

coords = rev_slice(coords1)

相关问题 更多 >

    热门问题