Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Quaternion FFT using 1D FFT in C++

I am currently implemented research paper that use a Quaternion FFT in c++. However, I could not find any C++ library that support Quaternion FFT. After some surveys, I have found someone on the Internet said that it is possible to convert Quaternion FFT process to several 1D complex-to-complex FFT. Does anyone know how to do that?

I try to use FFTW++ library which supports some basic FFT methods to implement it. I would be really appreciate if anyone could help.

like image 449
wlee Avatar asked Oct 18 '22 09:10

wlee


1 Answers

Thanks for Severin's help, I finally follow the instruction mentioned in this paper to seperate the quaternion FFT into two Complex-to-Complex 2D FFT and has successfully reproduce the result that shown in the paper.

Something like this: ( Please tell me if I am wrong :) )

#include <Array.h>
#include <fftw++.h>

using namespace std;
using namespace utils;
using namespace Array;
using namespace fftwpp;

void SaliencyMapHandler::quaternionFourierTransform(int dim1, int dim2, double* d1, double* d2, double* d3, double* d4) {
    // dim1 is the 1-st dimension of data, dim2 is the 2-nd dimension of data
    fftw::maxthreads = get_max_threads();
    size_t align = sizeof(Complex);

    array2<Complex> f1(dim1, dim2, align);
    array2<Complex> f2(dim1, dim2, align);

    fft2d forward_1(-1, f1);
    fft2d backward_1(1, f1);
    fft2d forward_2(-1, f2);
    fft2d backward_2(1, f2);

    for(int j=0; j<dim1; j++) {
        for(int i=0; i<dim2; i++) {
            f1(i,j) = Complex(d1[j*dim2 + i], d2[j*dim2 + i]);
            f2(i,j) = Complex(d3[j*dim2 + i], d4[j*dim2 + i]);
        }
    }

    forward_1.fft(f1);
    forward_2.fft(f2);

    // Do something on frequency domain

    backward_1.fftNormalized(f1);
    backward_2.fftNormalized(f2);

    for(int j=0; j<dim1; j++) {
        for(int i=0; i<dim2; i++) {
            double p1 = real(f1(i,j));
            double p2 = imag(f1(i,j));
            double p3 = real(f2(i,j));
            double p4 = imag(f2(i,j));

            // Do something after inverse transform
        }
    }
}
like image 179
wlee Avatar answered Nov 09 '22 16:11

wlee