I am trying to convert Sentinel-2 cloud masks from GML format to ESRI shapefiles. I have managed to use ogr2ogr in a python script as described in the answer to Using ogr2ogr to convert GML to shapefile in Python? to do this. I downloaded ogr2ogr from http://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/ogr2ogr.py and my Python call to that looks like this:
ogr2ogr.main(["", "-f", "ESRI Shapefile", 'output_shapefile.shp', 'cloud_mask.gml'])
However when I load the new shapefiles into ArcMap, it tells me that they are "missing spatial reference information". I've opened the GMLs in a text editor and I can see the EPSG code.
GDAL isn't creating a .prj file at all. If I try with '-lco SPATIAL_INDEX=YES' it creates a .qix file (which it doesn't otherwise), but the problem remains. The output of ogrinfo -al -so for the GML file is as follows:
INFO: Open of 'cloud_mask.gml' using driver `GML' successful.
Metadata: NAME=MSK_CLOUDS pixels mask from data-strip S2A_OPER_MSK_CLOUDS_SGS__20190403T 205739_A019743_T13RGP_B00_MSIL1C
Layer name: MaskFeature
Feature Count: 4
Extent: (806760.000000, 3335340.000000) - (809760.000000, 3381900.000000)
Layer SRS WKT: (unknown)
gml_id: String (0.0) NOT NULL
maskType: String (6.0)
Is there some way I can get ogr2ogr to project the shapefiles that it creates correctly, or another way I can project them afterwards using some extra Python code?