ADD: factors label writer in main pipeline

app/otbNDVITimeSeriesIndices.cxx

@@ -12,6 +12,9 @@
 // Functors
 #include "otbNDVITimeSeriesFunctor.h"
 #include "otbUnaryFunctorImageFilter.h"
+#include "itkBinaryFunctorImageFilter.h"
+#include "otbBinaryFunctorImageFilter.h"
+#include "otbTernaryFunctorImageFilter.h"
 
 // LayerStack
 #include "otbImageListToVectorImageFilter.h"
@@ -23,9 +26,6 @@
 #include "otbStreamingResampleImageFilter.h"
 #include "itkLinearInterpolateImageFunction.h"
 
-// Binary functor
-#include "itkBinaryFunctorImageFilter.h"
-
 namespace otb
 {
 
@@ -80,6 +80,14 @@ public:
   typedef otb::functor::RainfallEstimatedNDVIResiduesFunctor<FloatVectorImageType::PixelType, FloatVectorImageType::PixelType> ResiduesFunctorType;
   typedef itk::BinaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType, FloatVectorImageType, ResiduesFunctorType> ResiduesFilterType;
 
+  /** typedefs for pearson correlation */
+  typedef otb::functor::PearsonCorrelationCoefficientFunctor<FloatVectorImageType::PixelType, FloatVectorImageType::PixelType, FloatVectorImageType::PixelType> PearsonCorrelationFunctorType;
+  typedef otb::BinaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType, FloatVectorImageType, PearsonCorrelationFunctorType> PearsonCorrelationFilterType;
+
+  /** typedefs for final labeling */
+  typedef otb::functor::FactorsLabelingFunctor<FloatVectorImageType::PixelType, UInt8VectorImageType::PixelType> FactorsLabelingFunctorType;
+  typedef otb::TernaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType, FloatVectorImageType, FloatVectorImageType, FactorsLabelingFunctorType> FactorsLabelingFilterType;
+
 private:
 
   /*
@@ -117,7 +125,7 @@ private:
 
     // Documentation
     SetDocName("NDVITimeSeriesIndices");
-    SetDocLongDescription("This application performs NDVITimeSeriesIndices of images");
+    SetDocLongDescription("This application implements the method of Leroux et al. (https://doi.org/10.1016/j.rse.2017.01.014)");
     SetDocLimitations("None");
     SetDocAuthors("Remi Cresson");
     SetDocSeeAlso(" ");
@@ -207,6 +215,9 @@ private:
     SetDefaultOutputPixelType("reslabel", ImagePixelType_uint8);
     MandatoryOff("reslabel");
 
+    // Output factors labels
+    AddParameter(ParameterType_OutputImage,  "factorslabel",   "Output image for factors labels");
+    SetDefaultOutputPixelType("factorslabel", ImagePixelType_uint8);
 
     AddRAMParameter();
 
@@ -309,6 +320,7 @@ private:
           << "and " << ndviDates.size() << " dates.");
 
     // Get the detection method
+    FloatVectorImageType::Pointer reducedNDVI;
     if (GetParameterInt("ndvi.reduce") == CUMUL)
       {
       m_NDVICumulReduceFilter = TSCumulatedRangeReduceFilterType::New();
@@ -318,18 +330,21 @@ private:
       m_NDVICumulReduceFilter->GetFunctor().SetDay2(GetParameterInt("ndvi.reduce.cumul.day2"));
       m_NDVICumulReduceFilter->GetFunctor().SetMonth1(GetParameterInt("ndvi.reduce.cumul.month1"));
       m_NDVICumulReduceFilter->GetFunctor().SetMonth2(GetParameterInt("ndvi.reduce.cumul.month2"));
+      reducedNDVI = m_NDVICumulReduceFilter->GetOutput();
       }
     else if (GetParameterInt("ndvi.reduce") == MAX)
       {
       m_NDVIMaxReduceFilter = TSMaxReduceFilterType::New();
       m_NDVIMaxReduceFilter->SetInput(m_NDVIConcatener->GetOutput());
       m_NDVIMaxReduceFilter->GetFunctor().SetDates(ndviDates);
+      reducedNDVI = m_NDVIMaxReduceFilter->GetOutput();
       }
     else if (GetParameterInt("ndvi.reduce") == AMPLITUDE)
       {
       m_NDVIAmplitudeReduceFilter = TSAmplitudeReduceFilterType::New();
       m_NDVIAmplitudeReduceFilter->SetInput(m_NDVIConcatener->GetOutput());
       m_NDVIAmplitudeReduceFilter->GetFunctor().SetDates(ndviDates);
+      reducedNDVI = m_NDVIAmplitudeReduceFilter->GetOutput();
       }
     else
       {
@@ -338,7 +353,7 @@ private:
 
     // NDVI Trend
     m_NDVITimeSeriesRegressionFilter = SlopeAndPValueFilterType::New();
-    m_NDVITimeSeriesRegressionFilter->SetInput(m_NDVIAmplitudeReduceFilter->GetOutput());
+    m_NDVITimeSeriesRegressionFilter->SetInput(reducedNDVI);
 
     // Output NDVI trend image
     if (HasValue("ndvitrend"))
@@ -437,7 +452,7 @@ private:
 
     // Compute residues
     m_ResiduesFilter = ResiduesFilterType::New();
-    m_ResiduesFilter->SetInput1(m_NDVIAmplitudeReduceFilter->GetOutput());
+    m_ResiduesFilter->SetInput1(reducedNDVI);
     m_ResiduesFilter->SetInput2(m_RFCumulReduceFilter->GetOutput());
 
     // Output residues
@@ -464,6 +479,19 @@ private:
       SetParameterOutputImage("reslabel", m_NDVIResiduesLabelFilter->GetOutput());
       }
 
+    // Correlation between NDVI and Rainfall
+    m_CorrelationFilter = PearsonCorrelationFilterType::New();
+    m_CorrelationFilter->SetInput1(reducedNDVI);
+    m_CorrelationFilter->SetInput2(m_RFCumulReduceFilter->GetOutput());
+
+    // Factors labeling
+    m_FactorsLabelFilter = FactorsLabelingFilterType::New();
+    m_FactorsLabelFilter->SetInput1(m_NDVITimeSeriesRegressionFilter->GetOutput());
+    m_FactorsLabelFilter->SetInput2(m_NDVIResiduesRegressionFilter->GetOutput());
+    m_FactorsLabelFilter->SetInput3(m_CorrelationFilter->GetOutput());
+
+    // Write factors image
+    SetParameterOutputImage("factorslabel", m_FactorsLabelFilter->GetOutput());
 
   }   // DOExecute()
 
@@ -499,6 +527,12 @@ private:
 
   // Filter for residues
   ResiduesFilterType::Pointer               m_ResiduesFilter;
+
+  // Pearson correlation filter
+  PearsonCorrelationFilterType::Pointer     m_CorrelationFilter;
+
+  // Final labeling filter
+  FactorsLabelingFilterType::Pointer        m_FactorsLabelFilter;
 };
 }
 }

include/otbNDVITimeSeriesFunctor.h

@@ -855,9 +855,12 @@ public:
 
   PearsonCorrelationCoefficientFunctor(){
     m_CorrelationNoDataValue = itk::NumericTraits<CorrelationPixelValueType>::NonpositiveMin();
+    numberOfComponentsPerPixel = 1;
   }
   ~PearsonCorrelationCoefficientFunctor(){}
 
+  unsigned int GetOutputSize(){ return numberOfComponentsPerPixel; }
+
   // Purposely not implemented
   bool operator!=( const PearsonCorrelationCoefficientFunctor & ) const { return false; }
   bool operator==( const PearsonCorrelationCoefficientFunctor & other ) const { return !(*this != other); }
@@ -940,6 +943,8 @@ protected:
   NDVIPixelValueType       m_NDVINoDataValue;
   RainfallPixelValueType   m_RainfallNoDataValue;
   CorrelationPixelValueType m_CorrelationNoDataValue;
+
+  unsigned int numberOfComponentsPerPixel;
 }; // PearsonCorrelationCoefficientFunctor
 
 } // namespace functor