diff --git a/external/mdal/api/mdal.h b/external/mdal/api/mdal.h
index e30d5b2201d..865a764fffe 100644
--- a/external/mdal/api/mdal.h
+++ b/external/mdal/api/mdal.h
@@ -56,6 +56,7 @@ enum MDAL_Status
   Err_IncompatibleDatasetGroup,
   Err_MissingDriver,
   Err_MissingDriverCapability,
+  Err_FailToWriteToDisk,
   // Warnings
   Warn_UnsupportedElement,
   Warn_InvalidElements,
@@ -105,6 +106,9 @@ MDAL_EXPORT bool MDAL_DR_meshLoadCapability( DriverH driver );
 //! Returns whether driver has capability to write/edit dataset (groups)
 MDAL_EXPORT bool MDAL_DR_writeDatasetsCapability( DriverH driver );
 
+//! Returns whether driver has capability to save mesh
+MDAL_EXPORT bool MDAL_DR_SaveMeshCapability( DriverH driver );
+
 /**
  * Returns name of MDAL driver
  * not thread-safe and valid only till next call
@@ -134,9 +138,13 @@ MDAL_EXPORT const char *MDAL_DR_filters( DriverH driver );
  * Caller must free memory with MDAL_CloseMesh() afterwards
  */
 MDAL_EXPORT MeshH MDAL_LoadMesh( const char *meshFile );
+
 //! Closes mesh, frees the memory
 MDAL_EXPORT void MDAL_CloseMesh( MeshH mesh );
 
+//! Saves mesh (only mesh structure) on a file with the specified driver. On error see MDAL_LastStatus for error type.
+MDAL_EXPORT void MDAL_SaveMesh( MeshH mesh, const char *meshFile, const char *driver );
+
 /**
  * Returns mesh projection
  * not thread-safe and valid only till next call
@@ -348,6 +356,11 @@ MDAL_EXPORT bool MDAL_G_isInEditMode( DatasetGroupH group );
  */
 MDAL_EXPORT void MDAL_G_closeEditMode( DatasetGroupH group );
 
+/**
+ * Returns reference time for dataset group
+ * If returned value begins with word JULIAN, following number represents date in Julian format
+ */
+MDAL_EXPORT const char *MDAL_G_referenceTime( DatasetGroupH group );
 
 ///////////////////////////////////////////////////////////////////////////////////////
 /// DATASETS
diff --git a/external/mdal/frmts/mdal_2dm.cpp b/external/mdal/frmts/mdal_2dm.cpp
index 5f1402d45f5..fa57626f0d1 100644
--- a/external/mdal/frmts/mdal_2dm.cpp
+++ b/external/mdal/frmts/mdal_2dm.cpp
@@ -82,7 +82,7 @@ MDAL::Driver2dm::Driver2dm():
   Driver( DRIVER_NAME,
           "2DM Mesh File",
           "*.2dm",
-          Capability::ReadMesh
+          Capability::ReadMesh | Capability::SaveMesh
         )
 {
 }
@@ -271,7 +271,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::Driver2dm::load( const std::string &meshFile,
     new Mesh2dm(
       vertices.size(),
       faces.size(),
-      4, //maximum quads
+      MAX_VERTICES_PER_FACE_2DM,
       computeExtent( vertices ),
       mMeshFile,
       vertexIDtoIndex
@@ -286,3 +286,66 @@ std::unique_ptr<MDAL::Mesh> MDAL::Driver2dm::load( const std::string &meshFile,
 
   return std::unique_ptr<Mesh>( mesh.release() );
 }
+
+void MDAL::Driver2dm::save( const std::string &uri, MDAL::Mesh *mesh, MDAL_Status *status )
+{
+  if ( status ) *status = MDAL_Status::None;
+
+  std::ofstream file( uri, std::ofstream::out );
+
+  if ( !file.is_open() )
+  {
+    if ( status ) *status = MDAL_Status::Err_FailToWriteToDisk;
+  }
+
+  std::string line = "MESH2D";
+  file << line << std::endl;
+
+  //write vertices
+  std::unique_ptr<MDAL::MeshVertexIterator> vertexIterator = mesh->readVertices();
+  double vertex[3];
+  for ( size_t i = 0; i < mesh->verticesCount(); ++i )
+  {
+    vertexIterator->next( 1, vertex );
+    line = "ND ";
+    line.append( std::to_string( i + 1 ) );
+    for ( size_t j = 0; j < 2; ++j )
+    {
+      line.append( " " );
+      line.append( MDAL::coordinateToString( vertex[j] ) );
+    }
+    line.append( " " );
+    line.append( MDAL::doubleToString( vertex[2] ) );
+
+    file << line << std::endl;
+  }
+
+  //write faces
+  std::unique_ptr<MDAL::MeshFaceIterator> faceIterator = mesh->readFaces();
+  for ( size_t i = 0; i < mesh->facesCount(); ++i )
+  {
+    int faceOffsets[1];
+    int vertexIndices[4]; //max 4 vertices for a face
+    faceIterator->next( 1, faceOffsets, 4, vertexIndices );
+
+    if ( faceOffsets[0] > 2 && faceOffsets[0] < 5 )
+    {
+      if ( faceOffsets[0] == 3 )
+        line = "E3T ";
+      if ( faceOffsets[0] == 4 )
+        line = "E4Q ";
+
+      line.append( std::to_string( i + 1 ) );
+
+      for ( int j = 0; j < faceOffsets[0]; ++j )
+      {
+        line.append( " " );
+        line.append( std::to_string( vertexIndices[j] + 1 ) );
+      }
+    }
+
+    file << line << std::endl;
+
+  }
+  file.close();
+}
diff --git a/external/mdal/frmts/mdal_2dm.hpp b/external/mdal/frmts/mdal_2dm.hpp
index a30f24afc0a..4af561664a8 100644
--- a/external/mdal/frmts/mdal_2dm.hpp
+++ b/external/mdal/frmts/mdal_2dm.hpp
@@ -14,6 +14,9 @@
 #include "mdal.h"
 #include "mdal_driver.hpp"
 
+
+#define MAX_VERTICES_PER_FACE_2DM 4
+
 namespace MDAL
 {
   class Mesh2dm: public MemoryMesh
@@ -80,8 +83,12 @@ namespace MDAL
       ~Driver2dm() override;
       Driver2dm *create() override;
 
+      int faceVerticesMaximumCount() const override
+      {return MAX_VERTICES_PER_FACE_2DM;}
+
       bool canRead( const std::string &uri ) override;
       std::unique_ptr< Mesh > load( const std::string &meshFile, MDAL_Status *status ) override;
+      void save( const std::string &uri, Mesh *mesh, MDAL_Status *status ) override;
 
     private:
       std::string mMeshFile;
diff --git a/external/mdal/frmts/mdal_ascii_dat.cpp b/external/mdal/frmts/mdal_ascii_dat.cpp
index efb4f94335b..9983ecd03a0 100644
--- a/external/mdal/frmts/mdal_ascii_dat.cpp
+++ b/external/mdal/frmts/mdal_ascii_dat.cpp
@@ -142,7 +142,7 @@ void MDAL::DriverAsciiDat::loadNewFormat( std::ifstream &in,
   std::shared_ptr<DatasetGroup> group; // DAT outputs data
   std::string groupName( MDAL::baseName( mDatFile ) );
   std::string line;
-
+  std::string referenceTime;
   // see if it contains face-centered results - supported by BASEMENT
   bool faceCentered = false;
   if ( contains( groupName, "_els_" ) )
@@ -201,6 +201,7 @@ void MDAL::DriverAsciiDat::loadNewFormat( std::ifstream &in,
               );
       group->setIsScalar( !isVector );
       group->setIsOnVertices( !faceCentered );
+      group->setReferenceTime( referenceTime );
     }
     else if ( cardType == "ENDDS" )
     {
@@ -226,9 +227,24 @@ void MDAL::DriverAsciiDat::loadNewFormat( std::ifstream &in,
       if ( quoteIdx1 != std::string::npos && quoteIdx2 != std::string::npos )
         group->setName( line.substr( quoteIdx1 + 1, quoteIdx2 - quoteIdx1 - 1 ) );
     }
+    else if ( cardType == "RT_JULIAN" && items.size() >= 2 )
+    {
+      referenceTime = "JULIAN " + items[1];
+    }
+    else if ( cardType == "TIMEUNITS" && items.size() >= 2 )
+    {
+      if ( !group )
+      {
+        debug( "TIMEUNITS card for no active dataset!" );
+        EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
+      }
+
+      group->setMetadata( "TIMEUNITS", items[1] );
+    }
     else if ( cardType == "TS" && items.size() >= 3 )
     {
       double t = toDouble( items[2] );
+      t = convertTimeDataToHours( t, group->getMetadata( "TIMEUNITS" ) );
 
       if ( faceCentered )
       {
@@ -320,7 +336,7 @@ void MDAL::DriverAsciiDat::readVertexTimestep(
   size_t faceCount = mesh->facesCount();
 
   std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
-  dataset->setTime( t / 3600. ); // TODO read TIMEUNITS
+  dataset->setTime( t );
 
   int *active = dataset->active();
   // only for new format
@@ -392,7 +408,7 @@ void MDAL::DriverAsciiDat::readFaceTimestep(
   size_t faceCount = mesh->facesCount();
 
   std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
-  dataset->setTime( t / 3600. );
+  dataset->setTime( t );
   double *values = dataset->values();
   // TODO: hasStatus
   for ( size_t index = 0; index < faceCount; ++index )
@@ -427,3 +443,22 @@ void MDAL::DriverAsciiDat::readFaceTimestep(
   dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
   group->datasets.push_back( dataset );
 }
+
+
+double MDAL::DriverAsciiDat::convertTimeDataToHours( double time, const std::string &originalTimeDataUnit ) const
+{
+  if ( originalTimeDataUnit == "se" || originalTimeDataUnit == "2" || originalTimeDataUnit == "Seconds"
+       || originalTimeDataUnit.empty() )
+  {
+    time /= 3600.0;
+  }
+  else if ( originalTimeDataUnit == "mi" || originalTimeDataUnit == "1" || originalTimeDataUnit == "Minutes" )
+  {
+    time /= 60.0;
+  }
+  else if ( originalTimeDataUnit == "days" )
+  {
+    time *= 24;
+  }
+  return time;
+}
diff --git a/external/mdal/frmts/mdal_ascii_dat.hpp b/external/mdal/frmts/mdal_ascii_dat.hpp
index ab2af7ec8fa..2a2cdae61d5 100644
--- a/external/mdal/frmts/mdal_ascii_dat.hpp
+++ b/external/mdal/frmts/mdal_ascii_dat.hpp
@@ -84,6 +84,8 @@ namespace MDAL
         bool isVector,
         std::ifstream &stream ) const;
 
+      double convertTimeDataToHours( double time, const std::string &originalTimeDataUnit ) const;
+
       std::string mDatFile;
   };
 
diff --git a/external/mdal/frmts/mdal_binary_dat.cpp b/external/mdal/frmts/mdal_binary_dat.cpp
index 115248c70ef..c562b91cf9d 100644
--- a/external/mdal/frmts/mdal_binary_dat.cpp
+++ b/external/mdal/frmts/mdal_binary_dat.cpp
@@ -33,8 +33,8 @@ static const int CT_NUMCELLS  = 180;
 static const int CT_NAME      = 190;
 static const int CT_TS        = 200;
 static const int CT_ENDDS     = 210;
-//static const int CT_RT_JULIAN = 240;
-//static const int CT_TIMEUNITS = 250;
+static const int CT_RT_JULIAN = 240;
+static const int CT_TIMEUNITS = 250;
 static const int CT_2D_MESHES = 3;
 static const int CT_FLOAT_SIZE = 4;
 static const int CF_FLAG_SIZE = 1;
@@ -142,6 +142,9 @@ void MDAL::DriverBinaryDat::load( const std::string &datFile, MDAL::Mesh *mesh,
   int numdata;
   int numcells;
   char groupName[40];
+  double referenceTime;
+  int timeUnit = 0;
+  std::string timeUnitStr;
   char istat;
   float time;
 
@@ -245,6 +248,44 @@ void MDAL::DriverBinaryDat::load( const std::string &datFile, MDAL::Mesh *mesh,
         groupMax->setName( group->name() + "/Maximums" );
         break;
 
+      case CT_RT_JULIAN:
+        // Reference time
+        if ( readIStat( in, sflg, &istat ) )
+          EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
+
+        if ( read( in, reinterpret_cast< char * >( &time ), 8 ) )
+          EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
+
+        referenceTime = static_cast<double>( time );
+        group->setReferenceTime( "JULIAN " + std::to_string( referenceTime ) );
+        break;
+
+      case CT_TIMEUNITS:
+        // Time unit
+        if ( read( in, reinterpret_cast< char * >( &timeUnit ), 4 ) )
+          EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
+
+        switch ( timeUnit )
+        {
+          case 0:
+            timeUnitStr = "hours";
+            break;
+          case 1:
+            timeUnitStr = "minutes";
+            break;
+          case 2:
+            timeUnitStr = "seconds";
+            break;
+          case 4:
+            timeUnitStr = "days";
+            break;
+          default:
+            timeUnitStr = "unknown";
+            break;
+        }
+        group->setMetadata( "TIMEUNITS", timeUnitStr );
+        break;
+
       case CT_TS:
         // Time step!
         if ( readIStat( in, sflg, &istat ) )
@@ -254,6 +295,8 @@ void MDAL::DriverBinaryDat::load( const std::string &datFile, MDAL::Mesh *mesh,
           EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
 
         double t = static_cast<double>( time );
+        t = convertTimeDataToHours( t, timeUnit );
+
         if ( readVertexTimestep( mesh, group, groupMax, t, istat, sflg, in ) )
           EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
 
@@ -461,3 +504,23 @@ bool MDAL::DriverBinaryDat::persist( MDAL::DatasetGroup *group )
 
   return false;
 }
+
+double MDAL::DriverBinaryDat::convertTimeDataToHours( double time, int originalTimeDataUnit )
+{
+  switch ( originalTimeDataUnit )
+  {
+    case 1:
+      time /= 60.0;
+      break;
+    case 2:
+      time /= 3600.0;
+      break;
+    case 4:
+      time *= 24;
+      break;
+    case 0:
+    default:
+      break;
+  }
+  return time;
+}
diff --git a/external/mdal/frmts/mdal_binary_dat.hpp b/external/mdal/frmts/mdal_binary_dat.hpp
index fab0ba6732f..acb5255e66f 100644
--- a/external/mdal/frmts/mdal_binary_dat.hpp
+++ b/external/mdal/frmts/mdal_binary_dat.hpp
@@ -40,6 +40,7 @@ namespace MDAL
                                int sflg,
                                std::ifstream &in );
 
+      double convertTimeDataToHours( double time, int originalTimeDataUnit );
       std::string mDatFile;
   };
 
diff --git a/external/mdal/frmts/mdal_cf.cpp b/external/mdal/frmts/mdal_cf.cpp
index 9cf64e204a5..ba8b3bcb064 100644
--- a/external/mdal/frmts/mdal_cf.cpp
+++ b/external/mdal/frmts/mdal_cf.cpp
@@ -71,8 +71,24 @@ MDAL::cfdataset_info_map MDAL::DriverCF::parseDatasetGroupInfo()
       }
       else
       {
+        /*UGRID convention says "The order of dimensions on a data variable is arbitrary",
+        *So time dimension is not necessary the first dimension, event if the convention recommands it.
+        *And prevent that any of the the two dimensions is not time dimension
+        */
+        if ( mDimensions.type( dimids[0] ) == CFDimensions::Time )
+        {
+          dimid = dimids[1];
+        }
+        else if ( mDimensions.type( dimids[1] ) == CFDimensions::Time )
+        {
+          dimid = dimids[0];
+        }
+        else
+        {
+          continue;
+        }
+
         nTimesteps = mDimensions.size( CFDimensions::Time );
-        dimid = dimids[1];
       }
 
       if ( !mDimensions.isDatasetType( mDimensions.type( dimid ) ) )
@@ -169,6 +185,34 @@ std::shared_ptr<MDAL::Dataset> MDAL::DriverCF::createFace2DDataset( std::shared_
   return dataset;
 }
 
+std::shared_ptr<MDAL::Dataset> MDAL::DriverCF::createVertex2DDataset( std::shared_ptr<MDAL::DatasetGroup> group,
+    size_t ts, const MDAL::CFDatasetGroupInfo &dsi,
+    const std::vector<double> &vals_x,
+    const std::vector<double> &vals_y,
+    double fill_val_x, double fill_val_y )
+{
+  assert( dsi.outputType == CFDimensions::Vertex2D );
+  size_t nVertices2D = mDimensions.size( CFDimensions::Vertex2D );
+
+  std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared<MDAL::MemoryDataset>( group.get() );
+
+  for ( size_t i = 0; i < nVertices2D; ++i )
+  {
+    size_t idx = ts * nVertices2D + i;
+    populate_vals( dsi.is_vector,
+                   dataset->values(),
+                   i,
+                   vals_x,
+                   vals_y,
+                   idx,
+                   fill_val_x,
+                   fill_val_y );
+
+  }
+
+  return dataset;
+}
+
 void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<double> &times, const MDAL::cfdataset_info_map &dsinfo_map )
 {
   /* PHASE 2 - add dataset groups */
@@ -183,7 +227,16 @@ void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
           dsi.name
         );
     group->setIsScalar( !dsi.is_vector );
-    group->setIsOnVertices( false );
+
+    if ( dsi.outputType == CFDimensions::Vertex2D )
+      group->setIsOnVertices( true );
+    else if ( dsi.outputType == CFDimensions::Face2D )
+      group->setIsOnVertices( false );
+    else
+    {
+      // unsupported
+      continue;
+    }
 
     // read X data
     double fill_val_x = mNcFile.getFillValue( dsi.ncid_x );
@@ -209,6 +262,13 @@ void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
       if ( dsi.outputType == CFDimensions::Face2D )
       {
         dataset = createFace2DDataset( group, ts, dsi, vals_x, vals_y, fill_val_x, fill_val_y );
+      }
+      else     // Vertex2D
+      {
+        dataset = createVertex2DDataset( group, ts, dsi, vals_x, vals_y, fill_val_x, fill_val_y );
+      }
+      if ( dataset )
+      {
         dataset->setTime( time );
         dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
         group->datasets.push_back( dataset );
@@ -228,6 +288,13 @@ void MDAL::DriverCF::parseTime( std::vector<double> &times )
 {
 
   size_t nTimesteps = mDimensions.size( CFDimensions::Time );
+  if ( 0 == nTimesteps )
+  {
+    //if no time dimension is present creates only one time step to store the potential time-independent variable
+    nTimesteps = 1;
+    times = std::vector<double>( 1, 0 );
+    return;
+  }
   times = mNcFile.readDoubleArr( "time", nTimesteps );
   std::string units = mNcFile.getAttrStr( "time", "units" );
   double div_by = MDAL::parseTimeUnits( units );
diff --git a/external/mdal/frmts/mdal_cf.hpp b/external/mdal/frmts/mdal_cf.hpp
index c73f2bd3eb2..e1c29409fc9 100644
--- a/external/mdal/frmts/mdal_cf.hpp
+++ b/external/mdal/frmts/mdal_cf.hpp
@@ -94,6 +94,14 @@ namespace MDAL
         const std::vector<double> &vals_y,
         double fill_val_x, double fill_val_y );
 
+      std::shared_ptr<MDAL::Dataset> createVertex2DDataset(
+        std::shared_ptr<MDAL::DatasetGroup> group,
+        size_t ts,
+        const MDAL::CFDatasetGroupInfo &dsi,
+        const std::vector<double> &vals_x,
+        const std::vector<double> &vals_y,
+        double fill_val_x, double fill_val_y );
+
       void addDatasetGroups( Mesh *mesh,
                              const std::vector<double> &times,
                              const cfdataset_info_map &dsinfo_map );
diff --git a/external/mdal/frmts/mdal_driver.cpp b/external/mdal/frmts/mdal_driver.cpp
index d717ab9a6b0..84b892f7900 100644
--- a/external/mdal/frmts/mdal_driver.cpp
+++ b/external/mdal/frmts/mdal_driver.cpp
@@ -42,6 +42,11 @@ bool MDAL::Driver::hasCapability( MDAL::Capability capability ) const
   return capability == ( mCapabilityFlags & capability );
 }
 
+int MDAL::Driver::faceVerticesMaximumCount() const
+{
+  return -1;
+}
+
 std::unique_ptr< MDAL::Mesh > MDAL::Driver::load( const std::string &uri, MDAL_Status *status )
 {
   MDAL_UNUSED( uri );
@@ -57,6 +62,13 @@ void MDAL::Driver::load( const std::string &uri, Mesh *mesh, MDAL_Status *status
   return;
 }
 
+void MDAL::Driver::save( const std::string &uri, MDAL::Mesh *mesh, MDAL_Status *status )
+{
+  MDAL_UNUSED( uri );
+  MDAL_UNUSED( mesh );
+  MDAL_UNUSED( status );
+}
+
 void MDAL::Driver::createDatasetGroup( MDAL::Mesh *mesh, const std::string &groupName, bool isOnVertices, bool hasScalarData, const std::string &datasetGroupFile )
 {
   std::shared_ptr<MDAL::DatasetGroup> grp(
diff --git a/external/mdal/frmts/mdal_driver.hpp b/external/mdal/frmts/mdal_driver.hpp
index 564827e6b19..36a779cbeef 100644
--- a/external/mdal/frmts/mdal_driver.hpp
+++ b/external/mdal/frmts/mdal_driver.hpp
@@ -18,6 +18,7 @@ namespace MDAL
     ReadMesh      = 1 << 0, //! Can read mesh and all datasets stored in the mesh file
     ReadDatasets  = 1 << 1, //! Can read only datasets (groups) from existing mesh
     WriteDatasets = 1 << 2, //! Can write datasets (groups)
+    SaveMesh      = 1 << 3, //! Can save the mesh
   };
 
   class Driver
@@ -39,11 +40,17 @@ namespace MDAL
 
       virtual bool canRead( const std::string &uri ) = 0;
 
+      //! returns the maximum vertices per face
+      virtual int faceVerticesMaximumCount() const;
+
       // loads mesh
       virtual std::unique_ptr< Mesh > load( const std::string &uri, MDAL_Status *status );
       // loads datasets
       virtual void load( const std::string &uri, Mesh *mesh, MDAL_Status *status );
 
+      // save mesh
+      virtual void save( const std::string &uri, Mesh *mesh, MDAL_Status *status );
+
       // create new dataset group
       virtual void createDatasetGroup(
         Mesh *mesh,
diff --git a/external/mdal/frmts/mdal_esri_tin.cpp b/external/mdal/frmts/mdal_esri_tin.cpp
new file mode 100644
index 00000000000..04723267c6f
--- /dev/null
+++ b/external/mdal/frmts/mdal_esri_tin.cpp
@@ -0,0 +1,300 @@
+/*
+ MDAL - Mesh Data Abstraction Library (MIT License)
+ Copyright (C) 2019 Vincent Cloarec (vcloarec at gmail dot com)
+*/
+
+#include "mdal_esri_tin.hpp"
+
+MDAL::DriverEsriTin::DriverEsriTin(): Driver( "ESRI_TIN",
+      "Esri TIN",
+      "*.adf",
+      Capability::ReadMesh )
+{}
+
+MDAL::Driver *MDAL::DriverEsriTin::create()
+{
+  return new DriverEsriTin();
+}
+
+std::unique_ptr<MDAL::Mesh> MDAL::DriverEsriTin::load( const std::string &uri, MDAL_Status *status )
+{
+  if ( status ) *status = MDAL_Status::None;
+
+  try
+  {
+    std::list<int> superpointIndexes;
+    Vertices vertices;
+    Faces faces;
+
+    readSuperpoints( uri, superpointIndexes );
+    populateVertices( uri, vertices, superpointIndexes );
+    populateFaces( uri, faces, superpointIndexes );
+
+    std::unique_ptr< MemoryMesh > mesh(
+      new MemoryMesh(
+        name(),
+        vertices.size(),
+        faces.size(),
+        3,
+        computeExtent( vertices ),
+        uri
+      )
+    );
+
+    mesh->faces = std::move( faces );
+    mesh->vertices = std::move( vertices );
+
+    addBedElevationDatasetGroup( mesh.get(), mesh->vertices );
+    mesh->datasetGroups.back()->setName( "Altitude" );
+
+    std::string crs = getCrsWkt( uri );
+    if ( ! crs.empty() )
+      mesh->setSourceCrsFromWKT( crs );
+
+    return std::unique_ptr<Mesh>( mesh.release() );
+  }
+  catch ( MDAL_Status error )
+  {
+    if ( status ) *status = ( error );
+    return std::unique_ptr<Mesh>();
+  }
+}
+
+bool MDAL::DriverEsriTin::canRead( const std::string &uri )
+{
+
+  std::string zFileName = zFile( uri );
+  std::string faceFileName = faceFile( uri );
+
+  std::ifstream xyIn( xyFile( uri ), std::ifstream::in | std::ifstream::binary );
+  if ( ! xyIn.is_open() )
+    return false;
+
+  std::ifstream zIn( zFile( uri ), std::ifstream::in | std::ifstream::binary );
+  if ( ! zIn.is_open() )
+    return false;
+
+  std::ifstream faceIn( faceFile( uri ), std::ifstream::in | std::ifstream::binary );
+  if ( ! faceIn.is_open() )
+    return false;
+
+  std::ifstream hullIn( hullFile( uri ), std::ifstream::in | std::ifstream::binary );
+  if ( ! hullIn.is_open() )
+    return false;
+
+  return true;
+}
+
+std::string MDAL::DriverEsriTin::xyFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "tnxy.adf" );
+}
+
+std::string MDAL::DriverEsriTin::zFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "tnz.adf" );
+}
+
+std::string MDAL::DriverEsriTin::faceFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "tnod.adf" );
+}
+
+std::string MDAL::DriverEsriTin::mskFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "tmsk.adf" );
+}
+
+std::string MDAL::DriverEsriTin::msxFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "tmsx.adf" );
+}
+
+std::string MDAL::DriverEsriTin::hullFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "thul.adf" );
+}
+
+std::string MDAL::DriverEsriTin::crsFile( const std::string &uri ) const
+{
+  return pathJoin( dirName( uri ), "prj.adf" );
+}
+
+void MDAL::DriverEsriTin::readSuperpoints( const std::string &uri, std::list<int> &superpointsIndexes ) const
+{
+  superpointsIndexes.clear();
+  bool isNativeLittleEndian = MDAL::isNativeLittleEndian();
+  std::ifstream inHull( hullFile( uri ), std::ifstream::in | std::ifstream::binary );
+
+  int32_t index;
+  while ( readValue( index, inHull, isNativeLittleEndian ) && index != -1 )
+    superpointsIndexes.push_back( index );
+
+  superpointsIndexes.sort();
+}
+
+std::string MDAL::DriverEsriTin::getTinName( const std::string &uri ) const
+{
+  std::string tinName = uri;
+  size_t last_slash_idx = tinName.find_last_of( "\\/" );
+
+  if ( last_slash_idx == std::string::npos )
+    return "";
+
+  tinName.erase( last_slash_idx, tinName.size() - last_slash_idx );
+  last_slash_idx = tinName.find_last_of( "\\/" );
+
+  if ( last_slash_idx == std::string::npos )
+    return "";
+
+  tinName.erase( 0, last_slash_idx + 1 );
+
+  return tinName;
+}
+
+void MDAL::DriverEsriTin::populateVertices( const std::string &uri, MDAL::Vertices &vertices, const std::list<int> &superpointIndexes ) const
+{
+  //first, read the (x,y) coordinates
+
+  bool isNativeLittleEndian = MDAL::isNativeLittleEndian();
+  std::ifstream inXY( xyFile( uri ), std::ifstream::in | std::ifstream::binary );
+  std::ifstream inZ( zFile( uri ), std::ifstream::in | std::ifstream::binary );
+
+  if ( ! inXY.is_open() )
+    throw MDAL_Status::Err_FileNotFound;
+
+  if ( ! inZ.is_open() )
+    throw MDAL_Status::Err_FileNotFound;
+
+  int fileIndex = 1;
+  while ( true )
+  {
+    Vertex vert;
+
+    if ( !readValue( vert.x, inXY, isNativeLittleEndian ) )
+      break;
+
+    if ( !readValue( vert.y, inXY, isNativeLittleEndian ) )
+      throw MDAL_Status::Err_UnknownFormat;
+
+    float zValue;
+    if ( !readValue( zValue, inZ, isNativeLittleEndian ) )
+      throw MDAL_Status::Err_UnknownFormat;
+    vert.z = double( zValue );
+
+    if ( correctedIndex( fileIndex, superpointIndexes ) >= 0 ) //store the vertex only if it is not a superpoint
+      vertices.push_back( vert );
+    fileIndex++;
+  }
+
+  inXY.close();
+  inZ.close();
+}
+
+void MDAL::DriverEsriTin::populateFaces( const std::string &uri, MDAL::Faces &faces, const std::list<int> &superpointIndexes ) const
+{
+  bool isNativeLittleEndian = MDAL::isNativeLittleEndian();
+  std::ifstream inFaces( faceFile( uri ), std::ifstream::in | std::ifstream::binary );
+  std::ifstream inMsk( mskFile( uri ), std::ifstream::in | std::ifstream::binary );
+  std::ifstream inMsx( msxFile( uri ), std::ifstream::in | std::ifstream::binary );
+
+  if ( ! inFaces.is_open() )
+    throw MDAL_Status::Err_FileNotFound;
+  if ( ! inMsk.is_open() )
+    throw MDAL_Status::Err_FileNotFound;
+  if ( ! inMsx.is_open() )
+    throw MDAL_Status::Err_FileNotFound;
+
+  //Find the beginning of data in mskFile
+  inMsx.seekg( -4, std::ios::end );
+  int32_t mskBegin;
+  if ( ! readValue( mskBegin, inMsx, true ) )
+    throw  MDAL_Status::Err_UnknownFormat;
+
+  //read information in mskFile
+  inMsk.seekg( -mskBegin * 2, std::ios::end );
+  int32_t maskIntergerCount;
+  if ( ! readValue( maskIntergerCount, inMsk, true ) )
+    throw MDAL_Status::Err_UnknownFormat;
+  inMsk.ignore( 4 ); //unused 4 bytes
+  int32_t maskBitsCount;
+  if ( ! readValue( maskBitsCount, inMsk, true ) )
+    throw MDAL_Status::Err_UnknownFormat;
+
+  int c = 0;
+  int32_t maskInt = 0;
+  while ( true )
+  {
+    //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_Status::Err_UnknownFormat;
+
+    Face f;
+    for ( int i = 0; i < 3; ++i )
+    {
+      int32_t index;
+      if ( ! readValue( index, inFaces, isNativeLittleEndian ) )
+        break;
+
+      index = correctedIndex( index, superpointIndexes );
+
+      f.push_back( static_cast<size_t>( index ) );
+    }
+
+    if ( f.size() == 0 ) //that's mean this is the end of the file
+      break;
+
+    if ( f.size() < 3 ) //that's mean the face is not complete
+      throw MDAL_Status::Err_UnknownFormat;
+
+    //exclude masked face
+    if ( !( maskInt & 0x01 ) )
+      faces.push_back( f );
+
+    c++;
+    maskInt = maskInt >> 1;
+  }
+
+  inFaces.close();
+  inMsk.close();
+  inMsx.close();
+}
+
+std::string MDAL::DriverEsriTin::getCrsWkt( const std::string &uri ) const
+{
+  std::ifstream inCRS( crsFile( uri ), std::ifstream::in );
+  if ( ! inCRS.is_open() )
+    return std::string();
+
+  std::string crsWkt;
+  std::getline( inCRS, crsWkt );
+
+  if ( crsWkt == "{B286C06B-0879-11D2-AACA-00C04FA33C20}" ) //com class id of the esri UnknownCoordinateSystem class
+    crsWkt = "";
+
+  return crsWkt;
+}
+
+int MDAL::DriverEsriTin::correctedIndex( int i, const std::list<int> &superPointsIndexes ) const
+{
+  //usualy i is bigger than the biggest superpoint's index
+  if ( i > superPointsIndexes.back() )
+    return i - int( superPointsIndexes.size() ) - 1;
+
+  //correction of index depending of the position of the superpoints
+  int correctedIndex = i - 1;
+  for ( auto si : superPointsIndexes )
+  {
+    if ( i == si )
+      return -1; //return -1 if the index is associated with superpoint
+    else if ( i > si )
+    {
+      correctedIndex--;
+    }
+    else
+      break;
+  }
+
+  return correctedIndex;
+}
diff --git a/external/mdal/frmts/mdal_esri_tin.hpp b/external/mdal/frmts/mdal_esri_tin.hpp
new file mode 100644
index 00000000000..eb6c3b7d4e4
--- /dev/null
+++ b/external/mdal/frmts/mdal_esri_tin.hpp
@@ -0,0 +1,100 @@
+/*
+ MDAL - Mesh Data Abstraction Library (MIT License)
+ Copyright (C) 2019 Vincent Cloarec (vcloarec at gmail dot com)
+*/
+
+#ifndef MDAL_ESRI_TIN_H
+#define MDAL_ESRI_TIN_H
+
+#include <string>
+#include <vector>
+#include <list>
+#include <memory>
+#include <iosfwd>
+#include <iostream>
+#include <fstream>
+
+#include "mdal_data_model.hpp"
+#include "mdal.h"
+#include "mdal_driver.hpp"
+#include "mdal_utils.hpp"
+
+namespace MDAL
+{
+
+  /** *********************************************
+   *                   Structure of files
+   * https://en.wikipedia.org/wiki/Esri_TIN
+   ************************************************
+   * thul.adf
+   *
+   * file containing only int (big-endianness)
+   * first, list of superpoint indexes (infinity point to North, South East and West)
+   * used by ArcGIS but unused here
+   * This list finished with value -1
+   * Then indexes of the outer boundary indexes or inner boundary indexes (holes)
+   * Those list are separated with value 0
+   *
+   * ***************************************
+   * tnxy.adf
+   *
+   * vertices position X,Y : for each value 8 bytes : double with big-endianness
+   *
+   * ***************************************
+   * tnz.adf : z value of each vertices
+   *
+   * vertices Z value: for each value 4 bytes : float with big endianness
+   *
+   * ***************************************
+   * tnod.adf : faces (triangles) of the TIN
+   *
+   * indexes of the 3 vertices for each faces (int32 big-endianness)
+   *
+   * ***************************************
+   * tmsk.adf and tmsx.adf
+   *
+   * files used to store information about masked faces. Those faces are not displayed
+   *
+   *****************************************
+   * prj.adf
+   *
+   * text file containing the WKT CRS
+   *
+  */
+
+
+  class DriverEsriTin: public Driver
+  {
+    public:
+      DriverEsriTin();
+      ~DriverEsriTin() override {}
+
+      Driver *create() override;
+
+      virtual std::unique_ptr< Mesh > load( const std::string &uri, MDAL_Status *status ) override;
+
+      bool canRead( const std::string &uri ) override;
+
+    private:
+      std::string xyFile( const std::string &uri ) const;
+      std::string zFile( const std::string &uri ) const;
+      std::string faceFile( const std::string &uri ) const;
+      std::string mskFile( const std::string &uri ) const;
+      std::string msxFile( const std::string &uri ) const;
+      std::string hullFile( const std::string &uri ) const;
+      std::string crsFile( const std::string &uri ) const;
+
+      void readSuperpoints( const std::string &uri, std::list<int> &superpointsIndexes ) const;
+      void populateVertices( const std::string &uri, Vertices &vertices, const std::list<int> &superpointIndexes ) const;
+      void populateFaces( const std::string &uri, Faces &faces, const std::list<int> &superpointIndexes ) const;
+
+      std::string getCrsWkt( const std::string &uri ) const;
+      std::string getTinName( const std::string &uri ) const;
+
+      // correction of vertex index because indexing in ESRI TIN files begins with 1 and to take account superpoints position
+      // return -1 if i is the index of a superpoint
+      int correctedIndex( int i, const std::list<int> &superPointsIndexes ) const;
+  };
+}
+
+#endif // MDAL_ESRI_TIN_H
diff --git a/external/mdal/frmts/mdal_hdf5.cpp b/external/mdal/frmts/mdal_hdf5.cpp
index 719e18e14ff..1b0f08d0b16 100644
--- a/external/mdal/frmts/mdal_hdf5.cpp
+++ b/external/mdal/frmts/mdal_hdf5.cpp
@@ -4,6 +4,7 @@
 */
 
 #include "mdal_hdf5.hpp"
+#include <cstring>
 
 
 HdfFile::HdfFile( const std::string &path )
@@ -68,17 +69,21 @@ hid_t HdfAttribute::id() const { return d->id; }
 
 std::string HdfAttribute::readString() const
 {
-  char name[HDF_MAX_NAME];
-  hid_t datatype = H5Tcopy( H5T_C_S1 );
-  H5Tset_size( datatype, HDF_MAX_NAME );
+  hid_t datatype = H5Aget_type( id() );
+  char name[HDF_MAX_NAME + 1];
+  std::memset( name, '\0', HDF_MAX_NAME + 1 );
   herr_t status = H5Aread( d->id, datatype, name );
   if ( status < 0 )
   {
     //MDAL::debug("Failed to read data!");
+    H5Tclose( datatype );
     return std::string();
   }
+
   H5Tclose( datatype );
-  return std::string( name );
+  std::string res( name );
+  res = MDAL::trim( res );
+  return res;
 }
 
 HdfDataset::HdfDataset( hid_t file, const std::string &path )
@@ -189,7 +194,6 @@ std::string HdfDataset::readString() const
   return std::string( name );
 }
 
-
 HdfDataspace::HdfDataspace( const std::vector<hsize_t> &dims )
   : d( std::make_shared< Handle >( H5Screate_simple(
                                      static_cast<int>( dims.size() ),
diff --git a/external/mdal/frmts/mdal_hdf5.hpp b/external/mdal/frmts/mdal_hdf5.hpp
index f8a7609fca7..581dd11aaa8 100644
--- a/external/mdal/frmts/mdal_hdf5.hpp
+++ b/external/mdal/frmts/mdal_hdf5.hpp
@@ -97,7 +97,6 @@ class HdfGroup
     inline HdfDataset dataset( const std::string &dsName ) const;
     inline HdfAttribute attribute( const std::string &attr_name ) const;
     inline bool pathExists( const std::string &path ) const;
-
   protected:
     std::vector<std::string> objects( H5G_obj_t type ) const;
 
@@ -176,6 +175,8 @@ class HdfDataset
     std::vector<double> readArrayDouble( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
     std::vector<int> readArrayInt( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
 
+    inline bool hasAttribute( const std::string &attr_name ) const;
+    inline HdfAttribute attribute( const std::string &attr_name ) const;
 
     template <typename T> std::vector<T> readArray( hid_t mem_type_id ) const
     {
@@ -233,10 +234,18 @@ inline HdfGroup HdfGroup::group( const std::string &groupName ) const { return H
 
 inline HdfDataset HdfGroup::dataset( const std::string &dsName ) const { return HdfDataset( file_id(), childPath( dsName ) ); }
 
+inline bool HdfDataset::hasAttribute( const std::string &attr_name ) const
+{
+  htri_t res = H5Aexists( d->id, attr_name.c_str() );
+  return  res > 0 ;
+}
+
 inline HdfAttribute HdfFile::attribute( const std::string &attr_name ) const { return HdfAttribute( d->id, attr_name ); }
 
 inline HdfAttribute HdfGroup::attribute( const std::string &attr_name ) const { return HdfAttribute( d->id, attr_name ); }
 
+inline HdfAttribute HdfDataset::attribute( const std::string &attr_name ) const { return HdfAttribute( d->id, attr_name ); }
+
 inline bool HdfFile::pathExists( const std::string &path ) const { return H5Lexists( d->id, path.c_str(), H5P_DEFAULT ) > 0; }
 
 inline bool HdfGroup::pathExists( const std::string &path ) const { return H5Lexists( d->id, path.c_str(), H5P_DEFAULT ) > 0; }
diff --git a/external/mdal/frmts/mdal_hec2d.cpp b/external/mdal/frmts/mdal_hec2d.cpp
index 2f75f19553b..73614dd570a 100644
--- a/external/mdal/frmts/mdal_hec2d.cpp
+++ b/external/mdal/frmts/mdal_hec2d.cpp
@@ -43,7 +43,6 @@ static HdfDataset openHdfDataset( const HdfGroup &hdfGroup, const std::string &n
   return dsFileType;
 }
 
-
 static std::string openHdfAttribute( const HdfFile &hdfFile, const std::string &name )
 {
   HdfAttribute attr = hdfFile.attribute( name );
@@ -51,6 +50,13 @@ static std::string openHdfAttribute( const HdfFile &hdfFile, const std::string &
   return attr.readString();
 }
 
+static std::string openHdfAttribute( const HdfDataset &hdfDataset, const std::string &name )
+{
+  HdfAttribute attr = hdfDataset.attribute( name );
+  if ( !attr.isValid() ) throw MDAL_Status::Err_UnknownFormat;
+  return attr.readString();
+}
+
 static HdfGroup getBaseOutputGroup( const HdfFile &hdfFile )
 {
   HdfGroup gResults = openHdfGroup( hdfFile, "Results" );
@@ -69,12 +75,64 @@ static HdfGroup get2DFlowAreasGroup( const HdfFile &hdfFile, const std::string l
   return g2DFlowRes;
 }
 
+static std::string getDataTimeUnit( HdfDataset &dsTime )
+{
+  // Initially we expect data to be in hours
+  std::string dataTimeUnit = "Hours";
+
+  // First we look for the Time attribute
+  if ( dsTime.hasAttribute( "Time" ) )
+  {
+    dataTimeUnit = openHdfAttribute( dsTime, "Time" );
+    return dataTimeUnit;
+  }
+
+  // Some variants of HEC_RAS have time unit stored in Variables attribute
+  // With read values looking like "Time|Days       "
+  if ( dsTime.hasAttribute( "Variables" ) )
+  {
+    dataTimeUnit = openHdfAttribute( dsTime, "Variables" );
+    // Removing the "Time|" prefix
+    dataTimeUnit = MDAL::replace( dataTimeUnit, "Time|", "" );
+  }
+
+  return dataTimeUnit;
+}
+
+static void convertTimeDataToHours( std::vector<float> &times, const std::string &originalTimeDataUnit )
+{
+  if ( originalTimeDataUnit != "Hours" )
+  {
+    for ( size_t i = 0; i < times.size(); i++ )
+    {
+      if ( originalTimeDataUnit == "Seconds" ) { times[i] /= 3600.0f; }
+      else if ( originalTimeDataUnit == "Minutes" ) { times[i] /= 60.0f; }
+      else if ( originalTimeDataUnit == "Days" ) { times[i] *= 24; }
+    }
+  }
+}
+
+static std::string readReferenceTime( const HdfFile &hdfFile )
+{
+  std::string refTime;
+  HdfGroup gBaseO = getBaseOutputGroup( hdfFile );
+  HdfGroup gUnsteadTS = openHdfGroup( gBaseO, "Unsteady Time Series" );
+  HdfDataset dsTimeDateStamp = openHdfDataset( gUnsteadTS, "Time Date Stamp" );
+  std::vector<std::string> timeStamps = dsTimeDateStamp.readArrayString();
+
+  if ( timeStamps.size() > 0 )
+    refTime = timeStamps[0];
+  return refTime;
+}
+
 static std::vector<float> readTimes( const HdfFile &hdfFile )
 {
   HdfGroup gBaseO = getBaseOutputGroup( hdfFile );
   HdfGroup gUnsteadTS = openHdfGroup( gBaseO, "Unsteady Time Series" );
-  HdfDataset dsTimes = openHdfDataset( gUnsteadTS, "Time" );
-  std::vector<float> times = dsTimes.readArray();
+  HdfDataset dsTime = openHdfDataset( gUnsteadTS, "Time" );
+  std::string dataTimeUnits = getDataTimeUnit( dsTime );
+  std::vector<float> times = dsTime.readArray();
+  convertTimeDataToHours( times, dataTimeUnits );
   return times;
 }
 
@@ -99,7 +157,8 @@ void MDAL::DriverHec2D::readFaceOutput( const HdfFile &hdfFile,
                                         const std::vector<std::string> &flowAreaNames,
                                         const std::string rawDatasetName,
                                         const std::string datasetName,
-                                        const std::vector<float> &times )
+                                        const std::vector<float> &times,
+                                        const std::string &referenceTime )
 {
   double eps = std::numeric_limits<double>::min();
 
@@ -111,6 +170,7 @@ void MDAL::DriverHec2D::readFaceOutput( const HdfFile &hdfFile,
                                         );
   group->setIsOnVertices( false );
   group->setIsScalar( true );
+  group->setReferenceTime( referenceTime );
 
   std::vector<std::shared_ptr<MDAL::MemoryDataset>> datasets;
 
@@ -177,17 +237,17 @@ void MDAL::DriverHec2D::readFaceResults( const HdfFile &hdfFile,
   // UNSTEADY
   HdfGroup flowGroup = get2DFlowAreasGroup( hdfFile, "Unsteady Time Series" );
   std::vector<float> times = readTimes( hdfFile );
-
-  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Face Shear Stress", "Face Shear Stress", times );
-  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Face Velocity", "Face Velocity", times );
+  const std::string referenceTime = readReferenceTime( hdfFile );
+  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Face Shear Stress", "Face Shear Stress", times, referenceTime );
+  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Face Velocity", "Face Velocity", times, referenceTime );
 
   // SUMMARY
   flowGroup = get2DFlowAreasGroup( hdfFile, "Summary Output" );
   times.clear();
   times.push_back( 0.0f );
 
-  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Maximum Face Shear Stress", "Face Shear Stress/Maximums", times );
-  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Maximum Face Velocity", "Face Velocity/Maximums", times );
+  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Maximum Face Shear Stress", "Face Shear Stress/Maximums", times, referenceTime );
+  readFaceOutput( hdfFile, flowGroup, areaElemStartIndex, flowAreaNames, "Maximum Face Velocity", "Face Velocity/Maximums", times, referenceTime );
 }
 
 
@@ -197,7 +257,8 @@ std::shared_ptr<MDAL::MemoryDataset> MDAL::DriverHec2D::readElemOutput( const Hd
     const std::string rawDatasetName,
     const std::string datasetName,
     const std::vector<float> &times,
-    std::shared_ptr<MDAL::MemoryDataset> bed_elevation )
+    std::shared_ptr<MDAL::MemoryDataset> bed_elevation,
+    const std::string &referenceTime )
 {
   double eps = std::numeric_limits<double>::min();
 
@@ -209,6 +270,7 @@ std::shared_ptr<MDAL::MemoryDataset> MDAL::DriverHec2D::readElemOutput( const Hd
                                         );
   group->setIsOnVertices( false );
   group->setIsScalar( true );
+  group->setReferenceTime( referenceTime );
 
   std::vector<std::shared_ptr<MDAL::MemoryDataset>> datasets;
 
@@ -287,6 +349,8 @@ std::shared_ptr<MDAL::MemoryDataset> MDAL::DriverHec2D::readBedElevation(
   const std::vector<std::string> &flowAreaNames )
 {
   std::vector<float> times( 1, 0.0f );
+  std::string referenceTime;
+
   return readElemOutput(
            gGeom2DFlowAreas,
            areaElemStartIndex,
@@ -294,7 +358,8 @@ std::shared_ptr<MDAL::MemoryDataset> MDAL::DriverHec2D::readBedElevation(
            "Cells Minimum Elevation",
            "Bed Elevation",
            times,
-           std::shared_ptr<MDAL::MemoryDataset>()
+           std::shared_ptr<MDAL::MemoryDataset>(),
+           referenceTime
          );
 }
 
@@ -307,6 +372,7 @@ void MDAL::DriverHec2D::readElemResults(
   // UNSTEADY
   HdfGroup flowGroup = get2DFlowAreasGroup( hdfFile, "Unsteady Time Series" );
   std::vector<float> times = readTimes( hdfFile );
+  std::string referenceTime = readReferenceTime( hdfFile );
 
   readElemOutput(
     flowGroup,
@@ -315,7 +381,8 @@ void MDAL::DriverHec2D::readElemResults(
     "Water Surface",
     "Water Surface",
     times,
-    bed_elevation );
+    bed_elevation,
+    referenceTime );
   readElemOutput(
     flowGroup,
     areaElemStartIndex,
@@ -323,12 +390,14 @@ void MDAL::DriverHec2D::readElemResults(
     "Depth",
     "Depth",
     times,
-    bed_elevation );
+    bed_elevation,
+    referenceTime );
 
   // SUMMARY
   flowGroup = get2DFlowAreasGroup( hdfFile, "Summary Output" );
   times.clear();
   times.push_back( 0.0f );
+  referenceTime.clear();
 
   readElemOutput(
     flowGroup,
@@ -337,7 +406,8 @@ void MDAL::DriverHec2D::readElemResults(
     "Maximum Water Surface",
     "Water Surface/Maximums",
     times,
-    bed_elevation
+    bed_elevation,
+    referenceTime
   );
 }
 
diff --git a/external/mdal/frmts/mdal_hec2d.hpp b/external/mdal/frmts/mdal_hec2d.hpp
index 2abd8359efe..2f7b70e8ee9 100644
--- a/external/mdal/frmts/mdal_hec2d.hpp
+++ b/external/mdal/frmts/mdal_hec2d.hpp
@@ -26,6 +26,14 @@ namespace MDAL
    * There is a small change in the format in HEC-RAS 5.0.5+, where
    *    - Header File Type is different (HEC-RAS Results vs HEC-RAS Geometry)
    *    - Names or areas are stored in different place (names array vs attributes array)
+   *
+   * Time data unit should be present in Time dataset and Time or Variable attribute for given dataset root,
+   * Since MDAL API is reporting times in float hours, the original values need to be corrected
+   * based on value found in the Time attribute.
+   *
+   * All reference times can be found in Time Data Stamp dataset.
+   * First value in the dataset is reported by MDAL as reference time
+   *
    */
   class DriverHec2D: public Driver
   {
@@ -56,7 +64,8 @@ namespace MDAL
                            const std::vector<std::string> &flowAreaNames,
                            const std::string rawDatasetName,
                            const std::string datasetName,
-                           const std::vector<float> &times );
+                           const std::vector<float> &times,
+                           const std::string &referenceTime );
 
       void readFaceResults( const HdfFile &hdfFile,
                             const std::vector<size_t> &areaElemStartIndex,
@@ -69,7 +78,8 @@ namespace MDAL
         const std::string rawDatasetName,
         const std::string datasetName,
         const std::vector<float> &times,
-        std::shared_ptr<MDAL::MemoryDataset> bed_elevation );
+        std::shared_ptr<MDAL::MemoryDataset> bed_elevation,
+        const std::string &referenceTime );
 
       std::shared_ptr<MDAL::MemoryDataset> readBedElevation(
         const HdfGroup &gGeom2DFlowAreas,
diff --git a/external/mdal/frmts/mdal_netcdf.cpp b/external/mdal/frmts/mdal_netcdf.cpp
index 5675384020c..2707561dc49 100644
--- a/external/mdal/frmts/mdal_netcdf.cpp
+++ b/external/mdal/frmts/mdal_netcdf.cpp
@@ -82,6 +82,19 @@ std::vector<std::string> NetCDFFile::readArrNames() const
   return res;
 }
 
+bool NetCDFFile::hasAttrInt( const std::string &name, const std::string &attr_name ) const
+{
+  assert( mNcid != 0 );
+
+  int arr_id;
+  if ( nc_inq_varid( mNcid, name.c_str(), &arr_id ) != NC_NOERR ) return false;
+
+  int val;
+  if ( nc_get_att_int( mNcid, arr_id, attr_name.c_str(), &val ) ) return false;
+
+  return true;
+}
+
 int NetCDFFile::getAttrInt( const std::string &name, const std::string &attr_name ) const
 {
   assert( mNcid != 0 );
@@ -155,6 +168,27 @@ void NetCDFFile::getDimension( const std::string &name, size_t *val, int *ncid_v
   if ( nc_inq_dimlen( mNcid, *ncid_val, val ) != NC_NOERR ) throw MDAL_Status::Err_UnknownFormat;
 }
 
+void NetCDFFile::getDimensions( const std::string &variableName, std::vector<size_t> &dimensions, std::vector<int> &dimensionIds )
+{
+  assert( mNcid != 0 );
+
+  int n;
+  int varId;
+  if ( nc_inq_varid( mNcid, variableName.c_str(), &varId ) ) throw MDAL_Status::Err_UnknownFormat;
+  if ( nc_inq_varndims( mNcid, varId, &n ) ) throw MDAL_Status::Err_UnknownFormat;
+
+  dimensionIds.resize( size_t( n ) );
+  dimensions.resize( size_t( n ) );
+
+  if ( nc_inq_vardimid( mNcid, varId, dimensionIds.data() ) ) throw MDAL_Status::Err_UnknownFormat;
+
+  for ( int i = 0; i < n; ++i )
+  {
+    nc_inq_dimlen( mNcid, dimensionIds[size_t( i )], &dimensions[size_t( i )] );
+  }
+
+}
+
 bool NetCDFFile::hasDimension( const std::string &name ) const
 {
   int ncid_val;
diff --git a/external/mdal/frmts/mdal_netcdf.hpp b/external/mdal/frmts/mdal_netcdf.hpp
index 9eadb46f76e..4daf0f241ec 100644
--- a/external/mdal/frmts/mdal_netcdf.hpp
+++ b/external/mdal/frmts/mdal_netcdf.hpp
@@ -26,6 +26,7 @@ class NetCDFFile
     bool hasArr( const std::string &name ) const;
     std::vector<std::string> readArrNames() const;
 
+    bool hasAttrInt( const std::string &name, const std::string &attr_name ) const;
     int getAttrInt( const std::string &name, const std::string &attr_name ) const;
     double getAttrDouble( int varid, const std::string &attr_name ) const;
     /**
@@ -39,6 +40,7 @@ class NetCDFFile
     double getFillValue( int varid ) const;
     int getVarId( const std::string &name );
     void getDimension( const std::string &name, size_t *val, int *ncid_val ) const;
+    void getDimensions( const std::string &variableName, std::vector<size_t> &dimensionsId, std::vector<int> &dimensionIds );
     bool hasDimension( const std::string &name ) const;
 
   private:
diff --git a/external/mdal/frmts/mdal_selafin.cpp b/external/mdal/frmts/mdal_selafin.cpp
index fbf847f5c53..76c837d9bb6 100644
--- a/external/mdal/frmts/mdal_selafin.cpp
+++ b/external/mdal/frmts/mdal_selafin.cpp
@@ -145,35 +145,18 @@ std::string MDAL::SerafinStreamReader::read_string_without_length( size_t len )
 double MDAL::SerafinStreamReader::read_double( )
 {
   double ret;
-  unsigned char buffer[8];
 
   if ( mStreamInFloatPrecision )
   {
     float ret_f;
-    if ( mIn.read( reinterpret_cast< char * >( &buffer ), 4 ) )
-      if ( !mIn )
-        throw MDAL_Status::Err_UnknownFormat;
-    if ( mIsNativeLittleEndian )
-    {
-      std::reverse( reinterpret_cast< char * >( &buffer ), reinterpret_cast< char * >( &buffer ) + 4 );
-    }
-    memcpy( reinterpret_cast< char * >( &ret_f ),
-            reinterpret_cast< char * >( &buffer ),
-            4 );
+    if ( !readValue( ret_f, mIn, mIsNativeLittleEndian ) )
+      throw MDAL_Status::Err_UnknownFormat;
     ret = static_cast<double>( ret_f );
   }
   else
   {
-    if ( mIn.read( reinterpret_cast< char * >( &buffer ), 8 ) )
-      if ( !mIn )
-        throw MDAL_Status::Err_UnknownFormat;
-    if ( mIsNativeLittleEndian )
-    {
-      std::reverse( reinterpret_cast< char * >( &buffer ), reinterpret_cast< char * >( &buffer ) + 8 );
-    }
-    memcpy( reinterpret_cast< char * >( &ret ),
-            reinterpret_cast< char * >( &buffer ),
-            8 );
+    if ( !readValue( ret, mIn, mIsNativeLittleEndian ) )
+      throw MDAL_Status::Err_UnknownFormat;
   }
   return ret;
 }
@@ -453,6 +436,20 @@ void MDAL::DriverSelafin::addData( const std::vector<std::string> &var_names, co
       var_name =  MDAL::replace( var_name, "along y", "" );
     }
 
+    if ( MDAL::contains( var_name, "vitesse u" ) || MDAL::contains( var_name, "suivant x" ) )
+    {
+      is_vector = true;
+      var_name = MDAL::replace( var_name, "vitesse u", "vitesse" );
+      var_name = MDAL::replace( var_name, "suivant x", "" );
+    }
+    else if ( MDAL::contains( var_name, "vitesse v" ) || MDAL::contains( var_name, "suivant y" ) )
+    {
+      is_vector = true;
+      is_x =  false;
+      var_name =  MDAL::replace( var_name, "vitesse v", "vitesse" );
+      var_name =  MDAL::replace( var_name, "suivant y", "" );
+    }
+
     std::shared_ptr<DatasetGroup> group = mMesh->group( var_name );
     if ( !group )
     {
@@ -601,6 +598,5 @@ std::unique_ptr<MDAL::Mesh> MDAL::DriverSelafin::load( const std::string &meshFi
     if ( status ) *status = ( error );
     mMesh.reset();
   }
-
   return std::unique_ptr<Mesh>( mMesh.release() );
 }
diff --git a/external/mdal/frmts/mdal_sww.cpp b/external/mdal/frmts/mdal_sww.cpp
index d134293e057..7f101062f3c 100644
--- a/external/mdal/frmts/mdal_sww.cpp
+++ b/external/mdal/frmts/mdal_sww.cpp
@@ -190,6 +190,11 @@ void MDAL::DriverSWW::readDatasetGroups(
 
   for ( const std::string &name : names )
   {
+    // currently we do not support variables like elevation_c, friction_c, stage_c, xmomentum_c, ymomentum_c
+    // which contain values per volume instead of per vertex
+    if ( MDAL::endsWith( name, "_c" ) )
+      continue;
+
     // skip already parsed variables
     if ( parsedVariableNames.find( name ) == parsedVariableNames.cend() )
     {
diff --git a/external/mdal/frmts/mdal_ugrid.cpp b/external/mdal/frmts/mdal_ugrid.cpp
index 1efe1d0e691..1e3f07295c1 100644
--- a/external/mdal/frmts/mdal_ugrid.cpp
+++ b/external/mdal/frmts/mdal_ugrid.cpp
@@ -45,8 +45,21 @@ std::string MDAL::DriverUgrid::findMeshName( int dimension, bool optional ) cons
 
 std::string MDAL::DriverUgrid::nodeZVariableName() const
 {
-  // looks like mesh attributes does not have node_z array name
-  // reference
+  const std::vector<std::string> variables = mNcFile.readArrNames();
+  for ( const std::string &varName : variables )
+  {
+    const std::string stdName = mNcFile.getAttrStr( varName, "standard_name" );
+    const std::string meshName = mNcFile.getAttrStr( varName, "mesh" );
+    const std::string location = mNcFile.getAttrStr( varName, "location" );
+
+    if ( stdName == "altitude" && meshName == mMesh2dName && location == "node" )
+    {
+      return varName;
+    }
+  }
+
+  // not found, the file in non UGRID standard conforming,
+  // but lets try the common name
   return mMesh2dName + "_node_z";
 }
 
@@ -61,33 +74,91 @@ MDAL::CFDimensions MDAL::DriverUgrid::populateDimensions( )
 
   // 2D Mesh
 
-  // node_dimension is usually something like nMesh2D_node
-  // number of nodes/vertices in the mesh
-  const std::string mesh2dNode = mNcFile.getAttrStr( mMesh2dName, "node_dimension" );
-  mNcFile.getDimension( mesh2dNode, &count, &ncid );
-  dims.setDimension( CFDimensions::Vertex2D, count, ncid );
+  //node dimension location is retrieved from the node variable
 
-  // face_dimension is usually something like nMesh2D_face
-  // number of faces in the mesh
-  const std::string mesh2dFace = mNcFile.getAttrStr( mMesh2dName, "face_dimension" );
-  mNcFile.getDimension( mesh2dFace, &count, &ncid );
-  dims.setDimension( CFDimensions::Face2D, count, ncid );
+  std::vector<std::string> nodeVariablesName = MDAL::split( mNcFile.getAttrStr( mMesh2dName, "node_coordinates" ), ' ' );
+  if ( nodeVariablesName.size() < 2 )
+    throw MDAL_Status::Err_UnknownFormat;
 
-  // edge_dimension is usually something like nMesh2D_edge
-  // number of edges in the mesh
+  std::vector<size_t> nodeDimension;
+  std::vector<int> nodeDimensionId;
+  mNcFile.getDimensions( nodeVariablesName.at( 0 ), nodeDimension, nodeDimensionId );
+  if ( nodeDimension.size() != 1 )
+    throw MDAL_Status::Err_UnknownFormat;
+
+  dims.setDimension( CFDimensions::Vertex2D, nodeDimension.at( 0 ), nodeDimensionId.at( 0 ) );
+
+
+  //face dimension location is retrieved from the face_node_connectivity variable
+  //if face_dimension is defined as attribute, the dimension at this location help to desambiguate vertex per faces and number of faces
+
+  std::string faceConnectivityVariablesName = mNcFile.getAttrStr( mMesh2dName, "face_node_connectivity" );
+  std::string faceDimensionLocation = mNcFile.getAttrStr( mMesh2dName, "face_dimension" );
+  if ( faceConnectivityVariablesName == "" )
+    throw MDAL_Status::Err_UnknownFormat;
+
+  size_t facesCount;
+  size_t maxVerticesPerFace;
+
+  std::vector<size_t> faceDimension;
+  std::vector<int> faceDimensionId;
+  int facesIndexDimensionId;
+  int maxVerticesPerFaceDimensionId;
+  mNcFile.getDimensions( faceConnectivityVariablesName, faceDimension, faceDimensionId );
+  if ( faceDimension.size() != 2 )
+    throw MDAL_Status::Err_UnknownFormat;
+
+  if ( faceDimensionLocation != "" )
+  {
+    mNcFile.getDimension( faceDimensionLocation, &facesCount, &ncid );
+    if ( facesCount == faceDimension.at( 0 ) )
+    {
+      facesIndexDimensionId = faceDimensionId.at( 0 );
+      maxVerticesPerFaceDimensionId = faceDimensionId.at( 1 );
+      maxVerticesPerFace = faceDimension.at( 1 );
+    }
+    else
+    {
+      facesIndexDimensionId = faceDimensionId.at( 1 );
+      maxVerticesPerFaceDimensionId = faceDimensionId.at( 0 );
+      maxVerticesPerFace = faceDimension.at( 0 );
+    }
+  }
+  else
+  {
+    facesIndexDimensionId = faceDimensionId.at( 0 );
+    facesCount = faceDimension.at( 0 );
+    maxVerticesPerFaceDimensionId = faceDimensionId.at( 1 );
+    maxVerticesPerFace = faceDimension.at( 1 );
+  }
+
+  dims.setDimension( CFDimensions::Face2D, facesCount, facesIndexDimensionId );
+  dims.setDimension( CFDimensions::MaxVerticesInFace, maxVerticesPerFace, maxVerticesPerFaceDimensionId );
+
+
+
+  // number of edges in the mesh, not required for UGRID format
   const std::string mesh2dEdge = mNcFile.getAttrStr( mMesh2dName, "edge_dimension" );
-  mNcFile.getDimension( mesh2dEdge, &count, &ncid );
-  dims.setDimension( CFDimensions::Face2DEdge, count, ncid );
+  if ( mNcFile.hasDimension( mesh2dEdge ) )
+  {
+    mNcFile.getDimension( mesh2dEdge, &count, &ncid );
+    dims.setDimension( CFDimensions::Face2DEdge, count, ncid );
+  }
+  else
+  {
+    dims.setDimension( CFDimensions::Face2DEdge, 0, -1 );
+  }
 
-  // max_face_nodes_dimension is usually something like max_nMesh2D_face_nodes
-  // maximum number of vertices in faces
-  const std::string mesh2dMaxNodesInFace = mNcFile.getAttrStr( mMesh2dName, "max_face_nodes_dimension" );
-  mNcFile.getDimension( mesh2dMaxNodesInFace, &count, &ncid );
-  dims.setDimension( CFDimensions::MaxVerticesInFace, count, ncid );
-
-  // Time
-  mNcFile.getDimension( "time", &count, &ncid );
-  dims.setDimension( CFDimensions::Time, count, ncid );
+  // Time not required for UGRID format
+  if ( mNcFile.hasDimension( "time" ) )
+  {
+    mNcFile.getDimension( "time", &count, &ncid );
+    dims.setDimension( CFDimensions::Time, count, ncid );
+  }
+  else
+  {
+    dims.setDimension( CFDimensions::Time, 0, -1 );
+  }
 
   return dims;
 }
@@ -138,7 +209,9 @@ void MDAL::DriverUgrid::populateFaces( MDAL::Faces &faces )
   const std::string mesh2dFaceNodeConnectivity = mNcFile.getAttrStr( mMesh2dName, "face_node_connectivity" );
 
   size_t verticesInFace = mDimensions.size( CFDimensions::MaxVerticesInFace );
-  int fill_val = mNcFile.getAttrInt( mesh2dFaceNodeConnectivity, "_FillValue" );
+  int fill_val = -1;
+  if ( mNcFile.hasAttrInt( mesh2dFaceNodeConnectivity, "_FillValue" ) )
+    fill_val = mNcFile.getAttrInt( mesh2dFaceNodeConnectivity, "_FillValue" );
   int start_index = mNcFile.getAttrInt( mesh2dFaceNodeConnectivity, "start_index" );
   std::vector<int> face_nodes_conn = mNcFile.readIntArr( mesh2dFaceNodeConnectivity, faceCount * verticesInFace );
 
@@ -179,11 +252,16 @@ std::string MDAL::DriverUgrid::getCoordinateSystemVariableName()
   std::string coordinate_system_variable;
 
   // first try to get the coordinate system variable from grid definition
-  if ( mNcFile.hasArr( nodeZVariableName() ) )
+  std::vector<std::string> nodeVariablesName = MDAL::split( mNcFile.getAttrStr( mMesh2dName, "node_coordinates" ), ' ' );
+  if ( nodeVariablesName.size() > 1 )
   {
-    coordinate_system_variable = mNcFile.getAttrStr( nodeZVariableName(), "grid_mapping" );
+    if ( mNcFile.hasArr( nodeVariablesName[0] ) )
+    {
+      coordinate_system_variable = mNcFile.getAttrStr( nodeVariablesName[0], "grid_mapping" );
+    }
   }
 
+
   // if automatic discovery fails, try to check some hardcoded common variables that store projection
   if ( coordinate_system_variable.empty() )
   {
@@ -281,16 +359,18 @@ void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::strin
   }
   else
   {
-    if ( MDAL::contains( long_name, ", x-component" ) )
+    if ( MDAL::contains( long_name, ", x-component" ) || MDAL::contains( long_name, "u component of " ) )
     {
       *is_vector = true;
       name = MDAL::replace( long_name, ", x-component", "" );
+      name = MDAL::replace( name, "u component of ", "" );
     }
-    else if ( MDAL::contains( long_name, ", y-component" ) )
+    else if ( MDAL::contains( long_name, ", y-component" ) || MDAL::contains( long_name, "v component of " ) )
     {
       *is_vector = true;
       *is_x = false;
       name = MDAL::replace( long_name, ", y-component", "" );
+      name = MDAL::replace( name, "v component of ", "" );
     }
     else
     {
diff --git a/external/mdal/frmts/mdal_xmdf.cpp b/external/mdal/frmts/mdal_xmdf.cpp
index 97ea7000c7a..947e1cca094 100644
--- a/external/mdal/frmts/mdal_xmdf.cpp
+++ b/external/mdal/frmts/mdal_xmdf.cpp
@@ -247,7 +247,8 @@ std::shared_ptr<MDAL::DatasetGroup> MDAL::DriverXmdf::readXmdfGroupAsDatasetGrou
   bool isVector = dimValues.size() == 3;
 
   std::vector<double> times = dsTimes.readArrayDouble();
-
+  std::string timeUnit = rootGroup.attribute( "TimeUnits" ).readString();
+  convertTimeDataToHours( times, timeUnit );
   // all fine, set group and return
   group = std::make_shared<MDAL::DatasetGroup>(
             name(),
@@ -257,6 +258,7 @@ std::shared_ptr<MDAL::DatasetGroup> MDAL::DriverXmdf::readXmdfGroupAsDatasetGrou
           );
   group->setIsScalar( !isVector );
   group->setIsOnVertices( true );
+  group->setMetadata( "TIMEUNITS", timeUnit );
 
   // lazy loading of min and max of the dataset group
   std::vector<float> mins = dsMins.readArray();
@@ -279,3 +281,16 @@ std::shared_ptr<MDAL::DatasetGroup> MDAL::DriverXmdf::readXmdfGroupAsDatasetGrou
 
   return group;
 }
+
+void MDAL::DriverXmdf::convertTimeDataToHours( std::vector<double> &times, const std::string &originalTimeDataUnit )
+{
+  if ( originalTimeDataUnit != "Hours" )
+  {
+    for ( size_t i = 0; i < times.size(); i++ )
+    {
+      if ( originalTimeDataUnit == "Seconds" ) { times[i] /= 3600.0; }
+      else if ( originalTimeDataUnit == "Minutes" ) { times[i] /= 60.0; }
+      else if ( originalTimeDataUnit == "Days" ) { times[i] *= 24; }
+    }
+  }
+}
diff --git a/external/mdal/frmts/mdal_xmdf.hpp b/external/mdal/frmts/mdal_xmdf.hpp
index a1ea90c51e7..08a35e08754 100644
--- a/external/mdal/frmts/mdal_xmdf.hpp
+++ b/external/mdal/frmts/mdal_xmdf.hpp
@@ -100,6 +100,8 @@ namespace MDAL
         const std::string &nameSuffix,
         size_t vertexCount,
         size_t faceCount );
+
+      void convertTimeDataToHours( std::vector<double> &times, const std::string &originalTimeDataUnit );
   };
 
 } // namespace MDAL
diff --git a/external/mdal/mdal.cpp b/external/mdal/mdal.cpp
index 35410d1d856..ca4e80dac30 100644
--- a/external/mdal/mdal.cpp
+++ b/external/mdal/mdal.cpp
@@ -22,7 +22,7 @@ static MDAL_Status sLastStatus;
 
 const char *MDAL_Version()
 {
-  return "0.3.3";
+  return "0.4.0";
 }
 
 MDAL_Status MDAL_LastStatus()
@@ -92,6 +92,18 @@ bool MDAL_DR_writeDatasetsCapability( DriverH driver )
   return d->hasCapability( MDAL::Capability::WriteDatasets );
 }
 
+bool MDAL_DR_SaveMeshCapability( DriverH driver )
+{
+  if ( !driver )
+  {
+    sLastStatus = MDAL_Status::Err_MissingDriver;
+    return false;
+  }
+
+  MDAL::Driver *d = static_cast< MDAL::Driver * >( driver );
+  return d->hasCapability( MDAL::Capability::SaveMesh );
+}
+
 const char *MDAL_DR_longName( DriverH driver )
 {
   if ( !driver )
@@ -143,6 +155,39 @@ MeshH MDAL_LoadMesh( const char *meshFile )
   return static_cast< MeshH >( MDAL::DriverManager::instance().load( filename, &sLastStatus ).release() );
 }
 
+void MDAL_SaveMesh( MeshH mesh, const char *meshFile, const char *driver )
+{
+  if ( !meshFile )
+  {
+    sLastStatus = MDAL_Status::Err_FileNotFound;
+    return;
+  }
+
+  std::string driverName( driver );
+  auto d = MDAL::DriverManager::instance().driver( driver );
+
+  if ( !d )
+  {
+    sLastStatus = MDAL_Status::Err_MissingDriver;
+    return;
+  }
+
+  if ( !d->hasCapability( MDAL::Capability::SaveMesh ) )
+  {
+    sLastStatus = MDAL_Status::Err_MissingDriverCapability;
+    return;
+  }
+
+  if ( d->faceVerticesMaximumCount() < MDAL_M_faceVerticesMaximumCount( mesh ) )
+  {
+    sLastStatus = MDAL_Status::Err_IncompatibleMesh;
+    return;
+  }
+
+  std::string filename( meshFile );
+  MDAL::DriverManager::instance().save( static_cast< MDAL::Mesh * >( mesh ), filename, driverName, &sLastStatus );
+}
+
 
 void MDAL_CloseMesh( MeshH mesh )
 {
@@ -686,6 +731,16 @@ void MDAL_G_closeEditMode( DatasetGroupH group )
   }
 }
 
+const char *MDAL_G_referenceTime( DatasetGroupH group )
+{
+  if ( !group )
+  {
+    sLastStatus = MDAL_Status::Err_IncompatibleDataset;
+    return EMPTY_STR;
+  }
+  MDAL::DatasetGroup *g = static_cast< MDAL::DatasetGroup * >( group );
+  return _return_str( g->referenceTime() );
+}
 
 void MDAL_G_setMetadata( DatasetGroupH group, const char *key, const char *val )
 {
@@ -867,3 +922,4 @@ void MDAL_D_minimumMaximum( DatasetH dataset, double *min, double *max )
   *min = stats.minimum;
   *max = stats.maximum;
 }
+
diff --git a/external/mdal/mdal_data_model.cpp b/external/mdal/mdal_data_model.cpp
index 07b779a6753..150f114f814 100644
--- a/external/mdal/mdal_data_model.cpp
+++ b/external/mdal/mdal_data_model.cpp
@@ -155,6 +155,16 @@ void MDAL::DatasetGroup::setStatistics( const Statistics &statistics )
   mStatistics = statistics;
 }
 
+std::string MDAL::DatasetGroup::referenceTime() const
+{
+  return mReferenceTime;
+}
+
+void MDAL::DatasetGroup::setReferenceTime( const std::string &referenceTime )
+{
+  mReferenceTime = referenceTime;
+}
+
 MDAL::Mesh *MDAL::DatasetGroup::mesh() const
 {
   return mParent;
diff --git a/external/mdal/mdal_data_model.hpp b/external/mdal/mdal_data_model.hpp
index a9616a05445..72d83fb6c6d 100644
--- a/external/mdal/mdal_data_model.hpp
+++ b/external/mdal/mdal_data_model.hpp
@@ -109,6 +109,9 @@ namespace MDAL
       Statistics statistics() const;
       void setStatistics( const Statistics &statistics );
 
+      std::string referenceTime() const;
+      void setReferenceTime( const std::string &referenceTime );
+
       Mesh *mesh() const;
 
       bool isInEditMode() const;
@@ -124,6 +127,7 @@ namespace MDAL
       bool mIsOnVertices = true;
       std::string mUri; // file/uri from where it came
       Statistics mStatistics;
+      std::string mReferenceTime;
   };
 
   typedef std::vector<std::shared_ptr<DatasetGroup>> DatasetGroups;
diff --git a/external/mdal/mdal_driver_manager.cpp b/external/mdal/mdal_driver_manager.cpp
index e34b40d44fb..10f0bb99bab 100644
--- a/external/mdal/mdal_driver_manager.cpp
+++ b/external/mdal/mdal_driver_manager.cpp
@@ -9,6 +9,7 @@
 #include "frmts/mdal_ascii_dat.hpp"
 #include "frmts/mdal_binary_dat.hpp"
 #include "frmts/mdal_selafin.hpp"
+#include "frmts/mdal_esri_tin.hpp"
 #include "mdal_utils.hpp"
 
 #ifdef HAVE_HDF5
@@ -92,6 +93,15 @@ void MDAL::DriverManager::loadDatasets( Mesh *mesh, const std::string &datasetFi
     *status = MDAL_Status::Err_UnknownFormat;
 }
 
+void MDAL::DriverManager::save( MDAL::Mesh *mesh, const std::string &uri, const std::string &driverName, MDAL_Status *status ) const
+{
+  auto selectedDriver = driver( driverName );
+
+  std::unique_ptr<Driver> drv( selectedDriver->create() );
+
+  drv->save( uri, mesh, status );
+}
+
 size_t MDAL::DriverManager::driversCount() const
 {
   return mDrivers.size();
@@ -124,6 +134,7 @@ MDAL::DriverManager::DriverManager()
   // MESH DRIVERS
   mDrivers.push_back( std::make_shared<MDAL::Driver2dm>() );
   mDrivers.push_back( std::make_shared<MDAL::DriverSelafin>() );
+  mDrivers.push_back( std::make_shared<MDAL::DriverEsriTin>() );
 
 #ifdef HAVE_HDF5
   mDrivers.push_back( std::make_shared<MDAL::DriverFlo2D>() );
@@ -155,3 +166,4 @@ MDAL::DriverManager::DriverManager()
   mDrivers.push_back( std::make_shared<MDAL::DriverXdmf>() );
 #endif
 }
+
diff --git a/external/mdal/mdal_driver_manager.hpp b/external/mdal/mdal_driver_manager.hpp
index 7c340274aaf..fbe5ec24f0d 100644
--- a/external/mdal/mdal_driver_manager.hpp
+++ b/external/mdal/mdal_driver_manager.hpp
@@ -32,6 +32,8 @@ namespace MDAL
       std::unique_ptr< Mesh > load( const std::string &meshFile, MDAL_Status *status ) const;
       void loadDatasets( Mesh *mesh, const std::string &datasetFile, MDAL_Status *status ) const;
 
+      void save( Mesh *mesh, const std::string &uri, const std::string &driver, MDAL_Status *status ) const;
+
       size_t driversCount() const;
       std::shared_ptr<MDAL::Driver> driver( const std::string &driverName ) const;
       std::shared_ptr<MDAL::Driver> driver( size_t index ) const;
diff --git a/external/mdal/mdal_utils.cpp b/external/mdal/mdal_utils.cpp
index e1f0938c341..0a9a977e823 100644
--- a/external/mdal/mdal_utils.cpp
+++ b/external/mdal/mdal_utils.cpp
@@ -586,3 +586,43 @@ bool MDAL::isNativeLittleEndian()
   int n = 1;
   return ( *( char * )&n == 1 );
 }
+
+std::string MDAL::coordinateToString( double coordinate, int precision )
+{
+
+  std::ostringstream oss;
+  oss.setf( std::ios::fixed );
+  if ( fabs( coordinate ) > 180 )
+    oss.precision( precision ); //seems to not be a geographic coordinate, so 'precision' digits after the digital point
+  else
+    oss.precision( 6 + precision ); //could be a geographic coordinate, so 'precision'+6 digits after the digital point
+
+  oss << coordinate;
+
+  std::string returnString = oss.str();
+  returnString.back();
+
+  //remove unnecessary '0' or '.'
+  if ( returnString.size() > 0 )
+  {
+    while ( '0' == returnString.back() )
+    {
+      returnString.pop_back();
+    }
+
+    if ( '.' == returnString.back() )
+      returnString.pop_back();
+  }
+
+  return returnString;
+}
+
+std::string MDAL::doubleToString( double value, int precision )
+{
+  std::ostringstream oss;
+  oss.precision( precision );
+  oss << value;
+  return oss.str();
+}
+
+
diff --git a/external/mdal/mdal_utils.hpp b/external/mdal/mdal_utils.hpp
index 4ef88474333..b637b399ecd 100644
--- a/external/mdal/mdal_utils.hpp
+++ b/external/mdal/mdal_utils.hpp
@@ -14,6 +14,9 @@
 #include <vector>
 #include <stddef.h>
 #include <limits>
+#include <sstream>
+#include <fstream>
+#include <algorithm>
 
 #include "mdal_data_model.hpp"
 #include "mdal_memory_data_model.hpp"
@@ -66,6 +69,15 @@ namespace MDAL
   double toDouble( const std::string &str );
   bool toBool( const std::string &str );
 
+  //! Returns the string with a adapted format to coordinate
+  //! precision is the number of digits after the digital point if fabs(value)>180 (seems to not be a geographic coordinate)
+  //! precision+6 is the number of digits after the digital point if fabs(value)<=180 (could be a geographic coordinate)
+  std::string coordinateToString( double coordinate, int precision = 2 );
+
+  //! Returns a string with scientific format
+  //! precision is the number of signifiant digits
+  std::string doubleToString( double value, int precision = 6 );
+
   /**
    * Splits by deliminer and skips empty parts.
    * Faster than version with std::string
@@ -112,10 +124,25 @@ namespace MDAL
   // mesh & datasets
   //! Adds bed elevatiom dataset group to mesh
   void addBedElevationDatasetGroup( MDAL::Mesh *mesh, const Vertices &vertices );
-  //! Adds face scalar dataset group to mesh
+  //! Adds altitude dataset group to mesh
   void addFaceScalarDatasetGroup( MDAL::Mesh *mesh, const std::vector<double> &values, const std::string &name );
   //! Loop through all faces and activate those which has all 4 values on vertices valid
   void activateFaces( MDAL::MemoryMesh *mesh, std::shared_ptr<MemoryDataset> dataset );
 
+  //! function used to read all of type of value. Option to change the endianness is provided
+  template<typename T>
+  bool readValue( T &value, std::ifstream &in, bool changeEndianness = false )
+  {
+    char *const p = reinterpret_cast<char *>( &value );
+
+    if ( !in.read( p, sizeof( T ) ) )
+      return false;
+
+    if ( changeEndianness )
+      std::reverse( p, p + sizeof( T ) );
+
+    return true;
+  }
+
 } // namespace MDAL
 #endif //MDAL_UTILS_HPP
diff --git a/src/providers/mdal/CMakeLists.txt b/src/providers/mdal/CMakeLists.txt
index c59221cc57b..a0a1ae09fe4 100644
--- a/src/providers/mdal/CMakeLists.txt
+++ b/src/providers/mdal/CMakeLists.txt
@@ -44,6 +44,7 @@ IF (WITH_INTERNAL_MDAL)
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.cpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.cpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_selafin.cpp
+    ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_esri_tin.cpp
   )
 
   SET(MDAL_LIB_HDRS
@@ -57,6 +58,7 @@ IF (WITH_INTERNAL_MDAL)
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.hpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.hpp
     ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_selafin.hpp
+    ${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_esri_tin.hpp
   )
 
   IF(HDF5_FOUND)
diff --git a/tests/src/core/testqgsmeshlayer.cpp b/tests/src/core/testqgsmeshlayer.cpp
index 54c6e76409d..cee16b153c9 100644
--- a/tests/src/core/testqgsmeshlayer.cpp
+++ b/tests/src/core/testqgsmeshlayer.cpp
@@ -227,7 +227,7 @@ void TestQgsMeshLayer::test_read_vertex_scalar_dataset()
       QgsMeshDatasetIndex ds( 1, i );
 
       const QgsMeshDatasetGroupMetadata meta = dp->datasetGroupMetadata( ds );
-      if ( meta.extraOptions().count() == 2 )
+      if ( dp->name() == QStringLiteral( "mesh_memory" ) )
       {
         QCOMPARE( meta.extraOptions()["description"], QString( "Vertex Scalar Dataset" ) );
         QCOMPARE( meta.extraOptions()["meta2"], QString( "best dataset" ) );
@@ -270,7 +270,7 @@ void TestQgsMeshLayer::test_read_vertex_vector_dataset()
       QgsMeshDatasetIndex ds( 2, i );
 
       const QgsMeshDatasetGroupMetadata meta = dp->datasetGroupMetadata( ds );
-      if ( meta.extraOptions().count() == 2 )
+      if ( dp->name() == QStringLiteral( "mesh_memory" ) )
         QCOMPARE( meta.extraOptions()["description"], QString( "Vertex Vector Dataset" ) );
       QCOMPARE( meta.name(), QString( "VertexVectorDataset" ) );
       QVERIFY( !meta.isScalar() );
@@ -310,7 +310,7 @@ void TestQgsMeshLayer::test_read_face_scalar_dataset()
       QgsMeshDatasetIndex ds( 3, i );
 
       const QgsMeshDatasetGroupMetadata meta = dp->datasetGroupMetadata( ds );
-      if ( meta.extraOptions().count() == 2 )
+      if ( dp->name() == QStringLiteral( "mesh_memory" ) )
         QCOMPARE( meta.extraOptions()["description"], QString( "Face Scalar Dataset" ) );
       QCOMPARE( meta.name(), QString( "FaceScalarDataset" ) );
       QVERIFY( meta.isScalar() );
@@ -348,7 +348,7 @@ void TestQgsMeshLayer::test_read_face_vector_dataset()
       QgsMeshDatasetIndex ds( 4, i );
 
       const QgsMeshDatasetGroupMetadata meta = dp->datasetGroupMetadata( ds );
-      if ( meta.extraOptions().count() == 2 )
+      if ( dp->name() == QStringLiteral( "mesh_memory" ) )
         QCOMPARE( meta.extraOptions()["description"], QString( "Face Vector Dataset" ) );
       QCOMPARE( meta.name(), QString( "FaceVectorDataset" ) );
       QVERIFY( !meta.isScalar() );