I have been looking at the FFT (in python) of various simple functions. Everything was working as expected, except when I take the Fourier transform of the function:
$$ f(t) = e^{-\gamma t}\cos^2\omega_1t $$
where I obtain a peak at $\omega=0$ which I would like to remove. This peak shows up in $g(t)=\cos^2(\omega_1 t)$ as well, and can be removed by subtracting the mean from the input signal. However, for $f(t)$ the exponential decay broadens the width of the peak in the Fourier spectrum and subtracting the mean does not fully remove it.
How would I go about removing the broadened peak at $\omega=0$ in $f(t)$ ?
import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate
# number of points
samplingFrequency = 50
# spacing between points one second apart
samplingInterval = 1 / samplingFrequency
# begin and end time of signals
beginTime = 0
endTime = 100
# signal frequency (to find through fft)
signal1Frequency = 0.2
# time points (in terms of time step)
time = np.arange(beginTime, endTime, samplingInterval)
# create default functions
y1 = np.cos(2*np.pi*signal1Frequency*time)
# create subplot
figure, axis = plt.subplots(3, 1, figsize=(12,14), sharex=False)
plt.subplots_adjust(hspace=0.3)
# time domain representation for y1
axis[0].set_title('Cosine wave with a frequency of ' + str(signal1Frequency) + ' Hz')
axis[0].plot(time, y1)
axis[0].set_xlabel('Time')
axis[0].set_ylabel('Amplitude')
# create signal to fourier transform
y3 = y1**2*np.exp(-0.1*time)
mean = np.mean(y3)
y3 = y3 - np.mean(y3) # only removes 0Hz spike for cos squared without exp decay
# Time domain representation
axis[1].plot(time, y3)
axis[1].set_xlabel('Time')
axis[1].set_ylabel('Amplitude')
axis[1].set_xlim([0,10])
# Frequency domain representation
fourierTransform = np.fft.fft(y3)/len(y3) # Normalize amplitude
fourierTransform = fourierTransform[range(int(len(y3)/2))]
tpCount = len(y3)
values = np.arange(int(tpCount/2))
timePeriod = tpCount/samplingFrequency
frequencies = values/timePeriod
# Frequency domain plot
axis[2].set_title('Fourier transform depicting the frequency components')
axis[2].plot(frequencies, abs(fourierTransform))
axis[2].set_xlabel('Frequency')
axis[2].set_ylabel('Amplitude')
axis[2].set_xlim([0,1])
set_xlim([-.01, 1])
. $\omega=0$ is DC by definition. $\endgroup$