在geopandas 文档中它说
A
GeoDataFrame
还可能包含其他具有几何(形状)对象的列,但一次只能有一个列是活动几何。要更改哪一列是活动几何列,请使用该set_geometry
方法。
如果目标是将这些不同列中的几何数据灵活地重新投影到一个或多个其他坐标参考系统,我想知道如何使用这样的 GeoDataFrame 。这是我尝试过的。
第一次尝试
import geopandas as gpd
from shapely.geometry import Point
crs_lonlat = 'epsg:4326' #geometries entered in this crs (lon, lat in degrees)
crs_new = 'epsg:3395' #geometries needed in (among others) this crs
gdf = gpd.GeoDataFrame(crs=crs_lonlat)
gdf['geom1'] = [Point(9,53), Point(9,54)]
gdf['geom2'] = [Point(8,63), Point(8,64)]
#Working: setting geometry and reprojecting for first time.
gdf = gdf.set_geometry('geom1')
gdf = gdf.to_crs(crs_new) #geom1 is reprojected to crs_new, geom2 still in crs_lonlat
gdf
Out:
geom1 geom2
0 POINT (1001875.417 6948849.385) POINT (8 63)
1 POINT (1001875.417 7135562.568) POINT (8 64)
gdf.crs
Out: 'epsg:3395'
到目前为止,一切都很好。如果我想设置geom2
为几何列,事情就会出轨,并重新投影该列:
#Not working: setting geometry and reprojecting for second time.
gdf = gdf.set_geometry('geom2') #still in crs_lonlat...
gdf.crs #...but this still says crs_new...
Out: 'epsg:3395'
gdf = gdf.to_crs(crs_new) #...so this doesn't do anything! (geom2 unchanged)
gdf
Out:
geom1 geom2
0 POINT (1001875.417 6948849.385) POINT (8.00000 63.00000)
1 POINT (1001875.417 7135562.568) POINT (8.00000 64.00000)
好的,所以,显然,在更改用作几何图形的列时,.crs
属性并没有重置为其原始值 - 看起来,crs 没有为各个列存储。gdf
如果是这种情况,我看到对此数据框使用重新投影的唯一方法是回溯:开始->选择列作为几何图形->将gdf重新投影到crs_new->使用/可视化/...->将 gdf 重新投影回 crs_lonlat --> goto start。如果我想在一个图中可视化两列,这将不可用。
第二次尝试
我的第二次尝试是,crs
通过将上面脚本中的相应行更改为:
gdf = gpd.GeoDataFrame()
gdf['geom1'] = gpd.GeoSeries([Point(9,53), Point(9,54)], crs=crs_lonlat)
gdf['geom2'] = gpd.GeoSeries([Point(8,63), Point(8,64)], crs=crs_lonlat)
然而,很快就很清楚,虽然初始化为 a GeoSeries
,但这些列是 normal pandas
Series
,并且没有.crs
相同的属性GeoSeries
:
gdf['geom1'].crs
AttributeError: 'Series' object has no attribute 'crs'
s = gpd.GeoSeries([Point(9,53), Point(9,54)], crs=crs_lonlat)
s.crs
Out: 'epsg:4326'
我在这里缺少什么吗?
唯一的解决方案是事先决定“最终”crs - 并在添加列之前进行所有重新投影?像这样...
gdf = gpd.GeoDataFrame(crs=crs_new)
gdf['geom1'] = gpd.GeoSeries([Point(9,53), Point(9,54)], crs=crs_lonlat).to_crs(crs_new)
gdf['geom2'] = gpd.GeoSeries([Point(8,63), Point(8,64)], crs=crs_lonlat).to_crs(crs_new)
#no more reprojecting done/necessary/possible! :/
...然后,当需要另一个 crs 时,gdf
从头开始重建整个?这不可能是打算使用的方式。