Как получить проекцию из gml-файла

Question

Я пытаюсь получить проекцию из gml-файла. Это верхняя часть файла:

<?xml version="1.0" encoding="UTF-8"?> <eop:Mask xmlns:eop="http://www.opengis.net/eop/2.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xmlns:gml="http://www.opengis.net/gml/3.2" gml:id="S2A_OPER_MSK_CLOUDS_SGS__20160914T145755_A006426_T31UCT_B00_MSIL1C"> <gml:name>MSK_CLOUDS pixels mask from data-strip S2A_OPER_MSK_CLOUDS_SGS__20160914T145755_A006426_T31UCT_B00_MSIL1C</gml:name> <gml:boundedBy> 
    <gml:Envelope srsName="urn:ogc:def:crs:EPSG:8.7:32631"> 
     <gml:lowerCorner>300000 5690220</gml:lowerCorner> 
     <gml:upperCorner>368340 5777580</gml:upperCorner> 
    </gml:Envelope> </gml:boundedBy> <eop:maskMembers> 
    <eop:MaskFeature gml:id="OPAQUE.0"> 
     <eop:maskType codeSpace="urn:gs2:S2PDGS:maskType">OPAQUE</eop:maskType> 
     <eop:extentOf> 
     <gml:Polygon gml:id="OPAQUE.0_Polygon"> 
      <gml:exterior> 
      <gml:LinearRing> 
       <gml:posList srsDimension="2">320340 5776020 320520 5776020 320520 5775960 320700 5775960 320700 5775900 320760 5775900 320760 5775840 320820 5775840 320820 5775660 320760 5775660 320760 5775600 320700 5775600 320700 5775540 320340 5775540 320340 5775600 320280 5775600 320280 5775660 320220 5775660 320220 5775900 320280 5775900 320280 5775960 320340 5775960 320340 5776020</gml:posList> 
      </gml:LinearRing> 
      </gml:exterior> 
     </gml:Polygon> 
     </eop:extentOf> 
    </eop:MaskFeature> 
... 

Я попытался с помощью кода из https://pcjericks.github.io/py-gdalogr-cookbook/projection.html:

from osgeo import ogr, osr 
driver = ogr.GetDriverByName('ESRI Shapefile') 
dataset = driver.Open(r'c:\data\yourshpfile.shp') 

# from Layer 
layer = dataset.GetLayer() 
spatialRef = layer.GetSpatialRef() 
# from Geometry 
feature = layer.GetNextFeature() 
geom = feature.GetGeometryRef() 
spatialRef = geom.GetSpatialReference() 

Но обе версии spatialRef нет.

Из файла видно, что проекция указана в ограничивающей рамке (в самом конце первой строки она говорит, а затем конверт с кодом EPSG в строке 2.) (Это не так, t say 'crs' или 'EPSG' где-либо еще в файле).

Может ли кто-нибудь сказать мне, как я обращаюсь к проекционной информации?

Могу ли я как-нибудь добраться до ограничительной рамки, а затем получить проекцию?

Answer 1

Обычно вы можете найти информацию о проекции в файле JPEG 2000, предоставленном гранулами.

Использование GDAL это:

gdalinfo * .jp2

дает:

PROJCS["WGS 84/UTM zone 21N", 
GEOGCS["WGS 84", 
    DATUM["WGS_1984", 
     SPHEROID["WGS 84",6378137,298.257223563, 
      AUTHORITY["EPSG","7030"]], 
     AUTHORITY["EPSG","6326"]], 
    PRIMEM["Greenwich",0, 
     AUTHORITY["EPSG","8901"]], 
    UNIT["degree",0.0174532925199433, 
     AUTHORITY["EPSG","9122"]], 
    AXIS["Latitude",NORTH], 
    AXIS["Longitude",EAST], 
    AUTHORITY["EPSG","4326"]], 
PROJECTION["Transverse_Mercator"], 
PARAMETER["latitude_of_origin",0], 
PARAMETER["central_meridian",-57], 
PARAMETER["scale_factor",0.9996], 
PARAMETER["false_easting",500000], 
PARAMETER["false_northing",0], 
.... 

Это все у меня есть