Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
[FEATURE] raster_statistic expression function for retrieving
raster band stats from a loaded layer

Allows raster band stats (eg min, max, avg) to be used in
expressions
  • Loading branch information
nyalldawson committed Nov 22, 2016
1 parent 86ab302 commit 20dc7fb
Show file tree
Hide file tree
Showing 3 changed files with 112 additions and 0 deletions.
14 changes: 14 additions & 0 deletions resources/function_help/json/raster_statistic
@@ -0,0 +1,14 @@
{
"name": "raster_statistic",
"type": "function",
"description": "Returns statistics from a raster layer.",
"arguments": [
{"arg":"layer", "description":"a string, representing either a raster layer name or layer ID"},
{"arg":"band", "description":"integer representing the band number from the raster layer, starting at 1"},
{"arg":"property", "description":"a string corresponding to the property to return. Valid options are:<br /><ul><li>min: minimum value</li><li>max: maximum value</li><li>avg: average (mean) value</li><li>stdev: standard deviation of values</li><li>range: range of values (max - min)</li><li>sum: sum of all values from raster</li></ul>"}
],
"examples": [
{ "expression":"raster_statistic('lc',1,'avg')", "returns":"Average value from band 1 from 'lc' raster layer"},
{ "expression":"raster_statistic('ac2010',3,'min')", "returns":"Minimum value from band 3 from 'ac2010' raster layer"}
]
}
74 changes: 74 additions & 0 deletions src/core/qgsexpression.cpp
Expand Up @@ -54,6 +54,8 @@
#include "qgsmaptopixelgeometrysimplifier.h"
#include "qgsmessagelog.h"
#include "qgscsexception.h"
#include "qgsrasterlayer.h"
#include "qgsrasterdataprovider.h"

// from parser
extern QgsExpression::Node* parseExpression( const QString& str, QString& parserErrorMsg );
Expand Down Expand Up @@ -3423,6 +3425,75 @@ static QVariant fcnGetLayerProperty( const QVariantList& values, const QgsExpres
return QVariant();
}

static QVariant fcnGetRasterBandStat( const QVariantList& values, const QgsExpressionContext*, QgsExpression* parent )
{
QString layerIdOrName = getStringValue( values.at( 0 ), parent );

//try to find a matching layer by name
QgsMapLayer* layer = QgsMapLayerRegistry::instance()->mapLayer( layerIdOrName ); //search by id first
if ( !layer )
{
QList<QgsMapLayer *> layersByName = QgsMapLayerRegistry::instance()->mapLayersByName( layerIdOrName );
if ( !layersByName.isEmpty() )
{
layer = layersByName.at( 0 );
}
}

if ( !layer )
return QVariant();

QgsRasterLayer* rl = qobject_cast< QgsRasterLayer* >( layer );
if ( !rl )
return QVariant();

int band = getIntValue( values.at( 1 ), parent );
if ( band < 1 || band > rl->bandCount() )
{
parent->setEvalErrorString( QObject::tr( "Invalid band number %1 for layer %2" ).arg( band ).arg( layerIdOrName ) );
return QVariant();
}

QString layerProperty = getStringValue( values.at( 2 ), parent );
int stat = 0;

if ( QString::compare( layerProperty, QStringLiteral( "avg" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::Mean;
else if ( QString::compare( layerProperty, QStringLiteral( "stdev" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::StdDev;
else if ( QString::compare( layerProperty, QStringLiteral( "min" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::Min;
else if ( QString::compare( layerProperty, QStringLiteral( "max" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::Max;
else if ( QString::compare( layerProperty, QStringLiteral( "range" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::Range;
else if ( QString::compare( layerProperty, QStringLiteral( "sum" ), Qt::CaseInsensitive ) == 0 )
stat = QgsRasterBandStats::Sum;
else
{
parent->setEvalErrorString( QObject::tr( "Invalid raster statistic: '%1'" ).arg( layerProperty ) );
return QVariant();
}

QgsRasterBandStats stats = rl->dataProvider()->bandStatistics( band, stat );
switch ( stat )
{
case QgsRasterBandStats::Mean:
return stats.mean;
case QgsRasterBandStats::StdDev:
return stats.stdDev;
case QgsRasterBandStats::Min:
return stats.minimumValue;
case QgsRasterBandStats::Max:
return stats.maximumValue;
case QgsRasterBandStats::Range:
return stats.range;
case QgsRasterBandStats::Sum:
return stats.sum;
}
return QVariant();
}

static QVariant fcnArray( const QVariantList& values, const QgsExpressionContext*, QgsExpression* )
{
return values;
Expand Down Expand Up @@ -4007,6 +4078,9 @@ const QList<QgsExpression::Function*>& QgsExpression::Functions()
// **General** functions

<< new StaticFunction( QStringLiteral( "layer_property" ), 2, fcnGetLayerProperty, QStringLiteral( "General" ) )
<< new StaticFunction( QStringLiteral( "raster_statistic" ), ParameterList() << Parameter( QStringLiteral( "layer" ) )
<< Parameter( QStringLiteral( "band" ) )
<< Parameter( QStringLiteral( "statistic" ) ), fcnGetRasterBandStat, QStringLiteral( "General" ) )
<< new StaticFunction( QStringLiteral( "var" ), 1, fcnGetVariable, QStringLiteral( "General" ) )

//return all attributes string for referencedColumns - this is caught by
Expand Down
24 changes: 24 additions & 0 deletions tests/src/core/testqgsexpression.cpp
Expand Up @@ -32,6 +32,7 @@
#include "qgsmaplayerregistry.h"
#include "qgsvectordataprovider.h"
#include "qgsdistancearea.h"
#include "qgsrasterlayer.h"
#include "qgsproject.h"

static void _parseAndEvalExpr( int arg )
Expand All @@ -55,6 +56,7 @@ class TestQgsExpression: public QObject
, mMemoryLayer( nullptr )
, mAggregatesLayer( nullptr )
, mChildLayer( nullptr )
, mRasterLayer( nullptr )
{}

private:
Expand All @@ -63,6 +65,7 @@ class TestQgsExpression: public QObject
QgsVectorLayer* mMemoryLayer;
QgsVectorLayer* mAggregatesLayer;
QgsVectorLayer* mChildLayer;
QgsRasterLayer* mRasterLayer;

private slots:

Expand Down Expand Up @@ -94,6 +97,12 @@ class TestQgsExpression: public QObject
mPointsLayer->setMaximumScale( 500 );
mPointsLayer->setMinimumScale( 1000 );

QString rasterFileName = testDataDir + "tenbytenraster.asc";
QFileInfo rasterFileInfo( rasterFileName );
mRasterLayer = new QgsRasterLayer( rasterFileInfo.filePath(),
rasterFileInfo.completeBaseName() );
QgsMapLayerRegistry::instance()->addMapLayer( mRasterLayer );

// test memory layer for get_feature tests
mMemoryLayer = new QgsVectorLayer( QStringLiteral( "Point?field=col1:integer&field=col2:string" ), QStringLiteral( "test" ), QStringLiteral( "memory" ) );
QVERIFY( mMemoryLayer->isValid() );
Expand Down Expand Up @@ -1039,6 +1048,21 @@ class TestQgsExpression: public QObject
QTest::newRow( "layer_property storage_type" ) << QStringLiteral( "layer_property('%1','storage_type')" ).arg( mPointsLayer->name() ) << false << QVariant( "ESRI Shapefile" );
QTest::newRow( "layer_property geometry_type" ) << QStringLiteral( "layer_property('%1','geometry_type')" ).arg( mPointsLayer->name() ) << false << QVariant( "Point" );

// raster_statistic tests
QTest::newRow( "raster_statistic no layer" ) << "raster_statistic('',1,'min')" << false << QVariant();
QTest::newRow( "raster_statistic bad layer" ) << "raster_statistic('bad',1,'min')" << false << QVariant();
QTest::newRow( "raster_statistic bad band" ) << QStringLiteral( "raster_statistic('%1',0,'min')" ).arg( mRasterLayer->name() ) << true << QVariant();
QTest::newRow( "raster_statistic bad band 2" ) << QStringLiteral( "raster_statistic('%1',100,'min')" ).arg( mRasterLayer->name() ) << true << QVariant();
QTest::newRow( "raster_statistic no property" ) << QStringLiteral( "raster_statistic('%1',1,'')" ).arg( mRasterLayer->name() ) << true << QVariant();
QTest::newRow( "raster_statistic bad property" ) << QStringLiteral( "raster_statistic('%1',1,'bad')" ).arg( mRasterLayer->name() ) << true << QVariant();
QTest::newRow( "raster_statistic min by id" ) << QStringLiteral( "raster_statistic('%1',1,'min')" ).arg( mRasterLayer->id() ) << false << QVariant( 0.0 );
QTest::newRow( "raster_statistic min name" ) << QStringLiteral( "raster_statistic('%1',1,'min')" ).arg( mRasterLayer->name() ) << false << QVariant( 0.0 );
QTest::newRow( "raster_statistic max" ) << QStringLiteral( "raster_statistic('%1',1,'max')" ).arg( mRasterLayer->id() ) << false << QVariant( 9.0 );
QTest::newRow( "raster_statistic avg" ) << QStringLiteral( "round(10*raster_statistic('%1',1,'avg'))" ).arg( mRasterLayer->id() ) << false << QVariant( 45 );
QTest::newRow( "raster_statistic stdev" ) << QStringLiteral( "round(100*raster_statistic('%1',1,'stdev'))" ).arg( mRasterLayer->id() ) << false << QVariant( 287 );
QTest::newRow( "raster_statistic range" ) << QStringLiteral( "raster_statistic('%1',1,'range')" ).arg( mRasterLayer->id() ) << false << QVariant( 9.0 );
QTest::newRow( "raster_statistic sum" ) << QStringLiteral( "round(raster_statistic('%1',1,'sum'))" ).arg( mRasterLayer->id() ) << false << QVariant( 450 );

//test conversions to bool
QTest::newRow( "feature to bool false" ) << QStringLiteral( "case when get_feature('none','none',499) then true else false end" ) << false << QVariant( false );
QTest::newRow( "feature to bool true" ) << QStringLiteral( "case when get_feature('test','col1',10) then true else false end" ) << false << QVariant( true );
Expand Down

1 comment on commit 20dc7fb

@nirvn
Copy link
Contributor

@nirvn nirvn commented on 20dc7fb Nov 22, 2016

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's a really useful function, thanks. It'd be really cool for it to support an optional 4th parameter that would take a geometry (point, line, or polygon) to limit statistics to pixels overlapping the given geometry.

Please sign in to comment.