diff --git a/Modules/Core/Functor/include/otbFunctorImageFilter.h b/Modules/Core/Functor/include/otbFunctorImageFilter.h
index 7cb6096fcd666d39451853eb28697c1fb8535e56..7aeb4e7efc3f9ae35d00bb1eb2c3151f5f3d0640 100644
--- a/Modules/Core/Functor/include/otbFunctorImageFilter.h
+++ b/Modules/Core/Functor/include/otbFunctorImageFilter.h
@@ -141,8 +141,6 @@ template <class T> struct PixelTypeDeduction<T &>
   using PixelType = typename PixelTypeDeduction<T>::PixelType;
 };
 
-  
-
 /** 
  * \struct ImageTypeDeduction
  * \brief Helper struct to derive ImageType from template parameter
@@ -178,54 +176,92 @@ template <class T> struct ImageTypeDeduction<itk::VariableLengthVector<T>>
 
 template <typename T, typename TNameMap> struct FunctorFilterSuperclassHelper : public FunctorFilterSuperclassHelper<decltype(&T::operator()),TNameMap> {};
 
-/// Partial specialisation for R(*)(T...)
-template <typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<R(*)(T...),TNameMap>
+namespace functor_filter_details
 {
+template <typename R, typename TNameMap, typename...T> struct FunctorFilterSuperclassHelperImpl
+{
+  // OutputImageType is derived from return type R
   using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
+  // InputImageType is derived using pixel type deduction and image
+  // type deduction
   template <typename V> using InputImageType = typename ImageTypeDeduction<typename PixelTypeDeduction<V>::PixelType>::ImageType;
-  
+
+  // Filter type is either VariadicInputsImageFilter or
+  // VariadicNamedInputsImageFilter depending on if there is a
+  // TNameMap or not
   using FilterType = typename std::conditional<std::is_void<TNameMap>::value,
                                                VariadicInputsImageFilter<OutputImageType,InputImageType<T>...>,
                                                VariadicNamedInputsImageFilter<OutputImageType,TNameMap,InputImageType<T>...>>::type;
 
+  // InputHasNeighborhood is derived from IsNeighborhood
   using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
 };
+} // End namespace functor_filter_details
+
+/// Partial specialisation for R(*)(T...)
+template <typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<R(*)(T...),TNameMap>
+{
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
+};
 
 /// Partial specialisation for R(C::*)(T...) const
 template <typename C, typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<R(C::*)(T...) const, TNameMap>
 {
-  using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
-  template <typename V> using InputImageType = typename ImageTypeDeduction<typename PixelTypeDeduction<V>::PixelType>::ImageType;
-  
-  using FilterType = typename std::conditional<std::is_void<TNameMap>::value,
-                                               VariadicInputsImageFilter<OutputImageType,InputImageType<T>...>,
-                                               VariadicNamedInputsImageFilter<OutputImageType,TNameMap,InputImageType<T>...>>::type;
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
 
-  using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
 };
 
 /// Partial specialisation for R(C::*)(T...)
 template <typename C, typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<R(C::*)(T...), TNameMap>
-{ 
-  using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
-  template <typename V> using InputImageType = typename ImageTypeDeduction<typename PixelTypeDeduction<V>::PixelType>::ImageType;
-  
-  using FilterType = typename std::conditional<std::is_void<TNameMap>::value,
-                                               VariadicInputsImageFilter<OutputImageType,InputImageType<T>...>,
-                                               VariadicNamedInputsImageFilter<OutputImageType,TNameMap,InputImageType<T>...>>::type;
+{
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
 
-  using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
 };
 
+/// Partial specialisation for void(*)(R &,T...)
+template <typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<void(*)(R&, T...), TNameMap>
+{
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
+
+};
+
+/// Partial specialisation for void(C::*)(R&,T...) const
+template <typename C, typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<void(C::*)(R&,T...) const, TNameMap>
+{
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
+
+};
+
+/// Partial specialisation for void(C::*)(R&,T...)
+template <typename C, typename R, typename... T, typename TNameMap> struct FunctorFilterSuperclassHelper<void(C::*)(R&,T...), TNameMap>
+{ 
+  using OutputImageType      = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::OutputImageType;
+  using FilterType           = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::FilterType;
+  using InputHasNeighborhood = typename functor_filter_details::FunctorFilterSuperclassHelperImpl<R,TNameMap,T...>::InputHasNeighborhood;
+};
+
+
+
 /**
  * \brief This helper method builds a fully functional FunctorImageFilter from a functor instance
  * 
- * Functor can be any operator() that matches the following:
+ * Functor can be any operator() (const or non-const) that matches the following:
  * - Accepts any number of arguments of T,
  * (const) itk::VariableLengthVector<T> (&),(const)
  * itk::Neighborhood<T> (&), (const)
  * itk::Neighborhood<itk::VariableLengthVector<T>> (&) with T a scalar type
  * - returns T or itk::VariableLengthVector<T>, with T a scalar type
+ * or returns void and has first parameter as output (i.e. T& or itk::VariableLengthVector<T>&)
  *
  * The returned filter is ready to use. Inputs can be set through the
  * SetVariadicInputs() method (see VariadicInputsImageFilter class for
@@ -248,6 +284,16 @@ template <typename Functor, typename TNameMap = void> auto NewFunctorFilter(cons
 /** \class FunctorImageFilter
  * \brief A generic functor filter templated by its functor
  * 
+ * TFunction can be any operator() (const or non-const) that matches the following:
+ * - Accepts any number of arguments of T,
+ * (const) itk::VariableLengthVector<T> (&),(const)
+ * itk::Neighborhood<T> (&), (const)
+ * itk::Neighborhood<itk::VariableLengthVector<T>> (&) with T a scalar type
+ * - returns T or itk::VariableLengthVector<T>, with T a scalar type
+ * or returns void and has first parameter as output (i.e. T& or itk::VariableLengthVector<T>&)
+ *
+ * All image types will be deduced from the TFunction operator().
+ *
  * \sa VariadicInputsImageFilter
  * \sa NewFunctorFilter
  * 
@@ -272,8 +318,6 @@ public:
   using Superclass = typename SuperclassHelper::FilterType;
   using OutputImageType = typename Superclass::OutputImageType;
   using OutputImageRegionType = typename OutputImageType::RegionType;
-  
-  using ProcessObjectType = itk::ProcessObject;
 
   // A tuple of bool of the same size as the number of arguments in
   // the functor
@@ -283,7 +327,7 @@ public:
   using Superclass::NumberOfInputs;
   
   /** Run-time type information (and related methods). */
-  itkTypeMacro(FunctorImageFilter, ImageToImageFilter);
+  itkTypeMacro(FunctorImageFilter, VariadicInputsImageFilter);
   
   /** Get the functor object.
    * 
diff --git a/Modules/Core/Functor/include/otbFunctorImageFilter.hxx b/Modules/Core/Functor/include/otbFunctorImageFilter.hxx
index e1197d3c192f342ec03a3f10bcd4e82538c39da7..897839740023f9f3d9678b168239639db847d4ac 100644
--- a/Modules/Core/Functor/include/otbFunctorImageFilter.hxx
+++ b/Modules/Core/Functor/include/otbFunctorImageFilter.hxx
@@ -138,16 +138,67 @@ template <typename T> struct GetProxy<itk::ConstNeighborhoodIterator<T> >
 }
 };
 
+template <class Oper> struct OperProxy : public OperProxy<decltype(&Oper::operator())> {};
+
+template<class Out, class ... In> struct OperProxy<Out(*)(In...)>
+{
+  template <class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    out = oper(in...);
+  }
+};
+
+template<class C, class Out, class ... In> struct OperProxy<Out(C::*)(In...)>
+{
+  template<class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    out = oper(in...);
+  }
+};
+
+template<class C, class Out, class ... In> struct OperProxy<Out(C::*)(In...) const>
+{
+  template<class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    out = oper(in...);
+  }
+};
+
+template<class Out, class ... In> struct OperProxy<void(*)(Out&, In...)>
+{
+  template<class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    oper(out,in...);
+  }
+};
+
+template<class C, class Out, class ... In> struct OperProxy<void(C::*)(Out&, In...)>
+{
+  template<class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    oper(out,in...);
+  }
+};
+
+template<class C, class Out, class ... In> struct OperProxy<void(C::*)(Out&, In...) const>
+{
+  template template<class Oper> static void Compute(Oper& oper, Out& out, const In& ... in)
+  {
+    oper(out,in...);
+  }
+};
+
+
 // Will be easier to write in c++17 with std::apply and fold expressions
-template <class Tuple, class Oper, size_t...Is> auto CallOperatorImpl(Tuple& t, Oper & oper,std::index_sequence<Is...>)
+template <class Tuple, class Out, class Oper, size_t...Is> auto CallOperatorImpl(Tuple& t, Out & out, Oper & oper,std::index_sequence<Is...>)
 {
-  return oper(GetProxy<typename std::remove_reference<decltype(std::get<Is>(t))>::type>::Get(std::get<Is>(t))...);
+  OperProxy<Oper>::Compute(oper,out,GetProxy<typename std::remove_reference<decltype(std::get<Is>(t))>::type>::Get(std::get<Is>(t))...);
 }
 
 // Will be easier to write in c++17 with std::apply and fold expressions
-template <class Oper, typename ... Args> auto CallOperator(Oper& oper, std::tuple<Args...> & t)
+template <class Out, class Oper, typename ... Args> auto CallOperator(Out & out, Oper& oper, std::tuple<Args...> & t)
 {
-  return CallOperatorImpl(t,oper,std::make_index_sequence<sizeof...(Args)>{});
+  CallOperatorImpl(t,out,oper,std::make_index_sequence<sizeof...(Args)>{});
 }
 
 // Variadic move of iterators
@@ -228,13 +279,22 @@ FunctorImageFilter<TFunction, TNameMap>
   itk::ImageScanlineIterator<OutputImageType> outIt(this->GetOutput(),outputRegionForThread);
   itk::ProgressReporter p(this,threadId,outputRegionForThread.GetNumberOfPixels());
 
+  // This will build a tuple of iterators to be used
   auto inputIterators = functor_filter_details::MakeIterators(this->GetVariadicInputs(),outputRegionForThread, m_Radius,InputHasNeighborhood{});
+
+  // Build a default value
+  typename OutputImageType::PixelType outputValueHolder;
+  itk::NumericTraits<typename OutputImageType::PixelType>::SetLength(outputValueHolder,this->GetOutput()->GetNumberOfComponentsPerPixel());
   
   while(!outIt.IsAtEnd())
     {
+    // MoveIterartors will ++ all iterators in the tuple
     for(;!outIt.IsAtEndOfLine();++outIt,functor_filter_details::MoveIterators(inputIterators))
       {
-      outIt.Set(functor_filter_details::CallOperator(m_Functor,inputIterators));
+      // This will call the operator with inputIterators Get() results
+      // and fill outputValueHolder with the result.
+      functor_filter_details::CallOperator(outputValueHolder,m_Functor,inputIterators);
+      outIt.Set(outputValueHolder);
       // Update progress
       p.CompletedPixel();
       }
diff --git a/Modules/Core/Functor/include/otbVariadicAddFunctor.h b/Modules/Core/Functor/include/otbVariadicAddFunctor.h
index 61d687146c253d89ee6aaf12f641b2be9cdbfd81..cc0030d21513864a1d9aeff0f278bb7f2888fbe5 100644
--- a/Modules/Core/Functor/include/otbVariadicAddFunctor.h
+++ b/Modules/Core/Functor/include/otbVariadicAddFunctor.h
@@ -28,7 +28,10 @@ namespace otb
 
 namespace Functor
 {
-
+/**
+ * \class VariadicAdd
+ * \brief This functor adds any number of compile time scalar inputs 
+ */
 template <typename TOut, typename ...TIns> struct VariadicAdd
 {
   auto operator()(TIns... ins) const
diff --git a/Modules/Core/Functor/include/otbVariadicConcatenateFunctor.h b/Modules/Core/Functor/include/otbVariadicConcatenateFunctor.h
index d1bb0c99808b6a8120699849c53c3fa900085536..4cff7f504b61e950473f5ce1f328e6aac43b03fd 100644
--- a/Modules/Core/Functor/include/otbVariadicConcatenateFunctor.h
+++ b/Modules/Core/Functor/include/otbVariadicConcatenateFunctor.h
@@ -33,7 +33,7 @@ namespace Functor
 {
 
 
-namespace internal
+namespace variadic_concatenate_details
 {
 // helper function to implement next functor (convert a scalar value
 // to a VariableLengthVector)
@@ -70,17 +70,23 @@ template <typename v1, typename v2, typename ...vn> void concatenateVectors(v1 &
   concatenateVectors(a,b);
   concatenateVectors(a,z...);
 }
-} // end namespace internal
+} // end namespace variadic_concatenate_details
 
 // N  images (all types) -> vector image
 // This functor concatenates N images (N = variadic) of type
 // VectorImage and or Image, into a single VectorImage
+
+/** 
+ * \class VariadicConcatenate
+ * \brief This functor concatenates any number of input of scalar type
+ * or VariableLengthVector.
+ */
 template<typename TOut, typename ...TIns> struct VariadicConcatenate
 {
   auto operator()(const TIns &...  ins) const
   {
     itk::VariableLengthVector<TOut> out;
-    internal::concatenateVectors(out, internal::toVector(ins)...);
+    variadic_concatenate_details::concatenateVectors(out, internal::toVector(ins)...);
     
     return out;
   }
diff --git a/Modules/Core/Functor/test/otbFunctorImageFilter.cxx b/Modules/Core/Functor/test/otbFunctorImageFilter.cxx
index dd09fe05a53289dc0d2fa628803d153ddd8be6ea..67e418e845e0bfc3f3e292c23abba8e4e20934a4 100644
--- a/Modules/Core/Functor/test/otbFunctorImageFilter.cxx
+++ b/Modules/Core/Functor/test/otbFunctorImageFilter.cxx
@@ -91,6 +91,19 @@ template <typename T> struct TypesCheck
     }
   };
 
+template <typename TOut,typename TIn> struct TestOperatorVoidReturn
+  {
+    void operator()(TOut& out,const TIn&) const
+    {
+      out = TOut(OutputSize());
+    }
+    
+    constexpr size_t OutputSize(...) const
+    {
+      return 1;
+    }
+  };
+
   // Fake test operator non const
   template <typename TOut,typename TIn> struct TestOperatorNonConst
   {
@@ -106,7 +119,7 @@ template <typename T> struct TypesCheck
     }
   };
 
-    
+  
   template <typename TOut, typename TIn> void TestFilter()
   {
   // Deduce types
@@ -145,13 +158,23 @@ template <typename T> struct TypesCheck
   res = filter1->template GetVariadicNamedInput<tag>();
   filter1->Update();
   
+  // Test with void return
+  auto functorWithVoidReturn = TestOperatorVoidReturn<TOut,TIn>{};
+  auto filterWithVoidReturn = NewFunctorFilter(functorWithVoidReturn);
+  using FilterWithVoidReturnType = typename decltype(filter)::ObjectType;
+  static_assert(FilterWithVoidReturnType::NumberOfInputs == 1,"");
+  static_assert(std::is_same<typename FilterWithVoidReturnType::template InputImageType<0>, InputImageType>::value, "");
+  
   // Test with non const operator
   auto functorNonConstOperator = TestOperatorNonConst<TOut,TIn>{};
   auto filterWithNonConstOperator = NewFunctorFilter(functorNonConstOperator);
   filterWithNonConstOperator->SetInput1(in);
-  filterWithNonConstOperator->Update();
-  
-  
+  filterWithNonConstOperator->Update();  
+
+  filterWithVoidReturn->SetVariadicInputs(in);
+  filterWithVoidReturn->SetInput1(in);
+  filterWithVoidReturn->template SetVariadicInput<0>(in); // template keyword to avoid C++ parse ambiguity
+  filterWithVoidReturn->Update();
   
   // Test with simple lambda
   auto lambda = [] (const TIn &)