Skip to content

Commit

Permalink
use ogr to convert gml more robustly - my solution was a hack that di…
Browse files Browse the repository at this point in the history
…dn't work
  • Loading branch information
GondekNP committed Jan 11, 2024
1 parent 7609e05 commit 24d5993
Show file tree
Hide file tree
Showing 4 changed files with 945 additions and 5 deletions.
1 change: 1 addition & 0 deletions .deployment/prod_environment.yml
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@ dependencies:
- google-cloud-logging
- rioxarray
- stackstac
- ogr
- gdal
- numpy
- geopandas
Expand Down
1 change: 1 addition & 0 deletions .devcontainer/dev_environment.yml
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@ dependencies:
- ipykernel
- rioxarray
- stackstac
- ogr
- contextily
- rich
- rich-with-jupyter
Expand Down
921 changes: 921 additions & 0 deletions sdm_gml.gml

Large diffs are not rendered by default.

27 changes: 22 additions & 5 deletions src/lib/query_soil.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@
from shapely.ops import transform
import xml.etree.ElementTree as ET
import tempfile
import ogr

SDM_ENDPOINT_TABULAR = "https://SDMDataAccess.sc.egov.usda.gov/Tabular/post.rest"
SDM_ENDPOINT_SPATIAL = "https://SDMDataAccess.sc.egov.usda.gov/Spatial/SDMWGS84Geographic.wfs"
Expand Down Expand Up @@ -90,17 +91,33 @@ def sdm_get_esa_mapunitid_poly(geojson):
with tempfile.NamedTemporaryFile(suffix='.gml') as tmp:
tmp.write(response.content)
tmp.seek(0)
gdf = gpd.read_file(tmp.name)

# GML is in lon/lat, so we need to transform to lat/lon
gdf['geometry'] = gdf['geometry'].apply(lambda geom: transform(lambda x, y: (y, x), geom))
# gdf = gpd.read_file(tmp.name)
# # GML is in lon/lat, so we need to transform to lat/lon
# gdf['geometry'] = gdf['geometry'].apply(lambda geom: transform(lambda x, y: (y, x), geom))

gdf = gml_to_gpd(tmp.name)

mapunitpoly_geojson = gdf.to_json()
return mapunitpoly_geojson
else:
elif response.status_code == 400:
print("Error:", response.status_code)
return None

except Exception as e:
print("Error:", str(e))
return None
return None

def gml_to_gpd(gml_file):
gml = ogr.Open(gml_file)
layer = gml.GetLayer(0)

# Get the features and create a list to put geometries
features=[]
for i in range(0,layer.GetFeatureCount()):
feature = layer.GetFeature(i)
geom = feature.GetGeometryRef()
features.append(wkt.loads(geom.ExportToWkt()))

# convert this list into a geopandas dataframe
gdf = gpd.GeoDataFrame(geometry=features)

0 comments on commit 24d5993

Please sign in to comment.