使用proj4参数构建cartopy CRS

3

我下载了一个MODIS .hdf文件。使用xarray加载它后,它给了我一个属性ds.Proj4String == ' +a = 6378137.0 +b=6356752.3142451793 +no_defs +proj=latlong\n'

我该如何使用这个字符串将原始坐标转换为例如经纬度?

以下是数据的样子:

<xarray.Dataset>
Dimensions:                  (XDim:mod06: 1503, YDim:mod06: 833)
Dimensions without coordinates: XDim:mod06, YDim:mod06
Data variables:
    Cloud_Optical_Thickness  (YDim:mod06, XDim:mod06) float32 ...
Attributes:
    HDFEOSVersion:     HDFEOS_V2.19
    StructMetadata.0:  GROUP=SwathStructure\nEND_GROUP=SwathStructure\nGROUP=...
    CoreMetadata:      \nGROUP = INVENTORYMETADATA\n  GROUPTYPE = MASTERGROUP...
    ArchiveMetadata:   GROUP = ARCHIVEDMETADATA\n  GROUPTYPE = MASTERGROUP\n\...
    Proj4String:        +a=6378137.0 +b=6356752.3142451793 +no_defs +proj=lat...

这份文档告诉我需要传递proj4_params作为键值对,所以我这么做:

ss = ds.Proj4String.split()
proj = {}
for s in ss:
    k = s.split('=')
    if len(k)==2:
        proj[k[0][1:]] = k[1]
print(proj)

并且到达 {'a': '6378137.0', 'b': '6356752.3142451793', 'proj': 'latlong'},但是ccrs.CRS(proj)仍然报错,指出缺少属性globe

1个回答

6
创建CartoPy投影仅从proj.4字符串不受支持,参见this PR。您正在使用的CRS类本身是不足以创建完全功能的CartoPy投影的,而且实际上是为其他投影设计的基类。
就您遇到的错误而言,无论传递了什么参数,'globe'都是必需的参数。您可以使用以下方式创建'Globe'实例:
globe = ccrs.Globe(semimajor_axis=proj.pop('a'), semiminor_axis=proj.pop('b'))

谢谢,我能看出这个问题比我最初想象的要复杂。如果有可能的话,我会保持问题开放以便回答。 - doppler

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接