1 /*========================================================================= 2 * 3 * Copyright Insight Software Consortium 4 * 5 * Licensed under the Apache License, Version 2.0 (the "License"); 6 * you may not use this file except in compliance with the License. 7 * You may obtain a copy of the License at 8 * 9 * http://www.apache.org/licenses/LICENSE-2.0.txt 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 * 17 *=========================================================================*/ 18 #ifndef itkGaussianDerivativeOperator_h 19 #define itkGaussianDerivativeOperator_h 20 21 #include "itkGaussianOperator.h" 22 #include "itkDerivativeOperator.h" 23 24 #include <algorithm> 25 26 namespace itk 27 { 28 /** 29 * \class GaussianDerivativeOperator 30 * \brief A NeighborhoodOperator whose coefficients are a one dimensional, 31 * discrete derivative Gaussian kernel. 32 * 33 * GaussianDerivativeOperator can be used to calculate Gaussian derivatives 34 * by taking its inner product with to a Neighborhood 35 * (NeighborhooIterator) that is swept across an image region. 36 * It is a directional operator. N successive applications 37 * oriented along each dimensional direction will calculate separable, 38 * efficient, N-D Gaussian derivatives of an image region. 39 * 40 * GaussianDerivativeOperator takes three parameters: 41 * 42 * (1) The floating-point variance of the desired Gaussian function. 43 * 44 * (2) The order of the derivative to be calculated (zero order means 45 * it performs only smoothing as a standard itk::GaussianOperator) 46 * 47 * (3) The "maximum error" allowed in the discrete Gaussian 48 * function. "Maximum errror" is defined as the difference between the area 49 * under the discrete Gaussian curve and the area under the continuous 50 * Gaussian. Maximum error affects the Gaussian operator size. Care should 51 * be taken not to make this value too small relative to the variance 52 * lest the operator size become unreasonably large. 53 * 54 * References: 55 * The Gaussian kernel contained in this operator was described 56 * by Tony Lindeberg (Discrete Scale-Space Theory and the Scale-Space 57 * Primal Sketch. Dissertation. Royal Institute of Technology, Stockholm, 58 * Sweden. May 1991.). 59 * 60 * \author Ivan Macia, VICOMTech, Spain, http://www.vicomtech.es 61 * 62 * This implementation is derived from the Insight Journal paper: 63 * https://hdl.handle.net/1926/1290 64 * 65 * \sa GaussianOperator 66 * \sa NeighborhoodOperator 67 * \sa NeighborhoodIterator 68 * \sa Neighborhood 69 * 70 * \ingroup Operators 71 * \ingroup ITKCommon 72 * 73 * \wiki 74 * \wikiexample{Operators/GaussianDerivativeOperator,Create a Gaussian derivative kernel} 75 * \endwiki 76 */ 77 template< typename TPixel, unsigned int VDimension = 2, 78 typename TAllocator = NeighborhoodAllocator< TPixel > > 79 class ITK_TEMPLATE_EXPORT GaussianDerivativeOperator : 80 public NeighborhoodOperator< TPixel, VDimension, TAllocator > 81 { 82 public: 83 /** Standard class type aliases. */ 84 using Self = GaussianDerivativeOperator; 85 using Superclass = NeighborhoodOperator< TPixel, VDimension, TAllocator >; 86 87 /** Neighborhood operator types. */ 88 using GaussianOperatorType = GaussianOperator< TPixel, VDimension, TAllocator >; 89 using DerivativeOperatorType = DerivativeOperator< TPixel, VDimension, TAllocator >; 90 91 /** Constructor. */ 92 GaussianDerivativeOperator() = default; 93 94 /** Copy constructor */ 95 GaussianDerivativeOperator(const Self & other); 96 97 /** Assignment operator */ 98 Self & operator=(const Self & other); 99 100 101 /** Set/Get the flag for calculating scale-space normalized 102 * derivatives. 103 * 104 * Normalized derivatives are obtained multiplying by the scale 105 * parameter $t^1/order$. This use useful for scale-space selection 106 * algorithms such as blob detection. The scaling results in the 107 * value of the derivatives being independent of the size of an 108 * object. */ SetNormalizeAcrossScale(bool flag)109 void SetNormalizeAcrossScale(bool flag) { m_NormalizeAcrossScale = flag; } GetNormalizeAcrossScale()110 bool GetNormalizeAcrossScale() const { return m_NormalizeAcrossScale; } 111 itkBooleanMacro(NormalizeAcrossScale); 112 113 /** Set/Get the variance of the Gaussian kernel. 114 * 115 */ SetVariance(const double variance)116 void SetVariance(const double variance) { m_Variance = variance; } GetVariance()117 double GetVariance() const { return m_Variance; } 118 119 /** Set/Get the spacing for the direction of this kernel. */ SetSpacing(const double spacing)120 void SetSpacing(const double spacing) { m_Spacing = spacing; } GetSpacing()121 double GetSpacing() const { return m_Spacing; } 122 123 /** Set/Get the desired maximum error of the gaussian approximation. Maximum 124 * error is the difference between the area under the discrete Gaussian curve 125 * and the area under the continuous Gaussian. Maximum error affects the 126 * Gaussian operator size. The value is clamped between 0.00001 and 0.99999. */ SetMaximumError(const double maxerror)127 void SetMaximumError(const double maxerror) 128 { 129 constexpr double Min = 0.00001; 130 const double Max = 1.0 - Min; 131 132 m_MaximumError = std::max( Min, std::min( Max, maxerror ) ); 133 } GetMaximumError()134 double GetMaximumError() { return m_MaximumError; } 135 136 /** Sets/Get a limit for growth of the kernel. Small maximum error values with 137 * large variances will yield very large kernel sizes. This value can be 138 * used to truncate a kernel in such instances. A warning will be given on 139 * truncation of the kernel. */ SetMaximumKernelWidth(unsigned int n)140 void SetMaximumKernelWidth(unsigned int n) 141 { 142 m_MaximumKernelWidth = n; 143 } 144 145 /** Sets/Get the order of the derivative. */ SetOrder(const unsigned int order)146 void SetOrder(const unsigned int order) { m_Order = order;} GetOrder()147 unsigned int GetOrder() const { return m_Order; } 148 149 /** Prints member variables */ 150 void PrintSelf(std::ostream & os, Indent i) const override; 151 152 protected: 153 154 using CoefficientVector = typename Superclass::CoefficientVector; 155 156 /** Returns the value of the modified Bessel function I0(x) at a point x >= 0. 157 */ 158 static double ModifiedBesselI0(double); 159 160 /** Returns the value of the modified Bessel function I1(x) at a point x, 161 * x real. */ 162 static double ModifiedBesselI1(double); 163 164 /** Returns the value of the modified Bessel function Ik(x) at a point x>=0, 165 * where k>=2. */ 166 static double ModifiedBesselI(int, double); 167 168 /** Calculates operator coefficients. */ 169 CoefficientVector GenerateCoefficients() override; 170 171 /** Arranges coefficients spatially in the memory buffer. */ Fill(const CoefficientVector & coeff)172 void Fill(const CoefficientVector & coeff) override 173 { this->FillCenteredDirectional(coeff); } 174 175 private: 176 177 /* methods for generations of the coeeficients for a gaussian 178 * operator of 0-order respecting the remaining parameters */ 179 CoefficientVector GenerateGaussianCoefficients() const; 180 181 /** For compatibility with itkWarningMacro */ GetNameOfClass()182 const char * GetNameOfClass() const override 183 { 184 return "itkGaussianDerivativeOperator"; 185 } 186 187 /** Normalize derivatives across scale space */ 188 bool m_NormalizeAcrossScale{true}; 189 190 /** Desired variance of the discrete Gaussian function. */ 191 double m_Variance{1.0}; 192 193 /** Difference between the areas under the curves of the continuous and 194 * discrete Gaussian functions. */ 195 double m_MaximumError{0.005}; 196 197 /** Maximum kernel size allowed. This value is used to truncate a kernel 198 * that has grown too large. A warning is given when the specified maximum 199 * error causes the kernel to exceed this size. */ 200 unsigned int m_MaximumKernelWidth{30}; 201 202 /** Order of the derivative. */ 203 unsigned int m_Order{1}; 204 205 /** Spacing in the direction of this kernel. */ 206 double m_Spacing{1.0}; 207 }; 208 } // namespace itk 209 210 #ifndef ITK_MANUAL_INSTANTIATION 211 #include "itkGaussianDerivativeOperator.hxx" 212 #endif 213 214 #endif 215