ENH : Add a new application SARORthoInterferogram with two filters :...

ENH : Add a new application SARORthoInterferogram with two filters : SARCompensatedComplexImageFilter and SARGroupedByOrthoImageFilter

app/CMakeLists.txt

@@ -103,3 +103,8 @@ OTB_CREATE_APPLICATION(NAME SARMetadataCorrection
                        SOURCES otbSARMetadataCorrection.cxx
                        LINK_LIBRARIES ${${otb-module}_LIBRARIES}
 )
+
+OTB_CREATE_APPLICATION(NAME SAROrthoInterferogram
+                       SOURCES otbSAROrthoInterferogram.cxx
+                       LINK_LIBRARIES ${${otb-module}_LIBRARIES}
+)

app/otbSAROrthoInterferogram.cxx

+/*
+ * Copyright (C) 2005-2018 Centre National d'Etudes Spatiales (CNES)
+ *
+ * This file is part of Orfeo Toolbox
+ *
+ *     https://www.orfeo-toolbox.org/
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "otbWrapperApplication.h"
+#include "otbWrapperApplicationFactory.h"
+
+#include "otbSARCompensatedComplexImageFilter.h"
+#include "otbSARGroupedByOrthoImageFilter.h"
+
+#include "otbWrapperOutputImageParameter.h"
+#include "otbWrapperTypes.h"
+
+#include <iostream>
+#include <string>
+#include <fstream>
+
+namespace otb
+{
+
+  namespace Wrapper
+  {
+
+    class SAROrthoInterferogram : public Application
+    {
+    public:
+      typedef SAROrthoInterferogram Self;
+      typedef itk::SmartPointer<Self> Pointer; 
+
+      itkNewMacro(Self);
+      itkTypeMacro(SAROrthoInterferogram, otb::Wrapper::Application);
+
+      // Filters
+      typedef otb::SARCompensatedComplexImageFilter<ComplexFloatImageType, FloatVectorImageType, FloatVectorImageType> CompensatedComplexFilterType;
+      typedef otb::SARGroupedByOrthoImageFilter<FloatVectorImageType, FloatImageType, FloatVectorImageType> GroupedByOrthoFilterType;
+
+    private:
+      void DoInit() override
+      {
+	SetName("SAROrthoInterferogram");
+	SetDescription("Interferogram into ground geometry between two SAR images.");
+
+	SetDocLongDescription("This application builds the interferogram into ground geometry between two"
+        " SAR images.");
+
+	//Optional descriptors
+	SetDocLimitations("Only Sentinel 1 (IW and StripMap mode) and Cosmo products are supported for now.");
+	SetDocAuthors("OTB-Team");
+	SetDocSeeAlso(" ");
+	AddDocTag(Tags::SAR);
+	AddDocTag("DiapOTB");
+
+	//Parameter declarations
+	AddParameter(ParameterType_InputImage,  "insarslave",   "Input SAR Slave image (Coregistrated image)");
+	SetParameterDescription("insarslave", "Input SAR Slave image (Coregistrated image).");
+
+	AddParameter(ParameterType_InputImage,  "insarmaster",   "Input SAR Master image");
+	SetParameterDescription("insarmaster", "Input SAR Master image.");
+
+	AddParameter(ParameterType_InputImage,  "topographicphase",   "Input Topographic Phase (estimation with DEM projection)");
+	SetParameterDescription("topographicphase", "Input Topographic Phase (estimation with DEM projection).");
+	MandatoryOff("topographicphase");
+
+        AddParameter(ParameterType_InputImage,  "incartmeanmaster",   "Input Cartesian Mean Master image");
+	SetParameterDescription("incartmeanmaster", "Input Cartesian Mean Master image.");
+
+        AddParameter(ParameterType_InputImage,  "indem",   "Input DEM");
+        SetParameterDescription("indem", "GetGround geometry thanks to DEM (only metadata).");
+	
+	AddParameter(ParameterType_Float, "gain", "Gain to apply for amplitude estimation");
+	SetParameterDescription("gain","Gain to apply for amplitude estimation");
+	SetDefaultParameterFloat("gain", 0.1);
+	SetMinimumParameterFloatValue("gain", 0);
+	MandatoryOff("gain");
+
+	AddParameter(ParameterType_OutputImage, "out", "Interferogram");
+	SetParameterDescription("out","Output Vector Image : Interferogram.");
+
+	AddRAMParameter();
+
+	SetDocExampleParameterValue("insarslave","s1b-s4-slc-vv-20160929t014610-20160929t014634-002277-003d71-002.tiff");
+	SetDocExampleParameterValue("insarmaster","s1a-s4-slc-vv-20160818t014650-20160818t014715-012648-013db1-002_SLC.tiff");
+        SetDocExampleParameterValue("incartmeanmaster","CartesianMaster.tiff");
+        SetDocExampleParameterValue("indem","S21E055.hgt");
+	SetDocExampleParameterValue("gain","0.1");
+	SetDocExampleParameterValue("out","s1b-s1a-s4-interferogram.tiff");
+      }
+
+      void DoUpdateParameters() override
+      {
+	// Nothing to do here : all parameters are independent
+      }
+
+      void DoExecute() override
+      { 
+	// Get numeric parameters
+	double factor_gain = GetParameterFloat("gain");
+	
+	otbAppLogINFO(<<"Gain Factor : "<<factor_gain);
+
+	/////////////////////////////////// Compensated Complex Filter ////////////////////////////////////////
+	// Instanciate the first filter
+	CompensatedComplexFilterType::Pointer filterCompensatedComplex = CompensatedComplexFilterType::New();
+	m_Ref.push_back(filterCompensatedComplex.GetPointer());
+	
+	// Execute the Pipeline
+	ComplexFloatImageType::Pointer SARMasterPtr;
+	SARMasterPtr = GetParameterComplexFloatImage("insarmaster");
+
+	ComplexFloatImageType::Pointer SARSlavePtr;
+	SARSlavePtr = GetParameterComplexFloatImage("insarslave");
+
+	// Two Main inputs
+	filterCompensatedComplex->SetMasterInput(SARMasterPtr);
+	filterCompensatedComplex->SetSlaveInput(SARSlavePtr);
+
+	// One optionnal input
+	if (GetParameterByKey("topographicphase")->HasValue())
+	  {
+	    filterCompensatedComplex->SetTopographicPhaseInput(GetParameterImage("topographicphase"));
+	  }
+	
+
+	/////////////////////////////// GroupedBy Filter (for ground geometry) //////////////////////////////////
+	// Instanciate the second filter
+	FloatImageType::Pointer inputDEM =  GetParameterFloatImage("indem");
+	GroupedByOrthoFilterType::Pointer filterGroupedBy = GroupedByOrthoFilterType::New();
+	m_Ref.push_back(filterGroupedBy.GetPointer());
+	filterGroupedBy->SetSARImageKeyWorList(SARMasterPtr->GetImageKeywordlist());
+	filterGroupedBy->SetDEMImagePtr(inputDEM);
+
+	
+	// Execute the Pipeline
+	FloatVectorImageType::Pointer CartMeanPtr;
+	CartMeanPtr = GetParameterFloatVectorImage("incartmeanmaster");
+
+
+	filterGroupedBy->SetCartesianMeanInput(CartMeanPtr);
+	filterGroupedBy->SetCompensatedComplexInput(filterCompensatedComplex->GetOutput());
+	
+
+	// Main Output
+	SetParameterOutputImage("out", filterGroupedBy->GetOutput());
+
+      }
+      // Vector for filters 
+      std::vector<itk::ProcessObject::Pointer> m_Ref;
+
+    };
+
+  }
+
+}
+
+
+ 
+OTB_APPLICATION_EXPORT(otb::Wrapper::SAROrthoInterferogram)

include/otbSARCompensatedComplexImageFilter.h

+/*
+ * Copyright (C) 2005-2018 Centre National d'Etudes Spatiales (CNES)
+ *
+ * This file is part of Orfeo Toolbox
+ *
+ *     https://www.orfeo-toolbox.org/
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef otbSARCompensatedComplexImageFilter_h
+#define otbSARCompensatedComplexImageFilter_h
+
+#include "itkImageToImageFilter.h"
+#include "itkSmartPointer.h"
+#include "itkPoint.h"
+
+#include "itkSimpleFastMutexLock.h"
+
+#include "itkImageScanlineConstIterator.h"
+#include "itkImageScanlineIterator.h"
+
+#include "otbImageKeywordlist.h"
+
+
+namespace otb
+{
+/** \class SARCompensatedComplexImageFilter 
+ * \brief Creates an temporary image (to estimate conjugate multiplication) between two images. 
+ * 
+ * This filter built the conjugate multiplication between a master and a slave image.
+ * A topographic phase can be used to create this temporary image.
+ *
+ * The output is a vector image with real part, imag part and modulus of conjugate multiplication.
+ *
+ * \ingroup DiapOTBModule
+ */
+
+  template <typename TImageSAR, typename TImagePhase, typename TImageOut> 
+  class ITK_EXPORT SARCompensatedComplexImageFilter :
+    public itk::ImageToImageFilter<TImageSAR,TImageOut>
+{
+public:
+
+  // Standard class typedefs
+  typedef SARCompensatedComplexImageFilter                    Self;
+  typedef itk::ImageToImageFilter<TImageSAR,TImageOut>    Superclass;
+  typedef itk::SmartPointer<Self>                        Pointer;
+  typedef itk::SmartPointer<const Self>                  ConstPointer;
+
+  // Method for creation through object factory
+  itkNewMacro(Self);
+  // Run-time type information
+  itkTypeMacro(SARCompensatedComplexImageFilter,ImageToImageFilter);
+
+  /** Typedef to image input type (master and slave reech images) :  otb::Image (Complex)  */
+  typedef TImageSAR                                  ImageSARType;
+  /** Typedef to describe the inout image pointer type. */
+  typedef typename ImageSARType::Pointer             ImageSARPointer;
+  typedef typename ImageSARType::ConstPointer        ImageSARConstPointer;
+  /** Typedef to describe the inout image region type. */
+  typedef typename ImageSARType::RegionType          ImageSARRegionType;
+  /** Typedef to describe the type of pixel and point for inout image. */
+  typedef typename ImageSARType::PixelType           ImageSARPixelType;
+  typedef typename ImageSARType::PointType           ImageSARPointType;
+  /** Typedef to describe the image index, size types and spacing for inout image. */
+  typedef typename ImageSARType::IndexType           ImageSARIndexType;
+  typedef typename ImageSARType::IndexValueType      ImageSARIndexValueType;
+  typedef typename ImageSARType::SizeType            ImageSARSizeType;
+  typedef typename ImageSARType::SizeValueType       ImageSARSizeValueType;
+  typedef typename ImageSARType::SpacingType         ImageSARSpacingType;
+  typedef typename ImageSARType::SpacingValueType    ImageSARSpacingValueType;
+  
+  /** Typedef to image output type : otb::VectorImage  */
+  typedef TImageOut                                  ImageOutType;
+  /** Typedef to describe the output image pointer type. */
+  typedef typename ImageOutType::Pointer             ImageOutPointer;
+  typedef typename ImageOutType::ConstPointer        ImageOutConstPointer;
+  /** Typedef to describe the output image region type. */
+  typedef typename ImageOutType::RegionType          ImageOutRegionType;
+  /** Typedef to describe the type of pixel and point for output image. */
+  typedef typename ImageOutType::PixelType           ImageOutPixelType;
+  typedef typename ImageOutType::PointType           ImageOutPointType;
+  /** Typedef to describe the image index, size types and spacing for output image. */
+  typedef typename ImageOutType::IndexType           ImageOutIndexType;
+  typedef typename ImageOutType::IndexValueType      ImageOutIndexValueType;
+  typedef typename ImageOutType::SizeType            ImageOutSizeType;
+  typedef typename ImageOutType::SizeValueType       ImageOutSizeValueType;
+  typedef typename ImageOutType::SpacingType         ImageOutSpacingType;
+  typedef typename ImageOutType::SpacingValueType    ImageOutSpacingValueType;
+
+
+  // Define Point2DType and Point3DType
+  using Point2DType = itk::Point<double,2>;
+  using Point3DType = itk::Point<double,3>;
+
+  // Typedef for topographic Phase
+  typedef TImagePhase                                      ImagePhaseType;
+  typedef typename ImagePhaseType::Pointer                 ImagePhasePointer;
+  typedef typename ImagePhaseType::ConstPointer            ImagePhaseConstPointer;
+  typedef typename ImagePhaseType::RegionType              ImagePhaseRegionType;
+  typedef typename ImagePhaseType::IndexType               ImagePhaseIndexType;
+  typedef typename ImagePhaseType::IndexValueType          ImagePhaseIndexValueType;
+  typedef typename ImagePhaseType::SizeType                ImagePhaseSizeType;
+  typedef typename ImagePhaseType::SizeValueType           ImagePhaseSizeValueType;
+
+  // Typedef for iterators
+  typedef itk::ImageScanlineConstIterator< ImageSARType > InputSARIterator;
+  typedef itk::ImageScanlineConstIterator< ImagePhaseType > InputPhaseIterator;
+  typedef itk::ImageScanlineIterator< ImageOutType > OutputIterator;
+
+  typedef itk::SimpleFastMutexLock         MutexType;
+
+    
+  // Setter/Getter for inputs (SAR images and potentially topographic phase) 
+  /** Connect one of the operands for interferometry : Master */
+  void SetMasterInput( const ImageSARType* image);
+
+  /** Connect one of the operands for interferometry : Coregistrated Slave */
+  void SetSlaveInput(const ImageSARType* image);
+
+  /** Connect one of the operands for interferometry : Coregistrated Slave */
+  void SetTopographicPhaseInput(const ImagePhaseType* imagePhaseTopo);
+
+  /** Get the inputs */
+  const ImageSARType* GetMasterInput() const;
+  const ImageSARType * GetSlaveInput() const;
+  const ImagePhaseType * GetTopographicPhaseInput() const;
+
+
+protected:
+  // Constructor
+  SARCompensatedComplexImageFilter();
+
+  // Destructor
+  virtual ~SARCompensatedComplexImageFilter() ITK_OVERRIDE;
+
+  // Print
+  void PrintSelf(std::ostream & os, itk::Indent indent) const ITK_OVERRIDE;
+  
+  /** SARCompensatedComplexImageFilter produces an image which are into output geometry. The differences between 
+   * output and input images can be the size of images and the dimensions. 
+   * As such, SARCompensatedComplexImageFilter needs to provide an implementation for 
+   * GenerateOutputInformation() in order to inform the pipeline execution model. 
+   */ 
+  virtual void GenerateOutputInformation() ITK_OVERRIDE;
+
+  /** SARCompensatedComplexImageFilter needs a input requested region that corresponds to the margin and shifts 
+   * into the requested region of our output requested region. The output requested region needs to be modified
+   * to construct as wanted tiles with input size.
+   * As such, TilesAnalysesImageFilter needs to provide an implementation for 
+   * GenerateInputRequestedRegion() in order to inform the pipeline execution model. 
+   * \sa ProcessObject::GenerateInputRequestedRegion() */
+  virtual void GenerateInputRequestedRegion() ITK_OVERRIDE;
+
+
+  /** 
+   * SARCompensatedComplexImageFilterr can be implemented as a multithreaded filter.
+   * Therefore, this implementation provides a ThreadedGenerateData() routine
+   * which is called for each processing thread. The main output image data is
+   * allocated automatically by the superclass prior to calling
+   * ThreadedGenerateData().  ThreadedGenerateData can only write to the
+   * portion of the output image specified by the parameter
+   * "outputRegionForThread"
+   *
+   * \sa ImageToImageFilter::ThreadedGenerateData(),
+   *     ImageToImageFilter::GenerateData() */
+  virtual void ThreadedGenerateData(const ImageOutRegionType& outputRegionForThread, 
+				    itk::ThreadIdType threadId) ITK_OVERRIDE;
+  
+  virtual void ThreadedGenerateDataWithoutTopographicPhase(const ImageOutRegionType& outputRegionForThread, 
+							   itk::ThreadIdType threadId);
+  virtual void ThreadedGenerateDataWithTopographicPhase(const ImageOutRegionType& outputRegionForThread, 
+							itk::ThreadIdType threadId);
+
+
+ private:
+  SARCompensatedComplexImageFilter(const Self&); // purposely not implemented
+  void operator=(const Self &); // purposely not 
+
+  // SAR Dimensions
+  int m_nbLinesSAR;
+  int m_nbColSAR;
+
+  
+
+};
+
+} // End namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbSARCompensatedComplexImageFilter.txx"
+#endif
+
+
+
+#endif

include/otbSARCompensatedComplexImageFilter.txx

+/*
+ * Copyright (C) 2005-2018 Centre National d'Etudes Spatiales (CNES)
+ *
+ * This file is part of Orfeo Toolbox
+ *
+ *     https://www.orfeo-toolbox.org/
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef otbSARCompensatedComplexImageFilter_txx
+#define otbSARCompensatedComplexImageFilter_txx
+
+#include "otbSARCompensatedComplexImageFilter.h"
+
+#include "itkImageScanlineConstIterator.h"
+#include "itkImageScanlineIterator.h"
+#include "itkProgressReporter.h"
+#include "itkNumericTraitsPointPixel.h"
+#include "itkContinuousIndex.h"
+
+#include "ossim/ossimSarSensorModel.h"
+
+#include <cmath>
+#include <algorithm>
+
+
+namespace otb
+{
+  /** 
+   * Constructor with default initialization
+   */
+  template <class TImageSAR, class TImagePhase, class TImageOut> 
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >::SARCompensatedComplexImageFilter()
+    :  m_nbLinesSAR(0), m_nbColSAR(0)
+  {
+    // Inputs required and/or needed
+    this->SetNumberOfRequiredInputs(2);
+    this->SetNumberOfIndexedInputs(3);
+  }
+    
+  /** 
+   * Destructor
+   */
+  template <class TImageSAR, class TImagePhase, class TImageOut> 
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >::~SARCompensatedComplexImageFilter()
+  {
+  }
+
+  /**
+   * Print
+   */
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::PrintSelf(std::ostream & os, itk::Indent indent) const
+  {
+    Superclass::PrintSelf(os, indent);
+  }
+
+  /**
+   * Set Master Image
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+ SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::SetMasterInput(const ImageSARType* image )
+  {
+    // Process object is not const-correct so the const casting is required.
+    this->SetNthInput(0, const_cast<ImageSARType *>(image));
+  }
+
+  /**
+   * Set Coregistrated Slave Image
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+ SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::SetSlaveInput(const ImageSARType* image )
+  {
+    // Process object is not const-correct so the const casting is required.
+    this->SetNthInput(1, const_cast<ImageSARType *>(image));
+  }
+
+  /**
+   * Set Topographic phase
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::SetTopographicPhaseInput(const ImagePhaseType* phaseTopoImage)
+  {
+    // Process object is not const-correct so the const casting is required.
+    this->SetNthInput(2, const_cast<ImagePhaseType *>(phaseTopoImage));
+  }
+
+  /**
+   * Get Master Image
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  const typename SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >::ImageSARType *
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::GetMasterInput() const
+  {
+    if (this->GetNumberOfInputs()<1)
+      {
+	return 0;
+      }
+    return static_cast<const ImageSARType *>(this->itk::ProcessObject::GetInput(0));
+  }
+
+  /**
+   * Get Coregistrated Slave Image
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  const typename SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >::ImageSARType *
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::GetSlaveInput() const
+  {
+    if (this->GetNumberOfInputs()<2)
+      {
+	return 0;
+      }
+    return static_cast<const ImageSARType *>(this->itk::ProcessObject::GetInput(1));
+  }
+
+  /**
+   * Get Phase Topographic
+   */ 
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  const typename SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >::ImagePhaseType *
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::GetTopographicPhaseInput() const
+  {
+    if (this->GetNumberOfInputs()<3)
+      {
+	return 0;
+      }
+    return static_cast<const ImagePhaseType *>(this->itk::ProcessObject::GetInput(2));
+  }
+
+  /**
+   * Method GenerateOutputInformaton()
+   **/
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::GenerateOutputInformation()
+  {
+    // Call superclass implementation
+    Superclass::GenerateOutputInformation();
+
+    // Get pointers to the inputs and output
+    ImageSARConstPointer masterPtr = this->GetMasterInput();
+    ImageSARConstPointer slavePtr = this->GetSlaveInput();
+    ImageOutPointer outputPtr = this->GetOutput();
+
+    // KeyWordList
+    ImageKeywordlist masterKWL = masterPtr->GetImageKeywordlist();
+
+    // Master SAR Dimensions
+    m_nbLinesSAR = this->GetMasterInput()->GetLargestPossibleRegion().GetSize()[1];
+    m_nbColSAR = this->GetMasterInput()->GetLargestPossibleRegion().GetSize()[0];
+  
+    /////////////////////// Main Output : conjugate multiplication (into SAR geometry) ///////////////////////
+    // Vector Image  :
+    // At Least 3 Components :  
+    //                _ real part
+    //                _ imag part
+    //                _ modulus
+    outputPtr->SetNumberOfComponentsPerPixel(3);
+
+    // The output is defined with the Master SAR Image
+    // Origin, Spacing and Size (SAR master geometry)
+    ImageOutPointType outOrigin;
+    outOrigin = masterPtr->GetOrigin();
+    ImageOutSpacingType outSP;
+    outSP = masterPtr->GetSpacing();
+
+    // Define Output Largest Region
+    ImageOutRegionType outputLargestPossibleRegion = masterPtr->GetLargestPossibleRegion();
+
+    outputPtr->SetLargestPossibleRegion(outputLargestPossibleRegion);
+    outputPtr->SetOrigin(outOrigin);
+    outputPtr->SetSpacing(outSP);
+
+    // Add ML factors and bands meaning into keyWordList
+    ImageKeywordlist outputKWL = masterKWL;
+    outputKWL.AddKey("support_data.band.Amplitude", std::to_string(0));
+    outputKWL.AddKey("support_data.band.Phase", std::to_string(1));
+    outputKWL.AddKey("support_data.band.Coherency", std::to_string(2));
+
+    // Set new keyword list to output image
+    outputPtr->SetImageKeywordList(outputKWL);    
+
+    ///////// Checks (with input topographic phase keywordlists/metadata) /////////////
+    // Check ML Factors (must be 1 to have the same geometry than SAR images)
+     if (this->GetTopographicPhaseInput() != nullptr)
+       {
+	 // Get ImagKeyWordList
+	 ImageKeywordlist topoKWL = this->GetTopographicPhaseInput()->GetImageKeywordlist();
+	 
+	 if (topoKWL.HasKey("support_data.ml_ran") && topoKWL.HasKey("support_data.ml_azi"))
+	   {
+	     // Get Master ML Factors
+	     unsigned int topo_MLRan = atoi(topoKWL.GetMetadataByKey("support_data.ml_ran").c_str());
+	     unsigned int topo_MLAzi = atoi(topoKWL.GetMetadataByKey("support_data.ml_azi").c_str());
+	
+	     if ((topo_MLRan != 1) || (topo_MLAzi != 1))
+	       {
+		 itkExceptionMacro(<<"Error, ML Factor for topographic phase are different than 1.");
+	       }
+	   }
+       }
+  }
+
+
+
+  /** 
+   * Method GenerateInputRequestedRegion
+   */
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::GenerateInputRequestedRegion()
+  {
+    // call the superclass' implementation of this method
+    Superclass::GenerateInputRequestedRegion();
+    
+    // Get Output requested region
+    ImageOutRegionType outputRequestedRegion = this->GetOutput()->GetRequestedRegion();
+
+    // Input Region = Output region
+    ImageSARPointer  masterPtr = const_cast< ImageSARType * >( this->GetMasterInput() );
+    ImageSARPointer  slavePtr = const_cast< ImageSARType * >( this->GetSlaveInput() );
+    
+    masterPtr->SetRequestedRegion(outputRequestedRegion);
+    slavePtr->SetRequestedRegion(outputRequestedRegion);
+
+    ///////////// Find the region into topographic phase image (if needed) /////////////
+    if (this->GetTopographicPhaseInput() != nullptr)
+      {
+	ImagePhasePointer  phasePtr = const_cast< ImagePhaseType * >( this->GetTopographicPhaseInput() );
+	phasePtr->SetRequestedRegion(outputRequestedRegion);
+      }
+  }
+  
+
+  /**
+   * Method ThreadedGenerateData
+   */
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::ThreadedGenerateData(const ImageOutRegionType & outputRegionForThread,
+			 itk::ThreadIdType threadId)
+  {
+    if (this->GetTopographicPhaseInput() == nullptr)
+      {
+	// Process for line
+	this->ThreadedGenerateDataWithoutTopographicPhase(outputRegionForThread,threadId);
+      } 
+    else
+      {
+	// Process for column
+	this->ThreadedGenerateDataWithTopographicPhase(outputRegionForThread,threadId);
+      } 
+  }
+
+  /**
+   * Method ThreadedGenerateDataWithoutTopographicPhase
+   */
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::ThreadedGenerateDataWithoutTopographicPhase(const ImageOutRegionType & outputRegionForThread,
+						itk::ThreadIdType threadId)
+  {
+    //  Define/declare an iterator that will walk the input regions for this
+    // thread. 
+    InputSARIterator  inMasterIt(this->GetMasterInput(), outputRegionForThread);
+    InputSARIterator  inSlaveIt(this->GetSlaveInput(), outputRegionForThread);
+
+    // Define/declare an iterator that will walk the output region for this
+    // thread.
+    OutputIterator outIt(this->GetOutput(), outputRegionForThread);
+
+    // Support progress methods/callbacks
+    itk::ProgressReporter progress( this, threadId, outputRegionForThread.GetNumberOfPixels() );
+
+    inMasterIt.GoToBegin();
+    inSlaveIt.GoToBegin();
+    outIt.GoToBegin();
+    
+    // Output Pixel (VectorImage Pixel)
+    ImageOutPixelType outPixel;
+    outPixel.Reserve(3);
+
+    // For each Line
+    while ( !inMasterIt.IsAtEnd() && !inSlaveIt.IsAtEnd() && !outIt.IsAtEnd())
+      {
+	inMasterIt.GoToBeginOfLine();
+	inSlaveIt.GoToBeginOfLine();
+	outIt.GoToBeginOfLine();
+	
+	    // For each column
+	    while (!inMasterIt.IsAtEndOfLine() && !inSlaveIt.IsAtEndOfLine() && !outIt.IsAtEndOfLine())
+	      {
+		// Conjugate multiplication (master * conj(slave))
+		double real_interfero = (inMasterIt.Get().real()*inSlaveIt.Get().real() + 
+					 inMasterIt.Get().imag()*inSlaveIt.Get().imag());
+
+		double imag_interfero = (inMasterIt.Get().imag()*inSlaveIt.Get().real() - 
+					 inMasterIt.Get().real()*inSlaveIt.Get().imag());
+
+		///////////// Output assignations ///////////////
+		outPixel[0] = real_interfero;
+
+			
+		outPixel[1] = imag_interfero;
+			
+			
+		outPixel[2] = sqrt((real_interfero*real_interfero) + 
+							 (imag_interfero*imag_interfero));
+
+		
+		outIt.Set(outPixel);
+	    
+		// Next colunm inputs
+		++inMasterIt;
+		++inSlaveIt;
+		// Next colunm output
+		++outIt;
+	      }
+	  
+	    // Next line intputs
+	    inMasterIt.NextLine();
+	    inSlaveIt.NextLine();
+	    // Next line output
+	    outIt.NextLine();
+      }
+  }
+
+   /**
+   * Method ThreadedGenerateDataWithtTopographicPhase
+   */
+  template<class TImageSAR, class TImagePhase, class TImageOut>
+  void
+  SARCompensatedComplexImageFilter< TImageSAR, TImagePhase, TImageOut >
+  ::ThreadedGenerateDataWithTopographicPhase(const ImageOutRegionType & outputRegionForThread,
+						itk::ThreadIdType threadId)
+  {
+    //  Define/declare an iterator that will walk the input regions for this
+    // thread. 
+    InputSARIterator  inMasterIt(this->GetMasterInput(), outputRegionForThread);
+    InputSARIterator  inSlaveIt(this->GetSlaveInput(), outputRegionForThread);
+
+    // Define/declare an iterator that will walk the output region for this
+    // thread.
+    OutputIterator outIt(this->GetOutput(), outputRegionForThread);
+
+    // Topographic phase
+    InputPhaseIterator inTopoPhaseIt(this->GetTopographicPhaseInput(), outputRegionForThread);
+
+    // Support progress methods/callbacks
+    itk::ProgressReporter progress( this, threadId, outputRegionForThread.GetNumberOfPixels() );
+
+    inMasterIt.GoToBegin();
+    inSlaveIt.GoToBegin();
+    inTopoPhaseIt.GoToBegin();
+    outIt.GoToBegin();
+    
+    // Output Pixel (VectorImage Pixel)
+    ImageOutPixelType outPixel;
+    outPixel.Reserve(3);
+
+    // For each Line
+    while ( !inMasterIt.IsAtEnd() && !inSlaveIt.IsAtEnd() && !outIt.IsAtEnd())
+      {
+	inMasterIt.GoToBeginOfLine();
+	inSlaveIt.GoToBeginOfLine();
+	inTopoPhaseIt.GoToBeginOfLine();
+	outIt.GoToBeginOfLine();
+	
+	    // For each column
+	    while (!inMasterIt.IsAtEndOfLine() && !inSlaveIt.IsAtEndOfLine() && !outIt.IsAtEndOfLine())
+	      {
+		// Check isData (grom topographic phase)
+		if (inTopoPhaseIt.Get()[1] != 0)
+		  {			
+		    // Complex Raw interferogram (master * conj(slave))
+		    double real_RawInterfero = (inMasterIt.Get().real()*inSlaveIt.Get().real() + 
+						inMasterIt.Get().imag()*inSlaveIt.Get().imag());
+
+		    double imag_RawInterfero = (inMasterIt.Get().imag()*inSlaveIt.Get().real() - 
+						inMasterIt.Get().real()*inSlaveIt.Get().imag());
+
+		    ///////////// Topographoc phase as complex number ////////////////
+		    double topoPhase = inTopoPhaseIt.Get()[0];
+				
+		    double complexTopoPhase_Re = std::cos(topoPhase);
+		    double complexTopoPhase_Im = std::sin(topoPhase);
+
+	
+		    // Multiply the conj(complexTopoPhase) with complex raw interferogram
+		    double real_interfero = (real_RawInterfero * complexTopoPhase_Re + 
+					     imag_RawInterfero * complexTopoPhase_Im);
+		    double imag_interfero = (imag_RawInterfero * complexTopoPhase_Re - 
+					     real_RawInterfero * complexTopoPhase_Im);
+
+		    ///////////// Output assignations ///////////////
+		    outPixel[0] = real_interfero;
+
+			
+		    outPixel[1] = imag_interfero;
+			
+			
+		    outPixel[2] = sqrt((real_interfero*real_interfero) + 
+				       (imag_interfero*imag_interfero));
+		  }
+		else
+		  {
+		    outPixel[0] = 0;
+		    outPixel[1] = 0;
+		    outPixel[2] = 0;
+		    outPixel[3] = 0;
+		  }
+	
+		outIt.Set(outPixel);
+	    
+		// Next colunm inputs
+		++inMasterIt;
+		++inSlaveIt;
+		++inTopoPhaseIt;
+		// Next colunm output
+		++outIt;
+	      }
+	  
+	    // Next line intputs
+	    inMasterIt.NextLine();
+	    inSlaveIt.NextLine();
+	    inTopoPhaseIt.NextLine();
+	    // Next line output
+	    outIt.NextLine();
+      }
+  }
+
+
+  } /*namespace otb*/
+
+#endif

include/otbSARGroupedByOrthoImageFilter.h

+/*
+ * Copyright (C) 2005-2018 Centre National d'Etudes Spatiales (CNES)
+ *
+ * This file is part of Orfeo Toolbox
+ *
+ *     https://www.orfeo-toolbox.org/
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef otbSARGroupedByOrthoImageFilter_h
+#define otbSARGroupedByOrthoImageFilter_h
+
+#include "itkImageToImageFilter.h"
+#include "itkSmartPointer.h"
+#include "itkPoint.h"
+
+#include "itkImageScanlineConstIterator.h"
+#include "itkImageScanlineIterator.h"
+
+#include "otbImageKeywordlist.h"
+#include "otbSarSensorModelAdapter.h"
+
+
+#if defined(__GNUC__) || defined(__clang__)
+# pragma GCC diagnostic push
+#   pragma GCC diagnostic ignored "-Woverloaded-virtual"
+#if defined(__GNUC__) && (__GNUC__ > 5)
+#pragma GCC diagnostic ignored "-Wnonnull-compare"
+#endif
+#   include "ossim/base/ossimGeoidEgm96.h"
+# pragma GCC diagnostic pop
+#else
+#   include "ossim/base/ossimGeoidEgm96.h"
+#endif
+
+
+namespace otb
+{
+/** \class SARGroupedByOrthoImageFilter 
+ * \brief Creates an interferogram into ground geometry between two images. 
+ * 
+ * This filter built the Ortho interferogram between a master and a slave image.
+ *
+ * The output is a vector image with amplitude, phase and coherency.
+ *
+ * \ingroup DiapOTBModule
+ */
+
+  template <typename TImageVector, typename TImageDEM, typename TImageOut> 
+  class ITK_EXPORT SARGroupedByOrthoImageFilter :
+    public itk::ImageToImageFilter<TImageVector,TImageOut>
+{
+public:
+
+  // Standard class typedefs
+  typedef SARGroupedByOrthoImageFilter                    Self;
+  typedef itk::ImageToImageFilter<TImageVector,TImageOut>    Superclass;
+  typedef itk::SmartPointer<Self>                        Pointer;
+  typedef itk::SmartPointer<const Self>                  ConstPointer;
+
+  // Method for creation through object factory
+  itkNewMacro(Self);
+  // Run-time type information
+  itkTypeMacro(SARGroupedByOrthoImageFilter,ImageToImageFilter);
+
+  /** Typedef to vector image input type : otb::VectorImage  */
+  typedef TImageVector                                  ImageVectorType;
+  /** Typedef to describe the output image pointer type. */
+  typedef typename ImageVectorType::Pointer             ImageVectorPointer;
+  typedef typename ImageVectorType::ConstPointer        ImageVectorConstPointer;
+  /** Typedef to describe the output image region type. */
+  typedef typename ImageVectorType::RegionType          ImageVectorRegionType;
+  /** Typedef to describe the type of pixel and point for output image. */
+  typedef typename ImageVectorType::PixelType           ImageVectorPixelType;
+  typedef typename ImageVectorType::PointType           ImageVectorPointType;
+  /** Typedef to describe the image index, size types and spacing for output image. */
+  typedef typename ImageVectorType::IndexType           ImageVectorIndexType;
+  typedef typename ImageVectorType::IndexValueType      ImageVectorIndexValueType;
+  typedef typename ImageVectorType::SizeType            ImageVectorSizeType;
+  typedef typename ImageVectorType::SizeValueType       ImageVectorSizeValueType;
+  typedef typename ImageVectorType::SpacingType         ImageVectorSpacingType;
+  typedef typename ImageVectorType::SpacingValueType    ImageVectorSpacingValueType;
+  
+  /** Typedef to image output type : otb::VectorImage  */
+  typedef TImageOut                                  ImageOutType;
+  /** Typedef to describe the output image pointer type. */
+  typedef typename ImageOutType::Pointer             ImageOutPointer;
+  typedef typename ImageOutType::ConstPointer        ImageOutConstPointer;
+  /** Typedef to describe the output image region type. */
+  typedef typename ImageOutType::RegionType          ImageOutRegionType;
+  /** Typedef to describe the type of pixel and point for output image. */
+  typedef typename ImageOutType::PixelType           ImageOutPixelType;
+  typedef typename ImageOutType::PointType           ImageOutPointType;
+  /** Typedef to describe the image index, size types and spacing for output image. */
+  typedef typename ImageOutType::IndexType           ImageOutIndexType;
+  typedef typename ImageOutType::IndexValueType      ImageOutIndexValueType;
+  typedef typename ImageOutType::SizeType            ImageOutSizeType;
+  typedef typename ImageOutType::SizeValueType       ImageOutSizeValueType;
+  typedef typename ImageOutType::SpacingType         ImageOutSpacingType;
+  typedef typename ImageOutType::SpacingValueType    ImageOutSpacingValueType;
+
+  /** Typedef to dem type : otb::Image  */
+  typedef TImageDEM                                  ImageDEMType;
+  /** Typedef to describe the output image pointer type. */
+  typedef typename ImageDEMType::Pointer             ImageDEMPointer;
+  typedef typename ImageDEMType::ConstPointer        ImageDEMConstPointer;
+  typedef typename ImageDEMType::SpacingType         ImageDEMSpacingType;
+  typedef typename ImageDEMType::SpacingValueType    ImageDEMSpacingValueType;
+
+  // DEMHandler
+  typedef otb::DEMHandler DEMHandlerType;
+  typedef typename DEMHandlerType::Pointer DEMHandlerPointerType;
+
+  
+  // Define Point2DType and Point3DType
+  using Point2DType = itk::Point<double,2>;
+  using Point3DType = itk::Point<double,3>;
+  typedef double  ValueType;
+
+  // Typedef for iterators
+  typedef itk::ImageScanlineConstIterator< ImageVectorType > InputIterator;
+  typedef itk::ImageScanlineIterator< ImageOutType > OutputIterator;
+
+  
+    
+  // Setter/Getter for inputs (Cartesian mean and GroupedByOrtho) 
+  /** Connect one of the operands for interferometry : Master */
+  void SetCartesianMeanInput( const ImageVectorType* image);
+
+  /** Connect one of the operands for interferometry : Coregistrated Slave */
+  void SetCompensatedComplexInput(const ImageVectorType* image);
+
+  /** Get the inputs */
+  const ImageVectorType* GetCartesianMeanInput() const;
+  const ImageVectorType* GetCompensatedComplexInput() const;
+
+  
+  // Setter for metadata
+  void SetSARImageKeyWorList(ImageKeywordlist sarImageKWL);
+  void SetDEMImagePtr(ImageDEMPointer demPtr);
+
+protected:
+  // Constructor
+  SARGroupedByOrthoImageFilter();
+
+  // Destructor
+  virtual ~SARGroupedByOrthoImageFilter() ITK_OVERRIDE;
+
+  // Print
+  void PrintSelf(std::ostream & os, itk::Indent indent) const ITK_OVERRIDE;
+  
+  /** SARGroupedByOrthoImageFilter produces an image which are into output geometry. The differences between 
+   * output and input images can be the size of images and the dimensions. 
+   * As such, SARGroupedByOrthoImageFilter needs to provide an implementation for 
+   * GenerateOutputInformation() in order to inform the pipeline execution model. 
+   */ 
+  virtual void GenerateOutputInformation() ITK_OVERRIDE;
+
+  /** SARGroupedByOrthoImageFilter needs a input requested region that corresponds to the margin and shifts 
+   * into the requested region of our output requested region. The output requested region needs to be modified
+   * to construct as wanted tiles with input size.
+   * As such, TilesAnalysesImageFilter needs to provide an implementation for 
+   * GenerateInputRequestedRegion() in order to inform the pipeline execution model. 
+   * \sa ProcessObject::GenerateInputRequestedRegion() */
+  virtual void GenerateInputRequestedRegion() ITK_OVERRIDE;
+
+  /**
+   * OutputRegionToInputRegion returns the input SAR region 
+   */
+  ImageVectorRegionType OutputRegionToInputRegion(const ImageOutRegionType& outputRegion,
+						  bool & intoInputImage) const;
+
+
+  /** 
+   * SARGroupedByOrthoImageFilterr can be implemented as a multithreaded filter.
+   * Therefore, this implementation provides a ThreadedGenerateData() routine
+   * which is called for each processing thread. The main output image data is
+   * allocated automatically by the superclass prior to calling
+   * ThreadedGenerateData().  ThreadedGenerateData can only write to the
+   * portion of the output image specified by the parameter
+   * "outputRegionForThread"
+   *
+   * \sa ImageToImageFilter::ThreadedGenerateData(),
+   *     ImageToImageFilter::GenerateData() */
+  virtual void ThreadedGenerateData(const ImageOutRegionType& outputRegionForThread, 
+				    itk::ThreadIdType threadId) ITK_OVERRIDE;
+
+
+ private:
+  SARGroupedByOrthoImageFilter(const Self&); // purposely not implemented
+  void operator=(const Self &); // purposely not 
+
+  // Instance of SarSensorModelAdapter
+  SarSensorModelAdapter::Pointer m_SarSensorModelAdapter; 
+
+  // SAR Image KeyWorldList
+  ImageKeywordlist m_SarImageKwl;
+
+  ossimGeoidEgm96 * m_geoidEmg96;
+  
+  // DEM
+  ImageDEMConstPointer m_DEMPtr;
+
+  float m_Gain;
+
+  int m_nbLinesSAR;
+  int m_nbColSAR;
+
+};
+
+} // End namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbSARGroupedByOrthoImageFilter.txx"
+#endif
+
+
+
+#endif