Remedial Filtering or FFT Question

I am pretty sure this is remedial, but I have beat my head on it too long to not ask...

I am trying to simply low-pass filter a signal, but I seem to get out frequency content that shouldn't be there. Following is a minimal working example:

fl = 1.024e6;
fh = 65e6;
fs = 135e6;

t = 0:(1/fs):(0.01 - 1/fs);
x = cos(2*pi*fl*t) + cos(2*pi*fh*t);
X = fftshift(fft(x))/length(x);

figure
subplot(2,1,1)
plot(real(X))
subplot(2,1,2)
plot(imag(X))
subplot(2,1,1)
title('Frequency Domain: Original Signal')

b = firpm(200,[0 0.2 0.23 1],[1 1 0 0]);
y = filter(b,1,x);

% Y = fftshift(fft(y(101:end)))/length(y(101:end));
Y = fftshift(fft(y))/length(y);

figure
subplot(2,1,1)
plot(real(Y))
subplot(2,1,2)
plot(imag(Y))
subplot(2,1,1)
title('Frequency Domain: Filtered Signal')

The higher frequency is eliminated, but the problem is in the lower frequencies. In this simple example I should have delta functions at +/- f_L for the real part of the transform and basically nothing for the imaginary part. Instead I get something nontrivial in the imaginary part.

As an added bit of weirdness, if I lop off the transient elements of 'y', the frequency response looks noticeably different even though that's tossing out 100 samples from more than 1,000,000.

I appreciate your input.

Answer 1

[Star Strider]

You will of course have imaginary parts to the fft of both the unfiltered and filtered signals. There is phase information in the signal before and after you filter it.

The only situation in which you would not have a complex result of the Fourier transform of a signal would be that of a signal that was symmetric (and real) about t = 0. If you have the Symbolic Math Toolbox (or feel like doing this by hand), consider:

syms w t 

F1 = int(exp(-j*w*t), t, -1, 1)
F1 = simplify(F1)

F2 = int(exp(-j*w*t), t, 0, 2)
F2 = simplify(F2)

Both signals are essentially the same, except that the time-domain signal in F1 is symmetric about t = 0 (defined on (-1 1) and the time-domain signal in F2 is defined on (0 2). The Fourier transforms are significantly different, most notably in this instance in that F1 is real and F2 is complex.

[John]

This pretty much clears up my confusion, with a couple caveats. I can't see how to edit the original question, so I will provide my own answer below.

[Star Strider]

‘First, the original time domain signal in my question is supposed to be real and symmetric about t=0. ’

It clearly is not symmetrical about (t = 0):

t = 0:(1/fs):(0.01 - 1/fs);
x = cos(2*pi*fl*t) + cos(2*pi*fh*t);

The symbolic Fourier transforms are easy enough to do by hand. The Symbolic Math Toolbox results are:

F1 = 2*sin(w))/w
F2 = (2*sin(w)*(cos(w) - sin(w)*i))/w

(Note the imaginary component in F2.)

Answer 2

[John]

Thanks to Star Strider for setting me on the right track. I don't have the Symbolic Toolbox, just Signal Processing, so I couldn't try the example offered in his/her answer. Below are a few additions to his/her input that I think make for a more complete answer.

Bottom Line Up Front: Apparently the calculated frequency response is VERY sensitive to where on the continuous time sinusoid the samples fall.

First, the original time domain signal in my question is supposed to be real and symmetric about t=0. It is simply the sum of two cosines. I have since adjusted the time window to make sure the signal was not symmetric about the center sample (e.g. t = 0:(1/fs):(0.0091 - 1/fs);) In this case, the maximum of the imaginary part of its FFT is still approx 1e-10. See the first figure included in the question's example.

Next, observe that the commented calculation of "Y" in the question accounts for the filter's group delay. Without cropping out the first 100 samples, the filtered signal will certainly not be symmetric about t=0. However, taking a close look at what happens after cropping the filter transients it seems you will still see a significant (compared to the unfiltered signal) imaginary component in the FFT of the filtered signal.

For example, run the code in the attached file which includes the code from the original question as well as some additional plots. When you align the filtered signal and the low frequency part of the original signal, the difference between their amplitudes after transients is on the order of 1e-3 (i.e. 0.1%). If you crop out the first few cycles of the cosine to get past transients, BOTH the filtered and original signals have significant imaginary components in their FFTs.

Since what I really care about is in the time domain, I'm going to quit worrying about this question. But it has been educational for me. Any additional feedback is still very welcome.

Thanks.

[John]

It appears the attachment referenced above didn't come through. Trying again...

[Star Strider]

‘First, the original time domain signal in my question is supposed to be real and symmetric about t=0. ’

It clearly is not symmetrical about (t = 0):

t = 0:(1/fs):(0.01 - 1/fs);
x = cos(2*pi*fl*t) + cos(2*pi*fh*t);

The symbolic Fourier transforms are easy enough to do by hand. The Symbolic Math Toolbox results are:

F1 = 2*sin(w))/w
F2 = (2*sin(w)*(cos(w) - sin(w)*i))/w

(Note the imaginary component in F2.)