Skip to content

Commit 83f286a

Browse files
elpasoelpaso
authored
Merge pull request #9305 from elpaso/bugfix-21405-raster-calc-wrong-results
[opencl] Fix raster calculator operator precedence
2 parents 7d83263 + 88a9612 commit 83f286a

File tree

3 files changed

+106
-55
lines changed

3 files changed

+106
-55
lines changed

‎src/analysis/raster/qgsrastercalcnode.cpp

Lines changed: 16 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -210,22 +210,17 @@ QString QgsRasterCalcNode::toString( bool cStyle ) const
210210
if ( mRight )
211211
right = mRight->toString( cStyle );
212212

213-
auto floatCast = [ ]( const QString s ) -> QString
214-
{
215-
return QStringLiteral( "(float) ( %1 )" ).arg( s );
216-
};
217-
218213
switch ( mType )
219214
{
220215
case tOperator:
221216
switch ( mOperator )
222217
{
223218
case opPLUS:
224-
result = QStringLiteral( "%1 + %2" ).arg( left ).arg( right );
219+
result = QStringLiteral( "( %1 + %2 )" ).arg( left ).arg( right );
225220
break;
226221
case opMINUS:
227222
case opSIGN:
228-
result = QStringLiteral( "%1 - %2" ).arg( left ).arg( right );
223+
result = QStringLiteral( "( %1 - %2 )" ).arg( left ).arg( right );
229224
break;
230225
case opMUL:
231226
result = QStringLiteral( "%1 * %2" ).arg( left ).arg( right );
@@ -235,7 +230,7 @@ QString QgsRasterCalcNode::toString( bool cStyle ) const
235230
break;
236231
case opPOW:
237232
if ( cStyle )
238-
result = QStringLiteral( "pow( %1, %2 )" ).arg( floatCast( left ) ).arg( floatCast( right ) );
233+
result = QStringLiteral( "pow( %1, %2 )" ).arg( left ).arg( right );
239234
else
240235
result = QStringLiteral( "%1^%2" ).arg( left ).arg( right );
241236
break;
@@ -273,58 +268,31 @@ QString QgsRasterCalcNode::toString( bool cStyle ) const
273268
result = QStringLiteral( "%1 OR %2" ).arg( left ).arg( right );
274269
break;
275270
case opSQRT:
276-
if ( cStyle )
277-
result = QStringLiteral( "sqrt( %1 )" ).arg( floatCast( left ) );
278-
else
279-
result = QStringLiteral( "sqrt( %1 )" ).arg( left );
271+
result = QStringLiteral( "sqrt( %1 )" ).arg( left );
280272
break;
281273
case opSIN:
282-
if ( cStyle )
283-
result = QStringLiteral( "sin( %1 )" ).arg( floatCast( left ) );
284-
else
285-
result = QStringLiteral( "sin( %1 )" ).arg( left );
274+
result = QStringLiteral( "sin( %1 )" ).arg( left );
286275
break;
287276
case opCOS:
288-
if ( cStyle )
289-
result = QStringLiteral( "cos( %1 )" ).arg( floatCast( left ) );
290-
else
291-
result = QStringLiteral( "cos( %1 )" ).arg( left );
277+
result = QStringLiteral( "cos( %1 )" ).arg( left );
292278
break;
293279
case opTAN:
294-
if ( cStyle )
295-
result = QStringLiteral( "tan( %1 )" ).arg( floatCast( left ) );
296-
else
297-
result = QStringLiteral( "tan( %1 )" ).arg( left );
280+
result = QStringLiteral( "tan( %1 )" ).arg( left );
298281
break;
299282
case opASIN:
300-
if ( cStyle )
301-
result = QStringLiteral( "asin( %1 )" ).arg( floatCast( left ) );
302-
else
303-
result = QStringLiteral( "asin( %1 )" ).arg( left );
283+
result = QStringLiteral( "asin( %1 )" ).arg( left );
304284
break;
305285
case opACOS:
306-
if ( cStyle )
307-
result = QStringLiteral( "acos( %1 )" ).arg( floatCast( left ) );
308-
else
309-
result = QStringLiteral( "acos( %1 )" ).arg( left );
286+
result = QStringLiteral( "acos( %1 )" ).arg( left );
310287
break;
311288
case opATAN:
312-
if ( cStyle )
313-
result = QStringLiteral( "atan( %1 )" ).arg( floatCast( left ) );
314-
else
315-
result = QStringLiteral( "atan( %1 )" ).arg( left );
289+
result = QStringLiteral( "atan( %1 )" ).arg( left );
316290
break;
317291
case opLOG:
318-
if ( cStyle )
319-
result = QStringLiteral( "log( %1 )" ).arg( floatCast( left ) );
320-
else
321-
result = QStringLiteral( "log( %1 )" ).arg( left );
292+
result = QStringLiteral( "log( %1 )" ).arg( left );
322293
break;
323294
case opLOG10:
324-
if ( cStyle )
325-
result = QStringLiteral( "log10( %1 )" ).arg( floatCast( left ) );
326-
else
327-
result = QStringLiteral( "log10( %1 )" ).arg( left );
295+
result = QStringLiteral( "log10( %1 )" ).arg( left );
328296
break;
329297
case opNONE:
330298
break;
@@ -335,6 +303,10 @@ QString QgsRasterCalcNode::toString( bool cStyle ) const
335303
break;
336304
case tNumber:
337305
result = QString::number( mNumber );
306+
if ( cStyle )
307+
{
308+
result = QStringLiteral( "(float) ( %1 )" ).arg( result );
309+
}
338310
break;
339311
case tMatrix:
340312
break;

‎tests/src/analysis/testqgsrastercalculator.cpp

Lines changed: 90 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,10 @@ Email : nyall dot dawson at gmail dot com
1414
***************************************************************************/
1515
#include "qgstest.h"
1616

17+
#ifdef HAVE_OPENCL
18+
#include "qgsopenclutils.h"
19+
#endif
20+
1721
#include "qgsrastercalculator.h"
1822
#include "qgsrastercalcnode.h"
1923
#include "qgsrasterdataprovider.h"
@@ -59,22 +63,29 @@ class TestQgsRasterCalculator : public QObject
5963
void findNodes();
6064

6165
void testRasterEntries();
66+
void calcFormulasWithReprojectedLayers();
6267

6368
private:
6469

6570
QgsRasterLayer *mpLandsatRasterLayer = nullptr;
6671
QgsRasterLayer *mpLandsatRasterLayer4326 = nullptr;
6772
};
6873

74+
6975
void TestQgsRasterCalculator::initTestCase()
7076
{
7177
//
7278
// Runs once before any tests are run
7379
//
74-
// init QGIS's paths - true means that all path will be inited from prefix
80+
// Set up the QgsSettings environment
81+
QCoreApplication::setOrganizationName( QStringLiteral( "QGIS" ) );
82+
QCoreApplication::setOrganizationDomain( QStringLiteral( "qgis.org" ) );
83+
QCoreApplication::setApplicationName( QStringLiteral( "QGIS-TEST" ) );
84+
7585
QgsApplication::init();
7686
QgsApplication::initQgis();
7787

88+
7889
QString testDataDir = QStringLiteral( TEST_DATA_DIR ) + '/'; //defined in CmakeLists.txt
7990

8091
QString landsatFileName = testDataDir + "landsat.tif";
@@ -99,8 +110,13 @@ void TestQgsRasterCalculator::cleanupTestCase()
99110

100111
void TestQgsRasterCalculator::init()
101112
{
102-
113+
#ifdef HAVE_OPENCL
114+
QgsOpenClUtils::setEnabled( false );
115+
// Reset to default in case some tests mess it up
116+
QgsOpenClUtils::setSourcePath( QDir( QgsApplication::pkgDataPath() ).absoluteFilePath( QStringLiteral( "resources/opencl_programs" ) ) );
117+
#endif
103118
}
119+
104120
void TestQgsRasterCalculator::cleanup()
105121
{
106122

@@ -676,12 +692,78 @@ void TestQgsRasterCalculator::toString()
676692
return error;
677693
return calcNode->toString( cStyle );
678694
};
679-
QCOMPARE( _test( QStringLiteral( "\"raster@1\" + 2" ), false ), QString( "\"raster@1\" + 2" ) );
680-
QCOMPARE( _test( QStringLiteral( "\"raster@1\" + 2" ), true ), QString( "\"raster@1\" + 2" ) );
681-
QCOMPARE( _test( QStringLiteral( "\"raster@1\" ^ 3 + 2" ), false ), QString( "\"raster@1\"^3 + 2" ) );
682-
QCOMPARE( _test( QStringLiteral( "\"raster@1\" ^ 3 + 2" ), true ), QString( "pow( (float) ( \"raster@1\" ), (float) ( 3 ) ) + 2" ) );
683-
QCOMPARE( _test( QStringLiteral( "atan(\"raster@1\") * cos( 3 + 2 )" ), false ), QString( "atan( \"raster@1\" ) * cos( 3 + 2 )" ) );
684-
QCOMPARE( _test( QStringLiteral( "atan(\"raster@1\") * cos( 3 + 2 )" ), true ), QString( "atan( (float) ( \"raster@1\" ) ) * cos( (float) ( 3 + 2 ) )" ) );
695+
QCOMPARE( _test( QStringLiteral( "\"raster@1\" + 2" ), false ), QString( "( \"raster@1\" + 2 )" ) );
696+
QCOMPARE( _test( QStringLiteral( "\"raster@1\" + 2" ), true ), QString( "( \"raster@1\" + (float) ( 2 ) )" ) );
697+
QCOMPARE( _test( QStringLiteral( "\"raster@1\" ^ 3 + 2" ), false ), QString( "( \"raster@1\"^3 + 2 )" ) );
698+
QCOMPARE( _test( QStringLiteral( "\"raster@1\" ^ 3 + 2" ), true ), QString( "( pow( \"raster@1\", (float) ( 3 ) ) + (float) ( 2 ) )" ) );
699+
QCOMPARE( _test( QStringLiteral( "atan(\"raster@1\") * cos( 3 + 2 )" ), false ), QString( "atan( \"raster@1\" ) * cos( ( 3 + 2 ) )" ) );
700+
QCOMPARE( _test( QStringLiteral( "atan(\"raster@1\") * cos( 3 + 2 )" ), true ), QString( "atan( \"raster@1\" ) * cos( ( (float) ( 3 ) + (float) ( 2 ) ) )" ) );
701+
QCOMPARE( _test( QStringLiteral( "0.5 * ( 1.4 * (\"raster@1\" + 2) )" ), false ), QString( "0.5 * 1.4 * ( \"raster@1\" + 2 )" ) );
702+
QCOMPARE( _test( QStringLiteral( "0.5 * ( 1.4 * (\"raster@1\" + 2) )" ), true ), QString( "(float) ( 0.5 ) * (float) ( 1.4 ) * ( \"raster@1\" + (float) ( 2 ) )" ) );
703+
}
704+
705+
void TestQgsRasterCalculator::calcFormulasWithReprojectedLayers()
706+
{
707+
QgsRasterCalculatorEntry entry1;
708+
entry1.bandNumber = 1;
709+
entry1.raster = mpLandsatRasterLayer;
710+
entry1.ref = QStringLiteral( "landsat@1" );
711+
712+
QgsRasterCalculatorEntry entry2;
713+
entry2.bandNumber = 2;
714+
entry2.raster = mpLandsatRasterLayer4326;
715+
entry2.ref = QStringLiteral( "landsat_4326@2" );
716+
717+
QVector<QgsRasterCalculatorEntry> entries;
718+
entries << entry1 << entry2;
719+
720+
QgsCoordinateReferenceSystem crs;
721+
crs.createFromId( 32633, QgsCoordinateReferenceSystem::EpsgCrsId );
722+
QgsRectangle extent( 783235, 3348110, 783350, 3347960 );
723+
724+
725+
auto _chk = [ = ]( const QString & formula, const std::vector<float> &values, bool useOpenCL )
726+
{
727+
728+
#ifdef HAVE_OPENCL
729+
if ( ! QgsOpenClUtils::available() )
730+
return ;
731+
QgsOpenClUtils::setEnabled( useOpenCL );
732+
#endif
733+
734+
QTemporaryFile tmpFile;
735+
tmpFile.open(); // fileName is not available until open
736+
QString tmpName = tmpFile.fileName();
737+
tmpFile.close();
738+
QgsRasterCalculator rc( formula,
739+
tmpName,
740+
QStringLiteral( "GTiff" ),
741+
extent, crs, 2, 3, entries );
742+
QCOMPARE( static_cast< int >( rc.processCalculation() ), 0 );
743+
//open output file and check results
744+
QgsRasterLayer *result = new QgsRasterLayer( tmpName, QStringLiteral( "result" ) );
745+
QCOMPARE( result->width(), 2 );
746+
QCOMPARE( result->height(), 3 );
747+
QgsRasterBlock *block = result->dataProvider()->block( 1, extent, 2, 3 );
748+
qDebug() << block->value( 0, 0 ) << block->value( 0, 1 ) << block->value( 1, 0 ) << block->value( 1, 1 ) << block->value( 2, 0 ) << block->value( 2, 1 );
749+
const float epsilon { 0.0001f };
750+
QVERIFY( std::abs( ( static_cast<float>( block->value( 0, 0 ) ) - values[0] ) / values[0] ) < epsilon );
751+
QVERIFY( std::abs( ( static_cast<float>( block->value( 0, 1 ) ) - values[1] ) / values[1] ) < epsilon );
752+
QVERIFY( std::abs( ( static_cast<float>( block->value( 1, 0 ) ) - values[2] ) / values[2] ) < epsilon );
753+
QVERIFY( std::abs( ( static_cast<float>( block->value( 1, 1 ) ) - values[3] ) / values[3] ) < epsilon );
754+
QVERIFY( std::abs( ( static_cast<float>( block->value( 2, 0 ) ) - values[4] ) / values[4] ) < epsilon );
755+
QVERIFY( std::abs( ( static_cast<float>( block->value( 2, 1 ) ) - values[5] ) / values[5] ) < epsilon );
756+
delete result;
757+
delete block;
758+
};
759+
760+
_chk( QStringLiteral( "\"landsat@1\" + \"landsat_4326@2\"" ), {264.0, 263.0, 264.0, 264.0, 266.0, 261.0}, false );
761+
_chk( QStringLiteral( "\"landsat@1\" + \"landsat_4326@2\"" ), {264.0, 263.0, 264.0, 264.0, 266.0, 261.0}, true );
762+
_chk( QStringLiteral( "\"landsat@1\"^2 + 3 + \"landsat_4326@2\"" ), {15767, 15766, 15519, 15767, 15769, 15516}, false );
763+
_chk( QStringLiteral( "\"landsat@1\"^2 + 3 + \"landsat_4326@2\"" ), {15767, 15766, 15519, 15767, 15769, 15516}, true );
764+
_chk( QStringLiteral( "0.5*((2*\"landsat@1\"+1)-sqrt((2*\"landsat@1\"+1)^2-8*(\"landsat@1\"-\"landsat_4326@2\")))" ), {-0.111504, -0.103543, -0.128448, -0.111504, -0.127425, -0.104374}, false );
765+
_chk( QStringLiteral( "0.5*((2*\"landsat@1\"+1)-sqrt((2*\"landsat@1\"+1)^2-8*(\"landsat@1\"-\"landsat_4326@2\")))" ), {-0.111504, -0.103543, -0.128448, -0.111504, -0.127425, -0.104374}, true );
766+
685767
}
686768

687769

‎tests/src/core/testqgsopenclutils.cpp

Lines changed: 0 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -29,9 +29,6 @@
2929
class TestQgsOpenClUtils: public QObject
3030
{
3131
Q_OBJECT
32-
public:
33-
34-
//void testRunMakeProgram();
3532

3633
private slots:
3734
void initTestCase();// will be called before the first testfunction is executed.

0 commit comments

Comments
 (0)