Do not print excessive decimals when identifying value on a Float32 raster

Was reported in https://lists.osgeo.org/pipermail/qgis-user/2016-April/036045.html

Author: rouault

src/core/qgis.h

@@ -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 )
 {

src/core/raster/qgsrasterblock.cpp

@@ -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 );

src/core/raster/qgsrasterblock.h

@@ -286,6 +286,12 @@ 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*/
+    static QString printValue( float value );
+
     /** \brief Convert data to different type.
      *  @param destDataType dest data type
      *  @return true on success */

src/gui/qgsmaptoolidentify.cpp

@@ -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 );
       }

src/providers/gdal/qgsgdalprovider.cpp

@@ -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;
   }

tests/src/app/testqgsmaptoolidentifyaction.cpp

@@ -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"

tests/testdata/raster/test.asc

@@ -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