Loading [MathJax]/extensions/tex2jax.js
ITK 6.0.0
Insight Toolkit
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
itkMaskedFFTNormalizedCorrelationImageFilter.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
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 * https://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 itkMaskedFFTNormalizedCorrelationImageFilter_h
19#define itkMaskedFFTNormalizedCorrelationImageFilter_h
20
22#include "itkImage.h"
23
24namespace itk
25{
138template <typename TInputImage, typename TOutputImage, typename TMaskImage = TInputImage>
140 : public ImageToImageFilter<TInputImage, TOutputImage>
141{
142public:
143 ITK_DISALLOW_COPY_AND_MOVE(MaskedFFTNormalizedCorrelationImageFilter);
144
150
152 itkNewMacro(Self);
153
155 itkOverrideGetNameOfClassMacro(MaskedFFTNormalizedCorrelationImageFilter);
156
159 static constexpr unsigned int ImageDimension = TOutputImage::ImageDimension;
160
162 using InputImageType = TInputImage;
163 using InputRegionType = typename InputImageType::RegionType;
164 using InputImagePointer = typename InputImageType::Pointer;
165 using InputImageConstPointer = typename InputImageType::ConstPointer;
166 using InputSizeType = typename InputImageType::SizeType;
168
169 using OutputImageType = TOutputImage;
170 using OutputImagePointer = typename OutputImageType::Pointer;
171 using OutputPixelType = typename OutputImageType::PixelType;
172
180
181 using MaskImageType = TMaskImage;
182 using MaskImagePointer = typename MaskImageType::Pointer;
183
186
209 itkSetMacro(RequiredNumberOfOverlappingPixels, SizeValueType);
210 itkGetMacro(RequiredNumberOfOverlappingPixels, SizeValueType);
214 itkGetMacro(RequiredFractionOfOverlappingPixels, RealPixelType);
215 itkSetClampMacro(RequiredFractionOfOverlappingPixels, RealPixelType, 0.0f, 1.0f);
218 itkGetMacro(MaximumNumberOfOverlappingPixels, SizeValueType);
219
220 itkConceptMacro(OutputPixelTypeIsFloatingPointCheck, (Concept::IsFloatingPoint<OutputPixelType>));
221
222protected:
224 {
225 // #0 "FixedImage" required
226 Self::SetPrimaryInputName("FixedImage");
227
228 // #1 "MaskImage" required
229 Self::AddRequiredInputName("MovingImage", 1);
230
231 // #2 "FixedImageMask" optional
232 Self::AddOptionalInputName("FixedImageMask", 2);
233
234 // #3 "MaskImageMask" optional
235 Self::AddOptionalInputName("MovingImageMask", 3);
236
241 }
243 void
244 PrintSelf(std::ostream & os, Indent indent) const override;
245
247 void
248 VerifyInputInformation() const override;
249
251 void
252 GenerateData() override;
253
259 void
261
266 void
268
269 void
271
272 typename TMaskImage::Pointer
273 PreProcessMask(const InputImageType * inputImage, const MaskImageType * inputMask);
274
275 typename TInputImage::Pointer
276 PreProcessImage(const InputImageType * inputImage, const MaskImageType * inputMask);
277
278 template <typename LocalInputImageType>
279 typename LocalInputImageType::Pointer
280 RotateImage(LocalInputImageType * inputImage);
281
282 template <typename LocalInputImageType, typename LocalOutputImageType>
283 typename LocalOutputImageType::Pointer
284 CalculateForwardFFT(LocalInputImageType * inputImage, InputSizeType & FFTImageSize);
285
286 template <typename LocalInputImageType, typename LocalOutputImageType>
287 typename LocalOutputImageType::Pointer
288 CalculateInverseFFT(LocalInputImageType * inputImage, RealSizeType & combinedImageSize);
289
290 // Helper math methods.
291 template <typename LocalInputImageType, typename LocalOutputImageType>
292 typename LocalOutputImageType::Pointer
293 ElementProduct(LocalInputImageType * inputImage1, LocalInputImageType * inputImage2);
294
295 template <typename LocalInputImageType>
296 typename LocalInputImageType::Pointer
297 ElementQuotient(LocalInputImageType * inputImage1, LocalInputImageType * inputImage2);
298
299 template <typename LocalInputImageType>
300 typename LocalInputImageType::Pointer
301 ElementSubtraction(LocalInputImageType * inputImage1, LocalInputImageType * inputImage2);
302
303 template <typename LocalInputImageType>
304 typename LocalInputImageType::Pointer
305 ElementPositive(LocalInputImageType * inputImage);
306
307 template <typename LocalInputImageType, typename LocalOutputImageType>
308 typename LocalOutputImageType::Pointer
309 ElementRound(LocalInputImageType * inputImage);
310
311 // This function factorizes the image size uses factors of 2, 3, and
312 // 5. After this factorization, if there are any remaining values,
313 // the function returns this value.
314 int
316
317 // Find the closest valid dimension above the desired dimension. This
318 // will be a combination of 2s, 3s, and 5s.
319 int
321
322 template <typename LocalInputImageType>
323 double
324 CalculatePrecisionTolerance(LocalInputImageType * inputImage);
325
326private:
331
335
338
340 const unsigned int m_TotalForwardAndInverseFFTs{ 12 };
341
344};
345} // end namespace itk
346
347#ifndef ITK_MANUAL_INSTANTIATION
348# include "itkMaskedFFTNormalizedCorrelationImageFilter.hxx"
349#endif
350
351#endif
Base class for all data objects in ITK.
Templated n-dimensional image class.
Definition itkImage.h:89
Point< PointValueType, VImageDimension > PointType
ImageRegion< VImageDimension > RegionType
Control indentation during Print() invocation.
Definition itkIndent.h:50
Image< std::complex< RealPixelType >, ImageDimension > FFTImageType
LocalOutputImageType::Pointer ElementProduct(LocalInputImageType *inputImage1, LocalInputImageType *inputImage2)
itkGetInputMacro(FixedImage, InputImageType)
LocalInputImageType::Pointer ElementQuotient(LocalInputImageType *inputImage1, LocalInputImageType *inputImage2)
itkSetInputMacro(FixedImage, InputImageType)
void PrintSelf(std::ostream &os, Indent indent) const override
itkSetInputMacro(MovingImage, InputImageType)
LocalInputImageType::Pointer RotateImage(LocalInputImageType *inputImage)
void EnlargeOutputRequestedRegion(DataObject *output) override
itkSetInputMacro(FixedImageMask, MaskImageType)
itkSetInputMacro(MovingImageMask, MaskImageType)
double CalculatePrecisionTolerance(LocalInputImageType *inputImage)
LocalInputImageType::Pointer ElementPositive(LocalInputImageType *inputImage)
TMaskImage::Pointer PreProcessMask(const InputImageType *inputImage, const MaskImageType *inputMask)
TInputImage::Pointer PreProcessImage(const InputImageType *inputImage, const MaskImageType *inputMask)
LocalOutputImageType::Pointer CalculateForwardFFT(LocalInputImageType *inputImage, InputSizeType &FFTImageSize)
LocalInputImageType::Pointer ElementSubtraction(LocalInputImageType *inputImage1, LocalInputImageType *inputImage2)
LocalOutputImageType::Pointer CalculateInverseFFT(LocalInputImageType *inputImage, RealSizeType &combinedImageSize)
itkGetInputMacro(MovingImage, InputImageType)
LocalOutputImageType::Pointer ElementRound(LocalInputImageType *inputImage)
~MaskedFFTNormalizedCorrelationImageFilter() override=default
itkGetInputMacro(FixedImageMask, MaskImageType)
itkGetInputMacro(MovingImageMask, MaskImageType)
virtual void SetPrimaryInputName(const DataObjectIdentifierType &key)
void AddOptionalInputName(const DataObjectIdentifierType &)
Define a named input that is not required and optionally associate with a numbered index.
bool AddRequiredInputName(const DataObjectIdentifierType &)
Define a required named input and optionally associate it with a numbered index.
Implements transparent reference counting.
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
unsigned long SizeValueType
Definition itkIntTypes.h:86