Skip to content
Permalink
Browse files

Merge pull request #3142 from rouault/identify_result_float_precision

Do not print excessive decimals when identifying value on a Float32 raster
  • Loading branch information
rouault committed May 30, 2016
2 parents 4b130cb + 9fe2102 commit 86f6e7eb0da9b5c1b792910660b9829d01ddc542
@@ -355,6 +355,16 @@ inline bool qgsDoubleNear( double a, double b, double epsilon = 4 * DBL_EPSILON
return diff > -epsilon && diff <= epsilon;
}

//! Compare two floats (but allow some difference)
//! @param a first float
//! @param b second float
//! @param epsilon maximum difference allowable between floats
inline bool qgsFloatNear( float a, float b, float epsilon = 4 * FLT_EPSILON )
{
const float diff = a - b;
return diff > -epsilon && diff <= epsilon;
}

//! Compare two doubles using specified number of significant digits
inline bool qgsDoubleNearSig( double a, double b, int significantDigits = 10 )
{
@@ -831,6 +831,27 @@ QString QgsRasterBlock::printValue( double value )
return s;
}

QString QgsRasterBlock::printValue( float value )
{
/*
* IEEE 754 double has 6-9 significant digits. See printValue(double)
*/

QString s;

for ( int i = 6; i <= 9; i++ )
{
s.setNum( value, 'g', i );
if ( qgsFloatNear( s.toFloat(), value ) )
{
return s;
}
}
// Should not happen
QgsDebugMsg( "Cannot correctly parse printed value" );
return s;
}

void * QgsRasterBlock::convert( void *srcData, QGis::DataType srcDataType, QGis::DataType destDataType, qgssize size )
{
int destDataTypeSize = typeSize( destDataType );
@@ -286,6 +286,15 @@ class CORE_EXPORT QgsRasterBlock
* @return string representing the value*/
static QString printValue( double value );

/** \brief Print float value with all necessary significant digits.
* It is ensured that conversion back to float gives the same number.
* @param value the value to be printed
* @return string representing the value
* @note added in QGIS 2.16
* @note not available in python bindings
*/
static QString printValue( float value );

/** \brief Convert data to different type.
* @param destDataType dest data type
* @return true on success */
@@ -554,8 +554,18 @@ bool QgsMapToolIdentify::identifyRasterLayer( QList<IdentifyResult> *results, Qg
}
else
{
double value = values.value( bandNo ).toDouble();
valueString = QgsRasterBlock::printValue( value );
QVariant value( values.value( bandNo ) );
// The cast is legit. Quoting QT doc :
// "Although this function is declared as returning QVariant::Type,
// the return value should be interpreted as QMetaType::Type"
if ( static_cast<QMetaType::Type>( value.type() ) == QMetaType::Float )
{
valueString = QgsRasterBlock::printValue( value.toFloat() );
}
else
{
valueString = QgsRasterBlock::printValue( value.toDouble() );
}
}
attributes.insert( dprovider->generateBandName( bandNo ), valueString );
}
@@ -139,9 +139,12 @@ QgsGdalProvider::QgsGdalProvider( const QString &uri, bool update )

QgsGdalProviderBase::registerGdalDrivers();

// GDAL tends to open AAIGrid as Float32 which results in lost precision
// and confusing values shown to users, force Float64
CPLSetConfigOption( "AAIGRID_DATATYPE", "Float64" );
if ( !CPLGetConfigOption( "AAIGRID_DATATYPE", nullptr ) )
{
// GDAL tends to open AAIGrid as Float32 which results in lost precision
// and confusing values shown to users, force Float64
CPLSetConfigOption( "AAIGRID_DATATYPE", "Float64" );
}

// To get buildSupportedRasterFileFilter the provider is called with empty uri
if ( uri.isEmpty() )
@@ -1068,7 +1071,14 @@ QgsRasterIdentifyResult QgsGdalProvider::identify( const QgsPoint & thePoint, Qg
}
else
{
results.insert( i, value );
if ( srcDataType( i ) == QGis::Float32 )
{
// Insert a float QVariant so that QgsMapToolIdentify::identifyRasterLayer()
// can print a string without an excessive precision
results.insert( i, static_cast<float>( value ) );
}
else
results.insert( i, value );
}
delete myBlock;
}
@@ -16,6 +16,7 @@
#include <QtTest/QtTest>
#include "qgsapplication.h"
#include "qgsvectorlayer.h"
#include "qgsrasterlayer.h"
#include "qgsfeature.h"
#include "qgsgeometry.h"
#include "qgsvectordataprovider.h"
@@ -24,6 +25,8 @@
#include "qgsunittypes.h"
#include "qgsmaptoolidentifyaction.h"

#include "cpl_conv.h"

class TestQgsMapToolIdentifyAction : public QObject
{
Q_OBJECT
@@ -40,9 +43,13 @@ class TestQgsMapToolIdentifyAction : public QObject
void lengthCalculation(); //test calculation of derived length attributes
void perimeterCalculation(); //test calculation of derived perimeter attribute
void areaCalculation(); //test calculation of derived area attribute
void identifyRasterFloat32(); // test pixel identification and decimal precision
void identifyRasterFloat64(); // test pixel identification and decimal precision

private:
QgsMapCanvas* canvas;

QString testIdentifyRaster( QgsRasterLayer* layer, double xGeoref, double yGeoref );
};

void TestQgsMapToolIdentifyAction::initTestCase()
@@ -241,6 +248,73 @@ void TestQgsMapToolIdentifyAction::areaCalculation()
QVERIFY( qgsDoubleNear( area, 389.6117, 0.001 ) );
}

QString TestQgsMapToolIdentifyAction::testIdentifyRaster( QgsRasterLayer* layer, double xGeoref, double yGeoref )
{
QScopedPointer< QgsMapToolIdentifyAction > action( new QgsMapToolIdentifyAction( canvas ) );
QgsPoint mapPoint = canvas->getCoordinateTransform()->transform( xGeoref, yGeoref );
QList<QgsMapToolIdentify::IdentifyResult> result = action->identify( mapPoint.x(), mapPoint.y(), QList<QgsMapLayer*>() << layer );
if ( result.length() != 1 )
return "";
return result[0].mAttributes["Band 1"];
}

void TestQgsMapToolIdentifyAction::identifyRasterFloat32()
{
//create a temporary layer
QString raster = QString( TEST_DATA_DIR ) + "/raster/test.asc";

// By default the QgsRasterLayer forces AAIGRID_DATATYPE=Float64
CPLSetConfigOption( "AAIGRID_DATATYPE", "Float32" );
QScopedPointer< QgsRasterLayer> tempLayer( new QgsRasterLayer( raster ) );
CPLSetConfigOption( "AAIGRID_DATATYPE", nullptr );

QVERIFY( tempLayer->isValid() );

canvas->setExtent( QgsRectangle( 0, 0, 7, 1 ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 0.5, 0.5 ), QString( "-999.9" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 1.5, 0.5 ), QString( "-999.987" ) );

// More than 6 significant digits : precision loss in Float32
QCOMPARE( testIdentifyRaster( tempLayer.data(), 2.5, 0.5 ), QString( "1.234568" ) ); // in .asc file : 1.2345678

QCOMPARE( testIdentifyRaster( tempLayer.data(), 3.5, 0.5 ), QString( "123456" ) );

// More than 6 significant digits: no precision loss here for that particular value
QCOMPARE( testIdentifyRaster( tempLayer.data(), 4.5, 0.5 ), QString( "1234567" ) );

// More than 6 significant digits: no precision loss here for that particular value
QCOMPARE( testIdentifyRaster( tempLayer.data(), 5.5, 0.5 ), QString( "-999.9876" ) );

// More than 6 significant digits: no precision loss here
QCOMPARE( testIdentifyRaster( tempLayer.data(), 6.5, 0.5 ), QString( "1.234568" ) ); // in .asc file : 1.2345678901234
}

void TestQgsMapToolIdentifyAction::identifyRasterFloat64()
{
//create a temporary layer
QString raster = QString( TEST_DATA_DIR ) + "/raster/test.asc";
QScopedPointer< QgsRasterLayer> tempLayer( new QgsRasterLayer( raster ) );
QVERIFY( tempLayer->isValid() );

canvas->setExtent( QgsRectangle( 0, 0, 7, 1 ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 0.5, 0.5 ), QString( "-999.9" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 1.5, 0.5 ), QString( "-999.987" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 2.5, 0.5 ), QString( "1.2345678" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 3.5, 0.5 ), QString( "123456" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 4.5, 0.5 ), QString( "1234567" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 5.5, 0.5 ), QString( "-999.9876" ) );

QCOMPARE( testIdentifyRaster( tempLayer.data(), 6.5, 0.5 ), QString( "1.2345678901234" ) );
}


QTEST_MAIN( TestQgsMapToolIdentifyAction )
#include "testqgsmaptoolidentifyaction.moc"
@@ -0,0 +1,6 @@
ncols 7
nrows 1
xllcorner 0
yllcorner 0
cellsize 1
-999.9 -999.987 1.2345678 123456 1234567 -999.9876 1.2345678901234

0 comments on commit 86f6e7e

Please sign in to comment.
You can’t perform that action at this time.