lundi 27 mars 2017

Amplitude spectrum change after altering phases of Fourier coefficients

I'm trying to generate some random 1d sequences (r) by changing the phases of the Fourier coefficients of a given array (y). I assume the amplitude spectrum of these two sequences (y and r) won't change after doing this, and they indeed don't if I use numpy.fft.fft(), like this:

import numpy as np
import numpy.fft as fft

n=512
x=np.linspace(0,10*np.pi,n)
y=np.sin(x)+2*np.sin(2*x)+2*np.sin(5*x)

#-------------Get Fourier coefficients-------------
cf1=fft.fft(y)
amp1=np.abs(cf1)
theta1=np.angle(cf1)

#-Randomly alter phase keeping amplitude unchanged-
theta2=np.random.normal(0,1.,size=theta1.shape)+theta1
cf2=amp1*np.cos(theta2)+1j*amp1*np.sin(theta2)

#-----------------------IFFT-----------------------
y2=fft.ifft(cf2)

#------------Compare amplitude spectrum------------
cf3=fft.fft(y2)
amp2=np.abs(cf3)

import matplotlib.pyplot as plt
figure=plt.figure()
ax=figure.add_subplot(111)
ax.plot(amp1-amp2,'k-')

plt.show(block=False)

The resultant plot is a random sequence at the order of 1e-13. So effectively unchanged.

But I'm interested in generating random real data rather than complex. Using numpy.fft.rfft() and numpy.fft.ifft() instead, the amplitudes agree at all frequencies but the last one (amp1[-1] and amp2[-1]), and the difference is at the order of 0.1. According to the docs this correspondes to the Nyquist frequency, and the difference doesn't go away if the data size is an odd number. I don't understand what is causing the difference or how should I generate a real valued array with identical amplitude spectrum.

Thanks in advance.




Aucun commentaire:

Enregistrer un commentaire