[MESH] scalar color settings depending on classification (#36313)

* [MESH] [FEATURE] Sets meh color ramp classification from metadata read by MDAL driver.
Some mesh layer formats can provide values that can be compressed by categorizing values in consecutive intervals, each represent by an integer or byte. MDAL has the capabilities to recognize this dataset type and store the bounds of each class an the units in the metadata.
QGIS uses this metadata to setup adapted color ramp shader.

* [MDAL] update to pre-release 0.5.92

external/mdal/frmts/mdal_3di.cpp

@@ -194,10 +194,16 @@ std::set<std::string> MDAL::Driver3Di::ignoreNetCDFVariables()
   return ignore_variables;
 }
 
-void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *is_vector, bool *is_x )
+void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid,
+    std::string &variableName,
+    std::string &name,
+    bool *is_vector,
+    bool *isPolar,
+    bool *is_x )
 {
   *is_vector = false;
   *is_x = true;
+  *isPolar = false;
 
   std::string long_name = mNcFile->getAttrStr( "long_name", varid );
   if ( long_name.empty() )
@@ -209,6 +215,7 @@ void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid, const std::string
     }
     else
     {
+      variableName = standard_name;
       if ( MDAL::contains( standard_name, "_x_" ) )
       {
         *is_vector = true;
@@ -228,6 +235,7 @@ void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid, const std::string
   }
   else
   {
+    variableName = long_name;
     if ( MDAL::contains( long_name, " in x direction" ) )
     {
       *is_vector = true;
@@ -245,3 +253,9 @@ void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid, const std::string
     }
   }
 }
+
+std::vector<std::pair<double, double>> MDAL::Driver3Di::parseClassification( int varid ) const
+{
+  MDAL_UNUSED( varid );
+  return std::vector<std::pair<double, double>>();
+}

external/mdal/frmts/mdal_3di.hpp

@@ -50,8 +50,13 @@ namespace MDAL
       std::string getCoordinateSystemVariableName() override;
       std::string getTimeVariableName() const override;
       std::set<std::string> ignoreNetCDFVariables() override;
-      void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
-                                        std::string &name, bool *is_vector, bool *is_x ) override;
+      void parseNetCDFVariableMetadata( int varid,
+                                        std::string &variableName,
+                                        std::string &name,
+                                        bool *is_vector,
+                                        bool *isPolar,
+                                        bool *is_x ) override;
+      std::vector<std::pair<double, double>> parseClassification( int varid ) const override;
 
       //! Returns number of vertices
       size_t parse2DMesh();

external/mdal/frmts/mdal_cf.cpp

@@ -6,7 +6,7 @@
 #include <vector>
 #include <string>
 #include <netcdf.h>
-#include "math.h"
+#include <cmath>
 #include <stdlib.h>
 #include <assert.h>
 #include <cstring>
@@ -16,6 +16,28 @@
 #include "mdal_utils.hpp"
 #include "mdal_logger.hpp"
 
+static std::pair<std::string, std::string> metadataFromClassification( const MDAL::Classification &classes )
+{
+  std::pair<std::string, std::string> classificationMeta;
+  classificationMeta.first = "classification";
+  std::string classification;
+  for ( const auto boundValues : classes )
+  {
+    if ( boundValues.first != NC_FILL_DOUBLE )
+      classification.append( MDAL::doubleToString( boundValues.first ) );
+    if ( boundValues.second != NC_FILL_DOUBLE )
+    {
+      classification.append( "," );
+      classification.append( MDAL::doubleToString( boundValues.second ) );
+    }
+    if ( boundValues != classes.back() )
+      classification.append( ";;" );
+  }
+  classificationMeta.second = classification;
+
+  return classificationMeta;
+}
+
 MDAL::cfdataset_info_map MDAL::DriverCF::parseDatasetGroupInfo()
 {
   /*
@@ -101,63 +123,188 @@ MDAL::cfdataset_info_map MDAL::DriverCF::parseDatasetGroupInfo()
         continue;
 
       // Get name, if it is vector and if it is x or y
-      std::string name;
+      std::string vectorName;
       bool is_vector = true;
+      bool is_polar = false;
       bool is_x = false;
+      Classification classes = parseClassification( varid );
+      bool isClassified = !classes.empty();
+      parseNetCDFVariableMetadata( varid, variable_name, vectorName, &is_vector, &is_polar, &is_x );
 
-      parseNetCDFVariableMetadata( varid, variable_name, name, &is_vector, &is_x );
-
-      // Add it to the map
-      auto it = dsinfo_map.find( name );
-      if ( it != dsinfo_map.end() )
+      Metadata meta;
+      // check for units
+      std::string units;
+      try
       {
+        units = mNcFile->getAttrStr( "units", varid );
+        std::pair<std::string, std::string> unitMeta;
+        unitMeta.first = "units";
+        unitMeta.second = units;
+        meta.push_back( unitMeta );
+      }
+      catch ( MDAL::Error & )
+      {}
+
+      //construct classification metadata
+      if ( isClassified && !is_vector )
+      {
+
+        meta.push_back( metadataFromClassification( classes ) );
+      }
+
+      // Add dsinfo to the map
+      auto it = dsinfo_map.find( vectorName );
+      if ( it != dsinfo_map.end() && is_vector )
+      {
+        // this dataset is already existing and it is a vector dataset
+
         if ( is_x )
         {
           it->second.ncid_x = varid;
+          it->second.classification_x = classes;
         }
         else
         {
           it->second.ncid_y = varid;
+          it->second.classification_y = classes;
         }
+
+        // If it is classified, we want to keep each component as scalar
+        // So create two scalar dataset groups
+        if ( isClassified )
+        {
+          CFDatasetGroupInfo scalarDsInfoX;
+          scalarDsInfoX = it->second;
+          scalarDsInfoX.is_vector = false;
+          scalarDsInfoX.is_polar = false;
+          CFDatasetGroupInfo scalarDsInfoY;
+          scalarDsInfoY = it->second;
+          scalarDsInfoY.is_vector = false;
+          scalarDsInfoX.is_polar = false;
+          scalarDsInfoX.ncid_x = it->second.ncid_x;
+          scalarDsInfoY.ncid_x = it->second.ncid_y;
+
+          if ( is_x )
+          {
+            scalarDsInfoX.name = variable_name;
+            scalarDsInfoX.classification_x = classes;
+            scalarDsInfoX.metadata = meta;
+          }
+          else
+          {
+            scalarDsInfoY.name = variable_name;
+            scalarDsInfoY.classification_x = classes;
+            scalarDsInfoY.metadata = meta;
+          }
+
+          scalarDsInfoX.metadata.push_back( metadataFromClassification( scalarDsInfoX.classification_x ) );
+          scalarDsInfoY.metadata.push_back( metadataFromClassification( scalarDsInfoY.classification_y ) );
+
+          dsinfo_map[scalarDsInfoX.name] = scalarDsInfoX;
+          dsinfo_map[scalarDsInfoY.name] = scalarDsInfoY;
+        }
+
+        it->second.name = vectorName;
       }
-      else
+      else if ( it == dsinfo_map.end() || ( isClassified && is_vector ) )
       {
         CFDatasetGroupInfo dsInfo;
         dsInfo.nTimesteps = nTimesteps;
-        dsInfo.is_vector = is_vector;
+        dsInfo.ncid_x = -1;
+        dsInfo.ncid_y = -1;
         if ( is_x )
         {
           dsInfo.ncid_x = varid;
+          dsInfo.classification_x = classes;
         }
         else
         {
           dsInfo.ncid_y = varid;
+          dsInfo.classification_y = classes;
         }
+
         dsInfo.outputType = mDimensions.type( dimid );
-        dsInfo.name = name;
+        dsInfo.is_vector = is_vector;
+        dsInfo.is_polar = is_polar;
         dsInfo.nValues = mDimensions.size( mDimensions.type( dimid ) );
         dsInfo.timeLocation = timeLocation;
-        dsinfo_map[name] = dsInfo;
+        dsInfo.metadata = meta;
+        if ( is_vector && !isClassified )
+          dsInfo.name = vectorName;
+        else
+          dsInfo.name = variable_name;
+        dsinfo_map[vectorName] = dsInfo; //if is not vector, vectorName=variableName
       }
     }
   }
   while ( true );
 
+  // check the dsinfo if there dataset group defined as vector without valid ncid_y
+  // if ncid_y<0 set the datasetinfo to scalar
+  for ( auto &it : dsinfo_map )
+  {
+    if ( it.second.is_vector && it.second.ncid_y < 0 )
+      it.second.is_vector = false;
+  }
+
+
   return dsinfo_map;
 }
 
-static void populate_vals( bool is_vector, double *vals, size_t i,
-                           const std::vector<double> &vals_x, const std::vector<double> &vals_y,
-                           size_t idx, double fill_val_x, double fill_val_y )
+static void populate_vector_vals( double *vals, size_t i,
+                                  const std::vector<double> &vals_x, const std::vector<double> &vals_y,
+                                  size_t idx, double fill_val_x, double fill_val_y )
 {
-  if ( is_vector )
+  vals[2 * i] = MDAL::safeValue( vals_x[idx], fill_val_x );
+  vals[2 * i + 1] = MDAL::safeValue( vals_y[idx], fill_val_y );
+}
+
+static void populate_polar_vector_vals( double *vals, size_t i,
+                                        const std::vector<double> &vals_x, const std::vector<double> &vals_y,
+                                        size_t idx, double fill_val_x, double fill_val_y, std::pair<double, double> referenceAngles )
+{
+  double magnitude = MDAL::safeValue( vals_x[idx], fill_val_x );
+  double direction = MDAL::safeValue( vals_y[idx], fill_val_y );
+
+  direction = 2 * M_PI * ( ( direction - referenceAngles.second ) / referenceAngles.first );
+
+  vals[2 * i] = magnitude * cos( direction );
+  vals[2 * i + 1] = magnitude * sin( direction );
+}
+
+static void populate_scalar_vals( double *vals, size_t i,
+                                  const std::vector<double> &rawVals,
+                                  size_t idx,
+                                  double fill_val )
+{
+  vals[i] = MDAL::safeValue( rawVals[idx], fill_val );
+}
+
+static void fromClassificationToValue( const MDAL::Classification &classification, std::vector<double> &values, size_t classStartAt = 0 )
+{
+  for ( size_t i = 0; i < values.size(); ++i )
   {
-    vals[2 * i] = MDAL::safeValue( vals_x[idx], fill_val_x );
-    vals[2 * i + 1] = MDAL::safeValue( vals_y[idx], fill_val_y );
-  }
-  else
-  {
-    vals[i] = MDAL::safeValue( vals_x[idx], fill_val_x );
+    if ( std::isnan( values[i] ) )
+      continue;
+
+    size_t boundIndex = size_t( values[i] ) - classStartAt;
+    if ( boundIndex >= classification.size() )
+    {
+      values[i] = std::numeric_limits<double>::quiet_NaN();
+      continue;
+    }
+
+    std::pair<double, double> bounds = classification.at( boundIndex );
+    double bound1 = bounds.first;
+    double bound2 = bounds.second;
+    if ( bound1 == NC_FILL_DOUBLE )
+      bound1 = bound2;
+    if ( bound2 == NC_FILL_DOUBLE )
+      bound2 = bound1;
+    if ( bound1 == NC_FILL_DOUBLE || bound2 == NC_FILL_DOUBLE )
+      values[i] = std::numeric_limits<double>::quiet_NaN();
+    else
+      values[i] = ( bound1 + bound2 ) / 2;
   }
 }
 
@@ -175,6 +322,8 @@ void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<Relat
           dsi.name
         );
     group->setIsScalar( !dsi.is_vector );
+    group->setIsPolar( dsi.is_polar );
+    group->setMetadata( dsi.metadata );
 
     if ( dsi.outputType == CFDimensions::Vertex )
       group->setDataLocation( MDAL_DataLocation::DataOnVertices );
@@ -275,6 +424,8 @@ std::shared_ptr<MDAL::Dataset> MDAL::DriverCF::create2DDataset( std::shared_ptr<
         fill_val_y,
         dsi.ncid_x,
         dsi.ncid_y,
+        dsi.classification_x,
+        dsi.classification_y,
         dsi.timeLocation,
         dsi.nTimesteps,
         dsi.nValues,
@@ -491,14 +642,24 @@ bool MDAL::CFDimensions::isDatasetType( MDAL::CFDimensions::Type type ) const
 
 //////////////////////////////////////////////////////////////////////////////////////
 MDAL::CFDataset2D::CFDataset2D( MDAL::DatasetGroup *parent,
-                                double fill_val_x, double fill_val_y,
-                                int ncid_x, int ncid_y, CFDatasetGroupInfo::TimeLocation timeLocation,
-                                size_t timesteps, size_t values, size_t ts, std::shared_ptr<NetCDFFile> ncFile )
+                                double fill_val_x,
+                                double fill_val_y,
+                                int ncid_x,
+                                int ncid_y,
+                                Classification classification_x,
+                                Classification classification_y,
+                                CFDatasetGroupInfo::TimeLocation timeLocation,
+                                size_t timesteps,
+                                size_t values,
+                                size_t ts,
+                                std::shared_ptr<NetCDFFile> ncFile )
   : Dataset2D( parent )
   , mFillValX( fill_val_x )
   , mFillValY( fill_val_y )
   , mNcidX( ncid_x )
   , mNcidY( ncid_y )
+  , mClassificationX( classification_x )
+  , mClassificationY( classification_y )
   , mTimeLocation( timeLocation )
   , mTimesteps( timesteps )
   , mValues( values )
@@ -547,14 +708,11 @@ size_t MDAL::CFDataset2D::scalarData( size_t indexStart, size_t count, double *b
 
   for ( size_t i = 0; i < copyValues; ++i )
   {
-    populate_vals( false,
-                   buffer,
-                   i,
-                   values_x,
-                   std::vector<double>(),
-                   i,
-                   mFillValX,
-                   mFillValY );
+    populate_scalar_vals( buffer,
+                          i,
+                          values_x,
+                          i,
+                          mFillValX );
   }
   return copyValues;
 }
@@ -611,16 +769,36 @@ size_t MDAL::CFDataset2D::vectorData( size_t indexStart, size_t count, double *b
                );
   }
 
+  //if values component are classified convert from index to value
+  if ( !mClassificationX.empty() )
+  {
+    fromClassificationToValue( mClassificationX, values_x, 1 );
+  }
+
+  if ( !mClassificationY.empty() )
+  {
+    fromClassificationToValue( mClassificationY, values_y, 1 );
+  }
+
   for ( size_t i = 0; i < copyValues; ++i )
   {
-    populate_vals( true,
-                   buffer,
-                   i,
-                   values_x,
-                   values_y,
-                   i,
-                   mFillValX,
-                   mFillValY );
+    if ( group()->isPolar() )
+      populate_polar_vector_vals( buffer,
+                                  i,
+                                  values_x,
+                                  values_y,
+                                  i,
+                                  mFillValX,
+                                  mFillValY,
+                                  group()->referenceAngles() );
+    else
+      populate_vector_vals( buffer,
+                            i,
+                            values_x,
+                            values_y,
+                            i,
+                            mFillValX,
+                            mFillValY );
   }
 
   return copyValues;

external/mdal/frmts/mdal_cf.hpp

@@ -61,11 +61,15 @@ namespace MDAL
     std::string name; //!< Dataset group name
     CFDimensions::Type outputType;
     bool is_vector;
+    bool is_polar;
     TimeLocation timeLocation;
     size_t nTimesteps;
     size_t nValues;
     int ncid_x; //!< NetCDF variable id
     int ncid_y; //!< NetCDF variable id
+    Metadata metadata;
+    Classification classification_x;
+    Classification classification_y;
   };
   typedef std::map<std::string, CFDatasetGroupInfo> cfdataset_info_map; // name -> DatasetInfo
 
@@ -77,6 +81,8 @@ namespace MDAL
                    double fill_val_y,
                    int ncid_x,
                    int ncid_y,
+                   Classification classification_x,
+                   Classification classification_y,
                    CFDatasetGroupInfo::TimeLocation timeLocation,
                    size_t timesteps,
                    size_t values,
@@ -93,6 +99,8 @@ namespace MDAL
       double mFillValY;
       int mNcidX; //!< NetCDF variable id
       int mNcidY; //!< NetCDF variable id
+      Classification mClassificationX; //!< Classification, void if not classified
+      Classification mClassificationY; //!< Classification, void if not classified
       CFDatasetGroupInfo::TimeLocation mTimeLocation;
       size_t mTimesteps;
       size_t mValues;
@@ -120,8 +128,13 @@ namespace MDAL
       virtual void addBedElevation( MDAL::MemoryMesh *mesh ) = 0;
       virtual std::string getCoordinateSystemVariableName() = 0;
       virtual std::set<std::string> ignoreNetCDFVariables() = 0;
-      virtual void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
-          std::string &name, bool *is_vector, bool *is_x ) = 0;
+      virtual void parseNetCDFVariableMetadata( int varid,
+          std::string &variableName,
+          std::string &name,
+          bool *is_vector,
+          bool *isPolar,
+          bool *is_x ) = 0;
+      virtual std::vector<std::pair<double, double> > parseClassification( int varid ) const = 0;
       virtual std::string getTimeVariableName() const = 0;
       virtual std::shared_ptr<MDAL::Dataset> create2DDataset(
         std::shared_ptr<MDAL::DatasetGroup> group,

external/mdal/frmts/mdal_esri_tin.cpp

@@ -61,17 +61,17 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverEsriTin::load( const std::string &uri, c
     inMsx.seekg( -4, std::ios::end );
     int32_t mskBegin;
     if ( ! readValue( mskBegin, inMsx, true ) )
-      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to find the beginning of data in msk file" );
+      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to find the beginning of data in tmsx.adf file" );
 
     //read information in mskFile
     inMsk.seekg( -mskBegin * 2, std::ios::end );
     int32_t maskIntergerCount;
     if ( ! readValue( maskIntergerCount, inMsk, true ) )
-      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in msk file" );
+      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in tmsk.adf file" );
     inMsk.ignore( 4 ); //unused 4 bytes
     int32_t maskBitsCount;
     if ( ! readValue( maskBitsCount, inMsk, true ) )
-      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in msk file" );
+      throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in tmsk.adf file" );
 
     int c = 0;
     int32_t maskInt = 0;
@@ -80,7 +80,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverEsriTin::load( const std::string &uri, c
       //read mask file
       if ( c % 32 == 0 && c < maskBitsCount ) //first bit in the mask array have to be used-->read next maskInt
         if ( ! readValue( maskInt, inMsk, true ) )
-          throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in msk file" );
+          throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in tmsk.adf file" );
 
       Face f;
       for ( int i = 0; i < 3; ++i )
@@ -96,7 +96,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverEsriTin::load( const std::string &uri, c
         break;
 
       if ( f.size() < 3 ) //that's mean the face is not complete
-        throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in mask file, face is not complete" );
+        throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Unable to read information in tnod.adf file, face is not complete" );
 
       //exclude masked face
       if ( !( maskInt & 0x01 ) )

external/mdal/frmts/mdal_tuflowfv.cpp

@@ -54,6 +54,8 @@ MDAL::TuflowFVDataset2D::TuflowFVDataset2D(
   double fillValY,
   int ncidX,
   int ncidY,
+  Classification classificationX,
+  Classification classificationY,
   int ncidActive,
   CFDatasetGroupInfo::TimeLocation timeLocation,
   size_t timesteps,
@@ -67,6 +69,8 @@ MDAL::TuflowFVDataset2D::TuflowFVDataset2D(
       fillValY,
       ncidX,
       ncidY,
+      classificationX,
+      classificationY,
       timeLocation,
       timesteps,
       values,
@@ -456,10 +460,16 @@ std::set<std::string> MDAL::DriverTuflowFV::ignoreNetCDFVariables()
   return ignore_variables;
 }
 
-void MDAL::DriverTuflowFV::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *is_vector, bool *is_x )
+void MDAL::DriverTuflowFV::parseNetCDFVariableMetadata( int varid,
+    std::string &variableName,
+    std::string &name,
+    bool *is_vector,
+    bool *isPolar,
+    bool *is_x )
 {
   *is_vector = false;
   *is_x = true;
+  *isPolar = false;
 
   std::string long_name = mNcFile->getAttrStr( "long_name", varid );
   if ( long_name.empty() || ( long_name == "??????" ) )
@@ -480,6 +490,8 @@ void MDAL::DriverTuflowFV::parseNetCDFVariableMetadata( int varid, const std::st
     if ( MDAL::startsWith( long_name, "time at minimum value of " ) )
       long_name = MDAL::replace( long_name, "time at minimum value of ", "" ) + "/Time at Minimums";
 
+    variableName = long_name;
+
     if ( MDAL::startsWith( long_name, "x_" ) )
     {
       *is_vector = true;
@@ -515,6 +527,8 @@ std::shared_ptr<MDAL::Dataset> MDAL::DriverTuflowFV::create2DDataset(
         fill_val_y,
         dsi.ncid_x,
         dsi.ncid_y,
+        dsi.classification_x,
+        dsi.classification_y,
         mNcFile->arrId( "stat" ),
         dsi.timeLocation,
         dsi.nTimesteps,
@@ -551,3 +565,9 @@ std::shared_ptr<MDAL::Dataset> MDAL::DriverTuflowFV::create3DDataset( std::share
   dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
   return std::move( dataset );
 }
+
+std::vector<std::pair<double, double>> MDAL::DriverTuflowFV::parseClassification( int varid ) const
+{
+  MDAL_UNUSED( varid );
+  return std::vector<std::pair<double, double>>();
+}

external/mdal/frmts/mdal_tuflowfv.hpp

@@ -40,6 +40,8 @@ namespace MDAL
                          double fillValY,
                          int ncidX,
                          int ncidY,
+                         Classification classificationX,
+                         Classification classificationY,
                          int ncidActive,
                          CFDatasetGroupInfo::TimeLocation timeLocation,
                          size_t timesteps,
@@ -122,8 +124,13 @@ namespace MDAL
       void addBedElevation( MemoryMesh *mesh ) override;
       std::string getCoordinateSystemVariableName() override;
       std::set<std::string> ignoreNetCDFVariables() override;
-      void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
-                                        std::string &name, bool *is_vector, bool *is_x ) override;
+      void parseNetCDFVariableMetadata( int varid,
+                                        std::string &variableName,
+                                        std::string &name,
+                                        bool *is_vector,
+                                        bool *isPolar,
+                                        bool *is_x ) override;
+      std::vector<std::pair<double, double>> parseClassification( int varid ) const override;
       std::string getTimeVariableName() const override;
       std::shared_ptr<MDAL::Dataset> create2DDataset(
         std::shared_ptr<MDAL::DatasetGroup> group,

external/mdal/frmts/mdal_ugrid.cpp

@@ -491,10 +491,17 @@ void MDAL::DriverUgrid::ignore2DMeshVariables( const std::string &mesh, std::set
   ignoreVariables.insert( mNcFile->getAttrStr( mesh, "edge_face_connectivity" ) );
 }
 
-void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *isVector, bool *isX )
+void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid,
+    std::string &variableName,
+    std::string &name,
+    bool *isVector,
+    bool *isPolar,
+    bool *isX )
 {
+
   *isVector = false;
   *isX = true;
+  *isPolar = false;
 
   std::string longName = mNcFile->getAttrStr( "long_name", varid );
   if ( longName.empty() )
@@ -506,6 +513,7 @@ void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::strin
     }
     else
     {
+      variableName = standardName;
       if ( MDAL::contains( standardName, "_x_" ) )
       {
         *isVector = true;
@@ -525,6 +533,7 @@ void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::strin
   }
   else
   {
+    variableName = longName;
     if ( MDAL::contains( longName, ", x-component" ) || MDAL::contains( longName, "u component of " ) )
     {
       *isVector = true;
@@ -538,6 +547,20 @@ void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::strin
       name = MDAL::replace( longName, ", y-component", "" );
       name = MDAL::replace( name, "v component of ", "" );
     }
+    else if ( MDAL::contains( longName, "velocity magnitude" ) )
+    {
+      *isVector = true;
+      *isPolar = true;
+      *isX = true;
+      name = MDAL::replace( longName, " magnitude", "" );
+    }
+    else if ( MDAL::contains( longName, "velocity direction" ) )
+    {
+      *isVector = true;
+      *isPolar = true;
+      *isX = false;
+      name = MDAL::replace( longName, " direction", "" );
+    }
     else
     {
       name = longName;
@@ -807,3 +830,40 @@ void MDAL::DriverUgrid::writeGlobals()
   mNcFile->putAttrStr( NC_GLOBAL, "date_created", MDAL::getCurrentTimeStamp() );
   mNcFile->putAttrStr( NC_GLOBAL, "Conventions", "CF-1.6 UGRID-1.0" );
 }
+
+std::vector<std::pair<double, double>> MDAL::DriverUgrid::parseClassification( int varid ) const
+{
+  std::vector<std::pair<double, double>> classes;
+  std::string flagBoundVarName = mNcFile->getAttrStr( "flag_bounds", varid );
+  if ( !flagBoundVarName.empty() )
+  {
+    try
+    {
+      int boundsVarId = mNcFile->getVarId( flagBoundVarName );
+      std::vector<size_t> classDims;
+      std::vector<int> classDimIds;
+      mNcFile->getDimensions( flagBoundVarName, classDims, classDimIds );
+      std::vector<double> boundValues = mNcFile->readDoubleArr( boundsVarId, 0, 0, classDims[0], classDims[1] );
+
+      if ( classDims[1] != 2 || classDims[0] <= 0 )
+        throw MDAL::Error( MDAL_Status::Err_UnknownFormat, "Invalid classification dimension" );
+
+      std::pair<std::string, std::string> classificationMeta;
+      classificationMeta.first = "classification";
+      std::string classification;
+      for ( size_t i = 0; i < classDims[0]; ++i )
+      {
+        std::pair<double, double> classBound;
+        classBound.first = boundValues[i * 2];
+        classBound.second = boundValues[i * 2 + 1];
+        classes.push_back( classBound );
+      }
+    }
+    catch ( MDAL::Error &err )
+    {
+      MDAL::Log::warning( err.status, err.driver, "Error when parsing class bounds: " + err.mssg + ", classification ignored" );
+    }
+  }
+
+  return classes;
+}

external/mdal/frmts/mdal_ugrid.hpp

@@ -38,8 +38,12 @@ namespace MDAL
       void addBedElevation( MemoryMesh *mesh ) override;
       std::string getCoordinateSystemVariableName() override;
       std::set<std::string> ignoreNetCDFVariables() override;
-      void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
-                                        std::string &name, bool *is_vector, bool *is_x ) override;
+      void parseNetCDFVariableMetadata( int varid,
+                                        std::string &variableName,
+                                        std::string &name,
+                                        bool *is_vector, bool *isPolar,
+                                        bool *is_x ) override;
+      std::vector<std::pair<double, double>> parseClassification( int varid ) const override;
       std::string getTimeVariableName() const override;
 
       void parse2VariablesFromAttribute( const std::string &name, const std::string &attr_name,

external/mdal/frmts/mdal_xms_tin.cpp

@@ -105,7 +105,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverXmsTin::load( const std::string &meshFil
   // Read triangles
   if ( !std::getline( in, line ) )
   {
-    MDAL::Log::error( MDAL_Status::Err_IncompatibleMesh, name(), meshFile + " does not contain valid triangle definitions" );
+    MDAL::Log::error( MDAL_Status::Err_IncompatibleMesh, name(), meshFile + " does not contain valid triangle definition" );
     return nullptr;
   }
   chunks = split( line,  ' ' );
@@ -127,7 +127,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverXmsTin::load( const std::string &meshFil
     if ( chunks.size() != 3 )
     {
       // should have 3 indexes
-      MDAL::Log::error( MDAL_Status::Err_IncompatibleMesh, name(), meshFile + " does not contain valid triangle definitions" );
+      MDAL::Log::error( MDAL_Status::Err_IncompatibleMesh, name(), meshFile + " does not contain valid triangle defintion" );
       return nullptr;
     }
 

external/mdal/mdal.cpp

@@ -21,7 +21,7 @@ static const char *EMPTY_STR = "";
 
 const char *MDAL_Version()
 {
-  return "0.5.91";
+  return "0.5.92";
 }
 
 MDAL_Status MDAL_LastStatus()

external/mdal/mdal_data_model.cpp

@@ -187,6 +187,12 @@ void MDAL::DatasetGroup::setMetadata( const std::string &key, const std::string
     metadata.push_back( std::make_pair( key, val ) );
 }
 
+void MDAL::DatasetGroup::setMetadata( const MDAL::Metadata &metadata )
+{
+  for ( const auto &meta : metadata )
+    setMetadata( meta.first, meta.second );
+}
+
 std::string MDAL::DatasetGroup::name()
 {
   return getMetadata( "name" );
@@ -254,6 +260,26 @@ void MDAL::DatasetGroup::stopEditing()
   mInEditMode = false;
 }
 
+void MDAL::DatasetGroup::setReferenceAngles( const std::pair<double, double> &referenceAngle )
+{
+  mReferenceAngles = referenceAngle;
+}
+
+bool MDAL::DatasetGroup::isPolar() const
+{
+  return mIsPolar;
+}
+
+void MDAL::DatasetGroup::setIsPolar( bool isPolar )
+{
+  mIsPolar = isPolar;
+}
+
+std::pair<double, double> MDAL::DatasetGroup::referenceAngles() const
+{
+  return mReferenceAngles;
+}
+
 MDAL_DataLocation MDAL::DatasetGroup::dataLocation() const
 {
   return mDataLocation;

external/mdal/mdal_data_model.hpp

@@ -38,6 +38,7 @@ namespace MDAL
   } Statistics;
 
   typedef std::vector< std::pair< std::string, std::string > > Metadata;
+  typedef std::vector<std::pair<double, double>> Classification;
 
   class Dataset
   {
@@ -150,6 +151,7 @@ namespace MDAL
 
       std::string getMetadata( const std::string &key );
       void setMetadata( const std::string &key, const std::string &val );
+      void setMetadata( const Metadata &metadata );
 
       std::string name();
       void setName( const std::string &name );
@@ -179,12 +181,20 @@ namespace MDAL
       void startEditing();
       void stopEditing();
 
+      //! First value is the angle for full rotation and second value is the start angle
+      void setReferenceAngles( const std::pair<double, double> &referenceAngle );
+      std::pair<double, double> referenceAngles() const;
+
+      bool isPolar() const;
+      void setIsPolar( bool isPolar );
     private:
       bool mInEditMode = false;
 
       const std::string mDriverName;
       Mesh *mParent = nullptr;
       bool mIsScalar = true;
+      bool mIsPolar = true;
+      std::pair<double, double> mReferenceAngles = {360, 0};
       MDAL_DataLocation mDataLocation = MDAL_DataLocation::DataOnVertices;
       std::string mUri; // file/uri from where it came
       Statistics mStatistics;

external/mdal/mdal_utils.hpp

@@ -23,6 +23,10 @@
 #define MDAL_UNUSED(x) (void)x;
 #define MDAL_NAN std::numeric_limits<double>::quiet_NaN()
 
+#ifndef M_PI
+#define M_PI 3.14159265358979323846264338327
+#endif
+
 namespace MDAL
 {
   // endianness

src/app/mesh/qgsmeshrendereractivedatasetwidget.cpp

@@ -185,6 +185,11 @@ QString QgsMeshRendererActiveDatasetWidget::metadata( QgsMeshDatasetIndex datase
   const auto options = gmeta.extraOptions();
   for ( auto it = options.constBegin(); it != options.constEnd(); ++it )
   {
+    if ( it.key() == QStringLiteral( "classification" ) )
+    {
+      msg += QStringLiteral( "<tr><td>%1</td></tr>" ).arg( tr( "Classified values" ) );
+      continue;
+    }
     msg += QStringLiteral( "<tr><td>%1</td><td>%2</td></tr>" ).arg( it.key() ).arg( it.value() );
   }
 

src/app/mesh/qgsmeshrendererscalarsettingswidget.cpp

@@ -112,7 +112,6 @@ void QgsMeshRendererScalarSettingsWidget::syncToLayer( )
   const double min = settings.classificationMinimum();
   const double max = settings.classificationMaximum();
 
-
   whileBlocking( mScalarMinLineEdit )->setText( QString::number( min ) );
   whileBlocking( mScalarMaxLineEdit )->setText( QString::number( max ) );
   whileBlocking( mScalarColorRampShaderWidget )->setFromShader( shader );

src/app/mesh/qgsrenderermeshpropertieswidget.cpp

@@ -116,6 +116,11 @@ void QgsRendererMeshPropertiesWidget::apply()
   if ( activeVectorDatasetGroupIndex > -1 )
     settings.setVectorSettings( activeVectorDatasetGroupIndex, mMeshRendererVectorSettingsWidget->settings() );
 
+  QgsMeshDatasetIndex staticScalarDatasetIndex( activeScalarDatasetGroupIndex, mMeshLayer->staticScalarDatasetIndex().dataset() );
+  QgsMeshDatasetIndex staticVectorDatasetIndex( activeVectorDatasetGroupIndex, mMeshLayer->staticVectorDatasetIndex().dataset() );
+  mMeshLayer->setStaticScalarDatasetIndex( staticScalarDatasetIndex );
+  mMeshLayer->setStaticVectorDatasetIndex( staticVectorDatasetIndex );
+
   //set the blend mode for the layer
   mMeshLayer->setBlendMode( mBlendModeComboBox->blendMode() );
   //set the averaging method for the layer

src/core/mesh/qgsmeshlayer.cpp

@@ -108,6 +108,10 @@ void QgsMeshLayer::setDefaultRendererSettings()
       case QgsMeshDatasetGroupMetadata::DataOnEdges:
         break;
     }
+
+    //override color ramp if the values in the dataset group are classified
+    applyClassificationOnScalarSettings( meta, scalarSettings );
+
     mRendererSettings.setScalarSettings( i, scalarSettings );
   }
 
@@ -407,6 +411,94 @@ QgsMeshDatasetIndex QgsMeshLayer::datasetIndexAtTime( const QgsDateTimeRange &ti
     return QgsMeshDatasetIndex();
 }
 
+void QgsMeshLayer::applyClassificationOnScalarSettings( const QgsMeshDatasetGroupMetadata &meta, QgsMeshRendererScalarSettings &scalarSettings ) const
+{
+  if ( meta.extraOptions().contains( QStringLiteral( "classification" ) ) )
+  {
+    QgsColorRampShader colorRampShader = scalarSettings.colorRampShader();
+    QgsColorRamp *colorRamp = colorRampShader.sourceColorRamp();
+    QStringList classes = meta.extraOptions()[QStringLiteral( "classification" )].split( QStringLiteral( ";;" ) );
+
+    QString units;
+    if ( meta.extraOptions().contains( QStringLiteral( "units" ) ) )
+      units = meta.extraOptions()[ QStringLiteral( "units" )];
+
+    QVector<QVector<double>> bounds;
+    for ( const QString classe : classes )
+    {
+      QStringList boundsStr = classe.split( ',' );
+      QVector<double> bound;
+      for ( const QString boundStr : boundsStr )
+        bound.append( boundStr.toDouble() );
+      bounds.append( bound );
+    }
+
+    if ( ( bounds.count() == 1  && bounds.first().count() > 2 ) || // at least a class with two value
+         ( bounds.count() > 1 ) ) // or at least two classes
+    {
+      const QVector<double> firstClass = bounds.first();
+      const QVector<double> lastClass = bounds.last();
+      double minValue = firstClass.count() > 1 ? ( firstClass.first() + firstClass.last() ) / 2 : firstClass.first();
+      double maxValue = lastClass.count() > 1 ? ( lastClass.first() + lastClass.last() ) / 2 : lastClass.first();
+      double diff = maxValue - minValue;
+      QList<QgsColorRampShader::ColorRampItem> colorRampItemlist;
+      for ( int i = 0; i < bounds.count(); ++i )
+      {
+        const QVector<double> &boundClass = bounds.at( i );
+        QgsColorRampShader::ColorRampItem item;
+        item.value = i + 1;
+        if ( !boundClass.isEmpty() )
+        {
+          double scalarValue = ( boundClass.first() + boundClass.last() ) / 2;
+          item.color = colorRamp->color( ( scalarValue - minValue ) / diff );
+          if ( i != 0 && i < bounds.count() - 1 ) //The first and last labels are treated after
+          {
+            item.label = QString( ( "%1 - %2 %3" ) ).
+                         arg( QString::number( boundClass.first() ) ).
+                         arg( QString::number( boundClass.last() ) ).
+                         arg( units );
+          }
+        }
+        colorRampItemlist.append( item );
+      }
+      //treat first and last labels
+      if ( firstClass.count() == 1 )
+        colorRampItemlist.first().label = QObject::tr( "below %1 %2" ).
+                                          arg( QString::number( firstClass.first() ) ).
+                                          arg( units );
+      else
+      {
+        colorRampItemlist.first().label = QString( ( "%1 - %2 %3" ) ).
+                                          arg( QString::number( firstClass.first() ) ).
+                                          arg( QString::number( firstClass.last() ) ).
+                                          arg( units );
+      }
+
+      if ( lastClass.count() == 1 )
+        colorRampItemlist.last().label = QObject::tr( "above %1 %2" ).
+                                         arg( QString::number( lastClass.first() ) ).
+                                         arg( units );
+      else
+      {
+        colorRampItemlist.last().label = QString( ( "%1 - %2 %3" ) ).
+                                         arg( QString::number( lastClass.first() ) ).
+                                         arg( QString::number( lastClass.last() ) ).
+                                         arg( units );
+      }
+
+      colorRampShader.setMinimumValue( 0 );
+      colorRampShader.setMaximumValue( colorRampItemlist.count() - 1 );
+      scalarSettings.setClassificationMinimumMaximum( 0, colorRampItemlist.count() - 1 );
+      colorRampShader.setColorRampItemList( colorRampItemlist );
+      colorRampShader.setColorRampType( QgsColorRampShader::Exact );
+      colorRampShader.setClassificationMode( QgsColorRampShader::EqualInterval );
+    }
+
+    scalarSettings.setColorRampShader( colorRampShader );
+    scalarSettings.setDataResamplingMethod( QgsMeshRendererScalarSettings::None );
+  }
+}
+
 QgsMeshDatasetIndex QgsMeshLayer::activeScalarDatasetAtTime( const QgsDateTimeRange &timeRange ) const
 {
   if ( mTemporalProperties->isActive() )

src/core/mesh/qgsmeshlayer.h

@@ -461,6 +461,9 @@ class CORE_EXPORT QgsMeshLayer : public QgsMapLayer
 
     QgsMeshDatasetIndex datasetIndexAtTime( const QgsDateTimeRange &timeRange, int datasetGroupIndex ) const;
 
+    //! Changes scalar settings for classified scalar value (information about is in the metadata
+    void applyClassificationOnScalarSettings( const QgsMeshDatasetGroupMetadata &meta, QgsMeshRendererScalarSettings &scalarSettings ) const;
+
   private slots:
     void onDatasetGroupsAdded( int count );
 

src/core/mesh/qgsmeshlayerutils.cpp

@@ -587,5 +587,4 @@ QVector<QVector3D> QgsMeshLayerUtils::calculateNormals( const QgsTriangularMesh
   return normals;
 }
 
-
 ///@endcond

tests/src/core/testqgsmeshlayerrenderer.cpp

@@ -94,6 +94,7 @@ class TestQgsMeshRenderer : public QObject
     void test_vertex_vector_traces_colorRamp();
     void test_stacked_3d_mesh_single_level_averaging();
     void test_simplified_triangular_mesh_rendering();
+    void test_classified_values();
 
     void test_signals();
 };
@@ -681,7 +682,19 @@ void TestQgsMeshRenderer::test_simplified_triangular_mesh_rendering()
   QVERIFY( imageCheck( "simplified_triangular_mesh", mMdal3DLayer ) );
 }
 
-// TODO test edge mesh rendering!
+void TestQgsMeshRenderer::test_classified_values()
+{
+  QgsMeshLayer classifiedMesh( mDataDir + "/simplebox_clm.nc", "Mesh with classified values", "mdal" );
+  QVERIFY( classifiedMesh.isValid() );
+
+  QgsProject::instance()->addMapLayer( &classifiedMesh );
+  mMapSettings->setLayers( QList<QgsMapLayer *>() << &classifiedMesh );
+
+  classifiedMesh.temporalProperties()->setIsActive( false );
+  classifiedMesh.setStaticScalarDatasetIndex( QgsMeshDatasetIndex( 3, 4 ) );
+
+  QVERIFY( imageCheck( "classified_values", &classifiedMesh ) );
+}
 
 QGSTEST_MAIN( TestQgsMeshRenderer )
 #include "testqgsmeshlayerrenderer.moc"