``````

In [1]:

%matplotlib inline

import numpy as np
import matplotlib.pyplot as plt

``````
``````

In [2]:

N = 2**16 # FFT size
M = 2**5 # decimation factor

Nsamples = 100*N
w = np.cos(np.linspace(-np.pi/2,np.pi/2,N))**2 # window before decimation
ww = np.cos(np.linspace(-np.pi/2,np.pi/2,N//M))**2 # window after decimation
f = np.linspace(-0.5,0.5,Nsamples)
x = np.exp(1j*2*np.pi*np.cumsum(f)) # input signal: frequency sweep
x0 = x.reshape((Nsamples//N,N))*w # signal for FFT with overlap phase 0
x1 = x[N//2:-N//2].reshape((Nsamples//N-1,N))*w # signal for FFT with overlap phase 1/2
y0 = (ww*np.fft.ifft(np.fft.fftshift(np.fft.fft(x0)[:,:N//M], axes=1))).ravel() # decimated signal with overlap phase 0
y1 = (ww*np.fft.ifft(np.fft.fftshift(np.fft.fft(x1)[:,:N//M], axes=1))).ravel() # decimated signal with overlap phase 1
y = y0 + y1 # output signal

``````
``````

In [3]:

plt.plot(f[::M][:-1],np.diff(np.unwrap(np.angle(y)))/(2*np.pi))
plt.title('Frequency response')
plt.xlabel('Input natural frequency (cycles/sample)')
plt.ylabel('Output natural frequency (cycles/sample)')
plt.figure()
plt.title('Amplitude response')
plt.xlabel('Input natural frequency (cycles/sample)')
plt.ylabel('Output amplitude (dB)')
plt.plot(f[::M], 20*np.log10(np.abs(y)))
plt.ylim([-200,50])

``````
``````

Out[3]:

(-200, 50)

``````
``````

In [ ]:

``````