kfrlib / kfr

Fast, modern C++ DSP framework, FFT, Sample Rate Conversion, FIR/IIR/Biquad Filters (SSE, AVX, AVX-512, ARM NEON)
https://www.kfrlib.com
GNU General Public License v2.0
1.66k stars 253 forks source link

FFT over tensor with higher resolution #174

Closed Elukej closed 1 year ago

Elukej commented 1 year ago

Hello. I want to apply dft to one of the dimensions of a tensor and expand that dimension to have more elements as to represent higher resolution of the dft. For example, row length of tensor is 4 and I want to get output tensor to have row length of 8. If I would do this with input univector that has 4 elements, and output univector of 8 elements and a dft_plan that is of size 8, it would treat input like it was padded with 4 zeros and the result would be higher resolution dft like I desire. The problem arises if I extract the row of a tensor and create univector_ref from it.

    kfr::univector<kfr::complex<double>> sig{1, 2, 3, 4};
    kfr::univector<kfr::complex<double>, 8> freq;
    const kfr::dft_plan<double> plan(8);
    kfr::univector<kfr::u8> temp(plan.temp_size);
    plan.execute(freq, sig, temp);

    for (int i = 0; i < freq.size(); i++) std::cout << "real: " << freq[i].real() << ", imaginary: " << freq[i].imag() << std::endl;

    kfr::tensor<kfr::complex<double>, 3> tens {{{1,2,3,4}, {5,6,7,8}}, {{9,10,11,12}, {13,14,15,16}}};
    auto arr = tens(0, 0, kfr::trange());
    plan.execute(freq, kfr::make_univector(arr), temp);

    std::cout << std::endl;
    for (int i = 0; i < freq.size(); i++) std::cout << "real: " << freq[i].real() << ", imaginary: " << freq[i].imag() << std::endl;

When I try to apply the transform like in the snippet above, the univector_ref is taking the next 4 items from the tensor so it does the tranformation over {1,2,3,4,5,6,7,8} instead of just {1,2,3,4}(or {1,2,3,4,0,0,0,0}) like in the first case. The output that snippet produces which ilustrates the point is:

real: 10, imaginary: 5.58294e-322
real: -0.414214, imaginary: -7.24264
real: -2, imaginary: 2
real: 2.41421, imaginary: -1.24264
real: -2, imaginary: -5.58294e-322
real: 2.41421, imaginary: 1.24264
real: -2, imaginary: -2
real: -0.414214, imaginary: 7.24264

real: 36, imaginary: 0
real: -4, imaginary: 9.65685
real: -4, imaginary: 4
real: -4, imaginary: 1.65685
real: -4, imaginary: 0
real: -4, imaginary: -1.65685
real: -4, imaginary: -4
real: -4, imaginary: -9.65685

The same behavior happens if the output array is univector_ref to row of some output tensor that has expanded rows to accept higher resolution dft. Is it possible with kfr to achieve the first behavior(like in the first plan execution) without copying input tensor into a tensor that matches the wanted output tensor layout with zero padding on the added dimension first and only then applying dft?

Elukej commented 1 year ago

@dancazarin Please

dancazarin commented 1 year ago

Hello, FFT size is 8 in your case. But then you're passing data that is less than 8 elements. So the code reads the data past the end of the input array. Pass an array with the correct size to make it work. To pad with zeros:

univector<complex<double>> full = truncate(padded(make_univector(arr)), 8);

The data will be copied but it's the only way to pass an array padded with zeros considering the fact that FFT requires the full array.

Some bounds checks are omitted for performance reasons but maybe throwing error is a good solution for detecting such cases as passing data with incorrect size.

KFR FFT algorithm always requires contiguous input and output data. So only contiguous tensors (or tensor slices that are tensors as well) are valid for FFT and can be passed directly (if size matches the selected FFT size). The above code that will make a copy will work well with any tensors.

Elukej commented 1 year ago

@dancazarin Thank you! My issue is resolved. The copying of data is the way to go then for non contigious tensor slices.