I have a model function for a power spectral density (PSD). I want to simulate a time series with a power PSD equal to my model function. After some browsing about how to do this, I came up with the following method:
Discretize my continuous PSD
Convert the discrete PSD into amplitude points given by $\sqrt{N F_s S}$, where $N$ is the number of points in the discrete PSD, $F_s$ is the sampling rate (two times the max frequency of the discrete PSD), and $S$ is a point from the discrete PSD.
Give each amplitude point a randomly chosen phase. To do this, I multiplied each point by $\exp(-j \phi)$, where $\phi$ is a random value from $0$ to $2\pi$ (chosen from a uniform distribution).
I read some forums where it was advocated that one can now perform an inverse fast Fourier transform and arrive at the appropriate time series. However, the resulting time series is complex, so I added this step (as suggested in a forum) to deal with this:
I now have a vector of complex numbers as a function of positive frequency values, running from $0$ to, say, $f_{max}$. I then constructed a new vector with positive and negative frequencies (running from $-f_{max}$ to $f_{max}$) with the following structure. If $A$ is an amplitude point at frequency $f_1$, then the corresponding amplitude point at $-f_1$ is $A^{*}$ (the complex conjugate of $A$).
I now perform an inverse fast Fourier transform of this. I still get an imaginary part but its much smaller than the real part so I discard it. I then treat this as my time series. The length is $2 N - 1$ points.
To confirm this approach, I compute the PSD (i.e., the periodogram) of the time series and verify that it agrees perfectly with my original continuous PSD model.
I have some questions about this approach:
Is this a reasonable method? Can anyone confirm whether it makes sense or whether it is flawed?
Is there a reference for this that you can share? If you know of a better approach, would you kindly let me know?