Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Pythran does not recognise fft2 #1095

Open
dxm447 opened this issue Nov 15, 2018 · 9 comments
Open

Pythran does not recognise fft2 #1095

dxm447 opened this issue Nov 15, 2018 · 9 comments

Comments

@dxm447
Copy link

dxm447 commented Nov 15, 2018

%%pythran
#pythran export normalized_correlation(float64[:,:],float64[:,:],float64)
#pythran export normalized_correlation(float32[:,:],float32[:,:],float32)
import numpy as np
def normalized_correlation(im1,im2):
im1_fft = np.fft.fft2(im1)
im2_fft = np.fft.fft2(im2)
im_mult = np.multiply(im1_fft,np.conj(im2_fft))
return np.fft.ifft2(im_mult)

Gives the error

PythranSyntaxError: Attribute 'fft2' unknown

@jeanlaroche
Copy link
Contributor

Indeed, np.fft.ifft2 has not been implemented at this point... rfft and irfft are though...

@dxm447
Copy link
Author

dxm447 commented Nov 15, 2018

@jeanlaroche thanks. Can I use general complex fft for 1d arrays?

@jeanlaroche
Copy link
Contributor

So far, pythran only has np.fft.rfft and np.fft.irfft, i.e. 1D ffts with real inputs.

@serge-sans-paille
Copy link
Owner

@jeanlaroche I'll happily review your patch for complex number support and/or 2D support :-) You did a great job with irfft and rfft

@jeanlaroche
Copy link
Contributor

@jeanlaroche I'll happily review your patch for complex number support and/or 2D support :-) You did a great job with irfft and rfft

ah ah ah ah! Flattery will get you nowhere! :D

@dxm447
Copy link
Author

dxm447 commented Nov 21, 2018

Numpy RFFT is actually good enough for now :)
I am attaching two working code snippets for working with 1D FFT in Pythran, which I have been using in my own project in the past two-three days

#pythran export pythran_fft1D(float32[])
#pythran export pythran_fft1D(float64[])
#pythran export pythran_fft1D(complex64[])
#pythran export pythran_fft1D(complex128[])
import numpy as np
def pythran_fft1D(input_array):
    number_points = np.int64(len(input_array))
    number_rfft_even = np.int64((0.5*number_points) + 1)
    number_rfft_odd = np.int64((number_points + 1)/2)
    if (number_points % 2 == 0):
        number_rfft = number_rfft_even
    else:
        number_rfft = number_rfft_odd
    real_part = np.real(input_array)
    imag_part = np.imag(input_array)
    real_rfft = np.fft.rfft(real_part)
    imag_rfft = np.fft.rfft(imag_part)
    output_fft_real = np.zeros(number_points)
    output_fft_imag = np.zeros(number_points)
    output_fft_real[0:number_rfft] = real_rfft
    output_fft_imag[0:number_rfft] = imag_rfft
    conj_real_rfft = np.conjugate(real_rfft)
    conj_imag_rfft = np.conjugate(imag_rfft)
    for ii in range(number_rfft,number_points):
        output_fft_real[ii] = conj_real_rfft[(number_points - (ii))]
        output_fft_imag[ii] = conj_imag_rfft[(number_points - (ii))]
    output_fft = output_fft_real + ((1j)*output_fft_imag)
    return output_fft

and

#pythran export pythran_ifft1D(float32[])
#pythran export pythran_ifft1D(float64[])
#pythran export pythran_ifft1D(complex64[])
#pythran export pythran_ifft1D(complex128[])
import numpy as np
def pythran_ifft1D(input_array):
    number_points = np.int64(len(input_array))
    conj_array = np.conjugate(input_array)
    number_rfft_even = np.int64((0.5*number_points) + 1)
    number_rfft_odd = np.int64((number_points + 1)/2)
    if (number_points % 2 == 0):
        number_rfft = number_rfft_even
    else:
        number_rfft = number_rfft_odd
    real_part = np.real(conj_array)
    imag_part = np.imag(conj_array)
    real_rfft = np.fft.rfft(real_part)
    imag_rfft = np.fft.rfft(imag_part)
    output_fft_real = np.zeros(number_points)
    output_fft_imag = np.zeros(number_points)
    output_fft_real[0:number_rfft] = real_rfft
    output_fft_imag[0:number_rfft] = imag_rfft
    conj_real_rfft = np.conjugate(real_rfft)
    conj_imag_rfft = np.conjugate(imag_rfft)
    for ii in range(number_rfft,number_points):
        output_fft_real[ii] = conj_real_rfft[(number_points - (ii))]
        output_fft_imag[ii] = conj_imag_rfft[(number_points - (ii))]
    first_fft = output_fft_real + ((1j)*output_fft_imag)
    second_conj = np.conjugate(first_fft)
    output_ifft = second_conj/number_points
    return output_ifft

Both work for real and complex inputs with identical results and times to numpy.fft.fft. If you think it's a good idea I can add the 2D support too and submit it.

@serge-sans-paille
Copy link
Owner

serge-sans-paille commented Nov 21, 2018 via email

@serge-sans-paille
Copy link
Owner

This should be fixed with #1134, at least it now compiles o/

@jeffwitz
Copy link

jeffwitz commented Jun 3, 2019

for fft2, normally performing two FFT works well, but it is efficient on big data, that's why there are a lot of hacks in FFTW in order to deal differently with the different sizes of matrix or vectors.
FFT_image = fft(fft(Image).T).T

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants