itkRecursiveSeparableImageFilter.h

Go to the documentation of this file.

/*=========================================================================
 *
 *  Copyright NumFOCUS
 *
 *  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
 *
 *         https://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  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 itkRecursiveSeparableImageFilter_h
#define itkRecursiveSeparableImageFilter_h
 
#include "itkInPlaceImageFilter.h"
#include "itkNumericTraits.h"
#include "itkVariableLengthVector.h"
 
namespace itk
{
template <typename TInputImage, typename TOutputImage = TInputImage>
class ITK_TEMPLATE_EXPORT RecursiveSeparableImageFilter : public InPlaceImageFilter<TInputImage, TOutputImage>
{
public:
  ITK_DISALLOW_COPY_AND_MOVE(RecursiveSeparableImageFilter);

  using Self = RecursiveSeparableImageFilter;
  using Superclass = InPlaceImageFilter<TInputImage, TOutputImage>;
  using Pointer = SmartPointer<Self>;
  using ConstPointer = SmartPointer<const Self>;

  itkOverrideGetNameOfClassMacro(RecursiveSeparableImageFilter);

  using InputImagePointer = typename TInputImage::Pointer;
  using InputImageConstPointer = typename TInputImage::ConstPointer;

  using InputPixelType = typename TInputImage::PixelType;
  using RealType = typename NumericTraits<InputPixelType>::RealType;
  using ScalarRealType = typename NumericTraits<InputPixelType>::ScalarRealType;
 
  using OutputImageRegionType = typename TOutputImage::RegionType;

  using InputImageType = TInputImage;

  using OutputImageType = TOutputImage;

  itkGetConstMacro(Direction, unsigned int);

  itkSetMacro(Direction, unsigned int);

  void
  SetInputImage(const TInputImage *);

  const TInputImage *
  GetInputImage();
 
protected:
  RecursiveSeparableImageFilter();
  ~RecursiveSeparableImageFilter() override = default;
  void
  PrintSelf(std::ostream & os, Indent indent) const override;
 
  void
  BeforeThreadedGenerateData() override;
 
  void
  GenerateData() override;
 
  void
  DynamicThreadedGenerateData(const OutputImageRegionType &) override;

  void
  EnlargeOutputRequestedRegion(DataObject * output) override;

  virtual void
  SetUp(ScalarRealType spacing) = 0;

  void
  FilterDataArray(RealType * outs, const RealType * data, RealType * scratch, SizeValueType ln) const;
 
protected:
  ScalarRealType m_N0{};
  ScalarRealType m_N1{};
  ScalarRealType m_N2{};
  ScalarRealType m_N3{};

  ScalarRealType m_D1{};
  ScalarRealType m_D2{};
  ScalarRealType m_D3{};
  ScalarRealType m_D4{};

  ScalarRealType m_M1{};
  ScalarRealType m_M2{};
  ScalarRealType m_M3{};
  ScalarRealType m_M4{};

  ScalarRealType m_BN1{};
  ScalarRealType m_BN2{};
  ScalarRealType m_BN3{};
  ScalarRealType m_BN4{};
 
  ScalarRealType m_BM1{};
  ScalarRealType m_BM2{};
  ScalarRealType m_BM3{};
  ScalarRealType m_BM4{};
 
 
  template <typename T1, typename T2>
  static inline void
  MathEMAMAMAM(T1 &       out,
               const T1 & a1,
               const T2 & b1,
               const T1 & a2,
               const T2 & b2,
               const T1 & a3,
               const T2 & b3,
               const T1 & a4,
               const T2 & b4)
  {
    out = a1 * b1 + a2 * b2 + a3 * b3 + a4 * b4;
  }
 
 
  template <typename T1, typename T2>
  static inline void
  MathEMAMAMAM(VariableLengthVector<T1> &       out,
               const VariableLengthVector<T1> & a1,
               const T2 &                       b1,
               const VariableLengthVector<T1> & a2,
               const T2 &                       b2,
               const VariableLengthVector<T1> & a3,
               const T2 &                       b3,
               const VariableLengthVector<T1> & a4,
               const T2 &                       b4)
  {
    const unsigned int sz = a1.GetSize();
    if (sz != out.GetSize())
    {
      out.SetSize(sz);
    }
    for (unsigned int i = 0; i < sz; ++i)
    {
      out[i] = a1[i] * b1 + a2[i] * b2 + a3[i] * b3 + a4[i] * b4;
    }
  }
 
  template <typename T1, typename T2>
  static inline void
  MathSMAMAMAM(T1 &       out,
               const T1 & a1,
               const T2 & b1,
               const T1 & a2,
               const T2 & b2,
               const T1 & a3,
               const T2 & b3,
               const T1 & a4,
               const T2 & b4)
  {
    out -= a1 * b1 + a2 * b2 + a3 * b3 + a4 * b4;
  }
 
  template <typename T1, typename T2>
  static inline void
  MathSMAMAMAM(VariableLengthVector<T1> &       out,
               const VariableLengthVector<T1> & a1,
               const T2 &                       b1,
               const VariableLengthVector<T1> & a2,
               const T2 &                       b2,
               const VariableLengthVector<T1> & a3,
               const T2 &                       b3,
               const VariableLengthVector<T1> & a4,
               const T2 &                       b4)
  {
    const unsigned int sz = a1.GetSize();
    if (sz != out.GetSize())
    {
      out.SetSize(sz);
    }
    for (unsigned int i = 0; i < sz; ++i)
    {
      out[i] -= a1[i] * b1 + a2[i] * b2 + a3[i] * b3 + a4[i] * b4;
    }
  }
 
private:
  unsigned int m_Direction{ 0 };
};
} // end namespace itk
 
#ifndef ITK_MANUAL_INSTANTIATION
#  include "itkRecursiveSeparableImageFilter.hxx"
#endif
 
#endif