diff --git a/app/otbTimeSeriesIndexTrend.cxx b/app/otbTimeSeriesIndexTrend.cxx
index f9fa4cca0ba3cb1a25d2d7966281d4a7acf1860d..6ee67a908e3550aa5ac72b8b8266fc9cf60524db 100644
--- a/app/otbTimeSeriesIndexTrend.cxx
+++ b/app/otbTimeSeriesIndexTrend.cxx
@@ -26,10 +26,7 @@
 
 // Functors
 #include "otbNDVITimeSeriesFunctor.h"
-#include "otbUnaryFunctorImageFilter.h"
-#include "itkBinaryFunctorImageFilter.h"
-#include "otbBinaryFunctorImageFilter.h"
-#include "otbTernaryFunctorImageFilter.h"
+#include "otbFunctorImageFilter.h"
 
 // LayerStack
 #include "otbImageListToVectorImageFilter.h"
@@ -78,24 +75,19 @@ public:
       FloatVectorImageType::PixelType>                                                 TSCumulatedRangeReduceFunctorType;
   typedef otb::functor::TSMaxReduceFunctor<FloatVectorImageType::PixelType,
       FloatVectorImageType::PixelType>                                                 TSMaxReduceFunctorType;
-  typedef otb::UnaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      TSAmplitudeReduceFunctorType>                                                    TSAmplitudeReduceFilterType;
-  typedef otb::UnaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      TSCumulatedRangeReduceFunctorType>                                               TSCumulatedRangeReduceFilterType;
-  typedef otb::UnaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      TSMaxReduceFunctorType>                                                          TSMaxReduceFilterType;
+  typedef otb::FunctorImageFilter<TSAmplitudeReduceFunctorType>                        TSAmplitudeReduceFilterType;
+  typedef otb::FunctorImageFilter<TSCumulatedRangeReduceFunctorType>                   TSCumulatedRangeReduceFilterType;
+  typedef otb::FunctorImageFilter<TSMaxReduceFunctorType>                              TSMaxReduceFilterType;
 
   /** typedef for trend computing (slope, p-value)  */
   typedef otb::functor::SlopeAndPValueFunctor<FloatVectorImageType::PixelType,
       FloatVectorImageType::PixelType>                                                 SlopeAndPValueFunctorType;
-  typedef otb::UnaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      SlopeAndPValueFunctorType>                                                       SlopeAndPValueFilterType;
+  typedef otb::FunctorImageFilter<SlopeAndPValueFunctorType>                           SlopeAndPValueFilterType;
 
   /** typedefs for trend labeling */
   typedef otb::functor::SlopeAndPValueLabelingFunctor<FloatVectorImageType::PixelType,
       UInt8VectorImageType::PixelType>                                                 SlopeAndPValueLabelFunctorType;
-  typedef otb::UnaryFunctorImageFilter<FloatVectorImageType, UInt8VectorImageType,
-      SlopeAndPValueLabelFunctorType>                                                  SlopeAndPValueLabelFilterType;
+  typedef otb::FunctorImageFilter<SlopeAndPValueLabelFunctorType>                      SlopeAndPValueLabelFilterType;
 
   /* typedefs for dates*/
   typedef otb::dates::SingleDate                                                       DateType;
@@ -109,20 +101,17 @@ public:
   /* typedefs for residues processing */
   typedef otb::functor::RainfallEstimatedNDVIResiduesFunctor<FloatVectorImageType::PixelType,
       FloatVectorImageType::PixelType>                                                 ResiduesFunctorType;
-  typedef itk::BinaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      FloatVectorImageType, ResiduesFunctorType>                                       ResiduesFilterType;
+  typedef otb::FunctorImageFilter<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;
+  typedef otb::FunctorImageFilter<PearsonCorrelationFunctorType>                       PearsonCorrelationFilterType;
 
   /** typedefs for final labeling */
   typedef otb::functor::FactorsLabelingFunctor<FloatVectorImageType::PixelType,
       UInt8VectorImageType::PixelType>                                                 FactorsLabelingFunctorType;
-  typedef otb::TernaryFunctorImageFilter<FloatVectorImageType, FloatVectorImageType,
-      FloatVectorImageType, FloatVectorImageType, FactorsLabelingFunctorType>          FactorsLabelingFilterType;
+  typedef otb::FunctorImageFilter<FactorsLabelingFunctorType>                          FactorsLabelingFilterType;
 
 private:
 
@@ -311,7 +300,6 @@ private:
     m_NDVIConcatener = ListConcatenerFilterType::New();
     m_NDVIImageList = ImageListType::New();
     m_NDVIExtractorList = ExtractROIFilterListType::New();
-
     // Split each input vector image into image
     // and generate an mono channel image list
     inNDVIList->GetNthElement(0)->UpdateOutputInformation();
@@ -349,24 +337,24 @@ private:
       {
       m_NDVICumulReduceFilter = TSCumulatedRangeReduceFilterType::New();
       m_NDVICumulReduceFilter->SetInput(m_NDVIConcatener->GetOutput());
-      m_NDVICumulReduceFilter->GetFunctor().SetDates(ndviDates);
-      m_NDVICumulReduceFilter->GetFunctor().SetStartDay(GetParameterInt  ("ndvi.reduce.cumul.startday"));
-      m_NDVICumulReduceFilter->GetFunctor().SetStartMonth(GetParameterInt("ndvi.reduce.cumul.startmonth"));
-      m_NDVICumulReduceFilter->GetFunctor().SetNbOfDays(GetParameterInt  ("ndvi.reduce.cumul.nbdays"));
+      m_NDVICumulReduceFilter->GetModifiableFunctor().SetDates(ndviDates);
+      m_NDVICumulReduceFilter->GetModifiableFunctor().SetStartDay(GetParameterInt  ("ndvi.reduce.cumul.startday"));
+      m_NDVICumulReduceFilter->GetModifiableFunctor().SetStartMonth(GetParameterInt("ndvi.reduce.cumul.startmonth"));
+      m_NDVICumulReduceFilter->GetModifiableFunctor().SetNbOfDays(GetParameterInt  ("ndvi.reduce.cumul.nbdays"));
       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);
+      m_NDVIMaxReduceFilter->GetModifiableFunctor().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);
+      m_NDVIAmplitudeReduceFilter->GetModifiableFunctor().SetDates(ndviDates);
       reducedNDVI = m_NDVIAmplitudeReduceFilter->GetOutput();
       }
     else
@@ -377,7 +365,7 @@ private:
     // NDVI Trend
     m_NDVITimeSeriesRegressionFilter = SlopeAndPValueFilterType::New();
     m_NDVITimeSeriesRegressionFilter->SetInput(reducedNDVI);
-    m_NDVITimeSeriesRegressionFilter->GetFunctor().SetDates(ndviDates);
+    m_NDVITimeSeriesRegressionFilter->GetModifiableFunctor().SetDates(ndviDates);
 
     // Write reduced ndvi
     if (HasValue("ndvicumul"))
@@ -444,7 +432,6 @@ private:
           otbAppLogFATAL("Number of rainfall images and number of dates are different! There is "
               << (m_RFConcatener->GetOutput()->GetNumberOfComponentsPerPixel()) << " images "
               << "and " << rainDates.size() << " dates.");
-
         // Check dates consistency
         int firstYearForNDVI = ndviDates[0].year;
         int lastYearForNDVI = ndviDates[ndviDates.size()-1].year;
@@ -462,7 +449,6 @@ private:
                 << " NDVI last year " << lastYearForNDVI
                 << "and rainfall last year " << lastYearForRF << ".");
           }
-
         // Computes new dates of the reduced time series
         std::vector<DateType> reducedTimeSeriesDates;
         for (int year = firstYearForNDVI; year <= lastYearForNDVI; year++)
@@ -474,11 +460,10 @@ private:
         // Reduce rainfall
         m_RFCumulReduceFilter = TSCumulatedRangeReduceFilterType::New();
         m_RFCumulReduceFilter->SetInput(m_RFConcatener->GetOutput());
-        m_RFCumulReduceFilter->GetFunctor().SetDates(rainDates);
-        m_RFCumulReduceFilter->GetFunctor().SetStartDay(GetParameterInt  ("rain.reduce.cumul.startday"));
-        m_RFCumulReduceFilter->GetFunctor().SetStartMonth(GetParameterInt("rain.reduce.cumul.startmonth"));
-        m_RFCumulReduceFilter->GetFunctor().SetNbOfDays(GetParameterInt  ("rain.reduce.cumul.nbdays"));
-
+        m_RFCumulReduceFilter->GetModifiableFunctor().SetDates(rainDates);
+        m_RFCumulReduceFilter->GetModifiableFunctor().SetStartDay(GetParameterInt  ("rain.reduce.cumul.startday"));
+        m_RFCumulReduceFilter->GetModifiableFunctor().SetStartMonth(GetParameterInt("rain.reduce.cumul.startmonth"));
+        m_RFCumulReduceFilter->GetModifiableFunctor().SetNbOfDays(GetParameterInt  ("rain.reduce.cumul.nbdays"));
         // Resample rainfall
         LinearInterpolatorType::Pointer linInterpolator = LinearInterpolatorType::New();
         NNInterpolatorType::Pointer nnInterpolator = NNInterpolatorType::New();
@@ -505,7 +490,6 @@ private:
           {
           SetParameterOutputImage("rainfcumul", m_ResampleFilter->GetOutput());
           }
-
         // Compute residues
         m_ResiduesFilter = ResiduesFilterType::New();
         m_ResiduesFilter->SetInput1(reducedNDVI);
@@ -516,18 +500,16 @@ private:
           {
             SetParameterOutputImage("residues", m_ResiduesFilter->GetOutput());
           }
-
         // Regression (residues)
         m_NDVIResiduesRegressionFilter = SlopeAndPValueFilterType::New();
         m_NDVIResiduesRegressionFilter->SetInput(m_ResiduesFilter->GetOutput());
-        m_NDVIResiduesRegressionFilter->GetFunctor().SetDates(reducedTimeSeriesDates);
+        m_NDVIResiduesRegressionFilter->GetModifiableFunctor().SetDates(reducedTimeSeriesDates);
 
         // Output residues trend
         if (HasValue("restrend"))
           {
             SetParameterOutputImage("restrend", m_NDVIResiduesRegressionFilter->GetOutput());
           }
-
         // Output residues trend labels
         if (HasValue("reslabel"))
           {
@@ -535,7 +517,6 @@ private:
             m_NDVIResiduesLabelFilter->SetInput(m_NDVIResiduesRegressionFilter->GetOutput());
             SetParameterOutputImage("reslabel", m_NDVIResiduesLabelFilter->GetOutput());
           }
-
         // Correlation between NDVI and Rainfall
         m_CorrelationFilter = PearsonCorrelationFilterType::New();
         m_CorrelationFilter->SetInput1(reducedNDVI);
@@ -545,7 +526,6 @@ private:
           {
           SetParameterOutputImage("pearsoncoef", m_CorrelationFilter->GetOutput());
           }
-
         // Factors labeling
         m_FactorsLabelFilter = FactorsLabelingFilterType::New();
         m_FactorsLabelFilter->SetInput1(m_NDVITimeSeriesRegressionFilter->GetOutput());
@@ -557,7 +537,6 @@ private:
           {
             SetParameterOutputImage("factorslabel", m_FactorsLabelFilter->GetOutput());
           }
-
       }
     else
       {
diff --git a/include/otbNDVITimeSeriesFunctor.h b/include/otbNDVITimeSeriesFunctor.h
index e824619ae6d1ce3e8554013eea5d4bf84cd6db96..5066d1f01d123d9676ef79bb3a2c051666d432f3 100644
--- a/include/otbNDVITimeSeriesFunctor.h
+++ b/include/otbNDVITimeSeriesFunctor.h
@@ -20,6 +20,7 @@
 #include "vnl/vnl_vector.h"
 #include <vector>
 #include "otbDates.h"
+#include <numeric>
 
 namespace otb
 {
@@ -92,6 +93,10 @@ public:
 
   inline unsigned int GetOutputSize() const { return nbOfYears; }
 
+  constexpr std::size_t OutputSize(std::array<long unsigned int, 1ul>&) const
+  {
+    return nbOfYears;
+  }
   // Set / Get input no-data value
   void SetInputNoDataValue(PixelValueType value) { m_InputNoDataValue = value; }
   PixelValueType GetInputNoDataValue() const { return m_InputNoDataValue; }
@@ -362,6 +367,13 @@ public:
   // Get output size
   inline unsigned int GetOutputSize(){ return numberOfComponentsPerPixel; }
 
+  static constexpr std::size_t outputPixelSize{1};
+
+  std::size_t OutputSize(std::array<long unsigned int, 1ul>&) const
+  {
+    return numberOfComponentsPerPixel;
+  }
+
   // Compute output pixel
   inline TOutputPixel operator ()(const TInputPixel& input) const
   {
@@ -453,7 +465,6 @@ public:
 protected:
   InputPixelValueType m_InputNoDataValue;
   OutputPixelValueType m_OutputNoDataValue;
-
   unsigned int numberOfComponentsPerPixel;
 }; // SlopeAndPValueFunctor
 
@@ -510,6 +521,13 @@ public:
   // Get output size
   inline unsigned int GetOutputSize(){ return numberOfComponentsPerPixel; }
 
+  static constexpr std::size_t outputPixelSize{1};
+
+  std::size_t OutputSize(std::array<long unsigned int, 1ul>&) const
+  {
+    return numberOfComponentsPerPixel;
+  }
+
   // Compute output pixel
   inline TLabelPixel operator ()(const TInputPixel& input) const
   {
@@ -591,6 +609,12 @@ public:
   }
   ~RainfallEstimatedNDVIResiduesFunctor(){}
 
+  std::size_t OutputSize(std::array<long unsigned int, 2ul> Inputs) const
+  {
+    unsigned int nbComp = std::accumulate(Inputs.begin(),Inputs.end(),0)/2;
+    return nbComp;
+  }
+
   // Purposely not implemented
   bool operator!=( const RainfallEstimatedNDVIResiduesFunctor & ) const { return false; }
   bool operator==( const RainfallEstimatedNDVIResiduesFunctor & other ) const { return !(*this != other); }
@@ -707,7 +731,7 @@ public:
 
     return residues;
   }
-
+  // void SetStartDay(const TNDVIPixel& inputNDVIPixel) { numberOfComponentsPerPixel = inputNDVIPixel.Size(); }
 protected:
   NDVIPixelValueType       m_NDVINoDataValue;
   RainfallPixelValueType   m_RainfallNoDataValue;
@@ -774,6 +798,13 @@ public:
   // Get output size
   inline unsigned int GetOutputSize(){ return numberOfComponentsPerPixel; }
 
+  static constexpr std::size_t outputPixelSize{1};
+
+  std::size_t OutputSize(std::array<long unsigned int, 3ul>&) const
+  {
+    return numberOfComponentsPerPixel;
+  }
+
   // Compute output pixel
   inline TLabelPixel operator ()(const TInputPixel& ndviTrend, const TInputPixel& rainTrend, const TInputPixel& corr) const
   {
@@ -892,6 +923,13 @@ public:
 
   unsigned int GetOutputSize(){ return numberOfComponentsPerPixel; }
 
+  static constexpr std::size_t outputPixelSize{1};
+
+  std::size_t OutputSize(std::array<long unsigned int, 2ul>&) const
+  {
+    return numberOfComponentsPerPixel;
+  }
+
   // Purposely not implemented
   bool operator!=( const PearsonCorrelationCoefficientFunctor & ) const { return false; }
   bool operator==( const PearsonCorrelationCoefficientFunctor & other ) const { return !(*this != other); }