diff --git a/Code/IO/otbDEMReaderFilter.h b/Code/IO/otbDEMReaderFilter.h new file mode 100755 index 0000000000000000000000000000000000000000..e677d37bda3304673daf12f7b714d94c62423eff --- /dev/null +++ b/Code/IO/otbDEMReaderFilter.h @@ -0,0 +1,115 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ +#ifndef __otbDEMReaderFilter_h +#define __otbDEMReaderFilter_h + +#include "itkIndent.h" +#include "itkImageSource.h" +#include "otbImage.h" +#include <iostream> +#include <stdio.h> +#include "elevation/ossimElevManager.h" +#include "base/ossimFilename.h" +#include "itkImageRegionIteratorWithIndex.h" +#include "itkIndent.h" + +namespace otb +{ +/** \class DEMReaderFilter + * + * \brief Class for Reading a DEM data + * + * This class is based on ossimElevManager. It takes in input the UL and LR geographic coordinates and the spacing. + * Handle DTED and SRTM formats. + * \ingroup Images + * + */ +template <class TDEMImage> +class ITK_EXPORT DEMReaderFilter: +public itk::ImageSource<Image<typename TDEMImage::PixelType,2, 0> > +{ +public : +/** Standard class typedefs. */ + typedef itk::Indent Indent; + typedef TDEMImage DEMImageType; + typedef typename DEMImageType::Pointer DEMImagePointerType; + typedef typename DEMImageType::PixelType PixelType; + + typedef DEMReaderFilter Self; + typedef itk::ImageSource<Image<typename DEMImageType::PixelType,2, 0> > Superclass; + typedef itk::SmartPointer<Self> Pointer; + typedef itk::SmartPointer<const Self> ConstPointer; + typedef Image<PixelType,2> OutputImageType; + + typedef typename Superclass::Pointer OutputImagePointer; + typedef typename OutputImageType::SpacingType SpacingType; + typedef typename OutputImageType::SizeType SizeType; + typedef typename OutputImageType::PointType PointType; + typedef typename OutputImageType::IndexType IndexType; + typedef typename Superclass::OutputImageRegionType OutputImageRegionType; + typedef itk::ImageRegionIteratorWithIndex< Image<PixelType,2, 0> > ImageIteratorType; + + /** Method for creation through the object factory. */ + itkNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(DEMReader,ImageSource); + + /** Set the spacing. */ + itkSetMacro(Spacing,SpacingType); + itkGetConstReferenceMacro(Spacing,SpacingType); + + /** Set the Upper Left coordinates. */ + itkSetMacro(Ul,PointType); + itkGetConstReferenceMacro(Ul,PointType); + + /** Set the Lower Right coordinates. */ + itkSetMacro(Lr,PointType); + itkGetConstReferenceMacro(Lr,PointType); + + /** Set the spacing. */ + void SetSpacing(const double* spacing); + + /** Try to open the DEM directory. */ + bool OpenDEMDirectory(char* &DEMDirectory); + + /** Compute the height above MSL(Mean Sea Level) of the point. */ + virtual double GetHeightAboveMSL(const PointType& worldPoint); + +protected: + DEMReaderFilter(); + ~DEMReaderFilter(); + + void PrintSelf(std::ostream& os, Indent indent) const; + void GenerateData(); + virtual void GenerateOutputInformation(); + + ossimElevManager* m_ElevManager; + //DEMImagePointerType m_DEMImage; + SpacingType m_Spacing; + PointType m_Ul; + PointType m_Lr; +}; + +} // namespace otb + +#ifndef OTB_MANUAL_INSTANTIATION +#include "otbDEMReaderFilter.txx" +#endif + +#endif diff --git a/Code/IO/otbDEMReaderFilter.txx b/Code/IO/otbDEMReaderFilter.txx new file mode 100755 index 0000000000000000000000000000000000000000..144a68e5c4d2cc72b5061ddf2e4ba292daad4e7c --- /dev/null +++ b/Code/IO/otbDEMReaderFilter.txx @@ -0,0 +1,171 @@ +/*========================================================================= + + Program: ORFEO Toolbox + Language: C++ + Date: $Date$ + Version: $Revision$ + + + Copyright (c) Centre National d'Etudes Spatiales. All rights reserved. + See OTBCopyright.txt for details. + + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notices for more information. + +=========================================================================*/ +#ifndef __otbDEMReaderFilter_txx +#define __otbDEMReaderFilter_txx + +#include "otbDEMReaderFilter.h" +#include "otbMacro.h" + +namespace otb +{ + + +template<class TDEMImage> +DEMReaderFilter<TDEMImage> +::DEMReaderFilter() +{ + m_ElevManager=ossimElevManager::instance(); +//m_DEMImage = DEMImageType::New(); + m_Spacing[0]=0; + m_Spacing[1]=0; + m_Ul[0]=0; + m_Ul[1]=0; + m_Lr[0]=0; + m_Lr[1]=0; +} + +template<class TDEMImage> +DEMReaderFilter<TDEMImage> +::~DEMReaderFilter() +{ + delete m_ElevManager; +} + +///Methode pour specifier un dossier contenant des DEM +template<class TDEMImage> +bool +DEMReaderFilter<TDEMImage>:: +OpenDEMDirectory(char* &DEMDirectory) +{ + ossimFilename ossimDEMDir; + ossimDEMDir=ossimFilename(DEMDirectory); + bool result= false; + if (m_ElevManager->openDirectory(ossimDEMDir)) + { + result= true; + } + return result; +} + +///Methode pour calculer l'altitude d'un point geographique +template<class TDEMImage> +double +DEMReaderFilter<TDEMImage>:: +GetHeightAboveMSL(const PointType& worldPoint) +{ + float height; + ossimGpt ossimWorldPoint; + ossimWorldPoint.lat=worldPoint[0]; + ossimWorldPoint.lon=worldPoint[1]; + height=m_ElevManager->getHeightAboveMSL(ossimWorldPoint); + return height; +} + +///Methode SetSpacing +template<class TDEMImage> + void DEMReaderFilter<TDEMImage> + ::SetSpacing( const double* spacing) +{ + SpacingType s(spacing); + this->SetSpacing(s); +} + +///Methode GenerateOutputInformation +template <class TDEMImage> +void DEMReaderFilter<TDEMImage> +::GenerateOutputInformation() +{ + DEMImageType *output; + IndexType start; + start[0]=0; + start[1]=0; + + PointType origin; + origin[0]=m_Ul[0]; //latitude de l'origine. + origin[1]=m_Ul[1]; //longitude de l'origine. + output = this->GetOutput(0); + + // Calcul de la taille de l'image: + SizeType size; + size[0]=int (abs(((m_Lr[0]-m_Ul[0])/m_Spacing[0]))+1.5); + size[1]=int (abs(((m_Lr[1]-m_Ul[1])/m_Spacing[1]))+1.5); + + // On specifie les parametres de la region + OutputImageRegionType largestPossibleRegion; + largestPossibleRegion.SetSize( size ); + largestPossibleRegion.SetIndex( start ); + output->SetLargestPossibleRegion( largestPossibleRegion ); + + output->SetSpacing(m_Spacing); + output->SetOrigin(m_Ul); +} + +///Methode GenerateData +template <class TDEMImage> +void +DEMReaderFilter<TDEMImage> +::GenerateData() +{ + DEMImagePointerType m_DEMImage = this->GetOutput(); + + // allocate the output buffer + m_DEMImage->SetBufferedRegion( m_DEMImage->GetRequestedRegion() ); + m_DEMImage->Allocate(); + + // Create an iterator that will walk the output region + ImageIteratorType outIt = ImageIteratorType(m_DEMImage,m_DEMImage->GetRequestedRegion()); + + // Walk the output image, evaluating the height at each pixel + IndexType currentindex; + PointType phyPoint; + double height; + for (outIt.GoToBegin(); !outIt.IsAtEnd(); ++outIt) + { + currentindex=outIt.GetIndex(); + m_DEMImage->TransformIndexToPhysicalPoint(currentindex, phyPoint); + ossimGpt ossimWorldPoint; + ossimWorldPoint.lat=phyPoint[0]; + ossimWorldPoint.lon=phyPoint[1]; + height=m_ElevManager->getHeightAboveMSL(ossimWorldPoint); //Calcul de l'altitude + otbMsgDebugMacro(<<" HeightAboveMSL: "<<height); + if (height>-static_cast<double>(32768)) //On teste si les fichiers MNT recouvre la zone g�ographique demand�e (-32768 = theNullHeightValue) + { + m_DEMImage->SetPixel(currentindex, static_cast<PixelType>(height) ); + } //On remplit l'image + else + { + m_DEMImage->SetPixel(currentindex, static_cast<PixelType>(0) ); + } + } +} + +template <class TDEMImage> +void +DEMReaderFilter<TDEMImage> +::PrintSelf(std::ostream& os, Indent indent) const +{ + Superclass::PrintSelf(os,indent); + + os << indent << "Spacing:"<< m_Spacing[0] << ","<< m_Spacing[1] << std::endl; + os << indent << "Lr:"<< m_Lr[0] << ","<< m_Lr[1] << std::endl; + os << indent << "Ul:"<< m_Ul[0] << ","<< m_Ul[1] << std::endl; +} + +} // namespace otb + +#endif