WIP: A really, really generic functor filter

Summary

This MR provides THE ultimate functor filter. Let see why and how.

Note: this is a WIP, also the code is functional. It misses documentation and reviews. I am submitting this now because it already represents a huge amount of work and I want to get feedbacks on the idea before proceeding any further.

Rationale

We have a lot of base classes to define different cases of functor filters (with a VectorImage or Image output, with 2, 3 or more inputs, with neighborhood ...). Moreover, the very principle of functors is to be simple to use, but most often there is a lot of boiler plate code needed even for the simplest thing. This leads to BandMathFilter being preferred, which is a disaster in terms of maintenance and performances. So let see if we can do something really, really simple for the user.

The idea is that we could replace a LOT of code with this filter, only keeping the bare functor definitions, that will work with the current implementation.

Examples of use

This MR is rather complex, with a lot of meta-programming, sfinae and the rest, so I would like to present it from the user perspective. For the sake of the demonstration, image, image1, image2 ... will denote instances of otb::Image<T>, and vimage, vimage0, vimage1 ... will denote instances of otb::VectorImage<T>.

Lets start with a very simple example. Let say I have a scalar image, and want to multiply it by 10. I could of course use the well known MultiplyByScalarImageFilter (really?), but lets use a lambda:

double scale = 10.;
// dummy lambda that multiply by a scalar
auto lambda = [scale](double p)
{
  return scale*p;
};
// Create a filter that will accept an otb::Image<double> as input and produce an otb::Image<double>
auto filterLambda = NewFunctorFilter(lambda);
// Set the inputs using variadic SetVInputs
filterLambda->SetVInputs(image);
// All good, lets run it
filterLambda->Update();

So what happens under the hood? Well NewFunctorFilter calls the FunctorImageFilter method which analyses the operator() prototype and knows how to call it with the correct types. A lambda, 3 lines of code, and you have a fully functioning filter. Cool? Wait there is more.

// dummy lambda that sums 3 images with different types
auto lambda = [](double p, unsigned int p1, short p2)
{
  return static_cast<double>(p)+static_cast<double>(p1) + static_cast<double>(p2);
};
// Create a filter that will accept otb::Image<double>, otb::Image<unsigned int>,otb::Image<short>  as inputs and produce an otb::Image<double>
auto filterLambda = NewFunctorFilter(lambda);
// Set the inputs using variadic SetVinputs
filterLambda->SetVInputs(image, image0, image1); // image is otb::Image<double>, image0 is otb::Image<unsigned int> and image1 is otb::Image<short>
// All good, lets run it
filterLambda->Update();

So this is practically the same example as before, except that (you get it) it shows that you can use a lambda with any number of input arguments. It will build the corresponding filter, and you just have to set the inputs and run it. Pretty cool no? Wait there is more.

// lambda that returns the norm of the input vector pixel
auto lambda = [](const itk::VariableLengthVector<double> & inv)
{
  return inv.GetNorm();
};
// Create a filter that will accept otb::VectorImage<double>, as inputs and produce an otb::Image<double>
auto filterLambda = NewFunctorFilter(lambda);
// Set the inputs using variadic SetVinputs
filterLambda->SetVInputs(vimage); // vimage is otb::VectorImage<double>
// All good, lets run it
filterLambda->Update();

Yes, it supports VariableLengthVector as operator() arguments, and yes before you ask, you can combine any number of scalar and VariableLengthVector types as you like. But wait, there is more.

// lambda that performs neighborhood averaging
auto lambda = [](const itk::Neighborhood<short> & in)
{
   double out(0);
   for(auto it = in.Begin(); it!=in.End();++it)
      out+=static_cast<TOut>(*it);

   out/=in.Size(); 
   return out;
};
// Create a filter that will accept otb::Image<short>, as inputs and produce an otb::Image<double>
auto filterLambda = NewFunctorFilter(lambda,{{3,3}});
// Set the inputs using variadic SetVinputs
filterLambda->SetVInputs(image); // image is otb::VectorImage<short>
// All good, lets run it
filterLambda->Update();

Yes, you can use itk::Neighborhood<T>, and even itk::Neighborhood<itk::VariableLengthVector<T>>. In this case, you just have to provide the radius to NewFunctorFilter. But wait, there is more.

This works with lambda but also with any class that has an operator(). For instance, the following is a variadic functor that can add any (compile time defined) number of images.

template <typename O, typename ...T> struct VariadicAdd
{
  auto operator()(T... ins) const
  {
    std::vector<O> outVector{static_cast<O>(ins)...};
    return std::accumulate(outVector.begin(), outVector.end(),0);
  }
};

And here is how to use it:

using AddFunctorType = VariadicAdd<double, double, double>;
auto add = NewFunctorFilter(AddFunctorType{});
add->SetVInputs(image,image);
add->Update();

Now comes the tricky part. What if the functor produces a VariableLengthVector ? Well in that case we must provide a way to tell the filter how many output bands it has to allocate. Lets see how to do that:

template<typename O, typename T> struct BandExtraction
{
  BandExtraction(unsigned int indices...) : m_Indices({indices}){}

  auto operator()(const itk::VariableLengthVector<T> & in) const
  {
    itk::VariableLengthVector<O> out(m_Indices.size());

    size_t idx = 0;
    for(auto v: m_Indices)
      {
      out[idx] = static_cast<O>(in[v]);
      ++idx;
      }

    return out;
  }

  // This time OutputSize does not depend on input image size, hence
  // the ...
  size_t OutputSize(...) const
  {
    return m_Indices.size();
  }
  
  // set of band indices to extract
  std::set<unsigned int> m_Indices;
};

Note that there is now a size_t OuptutSize(...) const method that tells the filter how many bands the output will have.

And how to use it:

using ExtractFunctorType = BandExtraction<double,double>;
ExtractFunctorType extractFunctor{1,2};
auto extract = NewFunctorFilter(extractFunctor);
extract->SetVInputs(vimage);
extract->Update();

But what if the number of ouptut bands depends on the number of input bands? Well this is still possible:

// helper function to implement next functor (convert a scalar value
// to a VariableLengthVector)
template <typename T> itk::VariableLengthVector<T> toVector(const T & in)
{
  itk::VariableLengthVector<T> out;
  out.SetSize(1);
  out[0] = in;
  return out;
}

// helper function to implement next functor, VariableLengthVectorVersion (returns in)
template <typename  T> const itk::VariableLengthVector<T> & toVector(const itk::VariableLengthVector<T> & in)
{
  return in;
}

// helper function to implement next functor, Merge two VariableLengthVector in-place
template <typename v1, typename v2> void concatenateVectors(v1 & a, const v2 & b)
{
  const size_t previousSizeOfA = a.GetSize();
  
  a.SetSize(previousSizeOfA+b.GetSize());
  
  for(size_t it = 0; it<b.Size();++it)
    {
    a[previousSizeOfA+it] = static_cast<typename v1::ValueType>(b[it]);
    }
}

// helper function to implement next functor, Merge N VariableLengthVector in-place
template <typename v1, typename v2, typename ...vn> void concatenateVectors(v1 & a, const v2 & b, const vn&... z)
{
  concatenateVectors(a,b);
  concatenateVectors(a,z...);
}

// N  images (all types) -> vector image
// This functor concatenates N images (N = variadic) of type
// VectorImage and or Image, into a single VectorImage
template<typename O, typename ...T> struct VariadicConcatenate
{
  auto operator()(const T &...  ins) const
  {
    itk::VariableLengthVector<O> out;
    concatenateVectors(out, toVector(ins)...);
    
    return out;
  }

  // Must define OutputSize because output pixel is vector
  constexpr size_t OutputSize(const std::array<size_t, sizeof...(T)> inputsNbBands) const
  {
    return std::accumulate(inputsNbBands.begin(),inputsNbBands.end(),0);
  }
};

As you can see, the prototype of OuptutSize() is a bit more complex because it will recieve from the filter an array of the inputs number of bands, which can be used to derive the output number of bands. Note that this is an awesome functor that concatenates any number of inputs, either Image or VectorImage!

How to use it:

using ConcatFunctorType = VariadicConcatenate<double, double, itk::VariableLengthVector<double> >;
auto concatenate = NewFunctorFilter(ConcatFunctorType{});
concatenate->SetVInputs(image,vimage);
concatenate->Update();

Can we define a lambda that returns a VariableLengthVector ? We can not add an OutputSize() method in a lambda ... But there is a solution for that:

/ test FunctorImageFilter with a lambda that returns a
  // VariableLengthVector
  // Converts a neighborhood to a VariableLengthVector
  auto Lambda2 = [](const itk::Neighborhood<double> & in)
                 {
                   itk::VariableLengthVector<double> out(in.Size());
                   std::size_t idx{0};
                   for(auto it = in.Begin(); it!=in.End();++it,++idx)
                     {
                     out[idx]=*it;
                     }
                   return out;
                 };

  // In this case, we use the helper function which allows to specify
  // the number of outputs
  auto filterLambda2  = NewFunctorFilter(Lambda2,vimage->GetNumberOfComponentsPerPixel(),{{3,3}});
  filterLambda2->SetVInputs(image);
  filterLambda2->Update();

As you can see, in this case we need to specify the number of output bands to the NewFunctorFilter.

Implementation Details

More on that coming soon.

Copyright

Thanks to Jordi Inglada for the starting this together, and for many inputs.

The copyright owner is CNES and has signed the ORFEO ToolBox Contributor License Agreement.

---

Check before merging:

• All discussions are resolved
• At least 2 votes from core developers, no vote.
• The feature branch is (reasonably) up-to-date with the base branch
• Dashboard is green
• Copyright owner has signed the ORFEO ToolBox Contributor License Agreement

Modules/Core/Functor/include/otbFunctorImageFilter.h

+/*
+ * Copyright (C) 2005-2017 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 otbFunctorImageFilter_h
+#define otbFunctorImageFilter_h
+
+#include "otbVariadicInputsImageFilter.h"
+#include "otbImage.h"
+#include "otbVectorImage.h"
+#include "itkRGBPixel.h"
+#include "itkRGBAPixel.h"
+#include "itkFixedArray.h"
+#include "itkDefaultConvertPixelTraits.h"
+#include <type_traits>
+
+namespace otb
+{
+/**  
+ * \struct IsNeighborhood 
+ * \brief Struct testing if T is a neighborhood
+ * 
+ * Provides:
+ * - ValueType type set to false_type or true_type
+ * - value set to true or false
+ */
+template <class T> struct IsNeighborhood
+{
+  using ValueType = std::false_type;
+  static constexpr bool value = false;
+};
+
+/// Partial specialisation for itk::Neighborhood<T>
+template <class T> struct IsNeighborhood<itk::Neighborhood<T>>
+{
+  using ValueType = std::true_type;
+  static constexpr bool value = true;
+};
+
+
+/// Partial specialisation for const itk::Neighborhood<T> &
+template <class T> struct IsNeighborhood<const itk::Neighborhood<T>&>
+{
+  using ValueType = std::true_type;
+  static constexpr bool value = true;
+};
+
+/**
+ * \struct IsSuitableType
+ * \brief Helper struct to check if a type can be used as pixel type.
+ *
+ * ::value maps to true if type can be used and false otherwhise.
+ */
+template <class T> struct IsSuitableType
+{
+  static constexpr bool value = std::is_scalar<T>::value;
+};
+
+/// Unwrap complex
+template <class T> struct IsSuitableType<std::complex<T>> {
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Unwrap VariableLengthVector
+template <class T> struct IsSuitableType<itk::VariableLengthVector<T>> {
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Unwrap FixedArray
+template <class T, size_t N> struct IsSuitableType<itk::FixedArray<T,N>> {
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Unwrap RGBPixel
+template <class T> struct IsSuitableType<itk::RGBPixel<T>> {
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Unwrap RGBAPixel
+template <class T> struct IsSuitableType<itk::RGBAPixel<T>> {
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Discard const qualifier
+template <class T> struct IsSuitableType<const T>{
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+/// Discard reference
+template <class T> struct IsSuitableType<T &>{
+  static constexpr bool value = IsSuitableType<T>::value;
+};
+
+
+/**
+ * \struct PixelTypeDeduction
+ * \brief Helper struct to derive PixelType from template parameter.
+ * 
+ * T                           -> PixelType = T
+ * itk::Neighborhood<T>        -> PixelType = T
+ * const itk::Neighborhood<T>& -> PixelType = T
+*/
+template <class T> struct PixelTypeDeduction
+{
+  static_assert(IsSuitableType<T>::value,"T can not be used as a template parameter for Image or VectorImage classes.");
+  using PixelType = T;
+};
+
+/// Partial specialisation for itk::Neighborhood<T>
+template <class T> struct PixelTypeDeduction<itk::Neighborhood<T>>
+{
+  static_assert(IsSuitableType<T>::value,"T can not be used as a template parameter for Image or VectorImage classes.");
+  using PixelType = T;
+};
+
+/// Discard const qualifier
+template <class T> struct PixelTypeDeduction<const  T>
+{
+  using PixelType = typename PixelTypeDeduction<T>::PixelType;
+};
+
+/// Discard reference
+template <class T> struct PixelTypeDeduction<T &>
+{
+  using PixelType = typename PixelTypeDeduction<T>::PixelType;
+};
+
+  
+
+/** 
+ * \struct ImageTypeDeduction
+ * \brief Helper struct to derive ImageType from template parameter
+ * 
+ * T                            -> ImageType = otb::Image<T>
+ * itk::VariableLengthVector<T> -> ImageType = otb::VectorImage<T>
+ * const T &                    -> ImageType = ImageTypeDeduction<T>::ImageType
+ */
+template <class T> struct ImageTypeDeduction 
+{
+  using ImageType = otb::Image<T>;
+};
+
+/// Partial specialisation for itk::VariableLengthVector<T>
+template <class T> struct ImageTypeDeduction<itk::VariableLengthVector<T>>
+{
+  using ImageType = otb::VectorImage<T>;
+};
+
+/**
+* \struct FunctorFilterSuperclassHelper 
+* \brief Struct allowing to derive the superclass prototype for the
+*        FunctorImageFilter class
+*
+* Provides the following:
+* - OutputImageType : type of the output image
+* - FilterType : correct instanciation of VariadicInputsImageFilter from
+* - the operator() prototype
+* - InputHasNeighborhood a tuple of N false_type or true_type to denote
+* - if Ith arg of operator() expects a neighborhood.
+*/
+template <typename T> struct FunctorFilterSuperclassHelper : public FunctorFilterSuperclassHelper<decltype(&T::operator())> {};
+
+/// Partial specialisation for R(*)(T...)
+template <typename R, typename... T> struct FunctorFilterSuperclassHelper<R(*)(T...)>
+{
+  using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
+  using FilterType = VariadicInputsImageFilter<OutputImageType,typename ImageTypeDeduction<typename PixelTypeDeduction<T>::PixelType>::ImageType...>;
+  using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
+};
+
+/// Partial specialisation for R(C::*)(T...) const
+template <typename C, typename R, typename... T> struct FunctorFilterSuperclassHelper<R(C::*)(T...) const>
+{
+  using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
+  using FilterType = VariadicInputsImageFilter<OutputImageType,typename ImageTypeDeduction<typename PixelTypeDeduction<T>::PixelType>::ImageType...>;
+  using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
+};
+
+/// Partial specialisation for R(C::*)(T...)
+template <typename C, typename R, typename... T> struct FunctorFilterSuperclassHelper<R(C::*)(T...)>
+{ 
+  using OutputImageType = typename ImageTypeDeduction<R>::ImageType;
+  using FilterType = VariadicInputsImageFilter<OutputImageType,typename ImageTypeDeduction<typename PixelTypeDeduction<T>::PixelType>::ImageType...>;
+  using InputHasNeighborhood = std::tuple<typename IsNeighborhood<T>::ValueType...>;
+};
+
+/**
+ * \brief This helper method builds a fully functional FunctorImageFilter from a functor instance
+ * 
+ * Functor can be any operator() 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
+ *
+ * The returned filter is ready to use. Inputs can be set through the
+ * SetVariadicInputs() method (see VariadicInputsImageFilter class for
+ * details)
+ * 
+ * \param f the Functor to build the filter from
+ * \param radius The size of neighborhood to use, if there is any
+ * itk::Neighborhood<T> in the operator() arguments.
+ * \return A ready to use OTB filter, which accepts n input image of
+ * type derived from the operator() arguments, and producing an image
+ * correpsonding to the operator() return type.
+ *
+ * Note that this function also works with a lambda as Functor,
+ * provided it returns a scalar type. If your lambda returns a
+ * VariableLengthVector, see the other NewFunctorFilter implementation.
+ */
+template <typename Functor> auto NewFunctorFilter(const Functor& f, itk::Size<2> radius = {{0,0}});
+
+
+/** \class FunctorImageFilter
+ * \brief A generic functor filter templated by its functor
+ * 
+ * \sa VariadicInputsImageFilter
+ * \sa NewFunctorFilter
+ * 
+ * \ingroup IntensityImageFilters   Multithreaded Streamed
+ *
+ * \ingroup OTBFunctor
+*/
+template <class TFunction>
+    class ITK_EXPORT FunctorImageFilter
+  : public FunctorFilterSuperclassHelper<TFunction>::FilterType
+{
+
+public:
+  // Standard typedefs
+  using Self = FunctorImageFilter;
+  using FunctorType = TFunction;
+  using Pointer = itk::SmartPointer<Self>;
+  using ConstPointer = itk::SmartPointer<const Self>;
+
+  // Superclass through the helper struct
+  using Superclass = typename FunctorFilterSuperclassHelper<TFunction>::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
+  using InputHasNeighborhood = typename FunctorFilterSuperclassHelper<TFunction>::InputHasNeighborhood;
+  using InputTypesTupleType = typename Superclass::InputTypesTupleType;
+  template<size_t I> using InputImageType = typename Superclass::template InputImageType<I>;
+  using Superclass::NumberOfInputs;
+  
+  /** Run-time type information (and related methods). */
+  itkTypeMacro(FunctorImageFilter, ImageToImageFilter);
+  
+  /** Get the functor object.
+   * 
+   *  \return A non-const reference to the Functor in use. Note that
+   *  this call will call the Modified() function of the filter, which
+   *  will trigger a full output computation upon Update(). For
+   *  inspection only, prefer using GetFunctor()
+   */
+  FunctorType& GetModifiableFunctor()
+  {
+    this->Modified();
+    return m_Functor;
+  }
+
+  /** Get the functor object.  
+   * 
+   *  \return A const reference to the Functor in use. For a
+   *  modifiable version of the functor, see GetModifiableFunctor().
+   */
+  const FunctorType& GetFunctor() const
+  {
+    return m_Functor;
+  }
+
+protected:
+  /// Constructor of functor filter, will copy the functor
+  FunctorImageFilter(const FunctorType& f, itk::Size<2> radius) : m_Functor(f), m_Radius(radius) {};
+  FunctorImageFilter(const Self &) = delete;
+  void operator =(const Self&) = delete;
+  ~FunctorImageFilter() = default;
+
+private:
+  /// Actual creation of the filter is handled by this free function
+  friend auto NewFunctorFilter<TFunction>(const TFunction& f, itk::Size<2> radius);
+
+  /** Overload of ThreadedGenerateData  */
+  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId) override;
+
+  /**
+   * Pad the input requested region by radius
+   */
+  void GenerateInputRequestedRegion(void) override;
+
+  /**
+   * Will use the OutputSize() method of the functor if output is VectorImage
+   */
+  void GenerateOutputInformation() override;
+
+
+  // The functor member
+  FunctorType m_Functor;
+
+  // Radius if needed
+  itk::Size<2> m_Radius;
+};
+
+// Actual implementation of NewFunctorFilter free function
+template <typename Functor> auto NewFunctorFilter(const Functor& f, itk::Size<2> radius)
+{
+  using FilterType = FunctorImageFilter<Functor>;
+  using PointerType = typename FilterType::Pointer;
+
+  PointerType  p = new FilterType(f,radius);
+  p->UnRegister();
+  return p;
+}
+
+/** 
+ * \struct NumberOfOutputBandsDecorator
+ * \brief This struct allows to forward the operator of template
+ *        parameter, while adding number of ouptut components service.
+ * 
+ * Its purpose is to enable the use of lambda or functor witht
+ * Outputsize() method with FunctorImageFilter.
+ * 
+ * It is used internally in NewFunctorFilter version with
+ * numberOfOutputBands parameter.
+ */ 
+
+template <typename F> struct NumberOfOutputBandsDecorator : F
+{
+public:
+  constexpr NumberOfOutputBandsDecorator(const F t, unsigned int nbComp) : F(t), m_NumberOfOutputBands(nbComp) {}
+  using F::operator();
+
+  constexpr size_t OutputSize(...) const
+  {
+    return m_NumberOfOutputBands;
+  }
+
+private:
+  unsigned int m_NumberOfOutputBands;
+};
+
+/** 
+ * brief This helper method builds a fully functional
+ * FunctorImageFilter from a functor instance which does not provide
+ * the OutputSize() service, or a lambda, returing a VariableLengthVector
+ * 
+ * \param f the Functor to build the filter from
+ * \param numberOfOutputBands The number of output bands that
+ * this filter will return
+ * \param radius The size of neighborhood to use, if there is any
+ * itk::Neighborhood<T> in the operator() arguments.
+ * \return A ready to use OTB filter, which accepts n input image of
+ * type derived from the operator() arguments, and producing an image
+ * correpsonding to the operator() return type.
+ *
+ * Note that this function also works with a lambda as Functor,
+ * provided it returns a scalar type. If your lambda returns a
+ * VariableLengthVector, see the other NewFunctorFilter implementation.
+
+ */ 
+
+template <typename Functor> auto NewFunctorFilter(const Functor& f, unsigned int numberOfOutputBands, itk::Size<2> radius)
+{
+  using FunctorType = NumberOfOutputBandsDecorator<Functor>;
+  FunctorType decoratedF(f,numberOfOutputBands);
+  return  NewFunctorFilter(decoratedF,radius);
+}
+
+
+}// namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbFunctorImageFilter.hxx"
+#endif
+
+#endif