elmisterioso
elmisterioso

Reputation: 1

FFT of rectangular pulses in python

I am working on an assignment to convert a Matlab script to fft a rectangular pulse to python. However, the output in python is unexpectedly a straight horizontal line.

The original script performs fft on a rectangular pulse with different low-pass filters.

rect=zeros(100,1);
ffrect=zeros(100,1);
for j=45:55
rect(j,1)=1;
end
frect=fft(rect);
frect=fftshift(frect);
for j = 1:4
klimlow=50-j*10
klimhigh=50+j*10
for k=1:100
if (k < klimlow) | (k > klimhigh)
ffrect(k,1) = 0;
else
ffrect(k,1)=frect(k,1);
end 
end
rrect=ifft(ifftshift(ffrect));
plot(real(rrect))
ylim([-0.2 1.5])
foo=input('Press RETURN to proceed to the next rectangle pulse');
end

The correct and expected output of the first loop iteration looks like this: Output in Matlab

I tried converting the script to python:

import numpy as np
import matplotlib.pyplot as plt

rect = np.zeros((100, 1))
ffrect = np.zeros((100, 1))

for j in range(45, 56):
    rect[j, 0] = 1

frect = np.fft.fft(rect)
frect = np.fft.fftshift(frect)

for j in range(0, 5):
    klimlow = 50 - j * 10
    klimhigh = 50 + j * 10

    for k in range(100):
        if (k < klimlow) or (k > klimhigh):
            ffrect[k, 0] = 0
        else:
            ffrect[k, 0] = frect[k, 0]

    rrect = np.fft.ifft(np.fft.ifftshift(ffrect))
    
    plt.figure()
    plt.plot(np.real(rrect))
    plt.ylim([-0.2, 1.5])

plt.show()

But the output looks like this: Output in python

Any idea what I did wrong?

Upvotes: 0

Views: 149

Answers (1)

Tino D
Tino D

Reputation: 2636

When you define rect in python as np.zeros((n,m)), you create n arrays in each row. This comes down to how arrays are defined in python.

To solve your issue, there are some stuff that you need to take into consideration. First, you have to create rect like this, to make sure that rect is indeed 1 dimensional:

rect = np.zeros(100)

Alternatively, you can call it as you are right now, but you need to flatten it:

rect = np.zeros((100,1))
rect = rect.flatten()

Second, not really important, but nicer for the eyes and also for the interpreter: frect needs to be defined as a complex array, not a float.

ffrect = np.zeros(100, dtype = complex)

Also, why are you using a for loop to index elements 45 to 56? Just use this:

rect[45:56] = 1

After that, just do what you were doing:

frect = np.fft.fft(rect)
frect = np.fft.fftshift(frect)

for j in range(1, 5):
    klimlow = 50 - j * 10
    klimhigh = 50 + j * 10

    for k in range(100):
        if (k < klimlow) or (k > klimhigh):
            ffrect[k] = 0
        else:
            ffrect[k] = frect[k]

    rrect = np.fft.ifft(np.fft.ifftshift(ffrect))
    
    plt.plot(np.real(rrect), label = "FFT of j = "+str(j))
    fig.canvas.draw()

plt.show()
plt.grid()
plt.legend()

You will get the following: FFT of rects

Upvotes: 0

Related Questions