[gdal-dev] gdal-dev Digest, Vol 198, Issue 14
Rahkonen Jukka (MML)
jukka.rahkonen at maanmittauslaitos.fi
Tue Nov 10 13:18:30 PST 2020
Hi,
Sorry for the -spat that belongs to ogr2ogr. In gdal_translate it is -projwin.
The WMTS documentation is not relevant for the GDAL WCS driver. Most that they have in common is that they create a xml file for storing info about the capabilities of the service. The WCS documentation alone should be enough for using the WCS driver, even I have to admit that I do not understand it totally, for example when it comes to cache control commands. Fortunately documentation can be improved.
Using WCS with GDAL is meant to be easy but unfortunately it seems to me that something is broken with it with some WCS services. It should be this easy to use:
1. Read general capabilities, including the names of the coverages. Example
gdalinfo "WCS:http://demo.geo-solutions.it/geoserver/wcs?version=2.0.1"
The result contain the names of the subdatasets (coverages) as the WCS driver understands them
Subdatasets:
SUBDATASET_1_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=nurc__Arc_Sample
SUBDATASET_2_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=sde__GRAY_HR_SR_OB_DR
SUBDATASET_3_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=sde__HYP_HR_SR_OB_DR
SUBDATASET_4_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=sde__NE2_HR_LC_SR_W_DR
SUBDATASET_5_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=nurc__Img_Sample
SUBDATASET_6_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=geosolutions__black_marble
SUBDATASET_7_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=geosolutions__g30
SUBDATASET_8_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=nurc__mosaic
SUBDATASET_9_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=sf__sfdem
SUBDATASET_10_NAME=WCS:http://demo.geo-solutions.it:80/geoserver/wcs?version=2.0.1&coverage=geosolutions__truemarble
1. Select one coverage and read more metadata about that, again with gdalinfo. This step fails now with my GDAL 3.2. Here some debug info that shows that GDAL is generating a bbox/subset with lower limit = higher limit when it tries to read a small sample of the coverage. Geoserver does not accept the request and gives an error about an empty subset.
gdalinfo "WCS:http://demo.geo-solutions.it/geoserver/wcs?version=2.0.1&coverage=nurc__mosaic" --debug on
…
WCS: Requesting http://demo.geo-solutions.it/geoserver/wcs?SERVICE=WCS&REQUEST=DescribeCoverage&VERSION=2.0.1&COVERAGEID=nurc__mosaic&FORMAT=text/xml
HTTP: Fetch(http://demo.geo-solutions.it/geoserver/wcs?SERVICE=WCS&REQUEST=DescribeCoverage&VERSION=2.0.1&COVERAGEID=nurc__mosaic&FORMAT=text/xml)
HTTP: libcurl/7.70.0-DEV OpenSSL/1.1.1g zlib/1.2.3
WCS: Requesting http://demo.geo-solutions.it/geoserver/wcs?SERVICE=WCS&REQUEST=GetCoverage&VERSION=2.0.1&COVERAGEID=nurc__mosaic&SUBSET=Long%286.3953399152693953,6.3953399152693953%29&SUBSET=Lat%2846.541601968340089,46.541601968340089%29&Format=image/tiff
HTTP: Fetch(http://demo.geo-solutions.it/geoserver/wcs?SERVICE=WCS&REQUEST=GetCoverage&VERSION=2.0.1&COVERAGEID=nurc__mosaic&SUBSET=Long%286.3953399152693953,6.3953399152693953%29&SUBSET=Lat%2846.541601968340089,46.541601968340089%29&Format=image/tiff)
ERROR 1: HTTP error code : 404
ERROR 1: InvalidSubsetting: Empty intersection after subsetting
1. The third step would be to use gdal_translate with -projwin option and get coverage from the given area. Now gdal_translate is also making a probing request by min=max extents and it fails.
However, the example request on the WCS document page that reads data from the National Land Survey of Finland works.
gdal_translate -oo CACHE=wcs_cache -oo CLEAR_CACHE -oo INTERLEAVE=PIXEL -projwin 377418 6683393.87938218 377717.879386966 6683094 -oo Subset="time(1985-01-01T00:00:00.000Z)" -outsize 60 0 "WCS:https://beta-karttakuva.maanmittauslaitos.fi/wcs/service/ows?version=2.0.1&coverage=ortokuva__ortokuva" scaled.tiff
gdalinfo works as well
gdalinfo "WCS:https://beta-karttakuva.maanmittauslaitos.fi/wcs/service/ows?version=2.0.1&coverage=ortokuva__ortokuva"
-Jukka Rahkonen-
Lähettäjä: 1520 gis <juliermeopensourcedeveloper at gmail.com>
Lähetetty: tiistai 10. marraskuuta 2020 18.10
Vastaanottaja: gdal-dev at lists.osgeo.org; Rahkonen Jukka (MML) <jukka.rahkonen at maanmittauslaitos.fi>
Aihe: Re: gdal-dev Digest, Vol 198, Issue 14
Dear Jukka Rahkonen,
Thank you very much for your reply.
First of all, my apologies for my lack of knowledge about using GDAL_WCS driver. Let me give more details about the issues I am facing.
Extract of http://172.21.14.45:8181/geoserver/ows/wcs?request=GetCapabilities is shown in [1].
Regarding gdal_translate GDAL_WCS driver Image subset extract, looking at https://gdal.org/drivers/raster/wcs.html#examples, I tried the following:
Similar to GDAL_WMTS <gdal_translate "WMTS:http://maps.wien.gv.at/wmts/1.0.0/WMTSCapabilities.xml,layer=lb" wmts.xml -of WMTS>,
I tried running:
(1)
gdal_translate "WCS:http://172.21.14.45:8181/geoserver/wcs?version=2.0.1&coverage=deter-amazonia__CBERS-4_AWFI_170_111_21052020" CBERS-4_AWFI_170_111_21052020.xml
The response: I did not get the xml as wmts.xml showed above
Note: Including the -of option (WCS), the outcome fails
(2)
gdal_translate "WCS:http://172.21.14.45:8181/geoserver/ows/wcs?request=GetCapabilities&coverage=deter-amazonia__CBERS-4_AWFI_170_111_21052020" \
-spat -58 -11 -57.75 -10.7 \
CBERS-4_AWFI_170_111_21052020.xml
The response: There is no -spat option for gdal_translate
gdal_translate --help
Usage: gdal_translate [--help-general] [--long-usage]
[-ot {Byte/Int16/UInt16/UInt32/Int32/Float32/Float64/
CInt16/CInt32/CFloat32/CFloat64}] [-strict]
[-of format] [-b band] [-mask band] [-expand {gray|rgb|rgba}]
[-outsize xsize[%]|0 ysize[%]|0] [-tr xres yres]
[-r {nearest,bilinear,cubic,cubicspline,lanczos,average,mode}]
[-unscale] [-scale[_bn] [src_min src_max [dst_min dst_max]]]* [-exponent[_bn] exp_val]*
[-srcwin xoff yoff xsize ysize] [-epo] [-eco]
[-projwin ulx uly lrx lry] [-projwin_srs srs_def]
[-a_srs srs_def] [-a_ullr ulx uly lrx lry] [-a_nodata value]
[-a_scale value] [-a_offset value]
[-nogcp] [-gcp pixel line easting northing [elevation]]*
|-colorinterp{_bn} {red|green|blue|alpha|gray|undefined}]
|-colorinterp {red|green|blue|alpha|gray|undefined},...]
[-mo "META-TAG=VALUE"]* [-q] [-sds]
[-co "NAME=VALUE"]* [-stats] [-norat]
[-oo NAME=VALUE]*
src_dataset dst_dataset
[1] WCS GetCapabilities request
<wcs:CoverageSummary>
<ows:Title>CBERS-4_AWFI_170_111_21052020</ows:Title>
<ows:Abstract>Generated from GeoTIFF</ows:Abstract>
<ows:Keywords>
<ows:Keyword>deter-amazonia-CBERS-4_AWFI_170_111_21052020</ows:Keyword>
<ows:Keyword>WCS</ows:Keyword>
<ows:Keyword>GeoTIFF</ows:Keyword>
</ows:Keywords>
<wcs:CoverageId>deter-amazonia__CBERS-4_AWFI_170_111_21052020</wcs:CoverageId>
<wcs:CoverageSubtype>RectifiedGridCoverage</wcs:CoverageSubtype>
<ows:WGS84BoundingBox>
<ows:LowerCorner>-63.280290728152664 -13.794296142085976</ows:LowerCorner>
<ows:UpperCorner>-53.624004749255874 -5.636644661059637</ows:UpperCorner>
</ows:WGS84BoundingBox>
<ows:BoundingBox crs="http://www.opengis.net/def/crs/EPSG/0/EPSG:4326">
<ows:LowerCorner>-63.280290728152664 -13.794296142085976</ows:LowerCorner>
<ows:UpperCorner>-53.624004749255874 -5.636644661059637</ows:UpperCorner>
</ows:BoundingBox>
</wcs:CoverageSummary>
Any comment will be very appreciated.
GDAL website<https://gdal.org/drivers/raster/wmts.html> provides examples on how to run GDAL_WMTS driver in order to get a wmts xml file. I have looked in https://gdal.org/drivers/raster/wcs.html and could not find an example pointing to -spat option. I got an -spat option
to gdal_grid<https://gdal.org/programs/gdal_grid.html?highlight=spat#cmdoption-gdal_grid-spat> though.
gdal_translate "WMTS:http://maps.wien.gv.at/wmts/1.0.0/WMTSCapabilities.xml,layer=lb" wmts.xml -of WMTS
gdal_translate "WCS:http://www.dpi.inpe.br/fipcerrado-geoserver/ows/wcs?request=GetCapabilities" \
-oo CoverageName=deter-cerrado:CBERS-4_AWFI_161_099_13072018 CBERS-4_AWFI_161_099_13072018.xml -of WCS
On Mon, Nov 9, 2020 at 5:00 PM <gdal-dev-request at lists.osgeo.org<mailto:gdal-dev-request at lists.osgeo.org>> wrote:
Send gdal-dev mailing list submissions to
gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>
To subscribe or unsubscribe via the World Wide Web, visit
https://lists.osgeo.org/mailman/listinfo/gdal-dev
or, via email, send a message with subject or body 'help' to
gdal-dev-request at lists.osgeo.org<mailto:gdal-dev-request at lists.osgeo.org>
You can reach the person managing the list at
gdal-dev-owner at lists.osgeo.org<mailto:gdal-dev-owner at lists.osgeo.org>
When replying, please edit your Subject line so it is more specific
than "Re: Contents of gdal-dev digest..."
Today's Topics:
1. Re: Selecting just the geometry column with OGR SQL (Even Rouault)
2. Re: Selecting just the geometry column with OGR SQL
(Rahkonen Jukka (MML))
3. Re: WCS_GDAL Driver Image SubSet Request (jratike80)
----------------------------------------------------------------------
Message: 1
Date: Mon, 09 Nov 2020 16:09:45 +0100
From: Even Rouault <even.rouault at spatialys.com<mailto:even.rouault at spatialys.com>>
To: gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>
Cc: "Rahkonen Jukka (MML)" <jukka.rahkonen at maanmittauslaitos.fi<mailto:jukka.rahkonen at maanmittauslaitos.fi>>
Subject: Re: [gdal-dev] Selecting just the geometry column with OGR
SQL
Message-ID: <1659214.2tRo4lr4Oo at even-i700<mailto:1659214.2tRo4lr4Oo at even-i700>>
Content-Type: text/plain; charset="us-ascii"
Jukka,
> I believe that OGR SQL dialect adds the geometry column into SQL selection
> by default. Am I right with this? I volunteer to edit the documentation
> https://gdal.org/user/ogr_sql_dialect.html if this is the case. But what I
> can't understand is another side of the story, how to select just the
> geometry field with OGR SQL dialect from some data sources.
>
> This command returns both time column and geometry from GML
> ogrinfo -sql "select time from timetest" timetest.gml
> I can select just the geometry from geometry column named as
> "geometryProperty" by ogrinfo -sql "select geometryProperty from timetest"
> timetest.gml
>
> The behavior is the same with SpatiaLite when I use "-dialect OGRSQL".
> Geometry is selected automatically but it can be selected also by name. But
> I cannot discover any way to select just the geometry field from shapefile.
> Is it because geometry column in shapefile does not have any name as this
> Python test seems to prove?
> >>> from osgeo import ogr
> >>> shp_ds = ogr.Open('timetest.shp')
> >>> shp_lyr = shp_ds.GetLayer(0)
> >>> shp_lyr.GetGeometryColumn()
>
> ''
>
> OpenJUMP JML driver does not report a name for geometry column either. Is it
> rather an exception that drivers give names for the geometry columns? I
> know I can select only the geometry column with the SQLite dialect, but I
> am curious.
Indeed, in OGR SQL, if the geometry column has a non-empty name (typically for
database drivers), you can use it to select it. When the geometry column name
is empty, such as in the shapefile or OpenJUMP driver, you can use
_ogr_geometry_ (with leading and trainling underscore) to select it. This is
likely undocumented.
As far as I can see, this was introduced per
https://github.com/OSGeo/gdal/commit/1519dcb00bd7802d2f701114f0e1cb9b524a547f
to fix an issue with comparison with empty string literals (this was before
OGR enforced the difference between single-quoted string literals and double-
quoted identifiers)
Even
--
Spatialys - Geospatial professional services
http://www.spatialys.com
------------------------------
Message: 2
Date: Mon, 9 Nov 2020 15:35:11 +0000
From: "Rahkonen Jukka (MML)" <jukka.rahkonen at maanmittauslaitos.fi<mailto:jukka.rahkonen at maanmittauslaitos.fi>>
To: Even Rouault <even.rouault at spatialys.com<mailto:even.rouault at spatialys.com>>,
"gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>" <gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>>
Subject: Re: [gdal-dev] Selecting just the geometry column with OGR
SQL
Message-ID: <d59795ddfd994dee9f5eca836361324a at maanmittauslaitos.fi<mailto:d59795ddfd994dee9f5eca836361324a at maanmittauslaitos.fi>>
Content-Type: text/plain; charset="us-ascii"
Even Rouault wrote:
maanantai 9. marraskuuta 2020 17.10
> Jukka,
>> I believe that OGR SQL dialect adds the geometry column into SQL
>> selection by default. Am I right with this? I volunteer to edit the
>> documentation https://gdal.org/user/ogr_sql_dialect.html if this is
>> the case. But what I can't understand is another side of the story,
>> how to select just the geometry field with OGR SQL dialect from some data sources.
>>
>> This command returns both time column and geometry from GML ogrinfo
>> -sql "select time from timetest" timetest.gml I can select just the
>> geometry from geometry column named as "geometryProperty" by ogrinfo
>> -sql "select geometryProperty from timetest"
>> timetest.gml
>>
>> The behavior is the same with SpatiaLite when I use "-dialect OGRSQL".
>> Geometry is selected automatically but it can be selected also by
>> name. But I cannot discover any way to select just the geometry field from shapefile.
>> Is it because geometry column in shapefile does not have any name as
>> this Python test seems to prove?
>> >>> from osgeo import ogr
>> >>> shp_ds = ogr.Open('timetest.shp')
>> >>> shp_lyr = shp_ds.GetLayer(0)
>> >>> shp_lyr.GetGeometryColumn()
>>
>> ''
>>
>> OpenJUMP JML driver does not report a name for geometry column either.
>> Is it rather an exception that drivers give names for the geometry
>> columns? I know I can select only the geometry column with the SQLite
>> dialect, but I am curious.
> Indeed, in OGR SQL, if the geometry column has a nonempty name (typically for
> database drivers), you can use it to select it. When the geometry column name
> is empty, such as in the shapefile or OpenJUMP driver, you can use
> _ogr_geometry_ (with leading and trainling underscore) to select it. This is likely undocumented.
> As far as I can see, this was introduced per
> https://github.com/OSGeo/gdal/commit/1519dcb00bd7802d2f701114f0e1cb9b524a547f
> to fix an issue with comparison with empty string literals (this was before OGR enforced the
> difference between single-quoted string literals and double- quoted identifiers)
I could not manage to guess the right syntax for Windows right ahead but solved the quiz easily
with your double-quote hint
First trial:
ogrinfo -sql "select _ogr_geometry_ from timetest" timetest.shp
INFO: Open of `timetest.shp'
using driver `ESRI Shapefile' successful.
ERROR 1: SQL Expression Parsing Error: syntax error, unexpected $undefined. Occurred around :
select _ogr_geometry_ from timetest
Second trial:
ogrinfo -sql "select \"_ogr_geometry_\" from timetest" timetest.shp
-Jukka-
> Even
--
Spatialys - Geospatial professional services http://www.spatialys.com
------------------------------
Message: 3
Date: Mon, 9 Nov 2020 08:47:05 -0700 (MST)
From: jratike80 <jukka.rahkonen at maanmittauslaitos.fi<mailto:jukka.rahkonen at maanmittauslaitos.fi>>
To: gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>
Subject: Re: [gdal-dev] WCS_GDAL Driver Image SubSet Request
Message-ID: <1604936825138-0.post at n6.nabble.com<mailto:1604936825138-0.post at n6.nabble.com>>
Content-Type: text/plain; charset=us-ascii
Hi,
Do you want to make is somehow impossible to request data beyound the
BoundingBox that you defined in the XML file? If you just want to get data
from the rectangle that you or your users define, gdal_translate with the
regular -spat option should work out-of-the-box.
-Jukka Rahkonen-
juliermeopensourcedeveloper wrote
> Hi all,
>
> Is it possible to define a bounding box in order to request a coverage
> subset using the GDAL WCS Driver?
> I have tried the code in [1]. I am trying to get information only for the
> bbox_interest=( -58,-11,-57.75,-10.75) as shown below, however the subset
> image or bbox_interest request has not been delivered. The West Longitude
> reaches -50.17 and The Lower Latitude reaches -18.92 which are limits out
> of the bbox_interest.
> The image extent is:
> <GetCoverageExtra>
> &BoundingBox=-63.280290728152664,-13.794296142085976,-53.624004749255874,-5.636644661059637
> </GetCoverageExtra>
> Any comment on this issue will be greatly appreciated.
> Thank you very much for your time in advance.
> *CODE [1]:*
> <WCS_GDAL>
> <ServiceURL>
> http://172.21.14.45:6060/geoserver/wcs?version=1.0.0&
> </
> ServiceURL>
> <CoverageName>
> cbers_deter_amazonia:CBERS-4_AWFI_170_111_21052020
> </
> CoverageName>
> <GetCoverageExtra>
> &BoundingBox=-58,-11,-57.75,-10.75
> </GetCoverageExtra>
> <CoverageOffering>
> <description>
> Generated from GeoTIFF
> </description>
> <name>
> cbers_deter_amazonia:CBERS-4_AWFI_170_111_21052020
> </name>
> <label>
> CBERS-4_AWFI_170_111_21052020
> </label>
> <lonLatEnvelope srsName="urn:ogc:def:crs:OGC:1.3:CRS84">
> <pos>
> -58 -11
> </pos>
> <pos>
> -57.75 -10.75
> </pos>
> </lonLatEnvelope>
> <keywords>
> <keyword>
> WCS
> </keyword>
> <keyword>
> GeoTIFF
> </keyword>
> </keywords>
> <domainSet>
> <spatialDomain>
> <Envelope srsName="EPSG:4326">
> <pos>
> -58 -11
> </pos>
> <pos>
> -57.75 -10.75
> </pos>
> </Envelope>
> <RectifiedGrid dimension="2" srsName="EPSG:4326">
> <limits>
> <GridEnvelope>
> <low>
> 0 0
> </low>
> <high>
> 0 0
> </high>
> </GridEnvelope>
> </limits>
> <axisName>
> x
> </axisName>
> <axisName>
> y
> </axisName>
> <origin>
> <pos>
> -58 -10.75
> </pos>
> </origin>
> <offsetVector>
> 0.005364603321609328 0.0
> </offsetVector>
> <offsetVector>
> 0.0 -0.004532028600570188
> </offsetVector>
> </RectifiedGrid>
> </spatialDomain>
> </domainSet>
> <rangeSet>
> <RangeSet>
> <name>
> CBERS-4_AWFI_170_111_21052020
> </name>
> <label>
> CBERS-4_AWFI_170_111_21052020
> </label>
> <axisDescription>
> <AxisDescription>
> <name>
> Band
> </name>
> <label>
> Band
> </label>
> <values>
> <interval>
> <min>
> 1
> </min>
> <max>
> 4
> </max>
> </interval>
> </values>
> </AxisDescription>
> </axisDescription>
> </RangeSet>
> </rangeSet>
> <supportedCRSs>
> <requestResponseCRSs>
> EPSG:4326
> </requestResponseCRSs>
> </supportedCRSs>
> <supportedFormats nativeFormat="GeoTIFF">
> <formats>
> GeoTIFF
> </formats>
> <formats>
> GIF
> </formats>
> <formats>
> JPEG
> </formats>
> <formats>
> PNG
> </formats>
> <formats>
> TIFF
> </formats>
> </supportedFormats>
> <supportedInterpolations default="nearest neighbor">
> <interpolationMethod>
> nearest neighbor
> </interpolationMethod>
> <interpolationMethod>
> bilinear
> </interpolationMethod>
> <interpolationMethod>
> bicubic
> </interpolationMethod>
> </supportedInterpolations>
> </CoverageOffering>
> <PreferredFormat>
> GeoTIFF
> </PreferredFormat>
> <BandCount>
> 4
> </BandCount>
> <BandType>
> Byte
> </BandType>
> </WCS_GDAL>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at .osgeo<mailto:gdal-dev at .osgeo>
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
--
Sent from: http://osgeo-org.1560.x6.nabble.com/GDAL-Dev-f3742093.html
------------------------------
Subject: Digest Footer
_______________________________________________
gdal-dev mailing list
gdal-dev at lists.osgeo.org<mailto:gdal-dev at lists.osgeo.org>
https://lists.osgeo.org/mailman/listinfo/gdal-dev
------------------------------
End of gdal-dev Digest, Vol 198, Issue 14
*****************************************
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20201110/d22c6805/attachment-0001.html>
More information about the gdal-dev
mailing list