wcs null values, enabled statistics

Author: blazek

src/core/qgsrasterdataprovider.cpp

@@ -243,6 +243,7 @@ QByteArray QgsRasterDataProvider::noValueBytes( int theBandNo )
 
 QgsRasterBandStats QgsRasterDataProvider::bandStatistics( int theBandNo )
 {
+  QgsDebugMsg( QString( "theBandNo = %1" ).arg( theBandNo ) );
   double myNoDataValue = noDataValue();
   QgsRasterBandStats myRasterBandStats;
   myRasterBandStats.elementCount = 0; // because we'll be counting only VALID pixels later
@@ -255,7 +256,7 @@ QgsRasterBandStats QgsRasterDataProvider::bandStatistics( int theBandNo )
   myXBlockSize = xBlockSize();
   myYBlockSize = yBlockSize();
 
-  if ( myXBlockSize == 0 || myYBlockSize == 0 )
+  if ( !( capabilities() & QgsRasterDataProvider::Size ) || xSize() == 0 || ySize() == 0 || myXBlockSize == 0 || myYBlockSize == 0 )
   {
     return QgsRasterBandStats(); //invalid raster band stats
   }
@@ -295,7 +296,7 @@ QgsRasterBandStats QgsRasterDataProvider::bandStatistics( int theBandNo )
         for ( int iX = 0; iX < nXValid; iX++ )
         {
           double myValue = readValue( myData, myDataType, iX + ( iY * myXBlockSize ) );
-          //QgsDebugMsg ( QString ( "%1 %2 value %3" ).arg (iX).arg(iY).arg( myValue ) );
+          QgsDebugMsgLevel( QString( "%1 %2 value %3" ).arg( iX ).arg( iY ).arg( myValue ), 10 );
 
           if ( mValidNoDataValue && ( qAbs( myValue - myNoDataValue ) <= TINY_VALUE ) )
           {
@@ -381,14 +382,13 @@ QgsRasterBandStats QgsRasterDataProvider::bandStatistics( int theBandNo )
                              ( myRasterBandStats.elementCount - 1 ) ) );
 
 #ifdef QGISDEBUG
-  QgsLogger::debug( "************ STATS **************", 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "VALID NODATA", mValidNoDataValue, 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "NULL", noDataValue() , 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "MIN", myRasterBandStats.minimumValue, 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "MAX", myRasterBandStats.maximumValue, 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "RANGE", myRasterBandStats.range, 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "MEAN", myRasterBandStats.mean, 1, __FILE__, __FUNCTION__, __LINE__ );
-  QgsLogger::debug( "STDDEV", myRasterBandStats.stdDev, 1, __FILE__, __FUNCTION__, __LINE__ );
+  QgsDebugMsg( "************ STATS **************" );
+  QgsDebugMsg( QString( "VALID NODATA %1" ).arg( mValidNoDataValue ) );
+  QgsDebugMsg( QString( "MIN %1" ).arg( myRasterBandStats.minimumValue ) );
+  QgsDebugMsg( QString( "MAX %1" ).arg( myRasterBandStats.maximumValue ) );
+  QgsDebugMsg( QString( "RANGE %1" ).arg( myRasterBandStats.range ) );
+  QgsDebugMsg( QString( "MEAN %1" ).arg( myRasterBandStats.mean ) );
+  QgsDebugMsg( QString( "STDDEV %1" ).arg( myRasterBandStats.stdDev ) );
 #endif
 
   CPLFree( myData );

src/gui/qgsowssourceselect.cpp

@@ -221,14 +221,15 @@ void QgsOWSSourceSelect::populateFormats()
       QgsDebugMsg( QString( "format %1 not supported." ).arg( format ) );
       //btn->setEnabled( false );
       btn->setEnabled( true );
+      if ( firstEnabled < 0 ) { firstEnabled = i; }
       tip += " " + tr( "is not supported by GDAL" );
     }
     btn->setText( label );
     btn->setToolTip( tip );
   }
   // Set prefered
   // TODO: all enabled for now, see above
-  //prefered = prefered >= 0 ? prefered : firstEnabled;
+  prefered = prefered >= 0 ? prefered : firstEnabled;
   if ( prefered >= 0 )
   {
     mImageFormatGroup->button( prefered )->setChecked( true );

src/providers/wcs/qgswcsprovider.cpp

@@ -71,6 +71,8 @@ static QString WCS_DESCRIPTION = "OGC Web Coverage Service version 1.0/1.1 data
 
 static QString DEFAULT_LATLON_CRS = "CRS:84";
 
+// TODO: colortable - use comon baseclass with gdal, mapserver does not support http://trac.osgeo.org/mapserver/ticket/1671
+
 QgsWcsProvider::QgsWcsProvider( QString const &uri )
     : QgsRasterDataProvider( uri )
     , mHttpUri( QString::null )
@@ -88,6 +90,8 @@ QgsWcsProvider::QgsWcsProvider( QString const &uri )
     , mCachedMemFile( 0 )
     , mWidth( 0 )
     , mHeight( 0 )
+    , mXBlockSize( 0 )
+    , mYBlockSize( 0 )
     , mHasSize( false )
     , mBandCount( 0 )
 {
@@ -174,6 +178,12 @@ QgsWcsProvider::QgsWcsProvider( QString const &uri )
     mSrcGdalDataType.append( myGdalDataType );
     // TODO: This could be shared with GDAL provider
     int isValid = false;
+
+    // UMN Mapserver does not automaticaly set null value, METADATA wcs_rangeset_nullvalue must be used
+    // http://lists.osgeo.org/pipermail/mapserver-users/2010-April/065328.html
+
+    // TODO:
+
     double myNoDataValue = GDALGetRasterNoDataValue( gdalBand, &isValid );
     if ( isValid )
     {
@@ -216,6 +226,7 @@ QgsWcsProvider::QgsWcsProvider( QString const &uri )
       }
     }
     mNoDataValue.append( myNoDataValue );
+    mValidNoDataValue = true;
 
     QgsDebugMsg( QString( "mSrcGdalDataType[%1] = %2" ).arg( i - 1 ).arg( mSrcGdalDataType[i-1] ) );
     QgsDebugMsg( QString( "mGdalDataType[%1] = %2" ).arg( i - 1 ).arg( mGdalDataType[i-1] ) );
@@ -224,6 +235,18 @@ QgsWcsProvider::QgsWcsProvider( QString const &uri )
 
   clearCache();
 
+  // Block size is used for for statistics
+  // TODO: How to find maximum block size supported by server?
+  if ( mHasSize )
+  {
+    // This is taken from GDAL, how they come to these numbers?
+    if ( mWidth > 1800 ) mXBlockSize = 1024;
+    else mXBlockSize = mWidth;
+
+    if ( mHeight > 900 ) mYBlockSize = 512;
+    else mYBlockSize = mHeight;
+  }
+
   mValid = true;
   QgsDebugMsg( "Constructed ok, provider valid." );
 }
@@ -491,6 +514,32 @@ void QgsWcsProvider::getCache( int bandNo, QgsRectangle  const & viewExtent, int
 
 }
 
+// For stats only, maybe change QgsRasterDataProvider::bandStatistics() to
+// use standard readBlock with extent
+void QgsWcsProvider::readBlock( int theBandNo, int xBlock, int yBlock, void *block )
+{
+  QgsDebugMsg( "Entered" );
+
+  QgsDebugMsg( QString( "xBlock = %1 yBlock = %2" ).arg( xBlock ).arg( yBlock ) );
+
+  if ( !mHasSize ) return;
+
+  double xRes = mCoverageExtent.width() / mWidth;
+  double yRes = mCoverageExtent.height() / mHeight;
+
+  // blocks on edges may run out of extent, that should not be problem (at least for
+  // stats - there is a check for it)
+  double xMin = mCoverageExtent.xMinimum() + xRes * xBlock * mXBlockSize;
+  double xMax = xMin + xRes * mXBlockSize;
+  double yMax = mCoverageExtent.yMaximum() - yRes * yBlock * mYBlockSize;
+  double yMin = yMax - yRes * mYBlockSize;
+  //QgsDebugMsg( QString("yMin = %1 yMax = %2").arg(yMin).arg(yMax) );
+
+  QgsRectangle extent( xMin, yMin, xMax, yMax );
+
+  readBlock( theBandNo, extent, mXBlockSize, mYBlockSize, block );
+}
+
 void QgsWcsProvider::cacheReplyFinished()
 {
   if ( mCacheReply->error() == QNetworkReply::NoError )
@@ -685,6 +734,29 @@ int QgsWcsProvider::bandCount() const
   return mBandCount;
 }
 
+double  QgsWcsProvider::noDataValue() const
+{
+  if ( mNoDataValue.size() > 0 )
+  {
+    return mNoDataValue[0];
+  }
+  return std::numeric_limits<int>::max(); // should not happen or be used
+}
+
+// this is only called once when statistics are calculated
+// TODO
+int QgsWcsProvider::xBlockSize() const
+{
+  return mXBlockSize;
+}
+int QgsWcsProvider::yBlockSize() const
+{
+  return mYBlockSize;
+}
+
+int QgsWcsProvider::xSize() const { return mWidth; }
+int QgsWcsProvider::ySize() const { return mHeight; }
+
 void QgsWcsProvider::clearCache()
 {
   QgsDebugMsg( "Entered" );

src/providers/wcs/qgswcsprovider.h

@@ -118,6 +118,8 @@ class QgsWcsProvider : public QgsRasterDataProvider
 
     void readBlock( int bandNo, QgsRectangle  const & viewExtent, int width, int height, void *data );
 
+    void readBlock( int theBandNo, int xBlock, int yBlock, void *block );
+
     /** Download cache */
     void getCache( int bandNo, QgsRectangle  const & viewExtent, int width, int height );
 
@@ -141,6 +143,11 @@ class QgsWcsProvider : public QgsRasterDataProvider
     int dataType( int bandNo ) const;
     int srcDataType( int bandNo ) const;
     int bandCount() const;
+    double noDataValue() const;
+    int xBlockSize() const;
+    int yBlockSize() const;
+    int xSize() const;
+    int ySize() const;
     QString metadata();
     QString identifyAsHtml( const QgsPoint& point );
     QString identifyAsText( const QgsPoint& point );
@@ -266,6 +273,10 @@ class QgsWcsProvider : public QgsRasterDataProvider
     /** Coverage width, may be 0 if it could not be found in DescribeCoverage */
     int mHeight;
 
+    /** Block size */
+    int mXBlockSize;
+    int mYBlockSize;
+
     /** Flag if size was parsed successfully */
     bool mHasSize;