ITK
4.0.0
Insight Segmentation and Registration Toolkit
|
00001 /*========================================================================= 00002 * 00003 * Copyright Insight Software Consortium 00004 * 00005 * Licensed under the Apache License, Version 2.0 (the "License"); 00006 * you may not use this file except in compliance with the License. 00007 * You may obtain a copy of the License at 00008 * 00009 * http://www.apache.org/licenses/LICENSE-2.0.txt 00010 * 00011 * Unless required by applicable law or agreed to in writing, software 00012 * distributed under the License is distributed on an "AS IS" BASIS, 00013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 00014 * See the License for the specific language governing permissions and 00015 * limitations under the License. 00016 * 00017 *=========================================================================*/ 00018 #ifndef __itkBSplineControlPointImageFilter_h 00019 #define __itkBSplineControlPointImageFilter_h 00020 00021 #include "itkImageToImageFilter.h" 00022 #include "itkBSplineKernelFunction.h" 00023 #include "itkCoxDeBoorBSplineKernelFunction.h" 00024 #include "itkFixedArray.h" 00025 #include "itkPointSet.h" 00026 #include "itkVariableSizeMatrix.h" 00027 #include "itkVector.h" 00028 #include "itkVectorContainer.h" 00029 00030 #include "vnl/vnl_matrix.h" 00031 #include "vnl/vnl_vector.h" 00032 00033 namespace itk 00034 { 00035 00059 template <class TInputImage, class TOutputImage = TInputImage> 00060 class BSplineControlPointImageFilter 00061 : public ImageToImageFilter<TInputImage, TOutputImage> 00062 { 00063 public: 00064 typedef BSplineControlPointImageFilter Self; 00065 typedef ImageToImageFilter<TInputImage, TOutputImage> Superclass; 00066 typedef SmartPointer<Self> Pointer; 00067 typedef SmartPointer<const Self> ConstPointer; 00068 00070 itkNewMacro(Self); 00071 00073 itkStaticConstMacro( ImageDimension, unsigned int, 00074 TInputImage::ImageDimension ); 00075 00076 typedef TInputImage ControlPointLatticeType; 00077 typedef TOutputImage OutputImageType; 00078 00080 typedef typename OutputImageType::PixelType PixelType; 00081 typedef typename OutputImageType::RegionType RegionType; 00082 typedef typename OutputImageType::IndexType IndexType; 00083 typedef typename OutputImageType::PointType PointType; 00084 typedef typename OutputImageType::RegionType OutputImageRegionType; 00085 00086 typedef typename OutputImageType::SpacingType SpacingType; 00087 typedef typename OutputImageType::PointType OriginType; 00088 typedef typename OutputImageType::SizeType SizeType; 00089 typedef typename OutputImageType::DirectionType DirectionType; 00090 00092 typedef float RealType; 00093 typedef Image<RealType, 00094 itkGetStaticConstMacro( ImageDimension )> RealImageType; 00095 typedef typename RealImageType::Pointer RealImagePointer; 00096 00097 typedef FixedArray<unsigned, 00098 itkGetStaticConstMacro( ImageDimension )> ArrayType; 00099 00101 typedef PointSet<PixelType, 00102 itkGetStaticConstMacro( ImageDimension )> PointSetType; 00103 typedef typename PointSetType::PixelType PointDataType; 00104 typedef typename PointSetType::PointDataContainer PointDataContainerType; 00105 typedef Image<PointDataType, 00106 itkGetStaticConstMacro( ImageDimension )> PointDataImageType; 00107 typedef typename PointDataImageType::Pointer PointDataImagePointer; 00109 00111 typedef CoxDeBoorBSplineKernelFunction<3> KernelType; 00112 typedef BSplineKernelFunction<0> KernelOrder0Type; 00113 typedef BSplineKernelFunction<1> KernelOrder1Type; 00114 typedef BSplineKernelFunction<2> KernelOrder2Type; 00115 typedef BSplineKernelFunction<3> KernelOrder3Type; 00116 00121 void SetSplineOrder( unsigned int ); 00122 00127 void SetSplineOrder( ArrayType ); 00128 00132 itkGetConstReferenceMacro( SplineOrder, ArrayType ); 00133 00150 itkSetMacro( CloseDimension, ArrayType ); 00151 00155 itkGetConstReferenceMacro( CloseDimension, ArrayType ); 00156 00160 itkSetMacro( Spacing, SpacingType ); 00161 00165 itkGetConstMacro( Spacing, SpacingType ); 00166 00170 itkSetMacro( Origin, OriginType ); 00171 00175 itkGetConstMacro( Origin, OriginType ); 00176 00180 itkSetMacro( Size, SizeType ); 00181 00185 itkGetConstMacro( Size, SizeType ); 00186 00201 itkSetMacro( Direction, DirectionType ); 00202 00206 itkGetConstMacro( Direction, DirectionType ); 00207 00215 typename ControlPointLatticeType::Pointer 00216 RefineControlPointLattice( ArrayType ); 00217 00218 protected: 00219 BSplineControlPointImageFilter(); 00220 virtual ~BSplineControlPointImageFilter(); 00221 void PrintSelf( std::ostream& os, Indent indent ) const; 00222 00226 void ThreadedGenerateData( const OutputImageRegionType &, ThreadIdType ); 00227 00228 private: 00229 BSplineControlPointImageFilter( const Self& ); //purposely not implemented 00230 void operator=( const Self& ); //purposely not implemented 00231 00232 00237 void BeforeThreadedGenerateData(); 00238 00243 unsigned int SplitRequestedRegion( unsigned int, unsigned int, OutputImageRegionType & ); 00244 00249 void CollapsePhiLattice( PointDataImageType *, PointDataImageType *, 00250 const RealType, const unsigned int ); 00251 00255 void SetNumberOfLevels( ArrayType ); 00256 00258 SizeType m_Size; 00259 SpacingType m_Spacing; 00260 OriginType m_Origin; 00261 DirectionType m_Direction; 00262 00263 bool m_DoMultilevel; 00264 bool m_GenerateOutputImage; 00265 unsigned int m_MaximumNumberOfLevels; 00266 unsigned int m_CurrentLevel; 00267 ArrayType m_NumberOfControlPoints; 00268 ArrayType m_CloseDimension; 00269 ArrayType m_SplineOrder; 00270 ArrayType m_NumberOfLevels; 00271 00272 vnl_matrix<RealType> m_RefinedLatticeCoefficients[ImageDimension]; 00273 00274 typename KernelType::Pointer m_Kernel[ImageDimension]; 00275 typename KernelOrder0Type::Pointer m_KernelOrder0; 00276 typename KernelOrder1Type::Pointer m_KernelOrder1; 00277 typename KernelOrder2Type::Pointer m_KernelOrder2; 00278 typename KernelOrder3Type::Pointer m_KernelOrder3; 00279 00280 RealType m_BSplineEpsilon; 00281 00282 inline typename RealImageType::IndexType 00283 NumberToIndex( unsigned int number, typename RealImageType::SizeType size ) 00284 { 00285 typename RealImageType::IndexType k; 00286 k[0] = 1; 00287 00288 for ( unsigned int i = 1; i < ImageDimension; i++ ) 00289 { 00290 k[i] = size[ImageDimension-i-1]*k[i-1]; 00291 } 00292 typename RealImageType::IndexType index; 00293 for ( unsigned int i = 0; i < ImageDimension; i++ ) 00294 { 00295 index[ImageDimension-i-1] 00296 = static_cast<unsigned int>( number/k[ImageDimension-i-1] ); 00297 number %= k[ImageDimension-i-1]; 00298 } 00299 return index; 00300 } 00301 00302 }; 00303 00304 } // end namespace itk 00305 00306 #ifndef ITK_MANUAL_INSTANTIATION 00307 #include "itkBSplineControlPointImageFilter.hxx" 00308 #endif 00309 00310 #endif 00311