Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
Nicer API for raster sampling
  • Loading branch information
nyalldawson committed Jul 19, 2018
1 parent 54e5119 commit 74c2ed1
Show file tree
Hide file tree
Showing 6 changed files with 55 additions and 28 deletions.
12 changes: 7 additions & 5 deletions python/core/auto_generated/raster/qgsrasterdataprovider.sip.in
Expand Up @@ -277,17 +277,19 @@ if point is outside data source extent.
.. seealso:: :py:func:`sample`
%End

virtual QVariant sample( const QgsPointXY &point, int band, const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );
virtual double sample( const QgsPointXY &point, int band,
bool *ok /Out/ = 0,
const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );
%Docstring
Samples a raster value from the specified ``band`` found at the ``point`` position. The context
parameters ``boundingBox``, ``width`` and ``height`` are important to identify
on the same zoom level as a displayed map and to do effective
caching (WCS). If context params are not specified the highest
resolution is used. capabilities() may be used to test if format
is supported by provider.
resolution is used.

A null QVariant will be returned if the point is outside data source extent, or an invalid
band number was specified.
If ``ok`` is specified and the point is outside data source extent, or an invalid
band number was specified, then ``ok`` will be set to false. In this case the function will return
a NaN value.

.. seealso:: :py:func:`identify`

Expand Down
12 changes: 9 additions & 3 deletions src/core/raster/qgsrasterdataprovider.cpp
Expand Up @@ -300,12 +300,16 @@ QgsRasterIdentifyResult QgsRasterDataProvider::identify( const QgsPointXY &point
return QgsRasterIdentifyResult( QgsRaster::IdentifyFormatValue, results );
}

QVariant QgsRasterDataProvider::sample( const QgsPointXY &point, int band, const QgsRectangle &boundingBox, int width, int height, int )
double QgsRasterDataProvider::sample( const QgsPointXY &point, int band,
bool *ok, const QgsRectangle &boundingBox, int width, int height, int )
{
if ( ok )
*ok = false;

if ( !extent().contains( point ) )
{
// Outside the raster
return QVariant();
return std::numeric_limits<double>::quiet_NaN();
}

QgsRectangle finalExtent = boundingBox;
Expand Down Expand Up @@ -335,8 +339,10 @@ QVariant QgsRasterDataProvider::sample( const QgsPointXY &point, int band, const
QgsRectangle pixelExtent( xMin, yMin, xMax, yMax );

std::unique_ptr< QgsRasterBlock > bandBlock( block( band, pixelExtent, 1, 1 ) );
if ( bandBlock && ok )
*ok = true;

return bandBlock ? bandBlock->value( 0 ) : QVariant();
return bandBlock ? bandBlock->value( 0 ) : std::numeric_limits<double>::quiet_NaN();
}

QString QgsRasterDataProvider::lastErrorFormat()
Expand Down
12 changes: 7 additions & 5 deletions src/core/raster/qgsrasterdataprovider.h
Expand Up @@ -364,16 +364,18 @@ class CORE_EXPORT QgsRasterDataProvider : public QgsDataProvider, public QgsRast
* parameters \a boundingBox, \a width and \a height are important to identify
* on the same zoom level as a displayed map and to do effective
* caching (WCS). If context params are not specified the highest
* resolution is used. capabilities() may be used to test if format
* is supported by provider.
* resolution is used.
*
* A null QVariant will be returned if the point is outside data source extent, or an invalid
* band number was specified.
* If \a ok is specified and the point is outside data source extent, or an invalid
* band number was specified, then \a ok will be set to false. In this case the function will return
* a NaN value.
*
* \see identify(), which is much more flexible but considerably less efficient.
* \since QGIS 3.4
*/
virtual QVariant sample( const QgsPointXY &point, int band, const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );
virtual double sample( const QgsPointXY &point, int band,
bool *ok SIP_OUT = nullptr,
const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );

/**
* \brief Returns the caption error text for the last error in this provider
Expand Down
18 changes: 12 additions & 6 deletions src/providers/gdal/qgsgdalprovider.cpp
Expand Up @@ -1101,32 +1101,38 @@ bool QgsGdalProvider::worldToPixel( double x, double y, int &col, int &row ) con
return true;
};

QVariant QgsGdalProvider::sample( const QgsPointXY &point, int band, const QgsRectangle &, int, int, int )
double QgsGdalProvider::sample( const QgsPointXY &point, int band, bool *ok, const QgsRectangle &, int, int, int )
{
if ( ok )
*ok = false;

QMutexLocker locker( mpMutex );
if ( !initIfNeeded() )
return tr( "Cannot read data" );
return std::numeric_limits<double>::quiet_NaN();

if ( !extent().contains( point ) )
{
// Outside the raster
return QVariant();
return std::numeric_limits<double>::quiet_NaN();
}

GDALRasterBandH hBand = GDALGetRasterBand( mGdalDataset, band );
if ( !hBand )
return QVariant();
return std::numeric_limits<double>::quiet_NaN();

int row;
int col;
if ( !worldToPixel( point.x(), point.y(), col, row ) )
return QVariant();
return std::numeric_limits<double>::quiet_NaN();

float value = 0;
CPLErr err = GDALRasterIO( hBand, GF_Read, col, row, 1, 1,
&value, 1, 1, GDT_Float32, 0, 0 );
if ( err != CE_None )
return QVariant();
return std::numeric_limits<double>::quiet_NaN();

if ( ok )
*ok = true;

return static_cast< double >( value ) * bandScale( band ) + bandOffset( band );
}
Expand Down
2 changes: 1 addition & 1 deletion src/providers/gdal/qgsgdalprovider.h
Expand Up @@ -103,7 +103,7 @@ class QgsGdalProvider : public QgsRasterDataProvider, QgsGdalProviderBase
QgsRectangle extent() const override;
bool isValid() const override;
QgsRasterIdentifyResult identify( const QgsPointXY &point, QgsRaster::IdentifyFormat format, const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 ) override;
QVariant sample( const QgsPointXY &point, int band, const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );
double sample( const QgsPointXY &point, int band, bool *ok = nullptr, const QgsRectangle &boundingBox = QgsRectangle(), int width = 0, int height = 0, int dpi = 96 );
QString lastErrorTitle() override;
QString lastError() override;
int capabilities() const override;
Expand Down
27 changes: 19 additions & 8 deletions tests/src/core/testqgsrasterlayer.cpp
Expand Up @@ -729,21 +729,32 @@ void TestQgsRasterLayer::sample()
std::unique_ptr< QgsRasterLayer > rl = qgis::make_unique< QgsRasterLayer> ( rasterFileInfo.filePath(),
rasterFileInfo.completeBaseName() );
QVERIFY( rl->isValid() );
QVERIFY( !rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1 ).isValid() );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 1 ).toInt(), 125 );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1 ) ) );
bool ok = false;
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1, &ok ) ) );
QVERIFY( !ok );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 1 ), 125.0 );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 1, &ok ), 125.0 );
QVERIFY( ok );
// bad bands
QVERIFY( !rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 0 ).isValid() );
QVERIFY( !rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 10 ).isValid() );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 0 ) ) );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 0, &ok ) ) );
QVERIFY( !ok );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 788461, 3344957 ), 10, &ok ) ) );
QVERIFY( !ok );

fileName = mTestDataDir + "landsat_4326.tif";
rasterFileInfo = QFileInfo( fileName );
rl = qgis::make_unique< QgsRasterLayer> ( rasterFileInfo.filePath(),
rasterFileInfo.completeBaseName() );
QVERIFY( rl->isValid() );
QVERIFY( !rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1 ).isValid() );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 1 ).toInt(), 125 );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 2 ).toInt(), 139 );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 3 ).toInt(), 111 );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1 ) ) );
QVERIFY( std::isnan( rl->dataProvider()->sample( QgsPointXY( 0, 0 ), 1, &ok ) ) );
QVERIFY( !ok );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 1, &ok ), 125.0 );
QVERIFY( ok );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 2 ), 139.0 );
QCOMPARE( rl->dataProvider()->sample( QgsPointXY( 17.943731, 30.230791 ), 3 ), 111.0 );
}

QGSTEST_MAIN( TestQgsRasterLayer )
Expand Down

0 comments on commit 74c2ed1

Please sign in to comment.