Last active
November 17, 2016 12:07
-
-
Save mikaem/7c67630c582f16acd9877cd6ec8f3c4d to your computer and use it in GitHub Desktop.
Minor script for padding in 1D
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
from numpy import * | |
N = 8 | |
x = sin(4*pi*arange(N)/N)+sin(6*pi*arange(N)/N) | |
x_hat_c = fft.fft(x) # Use regular fft | |
x_hat_r = fft.rfft(x) # Use rfft. Should be same | |
M = 3*N//2 | |
x_hat_c_pad = zeros(M, complex) | |
x_hat_r_pad = zeros(M//2+1, complex) | |
c_indices = fft.fftfreq(N, 1./N).astype(int) | |
r_indices = fft.rfftfreq(N, 1./N).astype(int) | |
x_hat_c_pad[c_indices] = x_hat_c | |
x_hat_r_pad[r_indices] = x_hat_r | |
x_c_pad = fft.ifft(x_hat_c_pad) | |
x_r_pad = fft.irfft(x_hat_r_pad) | |
# Alternative | |
x_c_pad2 = fft.ifft(x_hat_c, n=M) # Not ok because of ordering (fftfreq) | |
x_r_pad2 = fft.irfft(x_hat_r, n=M) # Ok because of 0,1,2,...N/2 ordering (rfftfreq) | |
assert allclose(x_r_pad, x_r_pad2) | |
assert allclose(x_c_pad, x_c_pad2) == False | |
assert allclose(x_c_pad, x_r_pad) | |
x_hat_c2 = fft.fft(x_c_pad)[c_indices] | |
x_hat_r2 = fft.rfft(x_r_pad)[r_indices] | |
x_hat_r3 = fft.rfft(x_r_pad2, n=N) # Not ok because input is cropped, not output | |
assert allclose(x_hat_c2, x_hat_c) | |
assert allclose(x_hat_r2, x_hat_r) | |
assert allclose(x_hat_r2, x_hat_r3) == False |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment