Fast Fourier Transform for Harmonic Analysis

Question

I'm analyzing the harmonics present in the wave motion as a function of where along the string the pluck occurs. I hope to obtain a plot like those exhibited on this page: https://softwaredevelopmentperestroika.wordpress.com/2013/12/10/fast-fourier-transforms-with-python-the-noise-and-the-signal/. To do this, I've written code modelling an asymmetric triangle and implemented numpy's fft. The data being output, however, is not what I expect, and it's peaked about a frequency of 0 Hz. Here is my code and its output:

from numpy.fft import fft as npfft, fftfreq as npfftfreq
#triangular pulse
def triangular_pulse(x, xmean, sigma):
    for i in x:
        if x[i]<=xmean:
            y[i] = x[i]*(sigma/xmean)
        else :
            y[i] = sigma-(x[i]-xmean)*(sigma/(200-xmean))
    return y

N_masses = 200
T  = 0.0669264714
mu = .03937
cSq = T/mu
c  = np.sqrt(cSq)
dx = 1.0

dt = dx/c
print dt

#Initialize some arrays
x0  = np.arange(N_masses)*dx
y   = np.zeros(N_masses)
vy  = np.zeros(N_masses)
ay  = np.zeros(N_masses)

#Set Initial conditions (pluck)
# # half-pluck
# y = 30*gaussian_pulse(x0,x0[N_masses/2],2)

# quarter-pluck
y = triangular_pulse(x0,x0[N_masses/4],1)

rhat=npfft(y)
freaq=npfftfreq(len(y),dt)
plt.plot(freaq,np.abs(rhat)/len(rhat))
plt.show()

Please let me know if you spot the source of my problem. Thanks!

Update

Added y = triangular_pulse(x0,x0[N_masses/40],1) y-=np.mean(y) with the result of a broader non-zero band; peak is still centered around "0", however.

Answer 1

Just subtract the mean of the signal before running the frequency analysis, i.e. after calling triangular_pulse:

y-=y.mean()

and you will obtain the peak at a non-zero frequency. This is because the signal has a mean component that is not zero, that will show up as the component at zero frequency.

EDIT: as a comment, you can rewrite the triangular pulse function using numpy where:

def triangular_pulse2(x,xmean,sigma):
    return where(x<xmean,x*sigma/xmean,sigma-(x-xmean)*(sigma/(200-xmean)))