RTK  1.4.0
Reconstruction Toolkit
rtkForwardDifferenceGradientImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright RTK 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 
19 #ifndef rtkForwardDifferenceGradientImageFilter_h
20 #define rtkForwardDifferenceGradientImageFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkCovariantVector.h>
24 #include <itkImageRegionIterator.h>
25 
26 #include "rtkMacro.h"
27 
28 namespace rtk
29 {
50 template< typename TInputImage,
51  typename TOperatorValueType = float,
52  typename TOuputValue = float,
53  typename TOuputImage = itk::Image< itk::CovariantVector
54  < TOuputValue, TInputImage::ImageDimension >, TInputImage::ImageDimension > >
56  public itk::ImageToImageFilter< TInputImage, TOuputImage >
57 {
58 public:
60  itkStaticConstMacro(InputImageDimension, unsigned int,
61  TInputImage::ImageDimension);
62  itkStaticConstMacro(OutputImageDimension, unsigned int,
63  TOuputImage::ImageDimension);
65 
68 
70  typedef TInputImage InputImageType;
71  typedef typename InputImageType::Pointer InputImagePointer;
72  typedef TOuputImage OutputImageType;
73  typedef typename OutputImageType::Pointer OutputImagePointer;
74 
79 
81  itkNewMacro(Self);
82 
85 
87  typedef TOperatorValueType OperatorValueType;
88  typedef TOuputValue OutputValueType;
89  typedef typename InputImageType::PixelType InputPixelType;
90  typedef typename OutputImageType::PixelType CovariantVectorType;
91  typedef typename OutputImageType::RegionType OutputImageRegionType;
92 
99  void GenerateInputRequestedRegion() ITK_OVERRIDE;
100 
104  { this->SetUseImageSpacing(true); }
105 
109  { this->SetUseImageSpacing(false); }
110 
113  itkSetMacro(UseImageSpacing, bool);
114  itkGetConstMacro(UseImageSpacing, bool);
116 
120  void SetDimensionsProcessed(bool* DimensionsProcessed);
121 
124 
125 #ifdef ITK_USE_CONCEPT_CHECKING
126  // Begin concept checking
127  itkConceptMacro( InputConvertibleToOutputCheck,
129  itkConceptMacro( OutputHasNumericTraitsCheck,
131  // End concept checking
132 #endif
133 
144  itkSetMacro(UseImageDirection, bool);
145  itkGetConstMacro(UseImageDirection, bool);
146  itkBooleanMacro(UseImageDirection);
148 
149 protected:
151  virtual ~ForwardDifferenceGradientImageFilter() ITK_OVERRIDE;
152  void PrintSelf(std::ostream & os, itk::Indent indent) const ITK_OVERRIDE;
153 
164 #if ITK_VERSION_MAJOR<5
165  void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread,
166  itk::ThreadIdType threadId) ITK_OVERRIDE;
167 #else
168  void DynamicThreadedGenerateData(const OutputImageRegionType & outputRegionForThread) ITK_OVERRIDE;
169 #endif
170 
171 
172  void GenerateOutputInformation() ITK_OVERRIDE;
173 
174 private:
175  ForwardDifferenceGradientImageFilter(const Self &); //purposely not implemented
176  void operator=(const Self &); //purposely not implemented
177 
178 // // An overloaded method which may transform the gradient to a
179 // // physical vector and converts to the correct output pixel type.
180 // template <typename TValue>
181 // void SetOutputPixel( ImageRegionIterator< VectorImage<TValue,OutputImageDimension> > &it, CovariantVectorType &gradient )
182 // {
183 // if ( this->m_UseImageDirection )
184 // {
185 // CovariantVectorType physicalGradient;
186 // it.GetImage()->TransformLocalVectorToPhysicalVector( gradient, physicalGradient );
187 // it.Set( OutputPixelType( physicalGradient.GetDataPointer(), InputImageDimension, false ) );
188 // }
189 // else
190 // {
191 // it.Set( OutputPixelType( gradient.GetDataPointer(), InputImageDimension, false ) );
192 // }
193 // }
194 
195  template <typename T >
196  void SetOutputPixel( itk::ImageRegionIterator< T > &it, CovariantVectorType &gradient )
197  {
198  it.Value() = gradient;
199  }
200 
201 
203 
204  // flag to take or not the image direction into account
205  // when computing the derivatives.
207 
208  // list of the dimensions along which the gradient has
209  // to be computed. The components on other dimensions
210  // are set to zero
211  bool m_DimensionsProcessed[TInputImage::ImageDimension];
212 
215 };
216 } // end namespace itk
217 
218 #ifndef ITK_MANUAL_INSTANTIATION
219 #include "rtkForwardDifferenceGradientImageFilter.hxx"
220 #endif
221 
222 #endif
void OverrideBoundaryCondition(itk::ImageBoundaryCondition< TInputImage > *boundaryCondition)
virtual void SetUseImageSpacing(bool _arg)
itk::ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
itk::ImageToImageFilter< InputImageType, OutputImageType > Superclass
void PrintSelf(std::ostream &os, itk::Indent indent) const override
Computes the gradient of an image using forward difference.
void SetOutputPixel(itk::ImageRegionIterator< T > &it, CovariantVectorType &gradient)
virtual ~ForwardDifferenceGradientImageFilter() override
unsigned int ThreadIdType
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, itk::ThreadIdType threadId) override
void SetDimensionsProcessed(bool *DimensionsProcessed)
#define itkConceptMacro(name, concept)
#define itkSetMacro(name, type)