ITK
4.1.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 __itkMaskedFFTNormalizedCorrelationImageFilter_h 00019 #define __itkMaskedFFTNormalizedCorrelationImageFilter_h 00020 00021 #include "itkImageToImageFilter.h" 00022 #include "itkImage.h" 00023 00024 namespace itk 00025 { 00115 00116 template <class TInputImage, class TOutputImage > 00117 class ITK_EXPORT MaskedFFTNormalizedCorrelationImageFilter : 00118 public ImageToImageFilter< TInputImage, TOutputImage > 00119 { 00120 public: 00122 typedef MaskedFFTNormalizedCorrelationImageFilter Self; 00123 typedef ImageToImageFilter < TInputImage, TOutputImage > Superclass; 00124 typedef SmartPointer<Self> Pointer; 00125 typedef SmartPointer<const Self> ConstPointer; 00126 00128 itkNewMacro(Self); 00129 00131 itkTypeMacro(MaskedFFTNormalizedCorrelationImageFilter, MaskedFFTNormalizedCrossCorrelationImageFilter); 00132 00135 itkStaticConstMacro(ImageDimension, unsigned int, 00136 TOutputImage::ImageDimension); 00137 00139 typedef TInputImage InputImageType; 00140 typedef TOutputImage OutputImageType; 00141 typedef typename InputImageType::RegionType InputRegionType; 00142 typedef typename InputImageType::Pointer InputImagePointer; 00143 typedef typename InputImageType::ConstPointer InputImageConstPointer; 00144 typedef typename InputImageType::SizeType InputSizeType; 00145 typedef typename OutputImageType::Pointer OutputImagePointer; 00146 typedef typename OutputImageType::PixelType OutputPixelType; 00147 00148 typedef OutputPixelType RealPixelType; 00149 typedef Image< RealPixelType, ImageDimension> RealImageType; 00150 typedef typename RealImageType::Pointer RealImagePointer; 00151 typedef typename RealImageType::IndexType RealIndexType; 00152 typedef typename RealImageType::SizeType RealSizeType; 00153 typedef typename RealImageType::RegionType RealRegionType; 00154 00155 typedef Image< std::complex<RealPixelType>, ImageDimension > FFTImageType; 00156 typedef typename FFTImageType::Pointer FFTImagePointer; 00157 00159 void SetFixedImage(InputImageType *input) 00160 { 00161 this->SetNthInput(0, const_cast<InputImageType *>(input) ); 00162 } 00163 InputImageType * GetFixedImage() 00164 { 00165 return static_cast<InputImageType*>(const_cast<DataObject *>(this->ProcessObject::GetInput(0))); 00166 } 00168 00170 void SetMovingImage(InputImageType *input) 00171 { 00172 this->SetNthInput(1, const_cast<InputImageType *>(input) ); 00173 } 00174 InputImageType * GetMovingImage() 00175 { 00176 return static_cast<InputImageType*>(const_cast<DataObject *>(this->ProcessObject::GetInput(1))); 00177 } 00179 00181 void SetFixedImageMask(InputImageType *input) 00182 { 00183 this->SetNthInput(2, const_cast<InputImageType *>(input) ); 00184 } 00185 InputImageType * GetFixedImageMask() 00186 { 00187 return static_cast<InputImageType*>(const_cast<DataObject *>(this->ProcessObject::GetInput(2))); 00188 } 00190 00192 void SetMovingImageMask(InputImageType *input) 00193 { 00194 this->SetNthInput(3, const_cast<InputImageType *>(input) ); 00195 } 00196 InputImageType * GetMovingImageMask() 00197 { 00198 return static_cast<InputImageType*>(const_cast<DataObject *>(this->ProcessObject::GetInput(3))); 00199 } 00201 00203 itkSetMacro(RequiredNumberOfOverlappingPixels,SizeValueType); 00204 itkGetMacro(RequiredNumberOfOverlappingPixels,SizeValueType); 00206 00207 #ifdef ITK_USE_CONCEPT_CHECKING 00208 00209 itkConceptMacro( OutputPixelTypeIsFloatingPointCheck, 00210 ( Concept::IsFloatingPoint< OutputPixelType > ) ); 00211 00213 #endif 00214 00215 protected: 00216 MaskedFFTNormalizedCorrelationImageFilter() 00217 { 00218 this->SetNumberOfRequiredInputs(2); 00219 m_RequiredNumberOfOverlappingPixels = 0; 00220 } 00221 virtual ~MaskedFFTNormalizedCorrelationImageFilter() {} 00222 void PrintSelf(std::ostream& os, Indent indent) const; 00223 00225 void VerifyInputInformation(); 00226 00228 void GenerateData(); 00229 00235 virtual void GenerateInputRequestedRegion(); 00236 00241 void GenerateOutputInformation(); 00242 00243 template< class LocalInputImageType > 00244 typename LocalInputImageType::Pointer PreProcessMask( const LocalInputImageType * inputImage, const LocalInputImageType * inputMask ); 00245 00246 template< class LocalInputImageType > 00247 typename LocalInputImageType::Pointer PreProcessImage( const LocalInputImageType * inputImage, LocalInputImageType * inputMask ); 00248 00249 typename TInputImage::Pointer RotateImage( InputImageType * inputImage ); 00250 00251 template< class LocalInputImageType, class LocalOutputImageType > 00252 typename LocalOutputImageType::Pointer CalculateForwardFFT( LocalInputImageType * inputImage, InputSizeType & FFTImageSize ); 00253 00254 template< class LocalInputImageType, class LocalOutputImageType > 00255 typename LocalOutputImageType::Pointer CalculateInverseFFT( LocalInputImageType * inputImage, RealSizeType & combinedImageSize ); 00256 00257 // Helper math methods. 00258 template< class LocalInputImageType, class LocalOutputImageType > 00259 typename LocalOutputImageType::Pointer ElementProduct( LocalInputImageType * inputImage1, LocalInputImageType * inputImage2 ); 00260 00261 template< class LocalInputImageType > 00262 typename LocalInputImageType::Pointer ElementQuotient( LocalInputImageType * inputImage1, LocalInputImageType * inputImage2 ); 00263 00264 template< class LocalInputImageType > 00265 typename LocalInputImageType::Pointer ElementSubtraction( LocalInputImageType * inputImage1, LocalInputImageType * inputImage2 ); 00266 00267 template< class LocalInputImageType > 00268 typename LocalInputImageType::Pointer ElementPositive( LocalInputImageType * inputImage ); 00269 00270 template< class LocalInputImageType, class LocalOutputImageType > 00271 typename LocalOutputImageType::Pointer ElementRound( LocalInputImageType * inputImage ); 00272 00273 // This function factorizes the image size uses factors of 2, 3, and 00274 // 5. After this factorization, if there are any remaining values, 00275 // the function returns this value. 00276 int FactorizeNumber( int n ); 00277 00278 // Find the closest valid dimension above the desired dimension. This 00279 // will be a combination of 2s, 3s, and 5s. 00280 int FindClosestValidDimension( int n ); 00281 00282 template< class LocalInputImageType > 00283 double CalculatePrecisionTolerance( LocalInputImageType * inputImage ); 00284 00285 private: 00286 MaskedFFTNormalizedCorrelationImageFilter(const Self&); //purposely not implemented 00287 void operator=(const Self&); //purposely not implemented 00288 00292 SizeValueType m_RequiredNumberOfOverlappingPixels; 00293 }; 00294 } // end namespace itk 00295 00296 #ifndef ITK_MANUAL_INSTANTIATION 00297 #include "itkMaskedFFTNormalizedCorrelationImageFilter.hxx" 00298 #endif 00299 00300 #endif 00301