search

stactools-packages / threedep

stactools package for working with elevation data from the USGS 3DEP program (formerly known as NED)
Other
3 stars 4 forks source link

Possible discrepancy between STAC metadata and COG #9

Closed TomAugspurger closed 2 years ago

TomAugspurger commented 2 years ago

As reported by @aazuspan in https://github.com/gjoseph92/stackstac/issues/147, the proj information of these STAC items (sometimes?) doesn't align with the information reported by gdalinfo.

Looks like there's a pixel size discrepancy causing the issue: 0.00001 according to STAC vs. 0.000277 according to GDAL.

STAC metadata ``` {'type': 'Feature', 'stac_version': '1.0.0', 'id': 'n44w125-1', 'properties': {'gsd': 30, 'datetime': '2013-01-01T00:00:00Z', 'proj:epsg': 5498, 'proj:shape': [3612, 3612], 'end_datetime': '2013-11-01T00:00:00Z', 'proj:transform': [1e-05, 0.0, -125.0016666667, 0.0, -1e-05, 44.00166666666, 0.0, 0.0, 1.0], 'start_datetime': '1999-02-01T00:00:00Z', 'threedep:region': 'n40w130'}, 'geometry': {'type': 'Polygon', 'coordinates': [[[-123.99833333, 42.99833333], [-123.99833333, 44.00166667], [-125.00166667, 44.00166667], [-125.00166667, 42.99833333], [-123.99833333, 42.99833333]]]}, 'links': [{'rel': 'collection', 'href': 'https://planetarycomputer.microsoft.com/api/stac/v1/collections/3dep-seamless', 'type': 'application/json'}, {'rel': 'parent', 'href': 'https://planetarycomputer.microsoft.com/api/stac/v1/collections/3dep-seamless', 'type': 'application/json'}, {'rel': , 'href': 'https://planetarycomputer.microsoft.com/api/stac/v1', 'type': , 'title': 'Microsoft Planetary Computer STAC API'}, {'rel': 'self', 'href': 'https://planetarycomputer.microsoft.com/api/stac/v1/collections/3dep-seamless/items/n44w125-1', 'type': 'application/geo+json'}, {'rel': 'via', 'href': 'https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/USGS_1_n44w125.xml'}, {'rel': 'preview', 'href': 'https://planetarycomputer.microsoft.com/api/data/v1/item/map?collection=3dep-seamless&item=n44w125-1', 'type': 'text/html', 'title': 'Map of item'}], 'assets': {'data': {'href': 'https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/USGS_1_n44w125.tif?st=2022-04-29T23%3A13%3A27Z&se=2022-04-30T23%3A58%3A27Z&sp=rl&sv=2020-06-12&sr=c&skoid=c85c15d6-d1ae-42d4-af60-e2ca0f81359b&sktid=72f988bf-86f1-41af-91ab-2d7cd011db47&skt=2022-04-30T19%3A43%3A28Z&ske=2022-05-07T19%3A43%3A28Z&sks=b&skv=2020-06-12&sig=AgVxx/WNq/thlRHe9wO3ixVxXCWnk5ZxXtMwynWPaVY%3D', 'type': 'image/tiff; application=geotiff; profile=cloud-optimized', 'title': 'USGS 1 arc-second n44w125 1 x 1 degree', 'description': 'This tile of the 3D Elevation Program (3DEP) seamless products is 1 arc-second resolution. 3DEP data serve as the elevation layer of The National Map, and provide basic elevation information for Earth science studies and mapping applications in the United States. Scientists and resource managers use 3DEP data for global change research, hydrologic modeling, resource monitoring, mapping and visualization, and many other applications. 3DEP data compose an elevation dataset that consists of seamless layers and a high resolution layer. Each of these layers consists of the best available raster elevation data of the conterminous United States, Alaska, Hawaii, territorial islands, Mexico and Canada. 3DEP data are updated continually as new data become available. Seamless 3DEP data are derived from diverse source data that are processed to a common coordinate system and unit of vertical measure. These data are distributed in geographic coordinates in units of decimal degrees, and in conformance with the North American Datum of 1983 (NAD 83). All elevation values are in meters and, over the conterminous United States, are referenced to the North American Vertical Datum of 1988 (NAVD 88). The vertical reference will vary in other areas. Seamless 3DEP data are available nationally (except for Alaska) at resolutions of 1 arc-second (approximately 30 meters) and 1/3 arc-second (approximately 10 meters). In most of Alaska, only lower resolution source data are available. As a result, most seamless 3DEP data for Alaska are at 2 arc-second (approximately 60 meters) grid spacing. Part of Alaska is available at the 1 and 1/3 arc-second resolutions from interferometric synthetic aperture radar (ifsar) collections starting in 2010. Plans are in place for collection of statewide ifsar in Alaska. All 3DEP products are public domain.', 'roles': ['data']}, 'gpkg': {'href': 'https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/n44w125.gpkg?st=2022-04-29T23%3A13%3A27Z&se=2022-04-30T23%3A58%3A27Z&sp=rl&sv=2020-06-12&sr=c&skoid=c85c15d6-d1ae-42d4-af60-e2ca0f81359b&sktid=72f988bf-86f1-41af-91ab-2d7cd011db47&skt=2022-04-30T19%3A43%3A28Z&ske=2022-05-07T19%3A43%3A28Z&sks=b&skv=2020-06-12&sig=AgVxx/WNq/thlRHe9wO3ixVxXCWnk5ZxXtMwynWPaVY%3D', 'type': 'application/geopackage+sqlite3', 'roles': ['metadata']}, 'metadata': {'href': 'https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/USGS_1_n44w125.xml?st=2022-04-29T23%3A13%3A27Z&se=2022-04-30T23%3A58%3A27Z&sp=rl&sv=2020-06-12&sr=c&skoid=c85c15d6-d1ae-42d4-af60-e2ca0f81359b&sktid=72f988bf-86f1-41af-91ab-2d7cd011db47&skt=2022-04-30T19%3A43%3A28Z&ske=2022-05-07T19%3A43%3A28Z&sks=b&skv=2020-06-12&sig=AgVxx/WNq/thlRHe9wO3ixVxXCWnk5ZxXtMwynWPaVY%3D', 'type': 'application/xml', 'roles': ['metadata']}, 'thumbnail': {'href': 'https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/USGS_1_n44w125.jpg?st=2022-04-29T23%3A13%3A27Z&se=2022-04-30T23%3A58%3A27Z&sp=rl&sv=2020-06-12&sr=c&skoid=c85c15d6-d1ae-42d4-af60-e2ca0f81359b&sktid=72f988bf-86f1-41af-91ab-2d7cd011db47&skt=2022-04-30T19%3A43%3A28Z&ske=2022-05-07T19%3A43%3A28Z&sks=b&skv=2020-06-12&sig=AgVxx/WNq/thlRHe9wO3ixVxXCWnk5ZxXtMwynWPaVY%3D', 'type': 'image/jpeg', 'roles': ['thumbnail']}, 'tilejson': {'href': 'https://planetarycomputer.microsoft.com/api/data/v1/item/tilejson.json?collection=3dep-seamless&item=n44w125-1&assets=data&colormap_name=terrain&rescale=-1000,4000', 'type': 'application/json', 'title': 'TileJSON with default rendering', 'roles': ['tiles']}, 'rendered_preview': {'href': 'https://planetarycomputer.microsoft.com/api/data/v1/item/preview.png?collection=3dep-seamless&item=n44w125-1&assets=data&colormap_name=terrain&rescale=-1000,4000', 'type': 'image/png', 'title': 'Rendered preview', 'rel': 'preview', 'roles': ['overview']}}, 'bbox': [-125.0016666667, 42.99833333333, -123.9983333334, 44.00166666666], 'stac_extensions': ['https://stac-extensions.github.io/projection/v1.0.0/schema.json'], 'collection': '3dep-seamless'} ```
gdalinfo ``` Driver: GTiff/GeoTIFF Files: USGS_1_n44w125.tif Size is 3612, 3612 Coordinate System is: GEOGCRS["NAD83", DATUM["North American Datum 1983", ELLIPSOID["GRS 1980",6378137,298.257222101004, LENGTHUNIT["metre",1]]], PRIMEM["Greenwich",0, ANGLEUNIT["degree",0.0174532925199433]], CS[ellipsoidal,2], AXIS["geodetic latitude (Lat)",north, ORDER[1], ANGLEUNIT["degree",0.0174532925199433]], AXIS["geodetic longitude (Lon)",east, ORDER[2], ANGLEUNIT["degree",0.0174532925199433]], ID["EPSG",4269]] Data axis to CRS axis mapping: 2,1 Origin = (-125.001666666666665,44.001666666666665) Pixel Size = (0.000277777777778,-0.000277777777778) Metadata: AREA_OR_POINT=Area Image Structure Metadata: COMPRESSION=LZW INTERLEAVE=BAND PREDICTOR=3 Corner Coordinates: Upper Left (-125.0016667, 44.0016667) (125d 0' 6.00"W, 44d 0' 6.00"N) Lower Left (-125.0016667, 42.9983333) (125d 0' 6.00"W, 42d59'54.00"N) Upper Right (-123.9983333, 44.0016667) (123d59'54.00"W, 44d 0' 6.00"N) Lower Right (-123.9983333, 42.9983333) (123d59'54.00"W, 42d59'54.00"N) Center (-124.5000000, 43.5000000) (124d30' 0.00"W, 43d30' 0.00"N) Band 1 Block=256x256 Type=Float32, ColorInterp=Gray Description = Layer_1 Min=-5.715 Max=517.677 Minimum=-5.715, Maximum=517.677, Mean=27.307, StdDev=64.316 NoData Value=-999999 Overviews: 1806x1806, 903x903, 452x452, 226x226, 113x113 Metadata: LAYER_TYPE=athematic STATISTICS_MAXIMUM=517.67651367188 STATISTICS_MEAN=27.306801648412 STATISTICS_MINIMUM=-5.7150769233704 STATISTICS_STDDEV=64.31571761273 STATISTICS_VALID_PERCENT=99.83 ```

Looking at the metadata in https://ai4edataeuwest.blob.core.windows.net/3dep/Elevation/1/TIFF/n44w125/USGS_1_n44w125.xml (needs to be signed) I see the 0.00001.

<spref>
<horizsys>
<geograph>
<latres>0.00001</latres>
<longres>0.00001</longres>
<geogunit>Decimal degrees</geogunit>
</geograph>
<geodetic>
<horizdn>North American Datum of 1983</horizdn>
<ellips>Geodetic Reference System 80</ellips>
<semiaxis>6378137.000000</semiaxis>
<denflat>298.2572221</denflat>
</geodetic>
</horizsys>
<vertdef>
<altsys>
<altdatum>North American Vertical Datum of 1988</altdatum>
<altres>0.001</altres>
<altunits>meters</altunits>
<altenc>Implicit coordinate</altenc>
</altsys>
</vertdef>
</spref>
TomAugspurger commented 2 years ago

cc @gadomski, if you recall anything about this dataset. Do you know why the metadata might differ from the values in the COG?

gadomski commented 2 years ago

I've taken a quick initial look, and there's some discrepancy in the resolution within the same dataset:

I'll have to look over the full holdings and dig into what's up.

gadomski commented 2 years ago

After further review, there's an almost even split between correct and incorrect metadata. @TomAugspurger let's catch up offline (aka at our meeting later) to talk through strategies -- #10 should probably be addressed as well.