Commit 5b0b2bed authored by Andi Braimllari's avatar Andi Braimllari Committed by David Frank
Browse files

#105 Add (I)FFT shift operation to DataContainer

For further details check for example the numpy docs in:
https://numpy.org/doc/stable/reference/generated/numpy.fft.fftshift.html
parent 14fe7ad3
Pipeline #743587 passed with stages
in 22 minutes and 15 seconds
......@@ -604,6 +604,58 @@ namespace elsa
return concatenated;
}
template <typename data_t>
DataContainer<data_t> fftShift2D(DataContainer<data_t> dc)
{
assert(dc.getDataDescriptor().getNumberOfDimensions() == 2
&& "DataContainer::fftShift2D: currently only supporting 2D signals");
const DataDescriptor& dataDescriptor = dc.getDataDescriptor();
IndexVector_t numOfCoeffsPerDim = dataDescriptor.getNumberOfCoefficientsPerDimension();
index_t m = numOfCoeffsPerDim[0];
index_t n = numOfCoeffsPerDim[1];
index_t firstShift = m / 2;
index_t secondShift = n / 2;
DataContainer<data_t> copyDC(dataDescriptor);
for (index_t i = 0; i < m; ++i) {
for (index_t j = 0; j < n; ++j) {
copyDC((i + firstShift) % m, (j + secondShift) % n) = dc(i, j);
}
}
return copyDC;
}
template <typename data_t>
DataContainer<data_t> ifftShift2D(DataContainer<data_t> dc)
{
assert(dc.getDataDescriptor().getNumberOfDimensions() == 2
&& "DataContainer::ifftShift2D: currently only supporting 2D signals");
const DataDescriptor& dataDescriptor = dc.getDataDescriptor();
IndexVector_t numOfCoeffsPerDim = dataDescriptor.getNumberOfCoefficientsPerDimension();
index_t m = numOfCoeffsPerDim[0];
index_t n = numOfCoeffsPerDim[1];
index_t firstShift = -m / 2;
index_t secondShift = -n / 2;
DataContainer<data_t> copyDC(dataDescriptor);
for (index_t i = 0; i < m; ++i) {
for (index_t j = 0; j < n; ++j) {
index_t leftIndex = (((i + firstShift) % m) + m) % m;
index_t rightIndex = (((j + secondShift) % n) + n) % n;
copyDC(leftIndex, rightIndex) = dc(i, j);
}
}
return copyDC;
}
// ------------------------------------------
// explicit template instantiation
template class DataContainer<float>;
......@@ -623,4 +675,18 @@ namespace elsa
concatenate<std::complex<double>>(const DataContainer<std::complex<double>>&,
const DataContainer<std::complex<double>>&);
template DataContainer<float> fftShift2D<float>(DataContainer<float>);
template DataContainer<std::complex<float>>
fftShift2D<std::complex<float>>(DataContainer<std::complex<float>>);
template DataContainer<double> fftShift2D<double>(DataContainer<double>);
template DataContainer<std::complex<double>>
fftShift2D<std::complex<double>>(DataContainer<std::complex<double>>);
template DataContainer<float> ifftShift2D<float>(DataContainer<float>);
template DataContainer<std::complex<float>>
ifftShift2D<std::complex<float>>(DataContainer<std::complex<float>>);
template DataContainer<double> ifftShift2D<double>(DataContainer<double>);
template DataContainer<std::complex<double>>
ifftShift2D<std::complex<double>>(DataContainer<std::complex<double>>);
} // namespace elsa
......@@ -516,6 +516,18 @@ namespace elsa
DataContainer<data_t> concatenate(const DataContainer<data_t>& dc1,
const DataContainer<data_t>& dc2);
/// Perform the FFT shift operation to the provided signal. Refer to
/// https://numpy.org/doc/stable/reference/generated/numpy.fft.fftshift.html for further
/// details.
template <typename data_t>
DataContainer<data_t> fftShift2D(DataContainer<data_t> dc);
/// Perform the IFFT shift operation to the provided signal. Refer to
/// https://numpy.org/doc/stable/reference/generated/numpy.fft.ifftshift.html for further
/// details.
template <typename data_t>
DataContainer<data_t> ifftShift2D(DataContainer<data_t> dc);
/// User-defined template argument deduction guide for the expression based constructor
template <typename Source>
DataContainer(Source const& source) -> DataContainer<typename Source::data_t>;
......
......@@ -1303,6 +1303,205 @@ TEST_CASE_TEMPLATE("DataContainer: Slice a DataContainer", data_t, float, double
}
}
TEST_CASE_TEMPLATE("DataContainer: FFT shift and IFFT shift a DataContainer", data_t, float, double,
std::complex<float>, std::complex<double>)
{
GIVEN("a one-element 2D data container")
{
DataContainer<data_t> dc(VolumeDescriptor{{1, 1}});
dc[0] = 8;
WHEN("running the FFT shift operation to the container")
{
DataContainer<data_t> fftShiftedDC = fftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(dc.getDataDescriptor(), fftShiftedDC.getDataDescriptor());
}
THEN("the data containers match") { REQUIRE_UNARY(fftShiftedDC == dc); }
}
WHEN("running the IFFT shift operation to the container")
{
DataContainer<data_t> ifftShiftedDC = ifftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(dc.getDataDescriptor(), ifftShiftedDC.getDataDescriptor());
}
THEN("the data containers match") { REQUIRE_UNARY(ifftShiftedDC == dc); }
}
}
GIVEN("a 3x3 2D data container")
{
DataContainer<data_t> dc(VolumeDescriptor{{3, 3}});
dc(0, 0) = 0;
dc(0, 1) = 1;
dc(0, 2) = 2;
dc(1, 0) = 3;
dc(1, 1) = 4;
dc(1, 2) = -4;
dc(2, 0) = -3;
dc(2, 1) = -2;
dc(2, 2) = -1;
DataContainer<data_t> expectedFFTShiftDC(VolumeDescriptor{{3, 3}});
expectedFFTShiftDC(0, 0) = -1;
expectedFFTShiftDC(0, 1) = -3;
expectedFFTShiftDC(0, 2) = -2;
expectedFFTShiftDC(1, 0) = 2;
expectedFFTShiftDC(1, 1) = 0;
expectedFFTShiftDC(1, 2) = 1;
expectedFFTShiftDC(2, 0) = -4;
expectedFFTShiftDC(2, 1) = 3;
expectedFFTShiftDC(2, 2) = 4;
WHEN("running the FFT shift operation to the container")
{
DataContainer<data_t> fftShiftedDC = fftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(fftShiftedDC.getDataDescriptor(),
expectedFFTShiftDC.getDataDescriptor());
}
THEN("the data containers match") { REQUIRE_UNARY(fftShiftedDC == expectedFFTShiftDC); }
}
DataContainer<data_t> expectedIFFTShiftDC(VolumeDescriptor{{3, 3}});
expectedIFFTShiftDC(0, 0) = 4;
expectedIFFTShiftDC(0, 1) = -4;
expectedIFFTShiftDC(0, 2) = 3;
expectedIFFTShiftDC(1, 0) = -2;
expectedIFFTShiftDC(1, 1) = -1;
expectedIFFTShiftDC(1, 2) = -3;
expectedIFFTShiftDC(2, 0) = 1;
expectedIFFTShiftDC(2, 1) = 2;
expectedIFFTShiftDC(2, 2) = 0;
WHEN("running the IFFT shift operation to the container")
{
DataContainer<data_t> ifftShiftedDC = ifftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(ifftShiftedDC.getDataDescriptor(),
expectedIFFTShiftDC.getDataDescriptor());
}
THEN("the data containers match")
{
REQUIRE_UNARY(ifftShiftedDC == expectedIFFTShiftDC);
}
}
}
GIVEN("a 5x5 2D data container")
{
DataContainer<data_t> dc(VolumeDescriptor{{5, 5}});
dc(0, 0) = 28;
dc(0, 1) = 1;
dc(0, 2) = 5;
dc(0, 3) = -18;
dc(0, 4) = 8;
dc(1, 0) = 5;
dc(1, 1) = 6;
dc(1, 2) = 50;
dc(1, 3) = -8;
dc(1, 4) = 9;
dc(2, 0) = 8;
dc(2, 1) = 9;
dc(2, 2) = 10;
dc(2, 3) = 11;
dc(2, 4) = 12;
dc(3, 0) = -12;
dc(3, 1) = -41;
dc(3, 2) = -10;
dc(3, 3) = -9;
dc(3, 4) = -8;
dc(4, 0) = -70;
dc(4, 1) = -6;
dc(4, 2) = 22;
dc(4, 3) = -10;
dc(4, 4) = -3;
DataContainer<data_t> expectedFFTShiftDC(VolumeDescriptor{{5, 5}});
expectedFFTShiftDC(0, 0) = -9;
expectedFFTShiftDC(0, 1) = -8;
expectedFFTShiftDC(0, 2) = -12;
expectedFFTShiftDC(0, 3) = -41;
expectedFFTShiftDC(0, 4) = -10;
expectedFFTShiftDC(1, 0) = -10;
expectedFFTShiftDC(1, 1) = -3;
expectedFFTShiftDC(1, 2) = -70;
expectedFFTShiftDC(1, 3) = -6;
expectedFFTShiftDC(1, 4) = 22;
expectedFFTShiftDC(2, 0) = -18;
expectedFFTShiftDC(2, 1) = 8;
expectedFFTShiftDC(2, 2) = 28;
expectedFFTShiftDC(2, 3) = 1;
expectedFFTShiftDC(2, 4) = 5;
expectedFFTShiftDC(3, 0) = -8;
expectedFFTShiftDC(3, 1) = 9;
expectedFFTShiftDC(3, 2) = 5;
expectedFFTShiftDC(3, 3) = 6;
expectedFFTShiftDC(3, 4) = 50;
expectedFFTShiftDC(4, 0) = 11;
expectedFFTShiftDC(4, 1) = 12;
expectedFFTShiftDC(4, 2) = 8;
expectedFFTShiftDC(4, 3) = 9;
expectedFFTShiftDC(4, 4) = 10;
WHEN("running the FFT shift operation to the container")
{
DataContainer<data_t> fftShiftedDC = fftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(fftShiftedDC.getDataDescriptor(),
expectedFFTShiftDC.getDataDescriptor());
}
THEN("the data containers match") { REQUIRE_UNARY(fftShiftedDC == expectedFFTShiftDC); }
}
DataContainer<data_t> expectedIFFTShiftDC(VolumeDescriptor{{5, 5}});
expectedIFFTShiftDC(0, 0) = 10;
expectedIFFTShiftDC(0, 1) = 11;
expectedIFFTShiftDC(0, 2) = 12;
expectedIFFTShiftDC(0, 3) = 8;
expectedIFFTShiftDC(0, 4) = 9;
expectedIFFTShiftDC(1, 0) = -10;
expectedIFFTShiftDC(1, 1) = -9;
expectedIFFTShiftDC(1, 2) = -8;
expectedIFFTShiftDC(1, 3) = -12;
expectedIFFTShiftDC(1, 4) = -41;
expectedIFFTShiftDC(2, 0) = 22;
expectedIFFTShiftDC(2, 1) = -10;
expectedIFFTShiftDC(2, 2) = -3;
expectedIFFTShiftDC(2, 3) = -70;
expectedIFFTShiftDC(2, 4) = -6;
expectedIFFTShiftDC(3, 0) = 5;
expectedIFFTShiftDC(3, 1) = -18;
expectedIFFTShiftDC(3, 2) = 8;
expectedIFFTShiftDC(3, 3) = 28;
expectedIFFTShiftDC(3, 4) = 1;
expectedIFFTShiftDC(4, 0) = 50;
expectedIFFTShiftDC(4, 1) = -8;
expectedIFFTShiftDC(4, 2) = 9;
expectedIFFTShiftDC(4, 3) = 5;
expectedIFFTShiftDC(4, 4) = 6;
WHEN("running the IFFT shift operation to the container")
{
DataContainer<data_t> ifftShiftedDC = ifftShift2D(dc);
THEN("the data descriptors match")
{
REQUIRE_EQ(ifftShiftedDC.getDataDescriptor(),
expectedIFFTShiftDC.getDataDescriptor());
}
THEN("the data containers match")
{
REQUIRE_UNARY(ifftShiftedDC == expectedIFFTShiftDC);
}
}
}
}
// "instantiate" the test templates for CPU types
TEST_CASE_TEMPLATE_APPLY(datacontainer_construction, CPUTypeTuple);
TEST_CASE_TEMPLATE_APPLY(datacontainer_reduction, CPUTypeTuple);
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment