ITK  4.6.0
Insight Segmentation and Registration Toolkit
itkDiffusionTensor3DReconstructionImageFilter.h
Go to the documentation of this file.
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 __itkDiffusionTensor3DReconstructionImageFilter_h
19 #define __itkDiffusionTensor3DReconstructionImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include "itkSpatialObject.h"
23 #include "itkDiffusionTensor3D.h"
24 #include "vnl/vnl_matrix.h"
25 #include "vnl/vnl_vector_fixed.h"
26 #include "vnl/vnl_matrix_fixed.h"
27 #include "vnl/algo/vnl_svd.h"
28 #include "itkVectorContainer.h"
29 #include "itkVectorImage.h"
30 
31 namespace itk
32 {
123 template< typename TReferenceImagePixelType,
124  typename TGradientImagePixelType = TReferenceImagePixelType,
125  typename TTensorPixelType = double,
126  typename TMaskImageType = Image<unsigned char, 3> >
128  public ImageToImageFilter< Image< TReferenceImagePixelType, 3 >,
129  Image< DiffusionTensor3D< TTensorPixelType >, 3 > >
130 {
131 public:
132 
139 
141  itkNewMacro(Self);
142 
146 
147  typedef TReferenceImagePixelType ReferencePixelType;
148 
149  typedef TGradientImagePixelType GradientPixelType;
150 
152 
156 
158 
160 
161  typedef typename Superclass::OutputImageRegionType
163 
166 
172 
175 
177  typedef TMaskImageType MaskImageType;
178 
180  typedef vnl_matrix_fixed< double, 6, 6 > TensorBasisMatrixType;
181 
182  typedef vnl_matrix< double > CoefficientMatrixType;
183 
185  typedef vnl_vector_fixed< double, 3 > GradientDirectionType;
186 
188  typedef VectorContainer< unsigned int,
190 
192  void AddGradientImage(const GradientDirectionType &, const GradientImageType *image);
193  const GradientImageType *GetGradientImage(unsigned index) const;
195 
203  const GradientImagesType *image);
204 
206  void SetReferenceImage(ReferenceImageType *referenceImage)
207  {
209  {
210  itkExceptionMacro(<< "Cannot call both methods:"
211  << "AddGradientImage and SetGradientImage. Please call only one of them.");
212  }
213 
214  this->ProcessObject::SetNthInput(0, referenceImage);
215 
217  }
218 
221  { return ( static_cast< ReferenceImageType * >( this->ProcessObject::GetInput(0) ) ); }
222 
224  virtual GradientDirectionType GetGradientDirection(unsigned int idx) const
225  {
226  if ( idx >= m_NumberOfGradientDirections )
227  {
228  itkExceptionMacro(<< "Gradient direction " << idx << "does not exist");
229  }
230  return m_GradientDirectionContainer->ElementAt(idx + 1);
231  }
233 
235  void SetMaskImage(MaskImageType *maskImage);
236 
238  void SetMaskSpatialObject(MaskSpatialObjectType *maskSpatialObject);
239 
240 
244  itkSetMacro(Threshold, ReferencePixelType);
245  itkGetConstMacro(Threshold, ReferencePixelType);
247 
254  itkSetMacro(BValue, TTensorPixelType);
255 #ifdef GetBValue
256 #undef GetBValue
257 #endif
258  itkGetConstReferenceMacro(BValue, TTensorPixelType);
260 
261 #ifdef ITK_USE_CONCEPT_CHECKING
262  // Begin concept checking
263  itkConceptMacro( ReferenceEqualityComparableCheck,
265  itkConceptMacro( TensorEqualityComparableCheck,
267  itkConceptMacro( GradientConvertibleToDoubleCheck,
269  itkConceptMacro( DoubleConvertibleToTensorCheck,
271  itkConceptMacro( GradientReferenceAdditiveOperatorsCheck,
273  ReferencePixelType > ) );
274  itkConceptMacro( GradientReferenceAdditiveAndAssignOperatorsCheck,
276  ReferencePixelType > ) );
277 
278  itkConceptMacro( ReferenceOStreamWritableCheck,
280  itkConceptMacro( TensorOStreamWritableCheck,
282  // End concept checking
283 #endif
284 
285 protected:
288  void PrintSelf(std::ostream & os, Indent indent) const;
289 
290  void ComputeTensorBasis();
291 
293 
294  void ThreadedGenerateData(const
295  OutputImageRegionType & outputRegionForThread, ThreadIdType);
296 
299  typedef enum {
304 
305 private:
306 
307  /* Tensor basis coeffs */
309 
311 
314 
317 
320 
323 
325  TTensorPixelType m_BValue;
326 
329 
332 };
333 }
334 
335 #ifndef ITK_MANUAL_INSTANTIATION
336 #include "itkDiffusionTensor3DReconstructionImageFilter.hxx"
337 #endif
338 
339 #endif
Light weight base class for most itk classes.
void PrintSelf(std::ostream &os, Indent indent) const
virtual GradientDirectionType GetGradientDirection(unsigned int idx) const
void AddGradientImage(const GradientDirectionType &, const GradientImageType *image)
Templated n-dimensional vector image class.
void SetMaskSpatialObject(MaskSpatialObjectType *maskSpatialObject)
const GradientImageType * GetGradientImage(unsigned index) const
DataObject * GetInput(const DataObjectIdentifierType &key)
void SetMaskImage(MaskImageType *maskImage)
This class takes as input one or more reference image (acquired in the absence of diffusion sensitizi...
Base class for filters that take an image as input and produce an image as output.
Define a front-end to the STL &quot;vector&quot; container that conforms to the IndexedContainerInterface.
void SetGradientImage(GradientDirectionContainerType *, const GradientImagesType *image)
Control indentation during Print() invocation.
Definition: itkIndent.h:49
virtual void SetNthInput(DataObjectPointerArraySizeType num, DataObject *input)
#define itkConceptMacro(name, concept)
Represent a diffusion tensor as used in DTI images.
ImageToImageFilter< Image< TReferenceImagePixelType, 3 >, Image< DiffusionTensor3D< TTensorPixelType >, 3 > > Superclass
Templated n-dimensional image class.
Definition: itkImage.h:75
VectorContainer< unsigned int, GradientDirectionType > GradientDirectionContainerType
unsigned int ThreadIdType
Definition: itkIntTypes.h:159
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, ThreadIdType)