This algorithm produces a sequence $y$ having a magnitude spectrum that resembles a normal distribution, peaking at frequency $\omega$:
$$b[n]=(1-\sigma)b[n-1]+\frac{\sigma}{\omega} C$$ $$p[n]=p[n-1]+\omega(1+b[n])$$ $$y[n]=cos(p)$$
where $b[n]$ is a leaky integrator having feedback $1-\sigma$ and taking white noise $C$ as input, $p[n]$ is a phase accumulator having frequency $\omega$ modulated by brownian noise $b[n]$, and y[n] is the final signal, cosine at phase $p[n]$.
Python function:
def signal(f, w, n):
"""
Arguments
f - peak frequency in radians/sample
w - bandwidth of frequency content
n - length
Returns
sig - list of samples
"""
br = 0
p = 0
sig = [0] * n
for i in range(n):
# br - leaky integrator having rate proportional to bandwidth,
# integrating white noise at magnitude inversely proportional
# to frequency.
br = (1 - w ) * br + w * random.uniform(-1, 1) / max([0.001, f])
# p - phase increments at frequency f modulated by brownian noise.
p = p + (1 + br) * f
sig[i] = math.cos(p)
return sig
To illustrate:
averaged magnitude spectrum at swept frequency from [0,π)
averaged magnitude spectrum with swept bandwidth from [0,0.5)
fourier transform of a sample sequence, hue indicates phase
audio sample at 1760hz with $\sigma=0.005$
averaged eight times
There are undesirable characteristics:
- The amplitude of the peak frequency has an unknown relationship to $\sigma$, the second animation shows the sharp curved response
- $\omega$ is not a perfect predictor of peak frequency, $\sigma$ allows noise to perturb the phase towards a range of $[0\omega,2\omega]$ from FM
Short of using optimisation and curve-fitting, how can the frequency content be predicted (e.g. by the gaussian function) and controlled?
Any additional insight appreciated