From 61c396fb710751fd4fc41d4e5193e61f9f8d972a Mon Sep 17 00:00:00 2001
From: Luc Hermitte <luc.hermitte@csgroup.eu>
Date: Wed, 31 May 2023 18:24:40 +0200
Subject: [PATCH] REFACT: Simplify geoid parameter in SARDEMProjection2

---
 app/otbSARDEMProjection2.cxx                | 55 ++++++++++++---------
 include/otbSARDEMProjectionImageFilter2.h   | 21 +++++++-
 include/otbSARDEMProjectionImageFilter2.hxx | 24 +++++----
 3 files changed, 62 insertions(+), 38 deletions(-)

diff --git a/app/otbSARDEMProjection2.cxx b/app/otbSARDEMProjection2.cxx
index a8767f6..2a88179 100644
--- a/app/otbSARDEMProjection2.cxx
+++ b/app/otbSARDEMProjection2.cxx
@@ -30,6 +30,7 @@
 #include "otbSARStreamingDEMCheckFilter.h"
 #include "otbWrapperElevationParametersHandler.h"
 #include "otbDEMHandler.h"
+#include "otbConfigurationManager.h"
 
 #include <iostream>
 #include <string>
@@ -119,7 +120,8 @@ private:
     SetParameterRole("directiontoscandeml", Role_Output);
     SetDefaultParameterInt("directiontoscandeml", 1);
 
-    AddParameter(ParameterType_Group, "elev", "Elevation parameters group");
+    AddParameter(ParameterType_Group, "elev", "Elevation management");
+    SetParameterDescription("elev", "This group of parameters allows managing elevation values.");
 
     AddParameter(ParameterType_InputFilename, "elev.geoid", "Input file path to Geoid.");
     SetParameterDescription("elev.geoid","Input file path to Geoid.");
@@ -134,32 +136,37 @@ private:
 
   void DoUpdateParameters() override
   {
-    // Elevation
-    if (IsParameterEnabled("elev.geoid") && HasValue("elev.geoid"))
-    {
-      ElevationParametersHandler::AddElevationParameters(this, "elev");
-
-      std::string geoid = GetParameterString("elev.geoid");
-      otbAppLogINFO(<< "elev.geoid : " << geoid);
-
-#if OTB_VERSION_MAJOR < 8
-      otb::DEMHandler::Instance()->ClearDEMs();
-#endif
-      otb::Wrapper::ElevationParametersHandler::SetupDEMHandlerFromElevationParameters(this, "elev");
-      // otbAppLogINFO(<< "Geoid used: " << otb::Wrapper::ElevationParametersHandler::IsGeoidUsed(this, "elev"));
-      // otbAppLogINFO(<< "Geoid active(rec): " << this->GetParameterByKey("elev.geoid")->GetActive(true));
-      // otbAppLogINFO(<< "Geoid active(no rec): " << this->GetParameterByKey("elev.geoid")->GetActive(false));
-      // otbAppLogINFO(<< "Geoid enabled: " << this->IsParameterEnabled("elev.geoid"));
-      // otbAppLogINFO(<< "Geoid has value: " << this->HasValue("elev.geoid"));
+    // ====[ Update Elevation information
+    // GEOID parameter is best taken, in priority order:
+    // 1. As -elev.geom
+    // 2. As $OTB_GEOID_FILE
+    // 3. As: ignored
+    auto get_geoid_filename = [&]() {
+      if (IsParameterEnabled("elev.geoid") && HasValue("elev.geoid"))
+        return GetParameterString("elev.geoid");
+      else
+        return ConfigurationManager::GetGeoidFile();
+    };
+    m_geoid_filename = get_geoid_filename();
+  }
 
+  void DoExecute() override
+  {
+    // Set elevation information
 #if OTB_VERSION_MAJOR >= 8
-      otbAppLogINFO(<< "configured geoid file: " << DEMHandler::GetInstance().GetGeoidFile());
+    auto& demHandler = DEMHandler::GetInstance();
+    demHandler.ClearElevationParameters(); // NO!!!!!
+#else
+    auto& demHandler = *DEMHandler::Instance();
+    demHandler.ClearDEMs();
 #endif
+    otbAppLogINFO(<< "Geoid file used : " << (m_geoid_filename.empty() ? "None" : m_geoid_filename));
+    demHandler.SetDefaultHeightAboveEllipsoid(0);
+    if (! m_geoid_filename.empty())
+    {
+      demHandler.OpenGeoidFile(m_geoid_filename);
     }
-  }
 
-  void DoExecute() override
-  {
     // Get No data value
     int noData = GetParameterInt("nodata");
     otbAppLogINFO(<<"No Data : "<<noData);
@@ -168,7 +175,7 @@ private:
     ComplexFloatImageType::Pointer SARPtr = GetParameterComplexFloatImage("insar");
 
     // DEMProjection Filter
-    DEMFilterType::Pointer filterDEMProjection = DEMFilterType::New();
+    DEMFilterType::Pointer filterDEMProjection = DEMFilterType::New(m_geoid_filename);
     m_Ref.push_back(filterDEMProjection.GetPointer());
 #if OTB_VERSION_MAJOR >= 8
     filterDEMProjection->SetSARImageMetadata(SARPtr->GetImageMetadata());
@@ -253,6 +260,8 @@ private:
 
   // Vector for filters
   std::vector<itk::ProcessObject::Pointer> m_Ref;
+
+  std::string m_geoid_filename;
 };
 
 }
diff --git a/include/otbSARDEMProjectionImageFilter2.h b/include/otbSARDEMProjectionImageFilter2.h
index a12f00b..b1e1b69 100644
--- a/include/otbSARDEMProjectionImageFilter2.h
+++ b/include/otbSARDEMProjectionImageFilter2.h
@@ -71,7 +71,15 @@ public:
     using ConstPointer = itk::SmartPointer<const Self>;
 
     // Method for creation through object factory
-    itkNewMacro(Self);
+    // itkNewMacro(Self);
+    /** Method for creation through the object factory. */
+    template <typename ... Args>
+    static Pointer New(Args&&... args)
+    {
+      Pointer smartPtr = new Self(std::forward<Args>(args)...);
+      smartPtr->UnRegister();
+      return smartPtr;
+    }
     // Run-time type information
     itkTypeMacro(SARDEMProjectionImageFilter2,ImageToImageFilter);
 
@@ -146,7 +154,16 @@ public:
 
 protected:
     // Constructor
-    SARDEMProjectionImageFilter2();
+    SARDEMProjectionImageFilter2() = delete;
+    /** Init constructor.
+     * \pre Expects `geoid_filename == DEMHandler.GetGeoidFile()`
+     *
+     * This constructor won't update the Geoid data on `DEMHandler` singleton.
+     * The singleton is expected to have been correctly configured beforehand.
+     * However with older OSSIM-based versions of OTB, it will configure the
+     * internal OSSIM based GEOID reader in the case of egm96 geoids.
+     */
+    SARDEMProjectionImageFilter2(std::string const& geoid_filename);
 
     // Destructor
     ~SARDEMProjectionImageFilter2() override = default;
diff --git a/include/otbSARDEMProjectionImageFilter2.hxx b/include/otbSARDEMProjectionImageFilter2.hxx
index 972316e..47f6fbe 100644
--- a/include/otbSARDEMProjectionImageFilter2.hxx
+++ b/include/otbSARDEMProjectionImageFilter2.hxx
@@ -53,6 +53,7 @@
 #include <boost/optional.hpp>
 #include <mutex>
 #include <thread>
+#include <cassert>
 
 namespace otb
 {
@@ -60,31 +61,26 @@ namespace otb
  * Constructor with default initialization
  */
 template <class TImageIn, class TImageOut>
-SARDEMProjectionImageFilter2<TImageIn ,TImageOut>::SARDEMProjectionImageFilter2()
+SARDEMProjectionImageFilter2<TImageIn ,TImageOut>::SARDEMProjectionImageFilter2(
+    std::string const& geoid_filename)
 {
 #if OTB_VERSION_MAJOR >= 8
     auto& demHandler = DEMHandler::GetInstance();
-    demHandler.ClearElevationParameters();
+    assert(demHandler.GetGeoidFile() == geoid_filename);
 #else
     auto& demHandler = *DEMHandler::Instance();
-#endif
 
-    if (!(ConfigurationManager::GetGeoidFile().empty()))
+    if (! geoid_filename.empty())
     {
-        // DEMHandler instance to specify the geoid file
-        demHandler.OpenGeoidFile(ConfigurationManager::GetGeoidFile());
-
-        // If Geoid by default (emg96) instanciate directly a ossimGeoidEgm96 (increase performance)
-#if OTB_VERSION_MAJOR < 8
         const std::string isEmg96 =  "egm96";
-        std::size_t found = ConfigurationManager::GetGeoidFile().find(isEmg96);
+        std::size_t found = geoid_filename.find(isEmg96);
 
         if (found!=std::string::npos)
         {
-            m_geoidEmg96 = std::make_unique<ossimGeoidEgm96>(ossimFilename(ConfigurationManager::GetGeoidFile()));
+            m_geoidEmg96 = std::make_unique<ossimGeoidEgm96>(ossimFilename(geoid_filename));
         }
-#endif
     }
+#endif
 }
 
 /**
@@ -311,6 +307,7 @@ SARDEMProjectionImageFilter2< TImageIn, TImageOut >
 #if OTB_VERSION_MAJOR >= 8
     auto& demHandler = DEMHandler::GetInstance();
     // demHandler.ClearElevationParameters(); // NO!!!!!
+    auto& demHandlerTLS = demHandler.GetHandlerForCurrentThread();
 #else
     auto& demHandler = *DEMHandler::Instance();
 #endif
@@ -336,7 +333,8 @@ SARDEMProjectionImageFilter2< TImageIn, TImageOut >
                 double h; // Elevation from earth geoid
 
 #if OTB_VERSION_MAJOR >= 8
-                h = demHandler.GetGeoidHeight(demLatLonPoint);
+                // h = demHandler.GetGeoidHeight(demLatLonPoint);
+                h = GetGeoidHeight(demHandlerTLS, demLatLonPoint);
 #else
                 if (m_geoidEmg96)
                 {
-- 
GitLab