Geopandas：缓冲区操作似乎忽略了CRS的度量单位

Question

我的目标是从现有数据框中的两列坐标创建一个地理数据框，获取1677个地理点，并在每个点周围添加一个缓冲圈，然后将生成的多边形合并成一个多多边形。我一直围绕着轴的地方是geopandas的.buffer（）部分，它似乎没有使用我选择的CRS的度量单位。

In  []: ven_coords

Out []:     VenLat      VenLon
       0    42.34768    -71.085359
       1    42.349014   -71.081096
       2    42.347627   -71.081685
       3    42.348718   -71.077984
       4    42.34896    -71.081467
     ...         ...           ...
    1672    42.308962   -71.073516
    1673    42.313169   -71.089027
    1674    42.309717   -71.08247
    1675    42.356336   -71.074386
    1676    42.313005   -71.089887
    1677 rows × 2 columns

In  []: ven_coords_gdf = geopandas.GeoDataFrame(ven_coords, 
                                        geometry=geopandas.points_from_xy(ven_coords.VenLon, ven_coords.VenLat))
        ven_coords_gdf

Out []: VenLat  VenLon  geometry
       0    42.34768    -71.085359  POINT (-71.08536 42.34768)
       1    42.349014   -71.081096  POINT (-71.08110 42.34901)
       2    42.347627   -71.081685  POINT (-71.08168 42.34763)
       3    42.348718   -71.077984  POINT (-71.07798 42.34872)
       4    42.34896    -71.081467  POINT (-71.08147 42.34896)
     ...         ...           ...                        ...
    1672    42.308962   -71.073516  POINT (-71.07352 42.30896)
    1673    42.313169   -71.089027  POINT (-71.08903 42.31317)
    1674    42.309717   -71.08247   POINT (-71.08247 42.30972)
    1675    42.356336   -71.074386  POINT (-71.07439 42.35634)
    1676    42.313005   -71.089887  POINT (-71.08989 42.31300)
    1677 rows × 3 columns

到目前为止还不错，让我们看看我得到了什么：

In  []: print('Type:', type(ven_coords_gdf), "/ current CRS is:",ven_coords_gdf.crs)

Out []: Type: <class 'geopandas.geodataframe.GeoDataFrame'> / current CRS is: None

它没有CRS，所以我给它分配了一个与我正在做的工作相关的：

In  []: ven_coords_gdf.crs = ("epsg:2249")
        print('Type:', type(ven_coords_gdf), "/ current CRS is:",ven_coords_gdf.crs)

Out []: Type: <class 'geopandas.geodataframe.GeoDataFrame'> / current CRS is: epsg:2249

它似乎已经“添加”了我添加的CRS，并且只是为了重复检查，让我们看一下CRS的细节问题：

In  []: CRS.from_epsg(2249)

Out []: <Projected CRS: EPSG:2249>
        Name: NAD83 / Massachusetts Mainland (ftUS)
        Axis Info [cartesian]:
        - X[east]: Easting (US survey foot)
        - Y[north]: Northing (US survey foot)
        Area of Use:
        - name: United States (USA) - Massachusetts onshore - counties of Barnstable; Berkshire; Bristol; Essex; Franklin; Hampden; Hampshire; Middlesex; Norfolk; Plymouth; Suffolk; Worcester.
        - bounds: (-73.5, 41.46, -69.86, 42.89)
        Coordinate Operation:
        - name: SPCS83 Massachusetts Mainland zone (US Survey feet)
        - method: Lambert Conic Conformal (2SP)
        Datum: North American Datum 1983
        - Ellipsoid: GRS 1980
        - Prime Meridian: Greenwich

2249使用美国调查英尺作为度量单位，所以我将我的缓冲区设置为1000，从我的数据中的每个点获得1000英尺的半径：

In  []: ven_coords_buffer = ven_coords_gdf.geometry.buffer(distance = 1000)
        ven_coords_buffer

Out []: 0       POLYGON ((928.915 42.348, 924.099 -55.669, 909...
        1       POLYGON ((928.919 42.349, 924.104 -55.668, 909...
        2       POLYGON ((928.918 42.348, 924.103 -55.670, 909...
        3       POLYGON ((928.922 42.349, 924.107 -55.668, 909...
        4       POLYGON ((928.919 42.349, 924.103 -55.668, 909...
                                     ...                        
        1672    POLYGON ((928.926 42.309, 924.111 -55.708, 909...
        1673    POLYGON ((928.911 42.313, 924.096 -55.704, 909...
        1674    POLYGON ((928.918 42.310, 924.102 -55.707, 909...
        1675    POLYGON ((928.926 42.356, 924.110 -55.661, 909...
        1676    POLYGON ((928.910 42.313, 924.095 -55.704, 909...
        Length: 1677, dtype: geometry

那些坐标只差一点点。显然，缓冲区本身是一个1000°，而不是1000英尺，形成了一个覆盖整个地球的1677个巨大重叠圆的球体。不完全是我要找的。显然我遗漏了什么，有什么建议吗？

和任何有趣的代码问题一样，我发誓它工作得更早，诚实。在我最终让它输出正确的东西之前，我花了一段时间摸索，然后我关闭它，去吃晚饭，回来重新运行它，然后得到上面的结果。显而易见的推论是，我在前面提到的Futuring around中所做的一些事情是让它工作的关键，一些重复使用的变量或其他什么，但我无法找出上面代码中缺少了什么。

GeoPandas 0.9.0，pyproj 3.0.1

《快乐时光》的截图，我把它放到了地图上

Answer 1

GeoPandas所做的正是人们所期望的。您必须将几何图形重新投影到目标CRS，只需指定它就不会有任何作用。

创建GeoDataFrame时，请确保指定数据所在的CRS。在这种情况下，它是EPSG:4326 aka地理投影，单位为度。

ven_coords_gdf = geopandas.GeoDataFrame(ven_coords, 
                                        geometry=geopandas.points_from_xy(ven_coords.VenLon, ven_coords.VenLat),
                                        crs=4326)

一旦正确设置，您必须使用to_crs将坐标重新投影（转换）到目标CRS。

ven_coords_gdf_projected = ven_coords_gdf.to_crs("epsg:2249")

现在你可以用脚来缓冲。如果要将结果再次存储在4326中，只需使用将其重新投影回crs（4326）。

老实说，我发誓它早些时候起作用了。

我很确定它没有：）。