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 itkLaplacianOperator_hxx
19 #define itkLaplacianOperator_hxx
20 #include "itkLaplacianOperator.h"
21 
22 namespace itk
23 {
24 template< typename TPixel, unsigned int VDimension, typename TAllocator >
25 void
26 LaplacianOperator< TPixel, VDimension, TAllocator >
SetDerivativeScalings(const double * s)27 ::SetDerivativeScalings(const double *s)
28 {
29   for ( unsigned int i = 0; i < VDimension; ++i )
30     {
31     m_DerivativeScalings[i] = s[i];
32     }
33 }
34 
35 //Create the operator
36 template< typename TPixel, unsigned int VDimension, typename TAllocator >
37 void
38 LaplacianOperator< TPixel, VDimension, TAllocator >
CreateOperator()39 ::CreateOperator()
40 {
41   CoefficientVector coefficients;
42 
43   coefficients = this->GenerateCoefficients();
44 
45   this->Fill(coefficients);
46 }
47 
48 //This function fills the coefficients into the corresponding neighborhodd.
49 template< typename TPixel, unsigned int VDimension, typename TAllocator >
50 void
51 LaplacianOperator< TPixel, VDimension, TAllocator >
Fill(const CoefficientVector & coeff)52 ::Fill(const CoefficientVector & coeff)
53 {
54   typename Superclass::CoefficientVector::const_iterator it;
55 
56   std::slice *temp_slice;
57   temp_slice = new std::slice(0, coeff.size(), 1);
58 
59   typename Self::SliceIteratorType data(this, *temp_slice);
60   delete temp_slice;
61 
62   it = coeff.begin();
63 
64   // Copy the coefficients into the neighborhood
65   for ( data = data.Begin(); data < data.End(); ++data, ++it )
66     {
67     *data = *it;
68     }
69 }
70 
71 template< typename TPixel, unsigned int VDimension, typename TAllocator >
72 typename LaplacianOperator< TPixel, VDimension, TAllocator >
73 ::CoefficientVector
74 LaplacianOperator< TPixel, VDimension, TAllocator >
GenerateCoefficients()75 ::GenerateCoefficients()
76 {
77   unsigned int i, w;
78 
79   // Here we set the radius to 1's, here the
80   // operator is 3x3 for 2D, 3x3x3 for 3D.
81   SizeType r;
82 
83   r.Fill(1);
84   this->SetRadius(r);
85 
86   // Create a vector of the correct size to hold the coefficients.
87   w = this->Size();
88   CoefficientVector coeffP(w);
89 
90   //Set the coefficients
91   double   sum = 0.0;
92   for ( i = 0; i < 2 * VDimension; i += 2 )
93     {
94     OffsetValueType stride = this->GetStride(i / 2);
95 
96     const double   hsq = m_DerivativeScalings[i / 2] * m_DerivativeScalings[i / 2];
97     coeffP[w / 2 - stride] =  coeffP[w / 2 + stride] = hsq;
98     sum += 2.0 * hsq;
99     }
100   coeffP[w / 2] = -sum;
101 
102   return coeffP;
103 }
104 } // namespace itk
105 
106 #endif
107