From bd08087c4bc89068e5b3306317d63eddf4e0cdc7 Mon Sep 17 00:00:00 2001
From: Emmanuel Christophe <emmanuel.christophe@orfeo-toolbox.org>
Date: Sat, 8 Aug 2009 14:59:24 +0800
Subject: [PATCH] ENH: Add geoid handling to the DEMHandler class
---
Code/IO/otbDEMHandler.cxx | 53 +++++++++++++++++++--------
Code/IO/otbDEMHandler.h | 23 +++++++-----
Code/IO/otbDEMToImageGenerator.h | 8 ++--
Code/IO/otbDEMToImageGenerator.txx | 11 +-----
Testing/Code/IO/CMakeLists.txt | 1 +
Testing/Code/IO/otbDEMHandlerTest.cxx | 26 +++++++++----
6 files changed, 76 insertions(+), 46 deletions(-)
diff --git a/Code/IO/otbDEMHandler.cxx b/Code/IO/otbDEMHandler.cxx
index 8fe5ace4e1..31ceed274d 100644
--- a/Code/IO/otbDEMHandler.cxx
+++ b/Code/IO/otbDEMHandler.cxx
@@ -18,26 +18,21 @@ PURPOSE. See the above copyright notices for more information.
#include "otbDEMHandler.h"
#include "otbMacro.h"
-namespace otb
-{
+#include "elevation/ossimElevManager.h"
+#include "base/ossimGeoidManager.h"
+#include "base/ossimFilename.h"
+#include "base/ossimDirectory.h"
+#include "base/ossimGeoidEgm96.h"
-
-DEMHandler
-::DEMHandler()
+namespace otb
{
- m_ElevManager=ossimElevManager::instance();
-}
-
DEMHandler
-::~DEMHandler()
+::DEMHandler():
+ m_ElevManager(ossimElevManager::instance())
{
- // not needed, m_ElevManager created with instance() method
- // delete m_ElevManager;
}
-
-
void
DEMHandler
::OpenDEMDirectory(const char* DEMDirectory)
@@ -52,14 +47,27 @@ DEMHandler
m_Mutex.Unlock();
itkExceptionMacro("Failed to open DEM Directory: "<<ossimDEMDir);
}
-
m_Mutex.Unlock();
}
+void
+DEMHandler
+::OpenGeoidFile(const char* geoidFile)
+{
+ ossimFilename geoid(geoidFile);
+ ossimGeoid* geoidPtr = new ossimGeoidEgm96(geoid);
+ if (geoidPtr->getErrorStatus() == ossimErrorCodes::OSSIM_OK)
+ {
+ m_Mutex.Lock();
+ ossimGeoidManager::instance()->addGeoid(geoidPtr);
+ m_Mutex.Unlock();
+ }
+}
+
double
DEMHandler
-::GetHeightAboveMSL(const PointType& geoPoint)
+::GetHeightAboveMSL(const PointType& geoPoint) const
{
double height;
ossimGpt ossimWorldPoint;
@@ -71,10 +79,23 @@ DEMHandler
return height;
}
+double
+DEMHandler
+::GetHeightAboveEllipsoid(const PointType& geoPoint) const
+{
+ double height;
+ ossimGpt ossimWorldPoint;
+ ossimWorldPoint.lon=geoPoint[0];
+ ossimWorldPoint.lat=geoPoint[1];
+ m_Mutex.Lock();
+ height=m_ElevManager->getHeightAboveEllipsoid(ossimWorldPoint);
+ m_Mutex.Unlock();
+ return height;
+}
void
DEMHandler
-::PrintSelf(std::ostream& os, Indent indent) const
+::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os,indent);
os << indent << "DEMHandler" << std::endl;
diff --git a/Code/IO/otbDEMHandler.h b/Code/IO/otbDEMHandler.h
index 6009bafa6d..954c8666af 100644
--- a/Code/IO/otbDEMHandler.h
+++ b/Code/IO/otbDEMHandler.h
@@ -23,13 +23,13 @@
#include "otbImage.h"
#include <iostream>
#include <stdio.h>
-#include "elevation/ossimElevManager.h"
-#include "base/ossimFilename.h"
-#include "base/ossimDirectory.h"
+
#include "itkImageRegionIteratorWithIndex.h"
#include "itkIndent.h"
#include "itkSimpleFastMutexLock.h"
+class ossimElevManager;
+
namespace otb
{
/** \class DEMHandler
@@ -47,13 +47,12 @@ class ITK_EXPORT DEMHandler: public itk::Object
{
public :
/** Standard class typedefs. */
- typedef itk::Indent Indent;
typedef DEMHandler Self;
- typedef itk::Object Superclass;
+ typedef itk::Object Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
- typedef itk::Point<double, 2> PointType;
+ typedef itk::Point<double, 2> PointType;
/** Method for creation through the object factory. */
itkNewMacro(Self);
@@ -64,14 +63,20 @@ public :
/** Try to open the DEM directory. */
virtual void OpenDEMDirectory(const char* DEMDirectory);
+ /** Open geoid file. */
+ virtual void OpenGeoidFile(const char* geoidFile);
+
/** Compute the height above MSL(Mean Sea Level) of a geographic point. */
- virtual double GetHeightAboveMSL(const PointType& geoPoint);
+ virtual double GetHeightAboveMSL(const PointType& geoPoint) const;
+
+ /** Compute the height above ellipsoid of a geographic point. */
+ virtual double GetHeightAboveEllipsoid(const PointType& geoPoint) const;
protected:
DEMHandler();
- ~DEMHandler();
+ ~DEMHandler(){};
- void PrintSelf(std::ostream& os, Indent indent) const;
+ void PrintSelf(std::ostream& os, itk::Indent indent) const;
ossimElevManager* m_ElevManager;
diff --git a/Code/IO/otbDEMToImageGenerator.h b/Code/IO/otbDEMToImageGenerator.h
index 6f7ae1db0e..9cf1965e4f 100644
--- a/Code/IO/otbDEMToImageGenerator.h
+++ b/Code/IO/otbDEMToImageGenerator.h
@@ -51,7 +51,6 @@ class ITK_EXPORT DEMToImageGenerator:
{
public :
/** Standard class typedefs. */
- typedef itk::Indent Indent;
typedef TDEMImage DEMImageType;
typedef typename DEMImageType::Pointer DEMImagePointerType;
typedef typename DEMImageType::PixelType PixelType;
@@ -60,7 +59,8 @@ public :
typedef itk::ImageSource<DEMImageType> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
- typedef Image<PixelType,2> OutputImageType;
+// typedef Image<PixelType,2> OutputImageType;
+ typedef DEMImageType OutputImageType;
typedef typename Superclass::Pointer OutputImagePointer;
typedef typename OutputImageType::SpacingType SpacingType;
@@ -100,9 +100,9 @@ public :
protected:
DEMToImageGenerator();
- ~DEMToImageGenerator();
+ ~DEMToImageGenerator(){};
- void PrintSelf(std::ostream& os, Indent indent) const;
+ void PrintSelf(std::ostream& os, itk::Indent indent) const;
void GenerateData();
virtual void GenerateOutputInformation();
diff --git a/Code/IO/otbDEMToImageGenerator.txx b/Code/IO/otbDEMToImageGenerator.txx
index 6b936449b4..550d700039 100644
--- a/Code/IO/otbDEMToImageGenerator.txx
+++ b/Code/IO/otbDEMToImageGenerator.txx
@@ -20,11 +20,11 @@
#include "otbDEMToImageGenerator.h"
#include "otbMacro.h"
+#include "base/ossimCommon.h"
namespace otb
{
-
template<class TDEMImage>
DEMToImageGenerator<TDEMImage>
::DEMToImageGenerator()
@@ -39,13 +39,6 @@ DEMToImageGenerator<TDEMImage>
m_DefaultUnknownValue = static_cast<PixelType>(-32768); // Value defined in the norm for points strm doesn't have information.
}
-template<class TDEMImage>
-DEMToImageGenerator<TDEMImage>
-::~DEMToImageGenerator()
-{
- // Nothing to be done...
-}
-
// DEM folder specification method
template<class TDEMImage>
void
@@ -124,7 +117,7 @@ DEMToImageGenerator<TDEMImage>
template <class TDEMImage>
void
DEMToImageGenerator<TDEMImage>
-::PrintSelf(std::ostream& os, Indent indent) const
+::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os,indent);
diff --git a/Testing/Code/IO/CMakeLists.txt b/Testing/Code/IO/CMakeLists.txt
index dda43aef1a..ef9dcadf75 100755
--- a/Testing/Code/IO/CMakeLists.txt
+++ b/Testing/Code/IO/CMakeLists.txt
@@ -1417,6 +1417,7 @@ ADD_TEST(ioTvDEMHandler ${IO_TESTS12}
${TEMP}/ioDEMGetHeightAboveMSL.txt
otbDEMHandlerTest
${INPUTDATA}/DEM/srtm_directory
+ ${INPUTDATA}/DEM/egm96
${TEMP}/ioDEMGetHeightAboveMSL.txt
3.6999 44.08
)
diff --git a/Testing/Code/IO/otbDEMHandlerTest.cxx b/Testing/Code/IO/otbDEMHandlerTest.cxx
index 6954d9f89b..139d93f07b 100644
--- a/Testing/Code/IO/otbDEMHandlerTest.cxx
+++ b/Testing/Code/IO/otbDEMHandlerTest.cxx
@@ -27,8 +27,10 @@ int otbDEMHandlerTest(int argc, char * argv[])
{
const unsigned int Dimension = 2;
char * srtm_directory(argv[1]);
- const char * outputfilename(argv[2]);
- double height(0.);
+ char * geoidFile(argv[2]);
+ const char * outputfilename(argv[3]);
+ double height = 0.0;
+ double height2 = 0.0;
typedef otb::Image<float,Dimension> ImageType;
typedef otb::DEMHandler DEMHandlerType;
@@ -36,23 +38,31 @@ int otbDEMHandlerTest(int argc, char * argv[])
// Instantiating object
DEMHandlerType::Pointer demHandler = DEMHandlerType::New();
demHandler->OpenDEMDirectory(srtm_directory);
+ demHandler->OpenGeoidFile(geoidFile);
typedef otb::UtmInverseProjection utmProjection;
typedef utmProjection::InputPointType InputPoint;
InputPoint geoPoint;
- geoPoint[0] = atof(argv[3]);//3.6999;
- geoPoint[1] = atof(argv[4]);//44.08;
+ geoPoint[0] = atof(argv[4]);//3.6999;
+ geoPoint[1] = atof(argv[5]);//44.08;
- height=demHandler->GetHeightAboveMSL(geoPoint);
+ height = demHandler->GetHeightAboveMSL(geoPoint);
+ height2 = demHandler->GetHeightAboveEllipsoid(geoPoint);
std::ofstream file;
file.open(outputfilename);
file << "--- HEIGHT ABOVE MSL TEST ---" << std::endl;
- file << " geoPoint: "<<geoPoint[1]<<" ; "<<geoPoint[0]<< std::endl;
- file << " -> Height: "<<height<< std::endl;
+ file << " geoPoint: " << geoPoint[1] << " ; " << geoPoint[0] << std::endl;
+ file << " -> Height above MSL: " << height << std::endl;
+ std::cout << "Height above MSL: " << height << std::endl;
+
+ file << " -> Height above Ellipsoid: " << height2 << std::endl;
+ std::cout << "Height above Ellipsoid: " << height2 << std::endl;
+
+
file.close();
- std::cout << "Height: "<<height<<std::endl;
+
return EXIT_SUCCESS;
--
GitLab