[gdal-dev] Retrieving 3DEP data via GDAL's WCS driver...

jratike80 jukka.rahkonen at maanmittauslaitos.fi
Sun Nov 1 23:36:59 PST 2020


Hi,

By adding "-- debug on" into your command it is possible to get details
about what GDAL is doing. Here are the main parts from the log

gdal_translate
"WCS:https://elevation.nationalmap.gov:443/arcgis/services/3DEPElevation/ImageServer/WCSServer?version=2.0.1&coverage=DEP3Elevation"
test.tif -projwin -107.03 37.28 -107 37.25 -projwin_srs EPSG:4326 -tr 100
100 --debug on

GDAL:
GDALOpen(WCS:https://elevation.nationalmap.gov:443/arcgis/services/3DEPElevation/ImageServer/WCSServer?version=2.0.1&coverage=DEP3Elevation,
this=000002944175BD50) succeeds as WCS.
GDAL: Using GTiff driver
Input file size is 40075015, 20498394
GDAL: Using GTiff driver
0GDAL: GDAL_CACHEMAX = 1632 MB
GDAL: GDALDatasetCopyWholeRaster(): 33*42 swaths, bInterleave=0
WCS: DirectRasterIO(8122982,14330794,3340,4196) -> (33,42) (1 bands)

WCS: Requesting
https://elevation.nationalmap.gov:443/arcgis/services/3DEPElevation/ImageServer/WCSServer?SERVICE=WCS&REQUEST=GetCoverage&VERSION=2.0.1&COVERAGEID=DEP3Elevation&SUBSET=x%28-20037507.0672000013,-11911185.14836644%29&SUBSET=y%284474012.1150791775,4478208.115037268%29&SCALESIZE=x%2833%29,y%2842%29&Format=image/tiff
HTTP:
Fetch(https://elevation.nationalmap.gov:443/arcgis/services/3DEPElevation/ImageServer/WCSServer?SERVICE=WCS&REQUEST=GetCoverage&VERSION=2.0.1&COVERAGEID=DEP3Elevation&SUBSET=x%28-20037507.0672000013,-11911185.14836644%29&SUBSET=y%284474012.1150791775,4478208.115037268%29&SCALESIZE=x%2833%29,y%2842%29&Format=image/tiff)
HTTP: libcurl/7.70.0-DEV OpenSSL/1.1.1g zlib/1.2.3
WCS: GDALOpenResult() on content-type: multipart/related; boundary="wcs";
start="GML-Part"; type="text/xml"
GDAL: GDALOpen(/vsimem/wcs/000002944175BD50/wcsresult.dat,
this=0000029441C0D9A0) succeeds as GTiff.
GDAL: GDALClose(/vsimem/wcs/000002944175BD50/wcsresult.dat,
this=0000029441C0D9A0)
...10...20...30...40...50...60...70...80...90...100 - done.
GDAL: Flushing dirty blocks:
0...10...20...30...40...50...60...70...80...90...100 - done.

Image test.tif looks somehow corrupted from the right side.
The only WCS request that GDAL sends can be found from the log. I removed
url-encoding for clarity

https://elevation.nationalmap.gov:443/arcgis/services/3DEPElevation/ImageServer/WCSServer?SERVICE=WCS&REQUEST=GetCoverage&VERSION=2.0.1&COVERAGEID=DEP3Elevation&SUBSET=x(-20037507.0672000013,-11911185.14836644)&SUBSET=y(4474012.1150791775,4478208.115037268)&SCALESIZE=x(33),y(42)&Format=image/tiff

Unfortunately the server returns the GeoTIFF within a multipart content with
some GML. I do not know how to extract the GeoTIFF from the multipart
package so I could compare the results obtained by using WCS service
directly with curl or a browser. But because GDAL does not do much for the
result, just extract the tiff part and write it into separate GeoTIFF, I
suppose that it is the WCS server that generates the corrupted output.
However, the request that GDAL sends feels odd to me. Compare the projwin in
gdal_translate command with subsets in the GetCoverage

 -projwin -107.03 37.28 -107 37.25
SUBSET=x(-20037507.0672000013,-11911185.14836644)
SUBSET=y(4474012.1150791775,4478208.115037268)

Subract the numbers and notice that subset x is 8126321.919 meters while
subset y is 4196 meters. 
I can see number 4196 also in this line of the log
WCS: DirectRasterIO(8122982,14330794,3340,4196) -> (33,42) (1 bands)

Perhaps the meaning is to create SUBSET=x() so that it is 3340 meters wide,
not 8126322? And maybe the ESRI WCS server creates a corrupted tiff when it
receives a somewhat lunatic request?

-Jukka Rahkonen-




  






--
Sent from: http://osgeo-org.1560.x6.nabble.com/GDAL-Dev-f3742093.html


More information about the gdal-dev mailing list