RTK  2.0.1
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:
59  ITK_DISALLOW_COPY_AND_ASSIGN(ForwardDifferenceGradientImageFilter);
60 
62  itkStaticConstMacro(InputImageDimension, unsigned int,
63  TInputImage::ImageDimension);
64  itkStaticConstMacro(OutputImageDimension, unsigned int,
65  TOuputImage::ImageDimension);
67 
70 
72  using InputImageType = TInputImage;
73  using InputImagePointer = typename InputImageType::Pointer;
74  using OutputImageType = TOuputImage;
75  using OutputImagePointer = typename OutputImageType::Pointer;
76 
81 
83  itkNewMacro(Self);
84 
87 
89  using OperatorValueType = TOperatorValueType;
90  using OutputValueType = TOuputValue;
91  using InputPixelType = typename InputImageType::PixelType;
92  using CovariantVectorType = typename OutputImageType::PixelType;
93  using OutputImageRegionType = typename OutputImageType::RegionType;
94 
101  void GenerateInputRequestedRegion() override;
102 
106  { this->SetUseImageSpacing(true); }
107 
111  { this->SetUseImageSpacing(false); }
112 
115  itkSetMacro(UseImageSpacing, bool);
116  itkGetConstMacro(UseImageSpacing, bool);
118 
122  void SetDimensionsProcessed(bool* DimensionsProcessed);
123 
126 
127 #ifdef ITK_USE_CONCEPT_CHECKING
128  // Begin concept checking
129  itkConceptMacro( InputConvertibleToOutputCheck,
131  itkConceptMacro( OutputHasNumericTraitsCheck,
133  // End concept checking
134 #endif
135 
146  itkSetMacro(UseImageDirection, bool);
147  itkGetConstMacro(UseImageDirection, bool);
148  itkBooleanMacro(UseImageDirection);
150 
151 protected:
154  void PrintSelf(std::ostream & os, itk::Indent indent) const override;
155 
166 #if ITK_VERSION_MAJOR<5
167  void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread,
168  itk::ThreadIdType threadId) override;
169 #else
170  void DynamicThreadedGenerateData(const OutputImageRegionType & outputRegionForThread) override;
171 #endif
172 
173 
174  void GenerateOutputInformation() override;
175 
176 private:
177 // // An overloaded method which may transform the gradient to a
178 // // physical vector and converts to the correct output pixel type.
179 // template <typename TValue>
180 // void SetOutputPixel( ImageRegionIterator< VectorImage<TValue,OutputImageDimension> > &it, CovariantVectorType &gradient )
181 // {
182 // if ( this->m_UseImageDirection )
183 // {
184 // CovariantVectorType physicalGradient;
185 // it.GetImage()->TransformLocalVectorToPhysicalVector( gradient, physicalGradient );
186 // it.Set( OutputPixelType( physicalGradient.GetDataPointer(), InputImageDimension, false ) );
187 // }
188 // else
189 // {
190 // it.Set( OutputPixelType( gradient.GetDataPointer(), InputImageDimension, false ) );
191 // }
192 // }
193 
194  template <typename T >
196  {
197  it.Value() = gradient;
198  }
199 
200 
202 
203  // flag to take or not the image direction into account
204  // when computing the derivatives.
206 
207  // list of the dimensions along which the gradient has
208  // to be computed. The components on other dimensions
209  // are set to zero
210  bool m_DimensionsProcessed[TInputImage::ImageDimension];
211 
214 };
215 } // end namespace itk
216 
217 #ifndef ITK_MANUAL_INSTANTIATION
218 #include "rtkForwardDifferenceGradientImageFilter.hxx"
219 #endif
220 
221 #endif
void OverrideBoundaryCondition(itk::ImageBoundaryCondition< TInputImage > *boundaryCondition)
virtual void SetUseImageSpacing(bool _arg)
itk::ImageBoundaryCondition< TInputImage, TInputImage > * m_BoundaryCondition
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)
unsigned int ThreadIdType
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, itk::ThreadIdType threadId) override
void SetDimensionsProcessed(bool *DimensionsProcessed)
#define itkConceptMacro(name, concept)
#define itkSetMacro(name, type)