COM: OTB-151/GeometriesProjectionFilter -> + projection (draft to test...

COM: OTB-151/GeometriesProjectionFilter -> + projection (draft to test compilation on all platforms)

Code/Projections/CMakeLists.txt

@@ -6,7 +6,7 @@ IF( NOT OTB_USE_CURL )
 ENDIF( NOT OTB_USE_CURL )
 
 ADD_LIBRARY(OTBProjections ${OTBProjections_SRCS})
-TARGET_LINK_LIBRARIES(OTBProjections OTBIO)
+TARGET_LINK_LIBRARIES(OTBProjections OTBIO OTBOGRAdapters)
 
 IF(OTB_LIBRARY_PROPERTIES)
   SET_TARGET_PROPERTIES(OTBProjections PROPERTIES ${OTB_LIBRARY_PROPERTIES})

Code/Projections/otbGeometriesProjectionFilter.cxx

+/*=========================================================================
+
+  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.
+
+=========================================================================*/
+
+/*===========================================================================*/
+/*===============================[ Includes ]================================*/
+/*===========================================================================*/
+#include "otbGeometriesProjectionFilter.h"
+#include <boost/type_traits/is_same.hpp>
+#include <boost/mpl/assert.hpp>
+#include "ogr_srs_api.h" // OCTDestroyCoordinateTransformation
+#include "ogr_spatialref.h" // OGRCoordinateTransformation
+#include "itkMacro.h"
+#include "otbGeometriesSet.h"
+#include "itkMetaDataObject.h"
+#include "otbOGRGeometryWrapper.h"
+#include "otbOGRGeometriesVisitor.h"
+
+
+/*===========================================================================*/
+/*===============================[ constants ]===============================*/
+/*===========================================================================*/
+namespace
+{ // Anonymous namespace
+const double k_averageElevation = -32768.0;
+} // Anonymous namespace
+
+/*===========================================================================*/
+/*======================[ Reprojection Transformation ]======================*/
+/*===========================================================================*/
+
+
+template <typename TGeometry>
+otb::ogr::UniqueGeometryPtr
+otb::ReprojectTransformationFunctor::ByCopy::operator()(TGeometry const* in) const
+{
+  boost::interprocess::unique_ptr<TGeometry, ogr::internal::GeometryDeleter>
+    out(in ? static_cast <TGeometry*>(in->clone()) : 0); // OGR clone doesn't use covariant return ...
+  if (out)
+    m_Reprojector.do_transform(*out);
+  ogr::UniqueGeometryPtr res(out.release());
+  return boost::move(res);
+}
+
+template <typename TGeometry>
+void
+otb::ReprojectTransformationFunctor::InPlace::operator()(TGeometry * inout) const
+{
+  if (inout)
+    m_Reprojector.do_transform(*inout);
+}
+
+#if 0
+void do_transform(OGRLinearRing & g) const
+{
+  assert(!"don't know how to reproject a ring yet");
+}
+#endif
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRPoint & g) const
+{
+  typedef InternalTransformType::InputPointType  InputPointType;
+  typedef InternalTransformType::OutputPointType OutputPointType;
+  // in our case:
+  BOOST_MPL_ASSERT((boost::is_same<InputPointType, OutputPointType>));
+  const double rawPoint[] = { g.getX(), g.getY() };
+  const InputPointType inPoint(rawPoint);
+  const OutputPointType outPoint = m_Transform->TransformPoint(inPoint);
+  g.setX(outPoint[0]);
+  g.setY(outPoint[1]);
+}
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRLineString & g) const
+{
+  OGRPoint point;
+  for (int i=0, N=g.getNumPoints(); i!=N; ++i)
+    {
+    g.getPoint(i, &point);
+    do_transform(point); // no need to go through the redispatch cycle
+    g.setPoint(i, &point);
+    }
+}
+
+#if 0
+void otb::ReprojectTransformationFunctor::do_transform(OGRCurve & g) const
+{
+  assert(!"OGRCurve is just an interface, they can't be reprojected yet");
+}
+#endif
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRPolygon & g) const
+{
+  // Note: OGRPolygon interface states that rings obtained through
+  // getInteriorRing() and getExteriorRing() shall not be modified.
+  // In this context, it seems to make sense as it gives a direct access to the
+  // ring (yes this is a de-encapsulation), and the polygon is a fresh one has
+  // just been cloned in the () operator.
+  // The fact that OGRPolygon::transform() directly works on each sub-ring seems
+  // to confirm that hyppothesis.
+  for (int i=0, N=g.getNumInteriorRings(); i!=N; ++i)
+    {
+    OGRLinearRing * ring = g.getInteriorRing(i);
+    assert(ring && "Rings in polygons are not expected to be nil");
+    if (ring)
+      {
+      do_transform(*ring); // no need to go through the redispatch cycle
+      // how to commit the change ? Are ring modifications automatically commited ?
+      }
+    }
+  OGRLinearRing * ring = g.getExteriorRing();
+  assert(ring && "Rings in polygons are not expected to be nil");
+  if (ring)
+    {
+    do_transform(*ring); // no need to go through the redispatch cycle
+    // how to commit the change ? Are ring modifications automatically commited ?
+    }
+}
+
+#if 0
+void otb::ReprojectTransformationFunctor::do_transform(OGRSurface & g) const
+{
+  assert(!"OGRSurface is just an interface, they can't be reprojected yet");
+}
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRMultiLineString & g) const
+{
+}
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRMultiPoint & g) const
+{
+}
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRMultiPolygon & g) const
+{
+}
+#endif
+
+void otb::ReprojectTransformationFunctor::do_transform(OGRGeometryCollection & col) const
+{
+  for (int i=0, N=col.getNumGeometries(); i!=N; ++i)
+    {
+    OGRGeometry * g = col.getGeometryRef(i);
+    assert(g && "geometries to transform shall not be null");
+    apply_inplace(g);
+    }
+}
+
+otb::ogr::UniqueGeometryPtr
+otb::ReprojectTransformationFunctor::operator()(OGRGeometry const* in) const
+// otb::ReprojectTransformationFunctor::apply(OGRGeometry const* in) const
+{
+  otb::ogr::UniqueGeometryPtr res
+    = ogr::apply<otb::ogr::UniqueGeometryPtr>(in, ByCopy(*this));
+  return boost::move(res);
+}
+
+void otb::ReprojectTransformationFunctor::apply_inplace(OGRGeometry * inout) const
+{
+  if (inout)
+    ogr::apply<void>(inout, InPlace(*this));
+}
+
+/*===========================================================================*/
+/*======================[ GeometriesProjectionFilter ]=======================*/
+/*===========================================================================*/
+
+otb::GeometriesProjectionFilter::GeometriesProjectionFilter()
+: InputImageReference(*this)
+, OutputImageReference(*this)
+{
+}
+
+/*virtual*/
+otb::GeometriesProjectionFilter::~GeometriesProjectionFilter()
+{
+}
+
+/*virtual*/
+void otb::GeometriesProjectionFilter::GenerateOutputInformation(void)
+{
+  Superclass::GenerateOutputInformation();
+
+  OutputGeometriesType::Pointer  output = this->GetOutput();
+  itk::MetaDataDictionary& dict = output->GetMetaDataDictionary();
+
+  itk::EncapsulateMetaData<std::string>(dict, MetaDataKey::ProjectionRefKey, m_OutputProjectionRef);
+}
+
+/*virtual*/
+void otb::GeometriesProjectionFilter::FinalizeInitialisation()
+{
+  m_Transform = InternalTransformType::New();
+
+  InputGeometriesType::ConstPointer   input      = this->GetInput();
+  const itk::MetaDataDictionary&      inputDict  = input->GetMetaDataDictionary();
+  OutputGeometriesType::Pointer       output     = this->GetOutput();
+  itk::MetaDataDictionary&            outputDict = output->GetMetaDataDictionary();
+
+  // m_Transform->SetInputDictionary(input->GetMetaDataDictionary());
+  m_Transform->SetInputDictionary(inputDict);
+  m_Transform->SetOutputDictionary(output->GetMetaDataDictionary());
+
+  m_Transform->SetInputProjectionRef(m_InputProjectionRef);
+  m_Transform->SetOutputProjectionRef(m_OutputProjectionRef);
+  m_Transform->SetInputKeywordList(m_InputKeywordList);
+  m_Transform->SetOutputKeywordList(m_OutputKeywordList);
+  m_Transform->SetDEMDirectory(m_DEMDirectory);
+  m_Transform->SetGeoidFile(m_GeoidFile);
+  m_Transform->SetAverageElevation(k_averageElevation);
+
+  m_Transform->SetInputSpacing(InputImageReference.GetSpacing());
+  m_Transform->SetInputOrigin(InputImageReference.GetOrigin());
+  m_Transform->SetOutputSpacing(OutputImageReference.GetSpacing());
+  m_Transform->SetOutputOrigin(OutputImageReference.GetOrigin());
+
+  m_Transform->InstanciateTransform();
+
+  m_TransformationFunctor->SetOnePointTransformation(m_Transform);
+
+  // retrieve the output projection ref
+  // if it is not specified and end up being geographic,
+  // only the m_Transform will know
+  m_OutputProjectionRef = m_Transform->GetOutputProjectionRef();
+
+  //If the projection information for the output is provided, propagate it
+
+  if (m_OutputKeywordList.GetSize() != 0)
+    {
+    itk::EncapsulateMetaData<ImageKeywordlist>(outputDict, MetaDataKey::OSSIMKeywordlistKey, m_OutputKeywordList);
+    }
+  if (!m_OutputProjectionRef.empty())
+    {
+    itk::EncapsulateMetaData<std::string>(outputDict, MetaDataKey::ProjectionRefKey, m_OutputProjectionRef);
+    }
+  output->SetImageReference(OutputImageReference);
+}
+
+/*virtual*/
+OGRSpatialReference* otb::GeometriesProjectionFilter::DoDefineNewLayerSpatialReference(
+  ogr::Layer const& source) const
+{
+  if (!m_OutputProjectionRef.empty())
+    {
+    char const* wkt_string = m_OutputProjectionRef.c_str();
+    OGRSpatialReference * res = static_cast <OGRSpatialReference*>(OSRNewSpatialReference(wkt_string));
+    return res;
+    }
+  else
+    {
+    return 0;
+    }
+}
+
+/*virtual*/
+void otb::GeometriesProjectionFilter::DoProcessLayer(ogr::Layer const& source, ogr::Layer & destination) const
+{
+  // std::cout << "GPF::DoProcessLayer: L("<<source.GetName()<<") -> L("<<destination.GetName()<<") ...\n";
+  if (source != destination)
+    {
+    m_TransformationFunctor(source, destination); // if TransformedElementType == layer
+    }
+  else
+    {
+    m_TransformationFunctor(destination); // if TransformedElementType == layer
+    }
+}

Code/Projections/otbGeometriesProjectionFilter.h

+/*=========================================================================
+
+  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 __otbGeometriesProjectionFilter_h
+#define __otbGeometriesProjectionFilter_h
+
+#include "otbGeometriesToGeometriesFilter.h"
+#include "otbOGRGeometryWrapper.h"
+#include "otbImageReference.h"
+#include "itkTransform.h"
+#include "otbGenericRSTransform.h"
+
+class OGRCoordinateTransformation;
+
+namespace otb
+{
+
+struct ReprojectTransformationFunctor
+{
+  typedef OGRGeometry TransformedElementType;
+
+  struct ByCopy
+    {
+    ByCopy(ReprojectTransformationFunctor const& reprojector) : m_Reprojector(reprojector){}
+    template <typename TGeometry>
+      ogr::UniqueGeometryPtr operator()(TGeometry const* in) const;
+  private:
+    ReprojectTransformationFunctor const& m_Reprojector;
+    };
+
+  struct InPlace
+    {
+    InPlace(ReprojectTransformationFunctor const& reprojector) : m_Reprojector(reprojector){}
+    template <typename TGeometry>
+      void operator()(TGeometry * inout) const;
+  private:
+    ReprojectTransformationFunctor const& m_Reprojector;
+    };
+
+
+  ogr::UniqueGeometryPtr operator()(OGRGeometry const* in) const;
+  // ogr::UniqueGeometryPtr apply(OGRGeometry const* in) const;
+  void apply_inplace          (OGRGeometry * inout) const;
+
+  typedef otb::GenericRSTransform<double, 2, 2>   InternalTransformType;
+  typedef InternalTransformType::Pointer          InternalTransformPointerType;
+  void SetOnePointTransformation(InternalTransformPointerType transform)
+    {
+    m_Transform = transform;
+    }
+private:
+  void do_transform(OGRPoint              & g) const;
+  // void do_transform(OGRLinearRing         & g) const;
+  void do_transform(OGRLineString         & g) const;
+  // void do_transform(OGRCurve              & g) const;
+  void do_transform(OGRPolygon            & g) const;
+#if 0
+  void do_transform(OGRSurface            & g) const;
+  void do_transform(OGRMultiLineString    & g) const;
+  void do_transform(OGRMultiPoint         & g) const;
+  void do_transform(OGRMultiPolygon       & g) const;
+#endif
+  void do_transform(OGRGeometryCollection & g) const;
+
+  InternalTransformPointerType m_Transform;
+};
+
+
+class ITK_EXPORT GeometriesProjectionFilter : public GeometriesToGeometriesFilter
+{
+public:
+  /**\name Standard ITK typedefs */
+  //@{
+  typedef GeometriesProjectionFilter    Self;
+  typedef GeometriesToGeometriesFilter  Superclass;
+  typedef itk::SmartPointer<Self>       Pointer;
+  typedef itk::SmartPointer<const Self> ConstPointer;
+  //@}
+
+  /**\name Standard macros */
+  //@{
+  /** Method for creation through the object factory. */
+  itkNewMacro(Self);
+
+  /** Run-time type information (and related methods). */
+  itkTypeMacro(GeometriesProjectionFilter, GeometriesToGeometriesFilter);
+  //@}
+
+  /**\name I/O typedefs */
+  //@{
+  typedef Superclass::InputGeometriesType        InputGeometriesType;
+  // typedef Superclass::InputGeometriesPointer  InputGeometriesPointer;
+  typedef Superclass::OutputGeometriesType       OutputGeometriesType;
+  // typedef Superclass::OutputGeometriesPointer OutputGeometriesPointer;
+  //@}
+
+  /**\name Class typedefs */
+  //@{
+  typedef ogr::ImageReference<double>                            ImageReference;
+  //@}
+
+private:
+  virtual OGRSpatialReference* DoDefineNewLayerSpatialReference(ogr::Layer const& source) const;
+  virtual void DoProcessLayer(ogr::Layer const& source, ogr::Layer & destination) const;
+  virtual void FinalizeInitialisation();
+
+protected:
+  /** Default constructor. */
+  GeometriesProjectionFilter();
+  /** Destructor. */
+  virtual ~GeometriesProjectionFilter();
+
+  virtual void GenerateOutputInformation(void);
+
+public:
+  ImageReference InputImageReference;
+  ImageReference OutputImageReference;
+
+  /**\name Keywords lists accessors and mutators */
+  //@{
+  itkGetMacro(InputKeywordList, ImageKeywordlist);
+  void SetInputKeywordList(const ImageKeywordlist& kwl)
+    {
+    this->m_InputKeywordList = kwl;
+    this->Modified();
+    }
+
+  itkGetMacro(OutputKeywordList, ImageKeywordlist);
+  void SetOutputKeywordList(const ImageKeywordlist& kwl)
+    {
+    this->m_OutputKeywordList = kwl;
+    this->Modified();
+    }
+  //@}
+
+  /**\name Projection references accessors and mutators */
+  //@{
+  itkSetStringMacro(InputProjectionRef);
+  itkGetStringMacro(InputProjectionRef);
+
+  itkSetStringMacro(OutputProjectionRef);
+  itkGetStringMacro(OutputProjectionRef);
+  //@}
+private:
+  /**\name Functor definition */
+  //@{
+  typedef ReprojectTransformationFunctor                          TransformationFunctorType;
+  typedef TransformationFunctorType::TransformedElementType       TransformedElementType;
+  typedef TransformationFunctorDispatcher<TransformationFunctorType, TransformedElementType>
+                                                                  TransformationFunctorDispatcherType;
+
+  TransformationFunctorDispatcherType             m_TransformationFunctor;
+  //@}
+
+  /**\name 1 Point Transformation definition */
+  //@{
+  typedef TransformationFunctorType::InternalTransformType        InternalTransformType;
+  typedef TransformationFunctorType::InternalTransformPointerType InternalTransformPointerType;
+
+  InternalTransformPointerType                    m_Transform;
+  //@}
+
+  std::string                                     m_InputProjectionRef; // in WKT format!
+  std::string                                     m_OutputProjectionRef; // in WKT format!
+  ImageKeywordlist                                m_InputKeywordList;
+  ImageKeywordlist                                m_OutputKeywordList;
+
+  std::string                                     m_DEMDirectory;
+  std::string                                     m_GeoidFile;
+};
+} // end namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+// #include "otbGeometriesProjectionFilter.txx"
+#endif
+
+#endif // __otbGeometriesProjectionFilter_h

Code/UtilitiesAdapters/OGRAdapters/otbGeometriesSet.cxx

@@ -58,17 +58,20 @@ otb::ogr::Layer & getRef(otb::ogr::Layer & layer)
 /*===========================================================================*/
 otb::GeometriesSet::GeometriesSet()
 : m_geometriesSet(otb::ogr::DataSource::New())
+, m_ImageReference(*this)
 {
 }
 
 otb::GeometriesSet::GeometriesSet(ogr::DataSource::Pointer datasource)
 : m_geometriesSet(datasource)
+, m_ImageReference(*this)
 {
   assert(datasource && "unexpected nil datasource");
 }
 
 otb::GeometriesSet::GeometriesSet(ogr::Layer layer)
 : m_geometriesSet(layer)
+, m_ImageReference(*this)
 {
 }
 
@@ -148,3 +151,28 @@ void otb::GeometriesSet::PrintSelf(std::ostream& os, itk::Indent indent) const
 {
   this->apply(Printer(os, indent));
 }
+
+/*===========================================================================*/
+/*=========================[ GetMetaDataDictionary ]=========================*/
+/*===========================================================================*/
+struct MetaDataDictionaryGetter : boost::static_visitor<itk::MetaDataDictionary&>
+{
+  itk::MetaDataDictionary& operator()(otb::ogr::Layer layer) const
+    {
+    return layer.GetMetaDataDictionary();
+    }
+  itk::MetaDataDictionary& operator()(otb::ogr::DataSource::Pointer datasource) const
+    {
+    return datasource->GetMetaDataDictionary();
+    }
+};
+
+itk::MetaDataDictionary & otb::GeometriesSet::GetMetaDataDictionary()
+{
+  return this->apply(MetaDataDictionaryGetter());
+}
+
+itk::MetaDataDictionary const& otb::GeometriesSet::GetMetaDataDictionary() const
+{
+  return const_cast <GeometriesSet *>(this)->GetMetaDataDictionary();
+}

Code/UtilitiesAdapters/OGRAdapters/otbGeometriesSet.h

@@ -22,6 +22,12 @@
 #include "itkDataObject.h"
 // #include "otbOGRLayerWrapper.h"
 #include "otbOGRDataSourceWrapper.h"
+#include "otbImageReference.h"
+
+namespace itk
+{
+class MetaDataDictionary;
+}
 
 namespace otb
 {
@@ -107,6 +113,23 @@ public:
    */
   bool IsSet() const;
 
+  /**\name Meta data dictionary */
+  //@{
+  itk::MetaDataDictionary      & GetMetaDataDictionary();
+  itk::MetaDataDictionary const& GetMetaDataDictionary() const;
+  //@}
+
+  /**\name Image reference (spacing + origin) */
+  //@{
+  typedef ogr::ImageReference<double>                            ImageReference;
+  /** imageReference setter. */
+  void SetImageReference(ImageReference const& imageReference)
+    { m_ImageReference = imageReference; }
+  /** imageReference getter. */
+  ImageReference const& GetImageReference() const
+    { return m_ImageReference; }
+  //@}
+
 protected:
   /** Default constructor.
    * This actual geometries set is an in-memory \c otb::ogr::DataSource.
@@ -130,6 +153,7 @@ protected:
 private:
   typedef boost::variant<ogr::DataSource::Pointer, ogr::Layer> AnyGeometriesSetType;
   AnyGeometriesSetType m_geometriesSet;
+  ImageReference       m_ImageReference;
   };
 
 } // end namespace otb

Code/UtilitiesAdapters/OGRAdapters/otbGeometriesToGeometriesFilter.cxx

@@ -70,15 +70,19 @@ struct ProcessVisitor : boost::static_visitor<>
     assert(destination && "can't filter to a nil datasource");
     // std::cout << "G2GF: Process Visitor: DS("<<source->ogr().GetName()<<") -> DS("<<destination->ogr().GetName()<<") ...\n";
     for (otb::ogr::DataSource::const_iterator b = source->begin(), e = source->end()
-; b != e
-; ++b
+      ; b != e
+      ; ++b
     )
       {
       otb::ogr::Layer const& sourceLayer = *b;
       assert(sourceLayer && "unexpected nil source layer");
       // std::cout << "source layer name:" << sourceLayer.GetName() << "\n";
       otb::ogr::Layer destLayer = destination->CreateLayer(
-        sourceLayer.GetName(), 0, sourceLayer.GetGeomType());
+        sourceLayer.GetName(),
+        m_filter.DoDefineNewLayerSpatialReference(sourceLayer),
+        m_filter.DoDefineNewLayerGeometryType(sourceLayer),
+        otb::ogr::StringListConverter(m_filter.DoDefineNewLayerOptions(sourceLayer)).to_ogr()
+        );
       // std::cout << "layer created!\n";