There is probably something trivial behind this, but I am missing something. I need to create a stationary random time series data v(t) which is the the sum of another time series u(t) and u(t) with a time shift $\tau$, i.e
$v(t) = u(t) + u(t - \tau)$
The series u(t) is itself a stationary colored Gaussian time-series with some PSD $S_u(f)$. No problem! I wrote a small function (using Python) to make u(t) for a given PSD. I calculate the PSD of the the generated time series using Welch's method and it checks out.
If we take the Fourier transform of u(t - \tau) it should just gain an extra phase compared to the Fourier transform of u(t)
$F[u(t - \tau)] = \int dt \, u(t -\tau) \, e^{2 \pi j f t } = \int dt \, u(t) \, e^{2 \pi j f t } e^{- 2 \pi j f \tau } = \tilde{u}(f) e^{- 2 \pi j f \tau }$
Since $v(t) = u(t) + u(t - \tau)$,
$\tilde{v}(f) = \tilde{u} (f) \, (1 + e^{- 2 \pi j f \tau }) $
And the PSD of V should then be (using Euler identity)
$S_v (f) = 4 S_u (f) \cos^2 (\pi f \tau )$
In particular $S_v$ should vanish when $f = \frac{2n + 1}{2 \tau}$ with the smallest zero being at $\frac{1}{2 \tau}$. However when I test my this with the u(t) series I generated above the PSD looks different. There are zeros low frequencies which shouldn't exist. For example here is the PSD of the v(t) series and the spectra I would expect for $\tau = 2$. (The first zero should be at 0.25 Hz)
Moreover the zeros in the data seems to change with the sampling frequency I am using. The v(t) series for the above plot has a sampling frequency of 0.25 Hz. Here is one with 0.5 Hz, keeping everything else the same.
I am kind of baffled by this. Any help would be greatly appreciated. Sorry about the long post, and thanks in advance!!
EDIT: After using a higher sample rate as suggested by Hilmar below, and fixing the error I described in the comment everything looks good. Spectrum for the fixed case below.