1 /***************************************************************************
2 qgswcsprovider.cpp - QGIS Data provider for
3 OGC Web Coverage Service layers
4 -------------------
5 begin : 2 July, 2012
6 copyright : (C) 2012 by Radim Blazek
7 email : radim dot blazek at gmail.com
8
9 Based on qgswmsprovider.cpp written by Brendan Morley.
10
11 ***************************************************************************/
12
13 /***************************************************************************
14 * *
15 * This program is free software; you can redistribute it and/or modify *
16 * it under the terms of the GNU General Public License as published by *
17 * the Free Software Foundation; either version 2 of the License, or *
18 * (at your option) any later version. *
19 * *
20 ***************************************************************************/
21
22 #include "qgslogger.h"
23 #include "qgswcsprovider.h"
24 #include "qgscoordinatetransform.h"
25 #include "qgsdatasourceuri.h"
26 #include "qgsrasteridentifyresult.h"
27 #include "qgsrasterlayer.h"
28 #include "qgsrectangle.h"
29 #include "qgscoordinatereferencesystem.h"
30 #include "qgsnetworkaccessmanager.h"
31 #include "qgsnetworkreplyparser.h"
32 #include "qgsmessagelog.h"
33 #include "qgsexception.h"
34 #include "qgswcsdataitems.h"
35
36 #include <QNetworkRequest>
37 #include <QNetworkReply>
38 #include <QNetworkProxy>
39
40 #include <QUrl>
41 #include <QEventLoop>
42 #include <QFile>
43 #include <QUrlQuery>
44
45 #ifdef QGISDEBUG
46 #include <QDir>
47 #endif
48
49 #include "gdalwarper.h"
50 #include "ogr_srs_api.h"
51 #include "cpl_conv.h"
52 #include "cpl_string.h"
53
54 #define ERR(message) QGS_ERROR_MESSAGE(message,"WCS provider")
55 #define SRVERR(message) QGS_ERROR_MESSAGE(message,"WCS server")
56 #define QGS_ERROR(message) QgsError(message,"WCS provider")
57
58 QString QgsWcsProvider::WCS_KEY = QStringLiteral( "wcs" );
59 QString QgsWcsProvider::WCS_DESCRIPTION = QStringLiteral( "OGC Web Coverage Service version 1.0/1.1 data provider" );
60
61 static QString DEFAULT_LATLON_CRS = QStringLiteral( "CRS:84" );
62
63 // TODO: colortable - use common baseclass with gdal, mapserver does not support http://trac.osgeo.org/mapserver/ticket/1671
64
QgsWcsProvider(const QString & uri,const ProviderOptions & options,QgsDataProvider::ReadFlags flags)65 QgsWcsProvider::QgsWcsProvider( const QString &uri, const ProviderOptions &options, QgsDataProvider::ReadFlags flags )
66 : QgsRasterDataProvider( uri, options, flags )
67 , mCachedViewExtent( 0 )
68 {
69 QgsDebugMsg( "constructing with uri '" + mHttpUri + "'." );
70
71 mValid = false;
72 mCachedMemFilename = QStringLiteral( "/vsimem/qgis/wcs/%0.dat" ).arg( reinterpret_cast<std::uintptr_t>( this ) );
73
74 if ( !parseUri( uri ) ) return;
75
76 // GetCapabilities and DescribeCoverage
77 // TODO(?): do only DescribeCoverage to avoid one request
78 // We need to get at least server version, which is not set in of URI (if not part of url)
79 // and probably also rangeSet
80
81 QgsDataSourceUri capabilitiesUri;
82 capabilitiesUri.setEncodedUri( uri );
83 // remove non relevant params
84 capabilitiesUri.removeParam( QStringLiteral( "identifier" ) );
85 capabilitiesUri.removeParam( QStringLiteral( "crs" ) );
86 capabilitiesUri.removeParam( QStringLiteral( "format" ) );
87 // TODO: check if successful (add return to capabilities)
88 mCapabilities.setUri( capabilitiesUri );
89
90 // 1.0 get additional coverage info
91 if ( !mCapabilities.describeCoverage( mIdentifier ) )
92 {
93 appendError( ERR( tr( "Cannot describe coverage" ) ) );
94 return;
95 }
96
97 mCoverageSummary = mCapabilities.coverage( mIdentifier );
98 if ( !mCoverageSummary.valid )
99 {
100 // Should not happen if describeCoverage() did not fail
101 appendError( ERR( tr( "Coverage not found" ) ) );
102 return;
103 }
104
105 // It may happen that format is empty (e.g. uri created in python script),
106 // in that casei select one from available formats
107 if ( mFormat.isEmpty() )
108 {
109 // TIFF is known by GDAL
110 mFormat = mCoverageSummary.supportedFormat.filter( QStringLiteral( "tif" ), Qt::CaseInsensitive ).value( 0 );
111 }
112 if ( mFormat.isEmpty() )
113 {
114 // Take the first if TIFF was not found
115 mFormat = mCoverageSummary.supportedFormat.value( 0 );
116 }
117
118 // We cannot continue without format, it is required
119 if ( mFormat.isEmpty() ) return;
120
121 // It could happen (usually not with current QgsWCSSourceSelect if at least
122 // one CRS is available) that crs is not set in uri, in that case we
123 // use the native, if available or WGS84 or the first supported
124 if ( mCoverageCrs.isEmpty() )
125 {
126 QgsDebugMsg( "nativeCrs = " + mCoverageSummary.nativeCrs );
127 QgsDebugMsg( "supportedCrs = " + mCoverageSummary.supportedCrs.join( "," ) );
128 if ( !mCoverageSummary.nativeCrs.isEmpty() )
129 {
130 setCoverageCrs( mCoverageSummary.nativeCrs );
131 }
132 else if ( mCoverageSummary.supportedCrs.contains( QStringLiteral( "EPSG:4326" ), Qt::CaseInsensitive ) )
133 {
134 setCoverageCrs( QStringLiteral( "EPSG:4326" ) );
135 }
136 else if ( !mCoverageSummary.supportedCrs.isEmpty() )
137 {
138 setCoverageCrs( mCoverageSummary.supportedCrs.value( 0 ) );
139 }
140 }
141 QgsDebugMsg( "mCoverageCrs = " + mCoverageCrs );
142
143 // It may happen that coverage CRS is not given or it is unknown
144 // in that case we continue without CRS and user is asked for it
145 //if ( mCoverageCrs.isEmpty() ) return;
146
147 // Native size
148 mWidth = mCoverageSummary.width;
149 mHeight = mCoverageSummary.height;
150 mHasSize = mCoverageSummary.hasSize;
151
152 QgsDebugMsg( QStringLiteral( "mWidth = %1 mHeight = %2" ).arg( mWidth ).arg( mHeight ) );
153
154 // TODO: Consider if/how to recalculate mWidth, mHeight if non native CRS is used
155
156 if ( !calculateExtent() )
157 {
158 appendError( ERR( tr( "Cannot calculate extent" ) ) );
159 return;
160 }
161
162 // Get small piece of coverage to find GDAL data type and number of bands
163 const int bandNo = 0; // All bands
164 int width;
165 int height;
166 QString crs;
167 QgsRectangle box; // box to use to calc extent
168 // Prefer native CRS
169 if ( !mCoverageSummary.nativeCrs.isEmpty() &&
170 !mCoverageSummary.nativeBoundingBox.isEmpty() &&
171 mCoverageSummary.supportedCrs.contains( mCoverageSummary.nativeCrs ) &&
172 mHasSize )
173 {
174 box = mCoverageSummary.nativeBoundingBox;
175 width = mWidth;
176 height = mHeight;
177 crs = mCoverageSummary.nativeCrs;
178 }
179 else
180 {
181 box = mCoverageExtent;
182 if ( mHasSize )
183 {
184 width = mWidth;
185 height = mHeight;
186 }
187 else
188 {
189 // just a number to get smaller piece of coverage
190 width = 1000;
191 height = 1000;
192 }
193 }
194 const double xRes = box.width() / width;
195 const double yRes = box.height() / height;
196 const QgsPointXY p = box.center();
197
198 // width and height different to recognize rotation
199 const int requestWidth = 6;
200 const int requestHeight = 3;
201
202 // extent to be used for test request
203 const double halfWidth = xRes * ( requestWidth / 2. );
204 const double halfHeight = yRes * ( requestHeight / 2. );
205
206 // Using non native CRS (if we don't know which is native) it could easily happen,
207 // that a small part of bbox in request CRS near margin falls outside
208 // coverage native bbox and server reports error => take a piece from center
209 const QgsRectangle extent = QgsRectangle( p.x() - halfWidth, p.y() - halfHeight, p.x() + halfWidth, p.y() + halfHeight );
210
211 getCache( bandNo, extent, requestWidth, requestHeight, crs );
212
213 if ( !mCachedGdalDataset )
214 {
215 setError( mCachedError );
216 appendError( ERR( tr( "Cannot get test dataset." ) ) );
217 return;
218 }
219
220 mBandCount = GDALGetRasterCount( mCachedGdalDataset.get() );
221 QgsDebugMsg( QStringLiteral( "mBandCount = %1" ).arg( mBandCount ) );
222
223 // Check for server particularities (bbox, rotation)
224 const int responseWidth = GDALGetRasterXSize( mCachedGdalDataset.get() );
225 const int responseHeight = GDALGetRasterYSize( mCachedGdalDataset.get() );
226
227 QgsDebugMsg( QStringLiteral( "requestWidth = %1 requestHeight = %2 responseWidth = %3 responseHeight = %4)" ).arg( requestWidth ).arg( requestHeight ).arg( responseWidth ).arg( responseHeight ) );
228 // GeoServer and ArcGIS are using for 1.1 box "pixel" edges
229 // Mapserver is using pixel centers according to 1.1. specification
230 if ( ( responseWidth == requestWidth - 1 && responseHeight == requestHeight - 1 ) ||
231 ( responseWidth == requestHeight - 1 && responseHeight == requestWidth - 1 ) )
232 {
233 mFixBox = true;
234 QgsDebugMsg( QStringLiteral( "Test response size is smaller by pixel, using mFixBox" ) );
235 }
236 // Geoserver is giving rotated raster for geographic CRS - switched axis,
237 // Geoserver developers argue that changed axis order applies also to
238 // returned raster, that is exaggerated IMO but we have to handle that.
239 if ( ( responseWidth == requestHeight && responseHeight == requestWidth ) ||
240 ( responseWidth == requestHeight - 1 && responseHeight == requestWidth - 1 ) )
241 {
242 mFixRotate = true;
243 QgsDebugMsg( QStringLiteral( "Test response is rotated, using mFixRotate" ) );
244 }
245
246 // Get types
247 // TODO: we are using the same data types like GDAL (not wider like GDAL provider)
248 // with expectation to replace 'no data' values by NaN
249 mSrcGdalDataType.reserve( mBandCount );
250 mGdalDataType.reserve( mBandCount );
251 for ( int i = 1; i <= mBandCount; i++ )
252 {
253 GDALRasterBandH gdalBand = GDALGetRasterBand( mCachedGdalDataset.get(), i );
254 const GDALDataType myGdalDataType = GDALGetRasterDataType( gdalBand );
255
256 QgsDebugMsg( QStringLiteral( "myGdalDataType[%1] = %2" ).arg( i - 1 ).arg( myGdalDataType ) );
257 mSrcGdalDataType.append( myGdalDataType );
258
259 // UMN Mapserver does not automatically set null value, METADATA wcs_rangeset_nullvalue must be used
260 // http://lists.osgeo.org/pipermail/mapserver-users/2010-April/065328.html
261
262 // TODO: This could be shared with GDAL provider
263 int isValid = false;
264 double myNoDataValue = GDALGetRasterNoDataValue( gdalBand, &isValid );
265 if ( isValid )
266 {
267 QgsDebugMsg( QStringLiteral( "GDALGetRasterNoDataValue = %1" ).arg( myNoDataValue ) );
268 myNoDataValue = QgsRaster::representableValue( myNoDataValue, dataTypeFromGdal( myGdalDataType ) );
269 mSrcNoDataValue.append( myNoDataValue );
270 mSrcHasNoDataValue.append( true );
271 mUseSrcNoDataValue.append( true );
272 }
273 else
274 {
275 mSrcNoDataValue.append( std::numeric_limits<double>::quiet_NaN() );
276 mSrcHasNoDataValue.append( false );
277 mUseSrcNoDataValue.append( false );
278 }
279 // It may happen that nodata value given by GDAL is wrong and it has to be
280 // disabled by user, in that case we need another value to be used for nodata
281 // (for reprojection for example) -> always internally represent as wider type
282 // with mInternalNoDataValue in reserve.
283 // No retyping, no internal values for now
284 #if 0
285 int myInternalGdalDataType = myGdalDataType;
286 double myInternalNoDataValue;
287 switch ( srcDataType( i ) )
288 {
289 case Qgis::DataType::Byte:
290 myInternalNoDataValue = -32768.0;
291 myInternalGdalDataType = GDT_Int16;
292 break;
293 case Qgis::DataType::Int16:
294 myInternalNoDataValue = -2147483648.0;
295 myInternalGdalDataType = GDT_Int32;
296 break;
297 case Qgis::DataType::UInt16:
298 myInternalNoDataValue = -2147483648.0;
299 myInternalGdalDataType = GDT_Int32;
300 break;
301 case Qgis::DataType::Int32:
302 // We believe that such values is no used in real data
303 myInternalNoDataValue = -2147483648.0;
304 break;
305 case Qgis::DataType::UInt32:
306 // We believe that such values is no used in real data
307 myInternalNoDataValue = 4294967295.0;
308 break;
309 default: // Float32, Float64
310 //myNoDataValue = std::numeric_limits<int>::max();
311 // NaN should work well
312 myInternalNoDataValue = std::numeric_limits<double>::quiet_NaN();
313 }
314 mGdalDataType.append( myInternalGdalDataType );
315 mInternalNoDataValue.append( myInternalNoDataValue );
316 QgsDebugMsg( QStringLiteral( "mInternalNoDataValue[%1] = %2" ).arg( i - 1 ).arg( mInternalNoDataValue[i - 1] ) );
317 #endif
318 mGdalDataType.append( myGdalDataType );
319
320 #if 0
321 // TODO: what to do if null values from DescribeCoverage differ?
322 if ( !mCoverageSummary.nullValues.contains( myNoDataValue ) )
323 {
324 QgsDebugMsg( QStringLiteral( "noDataValue %1 is missing in nullValues from CoverageDescription" ).arg( myNoDataValue ) );
325 }
326 #endif
327
328 QgsDebugMsg( QStringLiteral( "mSrcGdalDataType[%1] = %2" ).arg( i - 1 ).arg( mSrcGdalDataType.at( i - 1 ) ) );
329 QgsDebugMsg( QStringLiteral( "mGdalDataType[%1] = %2" ).arg( i - 1 ).arg( mGdalDataType.at( i - 1 ) ) );
330 QgsDebugMsg( QStringLiteral( "mSrcNoDataValue[%1] = %2" ).arg( i - 1 ).arg( mSrcNoDataValue.at( i - 1 ) ) );
331
332 // Create and store color table
333 // TODO: never tested because mapserver (6.0.3) does not support color tables
334 mColorTables.append( QgsGdalProviderBase::colorTable( mCachedGdalDataset.get(), i ) );
335 }
336
337 clearCache();
338
339 // Block size is used for for statistics
340 // TODO: How to find maximum block size supported by server?
341 if ( mHasSize )
342 {
343 // This is taken from GDAL, how they come to these numbers?
344 if ( mWidth > 1800 ) mXBlockSize = 1024;
345 else mXBlockSize = mWidth;
346
347 if ( mHeight > 900 ) mYBlockSize = 512;
348 else mYBlockSize = mHeight;
349 }
350
351 mValid = true;
352 QgsDebugMsg( QStringLiteral( "Constructed OK, provider valid." ) );
353 }
354
QgsWcsProvider(const QgsWcsProvider & other,const QgsDataProvider::ProviderOptions & providerOptions)355 QgsWcsProvider::QgsWcsProvider( const QgsWcsProvider &other, const QgsDataProvider::ProviderOptions &providerOptions )
356 : QgsRasterDataProvider( other.dataSourceUri(), providerOptions )
357 , mHttpUri( other.mHttpUri )
358 , mBaseUrl( other.mBaseUrl )
359 , mIdentifier( other.mIdentifier )
360 , mTime( other.mTime )
361 , mFormat( other.mFormat )
362 , mValid( other.mValid )
363 , mCapabilities( other.mCapabilities )
364 , mCoverageSummary( other.mCoverageSummary )
365 , mSrid( other.mSrid )
366 , mCoverageExtent( other.mCoverageExtent )
367 , mWidth( other.mWidth )
368 , mHeight( other.mHeight )
369 , mXBlockSize( other.mXBlockSize )
370 , mYBlockSize( other.mYBlockSize )
371 , mHasSize( other.mHasSize )
372 , mBandCount( other.mBandCount )
373 , mGdalDataType( other.mGdalDataType )
374 , mSrcGdalDataType( other.mSrcGdalDataType )
375 , mColorTables( other.mColorTables )
376 , mExtentForLayer( other.mExtentForLayer )
377 , mCrsForLayer( other.mCrsForLayer )
378 , mQueryableForLayer( other.mQueryableForLayer )
379 , mCoverageCrs( other.mCoverageCrs )
380 // intentionally omitted:
381 // - mCachedData
382 // - mCachedMemFilename
383 // - mCachedMemFile
384 // - mCachedGdalDataset
385 // - mCachedError
386 // - mCachedViewExtent
387 // - mCachedViewWidth
388 // - mCachedViewHeight
389 , mMaxWidth( other.mMaxWidth )
390 , mMaxHeight( other.mMaxHeight )
391 // intentionally omitted:
392 // - mErrorCaption
393 // - mError
394 // - mErrorFormat
395 , mCoordinateTransform( other.mCoordinateTransform )
396 , mExtentDirty( other.mExtentDirty )
397 , mGetFeatureInfoUrlBase( other.mGetFeatureInfoUrlBase )
398 , mServiceMetadataURL( other.mServiceMetadataURL )
399 , mAuth( other.mAuth )
400 , mIgnoreGetCoverageUrl( other.mIgnoreGetCoverageUrl )
401 , mIgnoreAxisOrientation( other.mIgnoreAxisOrientation )
402 , mInvertAxisOrientation( other.mInvertAxisOrientation )
403 , mCrs( other.mCrs )
404 , mFixBox( other.mFixBox )
405 , mFixRotate( other.mFixRotate )
406 , mCacheLoadControl( other.mCacheLoadControl )
407 {
408 mCachedMemFilename = QStringLiteral( "/vsimem/qgis/wcs/%0.dat" ).arg( reinterpret_cast<std::uintptr_t>( this ) );
409 }
410
411
parseUri(const QString & uriString)412 bool QgsWcsProvider::parseUri( const QString &uriString )
413 {
414
415 QgsDebugMsg( "uriString = " + uriString );
416 QgsDataSourceUri uri;
417 uri.setEncodedUri( uriString );
418
419 mMaxWidth = 0;
420 mMaxHeight = 0;
421
422 mHttpUri = uri.param( QStringLiteral( "url" ) );
423 mBaseUrl = prepareUri( mHttpUri );
424 QgsDebugMsg( "mBaseUrl = " + mBaseUrl );
425
426 mIgnoreGetCoverageUrl = uri.hasParam( QStringLiteral( "IgnoreGetMapUrl" ) );
427 mIgnoreAxisOrientation = uri.hasParam( QStringLiteral( "IgnoreAxisOrientation" ) ); // must be before parsing!
428 mInvertAxisOrientation = uri.hasParam( QStringLiteral( "InvertAxisOrientation" ) ); // must be before parsing!
429
430 mAuth.mUserName = uri.username();
431 QgsDebugMsg( "set username to " + mAuth.mUserName );
432
433 mAuth.mPassword = uri.password();
434 QgsDebugMsg( "set password to " + mAuth.mPassword );
435
436 if ( !uri.authConfigId().isEmpty() )
437 {
438 mAuth.mAuthCfg = uri.authConfigId();
439 }
440 QgsDebugMsg( "set authcfg to " + mAuth.mAuthCfg );
441
442 mIdentifier = uri.param( QStringLiteral( "identifier" ) );
443
444 mTime = uri.param( QStringLiteral( "time" ) );
445
446 setFormat( uri.param( QStringLiteral( "format" ) ) );
447
448 if ( !uri.param( QStringLiteral( "crs" ) ).isEmpty() )
449 {
450 setCoverageCrs( uri.param( QStringLiteral( "crs" ) ) );
451 }
452
453 const QString cache = uri.param( QStringLiteral( "cache" ) );
454 if ( !cache.isEmpty() )
455 {
456 mCacheLoadControl = QgsNetworkAccessManager::cacheLoadControlFromName( cache );
457 }
458 QgsDebugMsg( QStringLiteral( "mCacheLoadControl = %1" ).arg( mCacheLoadControl ) );
459
460 return true;
461 }
462
prepareUri(QString uri) const463 QString QgsWcsProvider::prepareUri( QString uri ) const
464 {
465 if ( !uri.contains( '?' ) )
466 {
467 uri.append( '?' );
468 }
469 else if ( uri.right( 1 ) != QLatin1String( "?" ) && uri.right( 1 ) != QLatin1String( "&" ) )
470 {
471 uri.append( '&' );
472 }
473
474 return uri;
475 }
476
~QgsWcsProvider()477 QgsWcsProvider::~QgsWcsProvider()
478 {
479 QgsDebugMsg( QStringLiteral( "deconstructing." ) );
480
481 // Dispose of any cached image as created by draw()
482 clearCache();
483 }
484
clone() const485 QgsWcsProvider *QgsWcsProvider::clone() const
486 {
487 QgsDataProvider::ProviderOptions providerOptions;
488 providerOptions.transformContext = transformContext();
489 QgsWcsProvider *provider = new QgsWcsProvider( *this, providerOptions );
490 provider->copyBaseSettings( *this );
491 return provider;
492 }
493
baseUrl() const494 QString QgsWcsProvider::baseUrl() const
495 {
496 return mBaseUrl;
497 }
498
format() const499 QString QgsWcsProvider::format() const
500 {
501 return mFormat;
502 }
503
setFormat(QString const & format)504 void QgsWcsProvider::setFormat( QString const &format )
505 {
506 QgsDebugMsg( "Setting format to " + format + '.' );
507 mFormat = format;
508 }
509
510
setCoverageCrs(QString const & crs)511 void QgsWcsProvider::setCoverageCrs( QString const &crs )
512 {
513 QgsDebugMsg( "Setting coverage CRS to " + crs + '.' );
514
515 if ( crs != mCoverageCrs && !crs.isEmpty() )
516 {
517 // delete old coordinate transform as it is no longer valid
518 mCoordinateTransform = QgsCoordinateTransform();
519
520 mExtentDirty = true;
521
522 mCoverageCrs = crs;
523
524 mCrs = QgsCoordinateReferenceSystem::fromOgcWmsCrs( mCoverageCrs );
525 }
526 }
527
setQueryItem(QUrl & url,const QString & item,const QString & value) const528 void QgsWcsProvider::setQueryItem( QUrl &url, const QString &item, const QString &value ) const
529 {
530 QUrlQuery query( url );
531 query.removeQueryItem( item );
532 query.addQueryItem( item, value );
533 url.setQuery( query );
534 }
535
readBlock(int bandNo,QgsRectangle const & viewExtent,int pixelWidth,int pixelHeight,void * block,QgsRasterBlockFeedback * feedback)536 bool QgsWcsProvider::readBlock( int bandNo, QgsRectangle const &viewExtent, int pixelWidth, int pixelHeight, void *block, QgsRasterBlockFeedback *feedback )
537 {
538 // TODO: set block to null values, move that to function and call only if fails
539 memset( block, 0, pixelWidth * pixelHeight * QgsRasterBlock::typeSize( dataType( bandNo ) ) );
540
541 // Requested extent must at least partially overlap coverage extent, otherwise
542 // server gives error. QGIS usually does not request blocks outside raster extent
543 // (higher level checks) but it is better to do check here as well
544 if ( !viewExtent.intersects( mCoverageExtent ) )
545 {
546 return false;
547 }
548
549 // Can we reuse the previously cached coverage?
550 if ( !mCachedGdalDataset ||
551 mCachedViewExtent != viewExtent ||
552 mCachedViewWidth != pixelWidth ||
553 mCachedViewHeight != pixelHeight )
554 {
555 getCache( bandNo, viewExtent, pixelWidth, pixelHeight, QString(), feedback );
556 }
557
558 if ( mCachedGdalDataset )
559 {
560 // It may happen (Geoserver) that if requested BBOX is larger than coverage
561 // extent, the returned data cover intersection of requested BBOX and coverage
562 // extent scaled to requested WIDTH/HEIGHT => check extent
563 // Unfortunately if the received raster does not have a CRS, the extent is the raster size
564 // and in that case it cannot be used to verify extent
565 QgsCoordinateReferenceSystem cacheCrs;
566 if ( !cacheCrs.createFromWkt( GDALGetProjectionRef( mCachedGdalDataset.get() ) ) &&
567 !cacheCrs.createFromWkt( GDALGetGCPProjection( mCachedGdalDataset.get() ) ) )
568 {
569 QgsDebugMsg( QStringLiteral( "Cached does not have CRS" ) );
570 }
571 QgsDebugMsg( "Cache CRS: " + cacheCrs.userFriendlyIdentifier() );
572
573 const QgsRectangle cacheExtent = QgsGdalProviderBase::extent( mCachedGdalDataset.get() );
574 QgsDebugMsg( "viewExtent = " + viewExtent.toString() );
575 QgsDebugMsg( "cacheExtent = " + cacheExtent.toString() );
576 // TODO: check also rotated
577 if ( cacheCrs.isValid() && !mFixRotate )
578 {
579 // using qgsDoubleNear is too precise, example acceptable difference:
580 // 179.9999999306699863 x 179.9999999306700431
581 if ( !qgsDoubleNearSig( cacheExtent.xMinimum(), viewExtent.xMinimum(), 10 ) ||
582 !qgsDoubleNearSig( cacheExtent.yMinimum(), viewExtent.yMinimum(), 10 ) ||
583 !qgsDoubleNearSig( cacheExtent.xMaximum(), viewExtent.xMaximum(), 10 ) ||
584 !qgsDoubleNearSig( cacheExtent.yMaximum(), viewExtent.yMaximum(), 10 ) )
585 {
586 // Just print a message so user is aware of a server side issue but don't left
587 // the tile blank so we can deal with eventually misconfigured WCS server
588 // https://github.com/qgis/QGIS/issues/33339
589
590 QgsDebugMsg( QStringLiteral( "cacheExtent and viewExtent differ" ) );
591 QgsMessageLog::logMessage( tr( "Received coverage has wrong extent %1 (expected %2)" ).arg( cacheExtent.toString(), viewExtent.toString() ), tr( "WCS" ) );
592 }
593 }
594
595 const int width = GDALGetRasterXSize( mCachedGdalDataset.get() );
596 const int height = GDALGetRasterYSize( mCachedGdalDataset.get() );
597 QgsDebugMsg( QStringLiteral( "cached data width = %1 height = %2 (expected %3 x %4)" ).arg( width ).arg( height ).arg( pixelWidth ).arg( pixelHeight ) );
598
599 GDALRasterBandH gdalBand = GDALGetRasterBand( mCachedGdalDataset.get(), bandNo );
600 // TODO: check type?, check band count?
601 if ( mFixRotate && width == pixelHeight && height == pixelWidth )
602 {
603 // Rotate counter clockwise
604 // If GridOffsets With GeoServer,
605 QgsDebugMsg( QStringLiteral( "Rotating raster" ) );
606 const int pixelSize = QgsRasterBlock::typeSize( dataType( bandNo ) );
607 QgsDebugMsg( QStringLiteral( "pixelSize = %1" ).arg( pixelSize ) );
608 const int size = width * height * pixelSize;
609 void *tmpData = malloc( size );
610 if ( ! tmpData )
611 {
612 QgsDebugMsg( QStringLiteral( "Couldn't allocate memory of %1 bytes" ).arg( size ) );
613 return false;
614 }
615 if ( GDALRasterIO( gdalBand, GF_Read, 0, 0, width, height, tmpData, width, height, ( GDALDataType ) mGdalDataType.at( bandNo - 1 ), 0, 0 ) != CE_None )
616 {
617 QgsDebugMsg( QStringLiteral( "Raster IO Error" ) );
618 }
619 for ( int i = 0; i < pixelHeight; i++ )
620 {
621 for ( int j = 0; j < pixelWidth; j++ )
622 {
623 const int destIndex = pixelSize * ( i * pixelWidth + j );
624 const int srcIndex = pixelSize * ( j * width + ( width - i - 1 ) );
625 memcpy( ( char * )block + destIndex, ( char * )tmpData + srcIndex, pixelSize );
626 }
627 }
628 free( tmpData );
629 }
630 else if ( width == pixelWidth && height == pixelHeight )
631 {
632 if ( GDALRasterIO( gdalBand, GF_Read, 0, 0, pixelWidth, pixelHeight, block, pixelWidth, pixelHeight, ( GDALDataType ) mGdalDataType.at( bandNo - 1 ), 0, 0 ) != CE_None )
633 {
634 QgsDebugMsg( QStringLiteral( "Raster IO Error" ) );
635 }
636 else
637 {
638 QgsDebugMsg( QStringLiteral( "Block read OK" ) );
639 }
640 }
641 else
642 {
643 // This should not happen, but it is better to give distorted result + warning
644 if ( GDALRasterIO( gdalBand, GF_Read, 0, 0, width, height, block, pixelWidth, pixelHeight, ( GDALDataType ) mGdalDataType.at( bandNo - 1 ), 0, 0 ) != CE_None )
645 {
646 QgsDebugMsg( QStringLiteral( "Raster IO Error" ) );
647 }
648 QgsMessageLog::logMessage( tr( "Received coverage has wrong size %1 x %2 (expected %3 x %4)" ).arg( width ).arg( height ).arg( pixelWidth ).arg( pixelHeight ), tr( "WCS" ) );
649 }
650 }
651 return true;
652 }
653
getCache(int bandNo,QgsRectangle const & viewExtent,int pixelWidth,int pixelHeight,QString crs,QgsRasterBlockFeedback * feedback) const654 void QgsWcsProvider::getCache( int bandNo, QgsRectangle const &viewExtent, int pixelWidth, int pixelHeight, QString crs, QgsRasterBlockFeedback *feedback ) const
655 {
656 Q_UNUSED( bandNo )
657 // delete cached data
658 clearCache();
659
660 if ( crs.isEmpty() )
661 {
662 crs = mCoverageCrs;
663 }
664
665 mCachedViewExtent = viewExtent;
666 mCachedViewWidth = pixelWidth;
667 mCachedViewHeight = pixelHeight;
668
669 // --------------- BOUNDING BOX --------------------------------
670 //according to the WCS spec for 1.1, some CRS have inverted axis
671 // box:
672 // 1.0.0: minx,miny,maxx,maxy
673 // 1.1.0, 1.1.2: OGC 07-067r5 (WCS 1.1.2) refers to OGC 06-121r3 which says
674 // "The number of axes included, and the order of these axes, shall be as specified
675 // by the referenced CRS." That means inverted for geographic.
676 bool changeXY = false;
677 if ( !mIgnoreAxisOrientation && ( mCapabilities.version().startsWith( QLatin1String( "1.1" ) ) ) )
678 {
679 //create CRS from string
680 const QgsCoordinateReferenceSystem srs = QgsCoordinateReferenceSystem::fromOgcWmsCrs( crs );
681 if ( srs.isValid() && srs.hasAxisInverted() )
682 {
683 changeXY = true;
684 }
685 }
686
687 if ( mInvertAxisOrientation ) changeXY = !changeXY;
688
689 const double xRes = viewExtent.width() / pixelWidth;
690 const double yRes = viewExtent.height() / pixelHeight;
691 QgsRectangle extent = viewExtent;
692 // WCS 1.1 grid is using grid points (surrounded by sample spaces) and
693 // "The spatial extent of a grid coverage extends only as far as the outermost
694 // grid points contained in the bounding box. It does NOT include any area
695 // (partial or whole grid cells or sample spaces) beyond those grid points."
696 // Mapserver and GDAL are using bbox defined by grid points, i.e. shrunk
697 // by 1 pixel, but Geoserver and ArcGIS are using full bbox including
698 // the space around edge grid points.
699 if ( mCapabilities.version().startsWith( QLatin1String( "1.1" ) ) && !mFixBox )
700 {
701 // shrink the extent to border cells centers by half cell size
702 extent = QgsRectangle( viewExtent.xMinimum() + xRes / 2., viewExtent.yMinimum() + yRes / 2., viewExtent.xMaximum() - xRes / 2., viewExtent.yMaximum() - yRes / 2. );
703 }
704
705 // Bounding box in WCS format (Warning: does not work with scientific notation)
706 QString bbox = QString( changeXY ? "%2,%1,%4,%3" : "%1,%2,%3,%4" )
707 .arg( qgsDoubleToString( extent.xMinimum() ),
708 qgsDoubleToString( extent.yMinimum() ),
709 qgsDoubleToString( extent.xMaximum() ),
710 qgsDoubleToString( extent.yMaximum() ) );
711
712 QUrl url( mIgnoreGetCoverageUrl ? mBaseUrl : mCapabilities.getCoverageUrl() );
713
714 // Version 1.0.0, 1.1.0, 1.1.2
715 setQueryItem( url, QStringLiteral( "SERVICE" ), QStringLiteral( "WCS" ) );
716 setQueryItem( url, QStringLiteral( "VERSION" ), mCapabilities.version() );
717 setQueryItem( url, QStringLiteral( "REQUEST" ), QStringLiteral( "GetCoverage" ) );
718 setQueryItem( url, QStringLiteral( "FORMAT" ), mFormat );
719
720 // Version 1.0.0
721 if ( mCapabilities.version().startsWith( QLatin1String( "1.0" ) ) )
722 {
723 setQueryItem( url, QStringLiteral( "COVERAGE" ), mIdentifier );
724 if ( !mTime.isEmpty() )
725 {
726 // It seems that Mmapserver (6.0.3) WCS 1.1 completely ignores
727 // TemporalDomain. Some code (copy-pasted from 1.0) is commented in
728 // msWCSDescribeCoverage_CoverageDescription11() and GetCoverage
729 // TimeSequence param is not supported at all. If a coverage is defined
730 // with timeposition in mapfile, the result of GetCoverage is empty
731 // raster (all values 0).
732 setQueryItem( url, QStringLiteral( "TIME" ), mTime );
733 }
734 setQueryItem( url, QStringLiteral( "BBOX" ), bbox );
735 setQueryItem( url, QStringLiteral( "CRS" ), crs ); // request BBOX CRS
736 setQueryItem( url, QStringLiteral( "RESPONSE_CRS" ), crs ); // response CRS
737 setQueryItem( url, QStringLiteral( "WIDTH" ), QString::number( pixelWidth ) );
738 setQueryItem( url, QStringLiteral( "HEIGHT" ), QString::number( pixelHeight ) );
739 }
740
741 // Version 1.1.0, 1.1.2
742 if ( mCapabilities.version().startsWith( QLatin1String( "1.1" ) ) )
743 {
744 setQueryItem( url, QStringLiteral( "IDENTIFIER" ), mIdentifier );
745 const QString crsUrn = QStringLiteral( "urn:ogc:def:crs:%1::%2" ).arg( crs.split( ':' ).value( 0 ), crs.split( ':' ).value( 1 ) );
746 bbox += ',' + crsUrn;
747
748 if ( !mTime.isEmpty() )
749 {
750 setQueryItem( url, QStringLiteral( "TIMESEQUENCE" ), mTime );
751 }
752
753 setQueryItem( url, QStringLiteral( "BOUNDINGBOX" ), bbox );
754
755 // Example:
756 // GridBaseCRS=urn:ogc:def:crs:SG:6.6:32618
757 // GridType=urn:ogc:def:method:CS:1.1:2dGridIn2dCrs
758 // GridCS=urn:ogc:def:cs:OGC:0Grid2dSquareCS
759 // GridOrigin=0,0
760 // GridOffsets=0.0707,-0.0707,0.1414,0.1414&
761
762 setQueryItem( url, QStringLiteral( "GRIDBASECRS" ), crsUrn ); // response CRS
763
764 setQueryItem( url, QStringLiteral( "GRIDCS" ), QStringLiteral( "urn:ogc:def:cs:OGC:0.0:Grid2dSquareCS" ) );
765
766 setQueryItem( url, QStringLiteral( "GRIDTYPE" ), QStringLiteral( "urn:ogc:def:method:WCS:1.1:2dSimpleGrid" ) );
767
768 // GridOrigin is BBOX minx, maxy
769 // Note: shifting origin to cell center (not really necessary nor making sense)
770 // does not work with Mapserver 6.0.3
771 // Mapserver 6.0.3 does not work with origin on yMinimum (lower left)
772 // Geoserver works OK with yMinimum (lower left)
773 const QString gridOrigin = QString( changeXY ? "%2,%1" : "%1,%2" )
774 .arg( qgsDoubleToString( extent.xMinimum() ),
775 qgsDoubleToString( extent.yMaximum() ) );
776 setQueryItem( url, QStringLiteral( "GRIDORIGIN" ), gridOrigin );
777
778 // GridOffsets WCS 1.1:
779 // GridType urn:ogc:def:method:WCS:1.1:2dSimpleGrid : 2 values
780 // GridType urn:ogc:def:method:WCS:1.1:2dGridIn2dCrs : 4 values, the center two of these offsets will be zero when the GridCRS is not rotated or skewed in the GridBaseCRS.
781 // The yOff must be negative because origin is in upper corner
782 // WCS 1.1.2: BaseX = origin(1) + offsets(1) * GridX
783 // BaseY = origin(2) + offsets(2) * GridY
784 // otherwise GeoServer gives mirrored result, however
785 // Mapserver works OK with both positive and negative
786 // OTOH, the yOff offset must not be negative with mFixRotate and Geoserver 2.1-SNAPSHOT
787 // but it must be negative with GeoServer 2.1.3 and mFixRotate. I am not sure
788 // at this moment 100% -> disabling positive yOff for now - TODO: try other servers
789 //double yOff = mFixRotate ? yRes : -yRes; // this was working with some servers I think
790 const double yOff = -yRes;
791 const QString gridOffsets = QString( changeXY ? "%2,%1" : "%1,%2" )
792 //setQueryItem( url, "GRIDTYPE", "urn:ogc:def:method:WCS:1.1:2dGridIn2dCrs" );
793 //QString gridOffsets = QString( changeXY ? "%2,0,0,%1" : "%1,0,0,%2" )
794 .arg( qgsDoubleToString( xRes ),
795 qgsDoubleToString( yOff ) );
796 setQueryItem( url, QStringLiteral( "GRIDOFFSETS" ), gridOffsets );
797 }
798
799 QgsDebugMsg( QStringLiteral( "GetCoverage: %1" ).arg( url.toString() ) );
800
801 // cache some details for if the user wants to do an identifyAsHtml() later
802
803 //mGetFeatureInfoUrlBase = mIgnoreGetFeatureInfoUrl ? mBaseUrl : getFeatureInfoUrl();
804
805 QgsWcsDownloadHandler handler( url, mAuth, mCacheLoadControl, mCachedData, mCapabilities.version(), mCachedError, feedback );
806 handler.blockingDownload();
807
808 QgsDebugMsg( QStringLiteral( "%1 bytes received" ).arg( mCachedData.size() ) );
809 if ( mCachedData.isEmpty() )
810 {
811 if ( !feedback || !feedback->isCanceled() )
812 {
813 QgsMessageLog::logMessage( tr( "No data received" ), tr( "WCS" ) );
814 }
815 clearCache();
816 return;
817 }
818
819 #if 0
820 QFile myFile( "/tmp/qgiswcscache.dat" );
821 if ( myFile.open( QIODevice::WriteOnly ) )
822 {
823 QDataStream myStream( &myFile );
824 myStream.writeRawData( mCachedData.data(), mCachedData.size() );
825 myFile.close();
826 }
827 #endif
828
829 mCachedMemFile = VSIFileFromMemBuffer( mCachedMemFilename.toUtf8().constData(),
830 ( GByte * )mCachedData.data(),
831 ( vsi_l_offset )mCachedData.size(),
832 FALSE );
833
834 if ( !mCachedMemFile )
835 {
836 QgsMessageLog::logMessage( tr( "Cannot create memory file" ), tr( "WCS" ) );
837 clearCache();
838 return;
839 }
840 QgsDebugMsg( QStringLiteral( "Memory file created" ) );
841
842 CPLErrorReset();
843 mCachedGdalDataset.reset( GDALOpen( mCachedMemFilename.toUtf8().constData(), GA_ReadOnly ) );
844 if ( !mCachedGdalDataset )
845 {
846 QgsMessageLog::logMessage( QString::fromUtf8( CPLGetLastErrorMsg() ), tr( "WCS" ) );
847 clearCache();
848 return;
849 }
850 QgsDebugMsg( QStringLiteral( "Dataset opened" ) );
851
852 }
853
854 // For stats only, maybe change QgsRasterDataProvider::bandStatistics() to
855 // use standard readBlock with extent
readBlock(int bandNo,int xBlock,int yBlock,void * block)856 bool QgsWcsProvider::readBlock( int bandNo, int xBlock, int yBlock, void *block )
857 {
858
859 QgsDebugMsg( QStringLiteral( "xBlock = %1 yBlock = %2" ).arg( xBlock ).arg( yBlock ) );
860
861 if ( !mHasSize )
862 return false;
863
864 const double xRes = mCoverageExtent.width() / mWidth;
865 const double yRes = mCoverageExtent.height() / mHeight;
866
867 // blocks on edges may run out of extent, that should not be problem (at least for
868 // stats - there is a check for it)
869 const double xMin = mCoverageExtent.xMinimum() + xRes * xBlock * mXBlockSize;
870 const double xMax = xMin + xRes * mXBlockSize;
871 const double yMax = mCoverageExtent.yMaximum() - yRes * yBlock * mYBlockSize;
872 const double yMin = yMax - yRes * mYBlockSize;
873 //QgsDebugMsg( QStringLiteral("yMin = %1 yMax = %2").arg(yMin).arg(yMax) );
874
875 const QgsRectangle extent( xMin, yMin, xMax, yMax );
876
877 return readBlock( bandNo, extent, mXBlockSize, mYBlockSize, block, nullptr );
878 }
879
880
881 // This could be shared with GDAL provider
sourceDataType(int bandNo) const882 Qgis::DataType QgsWcsProvider::sourceDataType( int bandNo ) const
883 {
884 if ( bandNo <= 0 || bandNo > mSrcGdalDataType.size() )
885 {
886 return Qgis::DataType::UnknownDataType;
887 }
888
889 return dataTypeFromGdal( mSrcGdalDataType[bandNo - 1] );
890 }
891
dataType(int bandNo) const892 Qgis::DataType QgsWcsProvider::dataType( int bandNo ) const
893 {
894 if ( bandNo <= 0 || bandNo > mGdalDataType.size() )
895 {
896 return Qgis::DataType::UnknownDataType;
897 }
898
899 return dataTypeFromGdal( mGdalDataType[bandNo - 1] );
900 }
901
bandCount() const902 int QgsWcsProvider::bandCount() const
903 {
904 return mBandCount;
905 }
906
907 // this is only called once when statistics are calculated
908 // TODO
xBlockSize() const909 int QgsWcsProvider::xBlockSize() const
910 {
911 return mXBlockSize;
912 }
yBlockSize() const913 int QgsWcsProvider::yBlockSize() const
914 {
915 return mYBlockSize;
916 }
917
xSize() const918 int QgsWcsProvider::xSize() const { return mWidth; }
ySize() const919 int QgsWcsProvider::ySize() const { return mHeight; }
920
clearCache() const921 void QgsWcsProvider::clearCache() const
922 {
923 if ( mCachedGdalDataset )
924 {
925 QgsDebugMsg( QStringLiteral( "Close mCachedGdalDataset" ) );
926 mCachedGdalDataset.reset();
927 QgsDebugMsg( QStringLiteral( "Closed" ) );
928 }
929 if ( mCachedMemFile )
930 {
931 QgsDebugMsg( QStringLiteral( "Close mCachedMemFile" ) );
932 VSIFCloseL( mCachedMemFile );
933 mCachedMemFile = nullptr;
934 QgsDebugMsg( QStringLiteral( "Closed" ) );
935 }
936 QgsDebugMsg( QStringLiteral( "Clear mCachedData" ) );
937 mCachedData.clear();
938 mCachedError.clear();
939 QgsDebugMsg( QStringLiteral( "Cleared" ) );
940 }
941
colorTable(int bandNumber) const942 QList<QgsColorRampShader::ColorRampItem> QgsWcsProvider::colorTable( int bandNumber )const
943 {
944 return mColorTables.value( bandNumber - 1 );
945 }
946
colorInterpretation(int bandNo) const947 int QgsWcsProvider::colorInterpretation( int bandNo ) const
948 {
949 GDALRasterBandH myGdalBand = GDALGetRasterBand( mCachedGdalDataset.get(), bandNo );
950 return colorInterpretationFromGdal( GDALGetRasterColorInterpretation( myGdalBand ) );
951 }
952
953
parseServiceExceptionReportDom(const QByteArray & xml,const QString & wcsVersion,QString & errorTitle,QString & errorText)954 bool QgsWcsProvider::parseServiceExceptionReportDom( const QByteArray &xml, const QString &wcsVersion, QString &errorTitle, QString &errorText )
955 {
956
957 #ifdef QGISDEBUG
958 //test the content of the QByteArray
959 const QString responsestring( xml );
960 QgsDebugMsg( "received the following data: " + responsestring );
961 #endif
962
963 // Convert completed document into a Dom
964 QDomDocument doc;
965 QString errorMsg;
966 int errorLine;
967 int errorColumn;
968 const bool contentSuccess = doc.setContent( xml, false, &errorMsg, &errorLine, &errorColumn );
969
970 if ( !contentSuccess )
971 {
972 errorTitle = tr( "Dom Exception" );
973 errorText = tr( "Could not get WCS Service Exception at %1 at line %2 column %3\n\nResponse was:\n\n%4" )
974 .arg( errorMsg )
975 .arg( errorLine )
976 .arg( errorColumn )
977 .arg( QString( xml ) );
978
979 QgsLogger::debug( "Dom Exception: " + errorText );
980
981 return false;
982 }
983
984 const QDomElement docElem = doc.documentElement();
985
986 // TODO: Assert the docElem.tagName() is
987 // ServiceExceptionReport // 1.0
988 // ows:ExceptionReport // 1.1
989
990 //QString version = docElem.attribute("version");
991
992 QDomElement e;
993 if ( wcsVersion.startsWith( QLatin1String( "1.0" ) ) )
994 {
995 e = QgsWcsCapabilities::domElement( docElem, QStringLiteral( "ServiceException" ) );
996 }
997 else // 1.1
998 {
999 e = QgsWcsCapabilities::domElement( docElem, QStringLiteral( "Exception" ) );
1000 }
1001 parseServiceException( e, wcsVersion, errorTitle, errorText );
1002
1003 QgsDebugMsg( QStringLiteral( "exiting." ) );
1004
1005 return true;
1006 }
1007
parseServiceException(QDomElement const & e,const QString & wcsVersion,QString & errorTitle,QString & errorText)1008 void QgsWcsProvider::parseServiceException( QDomElement const &e, const QString &wcsVersion, QString &errorTitle, QString &errorText )
1009 {
1010
1011 errorTitle = tr( "Service Exception" );
1012
1013 QMap<QString, QString> exceptions;
1014
1015 // Some codes are in both 1.0 and 1.1, in that case the 'meaning of code' is
1016 // taken from 1.0 specification
1017
1018 // set up friendly descriptions for the service exception
1019 // 1.0
1020 exceptions[QStringLiteral( "InvalidFormat" )] = tr( "Request contains a format not offered by the server." );
1021 exceptions[QStringLiteral( "CoverageNotDefined" )] = tr( "Request is for a Coverage not offered by the service instance." );
1022 exceptions[QStringLiteral( "CurrentUpdateSequence" )] = tr( "Value of (optional) UpdateSequence parameter in GetCapabilities request is equal to current value of service metadata update sequence number." );
1023 exceptions[QStringLiteral( "InvalidUpdateSequence" )] = tr( "Value of (optional) UpdateSequence parameter in GetCapabilities request is greater than current value of service metadata update sequence number." );
1024 // 1.0, 1.1
1025 exceptions[QStringLiteral( "MissingParameterValue" )] = tr( "Request does not include a parameter value, and the server instance did not declare a default value for that dimension." );
1026 exceptions[QStringLiteral( "InvalidParameterValue" )] = tr( "Request contains an invalid parameter value." );
1027 // 1.1
1028 exceptions[QStringLiteral( "NoApplicableCode" )] = tr( "No other exceptionCode specified by this service and server applies to this exception." );
1029 exceptions[QStringLiteral( "UnsupportedCombination" )] = tr( "Operation request contains an output CRS that can not be used within the output format." );
1030 exceptions[QStringLiteral( "NotEnoughStorage" )] = tr( "Operation request specifies to \"store\" the result, but not enough storage is available to do this." );
1031
1032 QString seCode;
1033 QString seText;
1034 if ( wcsVersion.startsWith( QLatin1String( "1.0" ) ) )
1035 {
1036 seCode = e.attribute( QStringLiteral( "code" ) );
1037 seText = e.text();
1038 }
1039 else
1040 {
1041 const QStringList codes;
1042 seCode = e.attribute( QStringLiteral( "exceptionCode" ) );
1043 // UMN Mapserver (6.0.3) has messed/switched 'locator' and 'exceptionCode'
1044 if ( ! exceptions.contains( seCode ) )
1045 {
1046 seCode = e.attribute( QStringLiteral( "locator" ) );
1047 if ( ! exceptions.contains( seCode ) )
1048 {
1049 seCode.clear();
1050 }
1051 }
1052 seText = QgsWcsCapabilities::firstChildText( e, QStringLiteral( "ExceptionText" ) );
1053 }
1054
1055 if ( seCode.isEmpty() )
1056 {
1057 errorText = tr( "(No error code was reported)" );
1058 }
1059 else if ( exceptions.contains( seCode ) )
1060 {
1061 errorText = exceptions.value( seCode );
1062 }
1063 else
1064 {
1065 errorText = seCode + ' ' + tr( "(Unknown error code)" );
1066 }
1067
1068 errorText += '\n' + tr( "The WCS vendor also reported: " );
1069 errorText += seText;
1070
1071 QgsMessageLog::logMessage( tr( "composed error message '%1'." ).arg( errorText ), tr( "WCS" ) );
1072 QgsDebugMsg( QStringLiteral( "exiting." ) );
1073 }
1074
1075
1076
extent() const1077 QgsRectangle QgsWcsProvider::extent() const
1078 {
1079 if ( mExtentDirty )
1080 {
1081 if ( calculateExtent() )
1082 {
1083 mExtentDirty = false;
1084 }
1085 }
1086
1087 return mCoverageExtent;
1088 }
1089
isValid() const1090 bool QgsWcsProvider::isValid() const
1091 {
1092 return mValid;
1093 }
1094
1095
wcsVersion()1096 QString QgsWcsProvider::wcsVersion()
1097 {
1098 return mCapabilities.version();
1099 }
1100
calculateExtent() const1101 bool QgsWcsProvider::calculateExtent() const
1102 {
1103
1104 // Make sure we know what extents are available
1105 if ( !mCoverageSummary.described )
1106 {
1107 return false;
1108 }
1109
1110 // Prefer to use extent from capabilities / coverage description because
1111 // transformation from WGS84 enlarges the extent
1112 mCoverageExtent = mCoverageSummary.boundingBoxes.value( mCoverageCrs );
1113 QgsDebugMsg( "mCoverageCrs = " + mCoverageCrs + " mCoverageExtent = " + mCoverageExtent.toString() );
1114 if ( !mCoverageExtent.isEmpty() && mCoverageExtent.isFinite() )
1115 {
1116 QgsDebugMsg( "mCoverageExtent = " + mCoverageExtent.toString() );
1117 //return true;
1118 }
1119 else
1120 {
1121 // Set up the coordinate transform from the WCS standard CRS:84 bounding
1122 // box to the user's selected CRS
1123 if ( !mCoordinateTransform.isValid() )
1124 {
1125 const QgsCoordinateReferenceSystem qgisSrsSource = QgsCoordinateReferenceSystem::fromOgcWmsCrs( QStringLiteral( "EPSG:4326" ) );
1126 const QgsCoordinateReferenceSystem qgisSrsDest = QgsCoordinateReferenceSystem::fromOgcWmsCrs( mCoverageCrs );
1127
1128 //QgsDebugMsg( "qgisSrsSource: " + qgisSrsSource.toWkt() );
1129 //QgsDebugMsg( "qgisSrsDest: " + qgisSrsDest.toWkt() );
1130
1131 mCoordinateTransform = QgsCoordinateTransform( qgisSrsSource, qgisSrsDest, transformContext() );
1132 }
1133
1134 QgsDebugMsg( "mCoverageSummary.wgs84BoundingBox= " + mCoverageSummary.wgs84BoundingBox.toString() );
1135
1136 // Convert to the user's CRS as required
1137 try
1138 {
1139 mCoverageExtent = mCoordinateTransform.transformBoundingBox( mCoverageSummary.wgs84BoundingBox, Qgis::TransformDirection::Forward );
1140 }
1141 catch ( QgsCsException &cse )
1142 {
1143 Q_UNUSED( cse )
1144 return false;
1145 }
1146
1147 //make sure extent does not contain 'inf' or 'nan'
1148 if ( !mCoverageExtent.isFinite() )
1149 {
1150 return false;
1151 }
1152 }
1153
1154 QgsDebugMsg( "mCoverageExtent = " + mCoverageExtent.toString() );
1155
1156 // It may happen (GeoServer) that extent reported in spatialDomain.Envelope is larger
1157 // than the coverage. Then if that larger BBOX is requested, the server returns
1158 // request BBOX intersected with coverage box scaled to requested WIDTH and HEIGHT.
1159 // GDAL WCS client does not suffer from this probably because it probably takes
1160 // extent from lonLatEnvelope (it probably does not calculate it from
1161 // spatialDomain.RectifiedGrid because calculated value is slightly different).
1162
1163 // To get the true extent (it can also be smaller than real if reported Envelope is
1164 // than real smaller, but smaller is safer because data cannot be shifted) we make
1165 // request of the whole extent cut the extent from spatialDomain.Envelope if
1166 // necessary
1167
1168 getCache( 1, mCoverageExtent, 10, 10 );
1169 if ( mCachedGdalDataset )
1170 {
1171 const QgsRectangle cacheExtent = QgsGdalProviderBase::extent( mCachedGdalDataset.get() );
1172 QgsDebugMsg( "mCoverageExtent = " + mCoverageExtent.toString() );
1173 QgsDebugMsg( "cacheExtent = " + cacheExtent.toString() );
1174 QgsCoordinateReferenceSystem cacheCrs;
1175 if ( !cacheCrs.createFromWkt( GDALGetProjectionRef( mCachedGdalDataset.get() ) ) &&
1176 !cacheCrs.createFromWkt( GDALGetGCPProjection( mCachedGdalDataset.get() ) ) )
1177 {
1178 QgsDebugMsg( QStringLiteral( "Cached does not have CRS" ) );
1179 }
1180 QgsDebugMsg( "Cache CRS: " + cacheCrs.userFriendlyIdentifier() );
1181
1182 // We can only verify extent if CRS is set
1183 // If dataset comes rotated, GDAL probably cuts latitude extend, disable
1184 // extent check for rotated, TODO: verify
1185 if ( cacheCrs.isValid() && !mFixRotate )
1186 {
1187 if ( !qgsDoubleNearSig( cacheExtent.xMinimum(), mCoverageExtent.xMinimum(), 10 ) ||
1188 !qgsDoubleNearSig( cacheExtent.yMinimum(), mCoverageExtent.yMinimum(), 10 ) ||
1189 !qgsDoubleNearSig( cacheExtent.xMaximum(), mCoverageExtent.xMaximum(), 10 ) ||
1190 !qgsDoubleNearSig( cacheExtent.yMaximum(), mCoverageExtent.yMaximum(), 10 ) )
1191 {
1192 QgsDebugMsg( QStringLiteral( "cacheExtent and mCoverageExtent differ, mCoverageExtent cut to cacheExtent" ) );
1193 mCoverageExtent = cacheExtent;
1194 }
1195 }
1196 }
1197 else
1198 {
1199 // Unfortunately it may also happen that a server (cubewerx.com) does not have
1200 // overviews and it is not able to respond for the whole extent within timeout.
1201 // It returns timeout error.
1202 // In that case (if request failed) we do not report error to allow working
1203 // with such servers on smaller portions of extent
1204 // (http://lists.osgeo.org/pipermail/qgis-developer/2013-January/024019.html)
1205
1206 // Unfortunately even if we get over this 10x10 check, QGIS also requests
1207 // 32x32 thumbnail where it is waiting for another timeout
1208
1209 QgsDebugMsg( QStringLiteral( "Cannot get cache to verify extent" ) );
1210 QgsMessageLog::logMessage( tr( "Cannot verify coverage full extent: %1" ).arg( mCachedError.message() ), tr( "WCS" ) );
1211 }
1212
1213 return true;
1214 }
1215
1216
capabilities() const1217 int QgsWcsProvider::capabilities() const
1218 {
1219 int capability = NoCapabilities;
1220 capability |= QgsRasterDataProvider::Identify;
1221 capability |= QgsRasterDataProvider::IdentifyValue;
1222 capability |= QgsRasterDataProvider::Prefetch;
1223
1224 if ( mHasSize )
1225 {
1226 capability |= QgsRasterDataProvider::Size;
1227 }
1228
1229 return capability;
1230 }
1231
coverageMetadata(const QgsWcsCoverageSummary & coverage)1232 QString QgsWcsProvider::coverageMetadata( const QgsWcsCoverageSummary &coverage )
1233 {
1234 QString metadata;
1235
1236 // Use a nested table
1237 metadata += QLatin1String( "<tr><td>" );
1238 metadata += QLatin1String( "<table width=\"100%\">" );
1239
1240 // Table header
1241 metadata += QLatin1String( "<tr><th class=\"strong\">" );
1242 metadata += tr( "Property" );
1243 metadata += QLatin1String( "</th>" );
1244 metadata += QLatin1String( "<th class=\"strong\">" );
1245 metadata += tr( "Value" );
1246 metadata += QLatin1String( "</th></tr>" );
1247
1248 metadata += htmlRow( tr( "Name (identifier)" ), coverage.identifier );
1249 metadata += htmlRow( tr( "Title" ), coverage.title );
1250 metadata += htmlRow( tr( "Abstract" ), coverage.abstract );
1251
1252 if ( !coverage.metadataLink.metadataType.isNull() &&
1253 !coverage.metadataLink.xlinkHref.isNull() )
1254 {
1255 metadata += htmlRow( tr( "Metadata Type" ), coverage.metadataLink.metadataType );
1256 metadata += htmlRow( tr( "Metadata Link" ), coverage.metadataLink.xlinkHref );
1257 }
1258
1259 #if 0
1260 // We don't have size, nativeCrs, nativeBoundingBox etc. until describe coverage which would be heavy for all coverages
1261 metadata += htmlRow( tr( "Fixed Width" ), QString::number( coverage.width ) );
1262 metadata += htmlRow( tr( "Fixed Height" ), QString::number( coverage.height ) );
1263 metadata += htmlRow( tr( "Native CRS" ), coverage.nativeCrs );
1264 metadata += htmlRow( tr( "Native Bounding Box" ), coverage.nativeBoundingBox.toString() );
1265 #endif
1266
1267 metadata += htmlRow( tr( "WGS 84 Bounding Box" ), coverage.wgs84BoundingBox.toString() );
1268
1269 // Layer Coordinate Reference Systems
1270 // TODO(?): supportedCrs and supportedFormat are not available in 1.0
1271 // until coverage is described - it would be confusing to show it only if available
1272 #if 0
1273 for ( int j = 0; j < std::min( coverage.supportedCrs.size(), 10 ); j++ )
1274 {
1275 metadata += htmlRow( tr( "Available in CRS" ), coverage.supportedCrs.value( j ) );
1276 }
1277
1278 if ( coverage.supportedCrs.size() > 10 )
1279 {
1280 metadata += htmlRow( tr( "Available in CRS" ), tr( "(and %n more)", "crs", coverage.supportedCrs.size() - 10 ) );
1281 }
1282
1283 for ( int j = 0; j < std::min( coverage.supportedFormat.size(), 10 ); j++ )
1284 {
1285 metadata += htmlRow( tr( "Available in format" ), coverage.supportedFormat.value( j ) );
1286 }
1287
1288 if ( coverage.supportedFormat.size() > 10 )
1289 {
1290 metadata += htmlRow( tr( "Available in format" ), tr( "(and %n more)", "crs", coverage.supportedFormat.size() - 10 ) );
1291 }
1292 #endif
1293
1294 // Close the nested table
1295 metadata += QLatin1String( "</table>" );
1296 metadata += QLatin1String( "</td></tr>" );
1297
1298 return metadata;
1299 }
1300
htmlMetadata()1301 QString QgsWcsProvider::htmlMetadata()
1302 {
1303 QString metadata;
1304 metadata += QStringLiteral( "<tr><td class=\"highlight\">" ) + tr( "WCS Info" ) + QStringLiteral( "</td><td><div>" );
1305
1306 metadata += QLatin1String( "</a> <a href=\"#coverages\">" );
1307 metadata += tr( "Coverages" );
1308 metadata += QLatin1String( "</a>" );
1309
1310 #if 0
1311 // TODO
1312 metadata += "<a href=\"#cachestats\">";
1313 metadata += tr( "Cache Stats" );
1314 metadata += "</a> ";
1315 #endif
1316
1317 metadata += QLatin1String( "<br /><table class=\"tabular-view\">" ); // Nested table 1
1318
1319 // Server Properties section
1320 metadata += QLatin1String( "<tr><th class=\"strong\"><a name=\"serverproperties\"></a>" );
1321 metadata += tr( "Server Properties" );
1322 metadata += QLatin1String( "</th></tr>" );
1323
1324 // Use a nested table
1325 metadata += QLatin1String( "<tr><td>" );
1326 metadata += QLatin1String( "<table width=\"100%\">" ); // Nested table 2
1327
1328 // Table header
1329 metadata += QLatin1String( "<tr><th class=\"strong\">" );
1330 metadata += tr( "Property" );
1331 metadata += QLatin1String( "</th>" );
1332 metadata += QLatin1String( "<th class=\"strong\">" );
1333 metadata += tr( "Value" );
1334 metadata += QLatin1String( "</th></tr>" );
1335
1336 metadata += htmlRow( ( "WCS Version" ), mCapabilities.version() );
1337 metadata += htmlRow( tr( "Title" ), mCapabilities.capabilities().title );
1338 metadata += htmlRow( tr( "Abstract" ), mCapabilities.capabilities().abstract );
1339 #if 0
1340 // TODO: probably apply stylesheet in QgsWcsCapabilities and save as html
1341 metadata += htmlRow( tr( "Keywords" ), mCapabilities.service.keywordList.join( "<br />" ) );
1342 metadata += htmlRow( tr( "Online Resource" ), "-" );
1343 metadata += htmlRow( tr( "Contact Person" ),
1344 mCapabilities.service.contactInformation.contactPersonPrimary.contactPerson
1345 + "<br />" + mCapabilities.service.contactInformation.contactPosition;
1346 + "<br />" + mCapabilities.service.contactInformation.contactPersonPrimary.contactOrganization );
1347 metadata += htmlRow( tr( "Fees" ), mCapabilities.service.fees );
1348 metadata += htmlRow( tr( "Access Constraints" ), mCapabilities.service.accessConstraints );
1349 metadata += htmlRow( tr( "Image Formats" ), mCapabilities.capability.request.getMap.format.join( "<br />" ) );
1350 metadata += htmlRow( tr( "GetCapabilitiesUrl" ), mBaseUrl );
1351 #endif
1352 metadata += htmlRow( tr( "Get Coverage Url" ), mCapabilities.getCoverageUrl() + ( mIgnoreGetCoverageUrl ? tr( " <font color=\"red\">(advertised but ignored)</font>" ) : QString() ) );
1353
1354 // Close the nested table
1355 metadata += QLatin1String( "</table>" ); // End nested table 2
1356 metadata += QLatin1String( "</td></tr>" );
1357
1358 // Coverage properties
1359 metadata += QLatin1String( "<tr><th class=\"strong\"><a name=\"coverages\"></a>" );
1360 metadata += tr( "Coverages" );
1361 metadata += QLatin1String( "</th></tr>" );
1362
1363 // Dialog takes too long to open if there are too many coverages (1000 for example)
1364 int count = 0;
1365 const auto constCoverages = mCapabilities.coverages();
1366 for ( const QgsWcsCoverageSummary &c : constCoverages )
1367 {
1368 metadata += coverageMetadata( c );
1369 count++;
1370 if ( count >= 100 ) break;
1371 }
1372 metadata += QLatin1String( "</table>" );
1373 if ( count < mCapabilities.coverages().size() )
1374 {
1375 metadata += tr( "And %1 more coverages" ).arg( mCapabilities.coverages().size() - count );
1376 }
1377
1378 metadata += QLatin1String( "</table></div></td></tr>\n" ); // End nested table 1
1379 return metadata;
1380 }
1381
htmlCell(const QString & text)1382 QString QgsWcsProvider::htmlCell( const QString &text )
1383 {
1384 return "<td>" + text + "</td>";
1385 }
1386
htmlRow(const QString & text1,const QString & text2)1387 QString QgsWcsProvider:: htmlRow( const QString &text1, const QString &text2 )
1388 {
1389 return "<tr>" + htmlCell( text1 ) + htmlCell( text2 ) + "</tr>";
1390 }
1391
identify(const QgsPointXY & point,QgsRaster::IdentifyFormat format,const QgsRectangle & boundingBox,int width,int height,int)1392 QgsRasterIdentifyResult QgsWcsProvider::identify( const QgsPointXY &point, QgsRaster::IdentifyFormat format, const QgsRectangle &boundingBox, int width, int height, int /*dpi*/ )
1393 {
1394 QgsDebugMsg( QStringLiteral( "thePoint = %1 %2" ).arg( point.x(), 0, 'g', 10 ).arg( point.y(), 0, 'g', 10 ) );
1395 QgsDebugMsg( QStringLiteral( "theWidth = %1 height = %2" ).arg( width ).arg( height ) );
1396 QgsDebugMsg( "boundingBox = " + boundingBox.toString() );
1397 QMap<int, QVariant> results;
1398
1399 if ( format != QgsRaster::IdentifyFormatValue )
1400 {
1401 return QgsRasterIdentifyResult( QGS_ERROR( tr( "Format not supported" ) ) );
1402 }
1403
1404 if ( !extent().contains( point ) )
1405 {
1406 // Outside the raster
1407 for ( int i = 1; i <= bandCount(); i++ )
1408 {
1409 results.insert( i, QVariant() );
1410 }
1411 return QgsRasterIdentifyResult( QgsRaster::IdentifyFormatValue, results );
1412 }
1413
1414 QgsRectangle finalExtent = boundingBox;
1415 const int maxSize = 2000;
1416 const int cacheSize = 1000; // tile cache size if context is not defined or small
1417 const double relResTol = 0.1; // relative resolution tolerance (10%)
1418
1419 // TODO: We are using cacheSize x cacheSize if context is not defined
1420 // or box is too small (in pixels). That is necessary to allow effective
1421 // map querying (valuetool for example). OTOH, it means reading of large
1422 // unnecessary amount of data for each identify() call if query points are too
1423 // distant (out of cache) for example when called from scripts
1424
1425 // if context size is to large we have to cut it, in that case caching big
1426 // big part does not make sense
1427 if ( finalExtent.isEmpty() || width == 0 || height == 0 ||
1428 width > maxSize || height > maxSize )
1429 {
1430 // context missing, use a small area around the point and highest resolution if known
1431
1432 double xRes, yRes;
1433 if ( mHasSize )
1434 {
1435 xRes = mCoverageExtent.width() / mWidth;
1436 yRes = mCoverageExtent.height() / mHeight;
1437 }
1438 else
1439 {
1440 // set resolution approximately to 1mm
1441 switch ( mCrs.mapUnits() )
1442 {
1443 case QgsUnitTypes::DistanceMeters:
1444 xRes = 0.001;
1445 break;
1446 case QgsUnitTypes::DistanceFeet:
1447 xRes = 0.003;
1448 break;
1449 case QgsUnitTypes::DistanceDegrees:
1450 // max length of degree of latitude on pole is 111694 m
1451 xRes = 1e-8;
1452 break;
1453 default:
1454 xRes = 0.001; // expecting meters
1455 }
1456 yRes = xRes;
1457 }
1458
1459 width = cacheSize;
1460 height = cacheSize;
1461
1462 finalExtent = QgsRectangle( point.x() - xRes * width / 2,
1463 point.y() - yRes * height / 2,
1464 point.x() + xRes * width / 2,
1465 point.y() + yRes * height / 2 );
1466
1467 const double xResDiff = std::fabs( mCachedViewExtent.width() / mCachedViewWidth - xRes );
1468 const double yResDiff = std::fabs( mCachedViewExtent.height() / mCachedViewHeight - yRes );
1469
1470 if ( !mCachedGdalDataset ||
1471 !mCachedViewExtent.contains( point ) ||
1472 mCachedViewWidth == 0 || mCachedViewHeight == 0 ||
1473 xResDiff / xRes > relResTol ||
1474 yResDiff / yRes > relResTol )
1475 {
1476 getCache( 1, finalExtent, width, height );
1477 }
1478 }
1479 else
1480 {
1481 // Use context -> effective caching (usually, if context is constant)
1482 QgsDebugMsg( QStringLiteral( "Using context extent and resolution" ) );
1483 // To use the cache it is sufficient to have point within cache and
1484 // similar resolution
1485 const double xRes = finalExtent.width() / width;
1486 const double yRes = finalExtent.height() / height;
1487 QgsDebugMsg( QStringLiteral( "width = %1 height = %2 xRes = %3 yRes = %4" ).arg( finalExtent.width() ).arg( finalExtent.height() ).arg( xRes ).arg( yRes ) );
1488 const double xResDiff = std::fabs( mCachedViewExtent.width() / mCachedViewWidth - xRes );
1489 const double yResDiff = std::fabs( mCachedViewExtent.height() / mCachedViewHeight - yRes );
1490 QgsDebugMsg( QStringLiteral( "xRes diff = %1 yRes diff = %2 relative xResDiff = %3 relative yResDiff = %4" ).arg( xResDiff ).arg( yResDiff ).arg( xResDiff / xRes ).arg( yResDiff / yRes ) );
1491 if ( !mCachedGdalDataset ||
1492 !mCachedViewExtent.contains( point ) ||
1493 xResDiff / xRes > relResTol ||
1494 yResDiff / yRes > relResTol )
1495 {
1496 // Identify map tool is now using fake context 1x1 pixel to get point/line
1497 // features from WMS. In such case we enlarge the extent to get data cached
1498 // for next queries around.
1499 // BTW: UMN Mapserver (6.0.3) seems to be buggy with 1x1 pixels request (returns 'no data' value
1500 if ( width < cacheSize )
1501 {
1502 const int buffer = ( cacheSize - width ) / 2;
1503 width += 2 * buffer;
1504 finalExtent.setXMinimum( finalExtent.xMinimum() - xRes * buffer );
1505 finalExtent.setXMaximum( finalExtent.xMaximum() + xRes * buffer );
1506 }
1507 if ( height < cacheSize )
1508 {
1509 const int buffer = ( cacheSize - height ) / 2;
1510 height += 2 * buffer;
1511 finalExtent.setYMinimum( finalExtent.yMinimum() - yRes * buffer );
1512 finalExtent.setYMaximum( finalExtent.yMaximum() + yRes * buffer );
1513 }
1514
1515 getCache( 1, finalExtent, width, height );
1516 }
1517 }
1518
1519 if ( !mCachedGdalDataset ||
1520 !mCachedViewExtent.contains( point ) )
1521 {
1522 return QgsRasterIdentifyResult( QGS_ERROR( tr( "Read data error" ) ) );
1523 }
1524
1525 const double x = point.x();
1526 const double y = point.y();
1527
1528 // Calculate the row / column where the point falls
1529 const double xRes = mCachedViewExtent.width() / mCachedViewWidth;
1530 const double yRes = mCachedViewExtent.height() / mCachedViewHeight;
1531
1532 // Offset, not the cell index -> flor
1533 const int col = ( int ) std::floor( ( x - mCachedViewExtent.xMinimum() ) / xRes );
1534 const int row = ( int ) std::floor( ( mCachedViewExtent.yMaximum() - y ) / yRes );
1535
1536 QgsDebugMsg( "row = " + QString::number( row ) + " col = " + QString::number( col ) );
1537
1538 for ( int i = 1; i <= GDALGetRasterCount( mCachedGdalDataset.get() ); i++ )
1539 {
1540 GDALRasterBandH gdalBand = GDALGetRasterBand( mCachedGdalDataset.get(), i );
1541
1542 double value;
1543 const CPLErr err = GDALRasterIO( gdalBand, GF_Read, col, row, 1, 1,
1544 &value, 1, 1, GDT_Float64, 0, 0 );
1545
1546 if ( err != CPLE_None )
1547 {
1548 QgsLogger::warning( "RasterIO error: " + QString::fromUtf8( CPLGetLastErrorMsg() ) );
1549 return QgsRasterIdentifyResult( QGS_ERROR( tr( "RasterIO error: " ) + QString::fromUtf8( CPLGetLastErrorMsg() ) ) );
1550 }
1551
1552 // Apply no data and user no data
1553 if ( ( sourceHasNoDataValue( i ) && useSourceNoDataValue( i ) &&
1554 ( std::isnan( value ) || qgsDoubleNear( value, sourceNoDataValue( i ) ) ) ) ||
1555 ( QgsRasterRange::contains( value, userNoDataValues( i ) ) ) )
1556 {
1557 results.insert( i, QVariant() );
1558 }
1559 else
1560 {
1561 results.insert( i, value );
1562 }
1563 }
1564
1565 return QgsRasterIdentifyResult( QgsRaster::IdentifyFormatValue, results );
1566 }
1567
crs() const1568 QgsCoordinateReferenceSystem QgsWcsProvider::crs() const
1569 {
1570 return mCrs;
1571 }
1572
lastErrorTitle()1573 QString QgsWcsProvider::lastErrorTitle()
1574 {
1575 return mErrorCaption;
1576 }
1577
lastError()1578 QString QgsWcsProvider::lastError()
1579 {
1580 QgsDebugMsg( "returning '" + mError + "'." );
1581 return mError;
1582 }
1583
lastErrorFormat()1584 QString QgsWcsProvider::lastErrorFormat()
1585 {
1586 return mErrorFormat;
1587 }
1588
name() const1589 QString QgsWcsProvider::name() const
1590 {
1591 return WCS_KEY;
1592 }
1593
providerKey()1594 QString QgsWcsProvider::providerKey()
1595 {
1596 return WCS_KEY;
1597 }
1598
description() const1599 QString QgsWcsProvider::description() const
1600 {
1601 return WCS_DESCRIPTION;
1602 }
1603
providerCapabilities() const1604 QgsRasterDataProvider::ProviderCapabilities QgsWcsProvider::providerCapabilities() const
1605 {
1606 return ProviderCapability::ReloadData;
1607 }
1608
reloadProviderData()1609 void QgsWcsProvider::reloadProviderData()
1610 {
1611 clearCache();
1612 }
1613
nodeAttribute(const QDomElement & e,const QString & name,const QString & defValue)1614 QString QgsWcsProvider::nodeAttribute( const QDomElement &e, const QString &name, const QString &defValue )
1615 {
1616 if ( e.hasAttribute( name ) )
1617 return e.attribute( name );
1618
1619 const QDomNamedNodeMap map( e.attributes() );
1620 for ( int i = 0; i < map.size(); i++ )
1621 {
1622 const QDomAttr attr( map.item( i ).toElement().toAttr() );
1623 if ( attr.name().compare( name, Qt::CaseInsensitive ) == 0 )
1624 return attr.value();
1625 }
1626
1627 return defValue;
1628 }
1629
supportedMimes()1630 QMap<QString, QString> QgsWcsProvider::supportedMimes()
1631 {
1632 QMap<QString, QString> mimes;
1633 GDALAllRegister();
1634
1635 QgsDebugMsg( QStringLiteral( "GDAL drivers cont %1" ).arg( GDALGetDriverCount() ) );
1636 for ( int i = 0; i < GDALGetDriverCount(); ++i )
1637 {
1638 GDALDriverH driver = GDALGetDriver( i );
1639 Q_CHECK_PTR( driver );
1640
1641 if ( !driver )
1642 {
1643 QgsLogger::warning( "unable to get driver " + QString::number( i ) );
1644 continue;
1645 }
1646
1647 QString desc = GDALGetDescription( driver );
1648
1649 const QString mimeType = GDALGetMetadataItem( driver, "DMD_MIMETYPE", "" );
1650
1651 if ( mimeType.isEmpty() ) continue;
1652
1653 desc = desc.isEmpty() ? mimeType : desc;
1654
1655 QgsDebugMsg( "add GDAL format " + mimeType + ' ' + desc );
1656
1657 mimes[mimeType] = desc;
1658 }
1659 return mimes;
1660 }
1661
createProvider(const QString & uri,const QgsDataProvider::ProviderOptions & options,QgsDataProvider::ReadFlags flags)1662 QgsWcsProvider *QgsWcsProviderMetadata::createProvider( const QString &uri, const QgsDataProvider::ProviderOptions &options, QgsDataProvider::ReadFlags flags )
1663 {
1664 return new QgsWcsProvider( uri, options, flags );
1665 }
1666
1667 // ----------
1668
1669 int QgsWcsDownloadHandler::sErrors = 0;
1670
QgsWcsDownloadHandler(const QUrl & url,QgsWcsAuthorization & auth,QNetworkRequest::CacheLoadControl cacheLoadControl,QByteArray & cachedData,const QString & wcsVersion,QgsError & cachedError,QgsRasterBlockFeedback * feedback)1671 QgsWcsDownloadHandler::QgsWcsDownloadHandler( const QUrl &url, QgsWcsAuthorization &auth, QNetworkRequest::CacheLoadControl cacheLoadControl, QByteArray &cachedData, const QString &wcsVersion, QgsError &cachedError, QgsRasterBlockFeedback *feedback )
1672 : mAuth( auth )
1673 , mEventLoop( new QEventLoop )
1674 , mCachedData( cachedData )
1675 , mWcsVersion( wcsVersion )
1676 , mCachedError( cachedError )
1677 , mFeedback( feedback )
1678 {
1679 if ( feedback )
1680 {
1681 connect( feedback, &QgsFeedback::canceled, this, &QgsWcsDownloadHandler::canceled, Qt::QueuedConnection );
1682
1683 // rendering could have been canceled before we started to listen to canceled() signal
1684 // so let's check before doing the download and maybe quit prematurely
1685 if ( feedback->isCanceled() )
1686 return;
1687 }
1688
1689 QNetworkRequest request( url );
1690 QgsSetRequestInitiatorClass( request, QStringLiteral( "QgsWcsDownloadHandler" ) );
1691 if ( !mAuth.setAuthorization( request ) )
1692 {
1693 QgsMessageLog::logMessage( tr( "Network request update failed for authentication config" ),
1694 tr( "WCS" ) );
1695 return;
1696 }
1697 request.setAttribute( QNetworkRequest::CacheSaveControlAttribute, true );
1698 request.setAttribute( QNetworkRequest::CacheLoadControlAttribute, cacheLoadControl );
1699
1700 mCacheReply = QgsNetworkAccessManager::instance()->get( request );
1701 if ( !mAuth.setAuthorizationReply( mCacheReply ) )
1702 {
1703 mCacheReply->deleteLater();
1704 mCacheReply = nullptr;
1705 QgsMessageLog::logMessage( tr( "Network reply update failed for authentication config" ),
1706 tr( "WCS" ) );
1707 finish();
1708 return;
1709 }
1710 connect( mCacheReply, &QNetworkReply::finished, this, &QgsWcsDownloadHandler::cacheReplyFinished );
1711 connect( mCacheReply, &QNetworkReply::downloadProgress, this, &QgsWcsDownloadHandler::cacheReplyProgress );
1712 }
1713
~QgsWcsDownloadHandler()1714 QgsWcsDownloadHandler::~QgsWcsDownloadHandler()
1715 {
1716 delete mEventLoop;
1717 }
1718
blockingDownload()1719 void QgsWcsDownloadHandler::blockingDownload()
1720 {
1721 if ( mFeedback && mFeedback->isCanceled() )
1722 return; // nothing to do
1723
1724 mEventLoop->exec( QEventLoop::ExcludeUserInputEvents );
1725
1726 Q_ASSERT( !mCacheReply );
1727 }
1728
cacheReplyFinished()1729 void QgsWcsDownloadHandler::cacheReplyFinished()
1730 {
1731 QgsDebugMsg( QStringLiteral( "mCacheReply->error() = %1" ).arg( mCacheReply->error() ) );
1732 if ( mCacheReply->error() == QNetworkReply::NoError )
1733 {
1734 const QVariant redirect = mCacheReply->attribute( QNetworkRequest::RedirectionTargetAttribute );
1735 if ( !redirect.isNull() )
1736 {
1737 mCacheReply->deleteLater();
1738
1739 QgsDebugMsg( QStringLiteral( "redirected getmap: %1" ).arg( redirect.toString() ) );
1740 QNetworkRequest request( redirect.toUrl() );
1741 QgsSetRequestInitiatorClass( request, QStringLiteral( "QgsWcsDownloadHandler" ) );
1742 if ( !mAuth.setAuthorization( request ) )
1743 {
1744 QgsMessageLog::logMessage( tr( "Network request update failed for authentication config" ),
1745 tr( "WCS" ) );
1746 return;
1747 }
1748 mCacheReply = QgsNetworkAccessManager::instance()->get( request );
1749 if ( !mAuth.setAuthorizationReply( mCacheReply ) )
1750 {
1751 mCacheReply->deleteLater();
1752 mCacheReply = nullptr;
1753 QgsMessageLog::logMessage( tr( "Network reply update failed for authentication config" ),
1754 tr( "WCS" ) );
1755 finish();
1756 return;
1757 }
1758 connect( mCacheReply, &QNetworkReply::finished, this, &QgsWcsDownloadHandler::cacheReplyFinished );
1759 connect( mCacheReply, &QNetworkReply::downloadProgress, this, &QgsWcsDownloadHandler::cacheReplyProgress );
1760
1761 return;
1762 }
1763
1764 const QVariant status = mCacheReply->attribute( QNetworkRequest::HttpStatusCodeAttribute );
1765 QgsDebugMsg( QStringLiteral( "status = %1" ).arg( status.toInt() ) );
1766 if ( !status.isNull() && status.toInt() >= 400 )
1767 {
1768 const QVariant phrase = mCacheReply->attribute( QNetworkRequest::HttpReasonPhraseAttribute );
1769
1770 QgsMessageLog::logMessage( tr( "Map request error (Status: %1; Reason phrase: %2; URL: %3)" )
1771 .arg( status.toInt() )
1772 .arg( phrase.toString(),
1773 mCacheReply->url().toString() ), tr( "WCS" ) );
1774
1775 mCacheReply->deleteLater();
1776 mCacheReply = nullptr;
1777
1778 finish();
1779 return;
1780 }
1781
1782 // Read response
1783
1784 const QString contentType = mCacheReply->header( QNetworkRequest::ContentTypeHeader ).toString();
1785 QgsDebugMsg( "contentType: " + contentType );
1786
1787 // Exception
1788 // Content type examples: text/xml
1789 // application/vnd.ogc.se_xml;charset=UTF-8
1790 // application/xml
1791 if ( contentType.startsWith( QLatin1String( "text/" ), Qt::CaseInsensitive ) ||
1792 contentType.compare( QLatin1String( "application/xml" ), Qt::CaseInsensitive ) == 0 ||
1793 contentType.startsWith( QLatin1String( "application/vnd.ogc.se_xml" ), Qt::CaseInsensitive ) )
1794 {
1795 QString errorTitle, errorText;
1796 const QByteArray text = mCacheReply->readAll();
1797 if ( ( contentType.compare( QLatin1String( "text/xml" ), Qt::CaseInsensitive ) == 0 ||
1798 contentType.compare( QLatin1String( "application/xml" ), Qt::CaseInsensitive ) == 0 ||
1799 contentType.startsWith( QLatin1String( "application/vnd.ogc.se_xml" ), Qt::CaseInsensitive ) )
1800 && QgsWcsProvider::parseServiceExceptionReportDom( text, mWcsVersion, errorTitle, errorText ) )
1801 {
1802 mCachedError.append( SRVERR( tr( "Map request error:<br>Title: %1<br>Error: %2<br>URL: <a href='%3'>%3</a>)" )
1803 .arg( errorTitle, errorText,
1804 mCacheReply->url().toString() ) ) );
1805 }
1806 else
1807 {
1808 QgsMessageLog::logMessage( tr( "Map request error (Status: %1; Response: %2; URL: %3)" )
1809 .arg( status.toInt() )
1810 .arg( QString::fromUtf8( text ),
1811 mCacheReply->url().toString() ), tr( "WCS" ) );
1812 }
1813
1814 mCacheReply->deleteLater();
1815 mCacheReply = nullptr;
1816
1817 finish();
1818 return;
1819 }
1820
1821 // WCS 1.1 gives response as multipart, e.g.
1822 if ( QgsNetworkReplyParser::isMultipart( mCacheReply ) )
1823 {
1824 QgsDebugMsg( QStringLiteral( "reply is multipart" ) );
1825 const QgsNetworkReplyParser parser( mCacheReply );
1826
1827 if ( !parser.isValid() )
1828 {
1829 QgsMessageLog::logMessage( tr( "Cannot parse multipart response: %1" ).arg( parser.error() ), tr( "WCS" ) );
1830 mCacheReply->deleteLater();
1831 mCacheReply = nullptr;
1832
1833 finish();
1834 return;
1835 }
1836
1837 if ( parser.parts() < 2 )
1838 {
1839 QgsMessageLog::logMessage( tr( "Expected 2 parts, %1 received" ).arg( parser.parts() ), tr( "WCS" ) );
1840 mCacheReply->deleteLater();
1841 mCacheReply = nullptr;
1842
1843 finish();
1844 return;
1845 }
1846 else if ( parser.parts() > 2 )
1847 {
1848 // We will try the second one
1849 QgsMessageLog::logMessage( tr( "More than 2 parts (%1) received" ).arg( parser.parts() ), tr( "WCS" ) );
1850 }
1851
1852 const QString transferEncoding = parser.rawHeader( 1, QStringLiteral( "Content-Transfer-Encoding" ).toLatin1() );
1853 QgsDebugMsg( "transferEncoding = " + transferEncoding );
1854
1855 // It may happen (GeoServer) that in part header is for example
1856 // Content-Type: image/tiff and Content-Transfer-Encoding: base64
1857 // but content is xml ExceptionReport which is not in base64
1858 const QByteArray body = parser.body( 1 );
1859 if ( body.startsWith( "<?xml" ) )
1860 {
1861 QString errorTitle, errorText;
1862 if ( QgsWcsProvider::parseServiceExceptionReportDom( body, mWcsVersion, errorTitle, errorText ) )
1863 {
1864 QgsMessageLog::logMessage( tr( "Map request error (Title: %1; Error: %2; URL: %3)" )
1865 .arg( errorTitle, errorText,
1866 mCacheReply->url().toString() ), tr( "WCS" ) );
1867 }
1868 else
1869 {
1870 QgsMessageLog::logMessage( tr( "Map request error (Response: %1; URL: %2)" )
1871 .arg( QString::fromUtf8( body ),
1872 mCacheReply->url().toString() ), tr( "WCS" ) );
1873 }
1874
1875 mCacheReply->deleteLater();
1876 mCacheReply = nullptr;
1877
1878 finish();
1879 return;
1880 }
1881
1882 if ( transferEncoding == QLatin1String( "binary" ) )
1883 {
1884 mCachedData = body;
1885 }
1886 else if ( transferEncoding == QLatin1String( "base64" ) )
1887 {
1888 mCachedData = QByteArray::fromBase64( body );
1889 }
1890 else
1891 {
1892 QgsMessageLog::logMessage( tr( "Content-Transfer-Encoding %1 not supported" ).arg( transferEncoding ), tr( "WCS" ) );
1893 }
1894 }
1895 else
1896 {
1897 // Mime types for WCS 1.0 response should be usually image/<format>
1898 // (image/tiff, image/png, image/jpeg, image/png; mode=8bit, etc.)
1899 // but other mime types (like application/*) may probably also appear
1900
1901 // Create memory file
1902 mCachedData = mCacheReply->readAll();
1903 }
1904
1905 mCacheReply->deleteLater();
1906 mCacheReply = nullptr;
1907
1908 finish();
1909 }
1910 else
1911 {
1912 // report any errors except for the one we have caused by canceling the request
1913 if ( mCacheReply->error() != QNetworkReply::OperationCanceledError )
1914 {
1915 // Resend request if AlwaysCache
1916 QNetworkRequest request = mCacheReply->request();
1917 QgsSetRequestInitiatorClass( request, QStringLiteral( "QgsWcsDownloadHandler" ) );
1918 if ( request.attribute( QNetworkRequest::CacheLoadControlAttribute ).toInt() == QNetworkRequest::AlwaysCache )
1919 {
1920 QgsDebugMsg( QStringLiteral( "Resend request with PreferCache" ) );
1921 request.setAttribute( QNetworkRequest::CacheLoadControlAttribute, QNetworkRequest::PreferCache );
1922
1923 mCacheReply->deleteLater();
1924
1925 mCacheReply = QgsNetworkAccessManager::instance()->get( request );
1926 if ( !mAuth.setAuthorizationReply( mCacheReply ) )
1927 {
1928 mCacheReply->deleteLater();
1929 mCacheReply = nullptr;
1930 QgsMessageLog::logMessage( tr( "Network reply update failed for authentication config" ),
1931 tr( "WCS" ) );
1932 finish();
1933 return;
1934 }
1935 connect( mCacheReply, &QNetworkReply::finished, this, &QgsWcsDownloadHandler::cacheReplyFinished, Qt::DirectConnection );
1936 connect( mCacheReply, &QNetworkReply::downloadProgress, this, &QgsWcsDownloadHandler::cacheReplyProgress, Qt::DirectConnection );
1937
1938 return;
1939 }
1940
1941 sErrors++;
1942 if ( sErrors < 100 )
1943 {
1944 QgsMessageLog::logMessage( tr( "Map request failed [error: %1 url: %2]" ).arg( mCacheReply->errorString(), mCacheReply->url().toString() ), tr( "WCS" ) );
1945 }
1946 else if ( sErrors == 100 )
1947 {
1948 QgsMessageLog::logMessage( tr( "Not logging more than 100 request errors." ), tr( "WCS" ) );
1949 }
1950 }
1951
1952 mCacheReply->deleteLater();
1953 mCacheReply = nullptr;
1954
1955 finish();
1956 }
1957 }
1958
cacheReplyProgress(qint64 bytesReceived,qint64 bytesTotal)1959 void QgsWcsDownloadHandler::cacheReplyProgress( qint64 bytesReceived, qint64 bytesTotal )
1960 {
1961 Q_UNUSED( bytesReceived )
1962 Q_UNUSED( bytesTotal )
1963 QgsDebugMsgLevel( QStringLiteral( "%1 of %2 bytes of map downloaded." ).arg( bytesReceived ).arg( bytesTotal < 0 ? QString( "unknown number of" ) : QString::number( bytesTotal ) ), 3 );
1964 }
1965
canceled()1966 void QgsWcsDownloadHandler::canceled()
1967 {
1968 QgsDebugMsg( QStringLiteral( "Caught canceled() signal" ) );
1969 if ( mCacheReply )
1970 {
1971 QgsDebugMsg( QStringLiteral( "Aborting WCS network request" ) );
1972 mCacheReply->abort();
1973 }
1974 }
1975
1976
dataItemProviders() const1977 QList< QgsDataItemProvider * > QgsWcsProviderMetadata::dataItemProviders() const
1978 {
1979 QList<QgsDataItemProvider *> providers;
1980 providers << new QgsWcsDataItemProvider;
1981 return providers;
1982 }
1983
QgsWcsProviderMetadata()1984 QgsWcsProviderMetadata::QgsWcsProviderMetadata():
1985 QgsProviderMetadata( QgsWcsProvider::WCS_KEY, QgsWcsProvider::WCS_DESCRIPTION ) {}
1986
1987
1988 #ifndef HAVE_STATIC_PROVIDERS
providerMetadataFactory()1989 QGISEXTERN QgsProviderMetadata *providerMetadataFactory()
1990 {
1991 return new QgsWcsProviderMetadata();
1992 }
1993 #endif
1994