/
itkVkInverseFFTImageFilter.h
144 lines (123 loc) · 4.9 KB
/
itkVkInverseFFTImageFilter.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
/*=========================================================================
*
* Copyright NumFOCUS
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
#ifndef itkVkInverseFFTImageFilter_h
#define itkVkInverseFFTImageFilter_h
#include "itkFFTImageFilterFactory.h"
#include "itkImage.h"
#include "itkInverseFFTImageFilter.h"
#include "itkVkCommon.h"
#include "itkVkGlobalConfiguration.h"
namespace itk
{
/**
*\class VkInverseFFTImageFilter
*
* \brief Vk-based inverse Fast Fourier Transform.
*
* This filter computes the inverse Fourier transform of an image. The
* implementation is based on the VkFFT library.
*
* This filter is multithreaded and supports input images with sizes which are
* divisible only by primes up to 13.
*
* \ingroup FourierTransform
* \ingroup MultiThreaded
* \ingroup ITKFFT
* \ingroup VkFFTBackend
*
* \sa VkGlobalConfiguration
* \sa InverseFFTImageFilter
*/
template <typename TInputImage,
typename TOutputImage = Image<typename TInputImage::PixelType::value_type, TInputImage::ImageDimension>>
class VkInverseFFTImageFilter : public InverseFFTImageFilter<TInputImage, TOutputImage>
{
public:
ITK_DISALLOW_COPY_AND_MOVE(VkInverseFFTImageFilter);
using InputImageType = TInputImage;
using OutputImageType = TOutputImage;
static_assert(std::is_same<typename TInputImage::PixelType, std::complex<float>>::value ||
std::is_same<typename TInputImage::PixelType, std::complex<double>>::value,
"Unsupported pixel type");
static_assert(std::is_same<typename TOutputImage::PixelType, float>::value ||
std::is_same<typename TOutputImage::PixelType, double>::value,
"Unsupported pixel type");
static_assert(TInputImage::ImageDimension >= 1 && TInputImage::ImageDimension <= 3, "Unsupported image dimension");
/** Standard class type aliases. */
using Self = VkInverseFFTImageFilter;
using Superclass = InverseFFTImageFilter<InputImageType, OutputImageType>;
using Pointer = SmartPointer<Self>;
using ConstPointer = SmartPointer<const Self>;
using InputPixelType = typename InputImageType::PixelType;
using OutputPixelType = typename OutputImageType::PixelType;
using ComplexType = InputPixelType;
using RealType = typename ComplexType::value_type;
using SizeType = typename InputImageType::SizeType;
using SizeValueType = typename InputImageType::SizeValueType;
using OutputImageRegionType = typename OutputImageType::RegionType;
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Run-time type information (and related methods). */
itkTypeMacro(VkInverseFFTImageFilter, InverseFFTImageFilter);
static constexpr unsigned int ImageDimension{ InputImageType::ImageDimension };
/** Determine whether local or global properties will be
* referenced for setting up GPU acceleration.
* Defaults to global so that the user can adjust default properties
* in filters constructed through the ITK object factory. */
itkSetMacro(UseVkGlobalConfiguration, bool);
itkGetMacro(UseVkGlobalConfiguration, bool);
/** Local platform identifier for accelerated backend.
* Ignored if `UseVkGlobalConfiguration` is true. */
itkSetMacro(DeviceID, uint64_t);
/** Return the enumerated GPU device to use for FFT
* according to current filter settings. */
uint64_t
GetDeviceID() const
{
return uint64_t{ m_UseVkGlobalConfiguration ? VkGlobalConfiguration::GetDeviceID() : m_DeviceID };
}
SizeValueType
GetSizeGreatestPrimeFactor() const override;
protected:
VkInverseFFTImageFilter();
~VkInverseFFTImageFilter() override = default;
void
GenerateData() override;
void
PrintSelf(std::ostream & os, Indent indent) const override;
private:
bool m_UseVkGlobalConfiguration{ true };
uint64_t m_DeviceID{ 0UL };
VkCommon m_VkCommon{};
};
// Describe whether input/output are real- or complex-valued
// for factory registration
template <>
struct FFTImageFilterTraits<VkInverseFFTImageFilter>
{
template <typename TUnderlying>
using InputPixelType = std::complex<TUnderlying>;
template <typename TUnderlying>
using OutputPixelType = TUnderlying;
using FilterDimensions = std::integer_sequence<unsigned int, 3, 2, 1>;
};
} // namespace itk
#ifndef ITK_MANUAL_INSTANTIATION
# include "itkVkInverseFFTImageFilter.hxx"
#endif
#endif // itkVkInverseFFTImageFilter_h