Reply by Martin Eisenberg●October 20, 20082008-10-20

Andor wrote:

> On 15 Okt., 18:16, Sampo Niskanen <spnis...@cc.hut.fi> wrote:

>> The simulation time scale, however, is quite short, so I don't
>> want a gust of wind to last the whole simulation duration (that
>> would effectively change the average wind speed).

> Ok, I see. I think I had a similar discussion with somebody here
> some time ago who wanted to use "pink noise" for audio testing
> purposes. He was also interested in generating series with PSD
> proportional to
>
> 1/(1+ K f),
>
> with a finite cutoff, going flat towards DC, and having the 1/f
> property for audible frequencies.

Another situation where a low-frequency cutoff would be desirable is
when the generating mechanism is intermittency. Some nonlinear
systems have an operating regime where the state grows very slowly,
until it crosses a threshold inducing bursty output, which eventually
returns the state to the slow region. If the duration of quiescent
intermissions is unbounded then the system can exhibit true low-
frequency divergence.
Now, feeding the intermittent signal to an allpass filter can turn
bursts into smoother oscillation suitable for modulating musical
quantities, like pitch. But we don't want that viola synth to sound
like an '80s Casio every now and then, so the intermission length
must not much exceed the effective length of the allpass impulse
response in this application.
I don't know if anyone has published on this in the musical context,
but here are two articles describing intermittent iteration
processes, one free and one not:
http://www.istia.univ-angers.fr/~chapeau/papers/psiplrc2.pdfhttp://prola.aps.org/abstract/PRA/v31/i3/p1830_1
Martin
--
When he had finally achieved a position that
allowed him to say everything he thought, he
only thought of his position anymore.
--Gabriel Laub

Reply by Sampo Niskanen●October 16, 20082008-10-16

Andor <andor.bariska@gmail.com> wrote:

> Ok, I see. I think I had a similar discussion with somebody here some
> time ago who wanted to use "pink noise" for audio testing purposes. He
> was also interested in generating series with PSD proportional to

> 1/(1+ K f),

> with a finite cutoff, going flat towards DC, and having the 1/f
> property for audible frequencies. Perhaps it would be a suitable
> nomenclature to call that kind (with finite cutoff frequency towards
> DC) of series "pink noise", as compared to 1/f^alpha noise.

That's true. However, as seen on the DSP generation of Pink Noise web
site[1] "pink noise" typically refers simply to noise with PSD
proportional to 1/f instead of 1/f^alpha. Therefore some distinction
must be made between these two. (Of course the correct way would be to
always talk about 1/f^alpha, with alpha=1 being a special case, but this
is currently uncommon.)
In any case, my code allows the limiting frequency to be made
arbitrarily low (though of course this slows down the generation a bit,
as generating a sample is O(n) for the number of poles).
[1] http://www.firstpr.com.au/dsp/pink-noise/
--
__________________________________________________
/____\ Sampo Niskanen <=> sampo.niskanen@iki.fi \
\ http://www.iki.fi/sampo.niskanen/ \
\ ________________________________________\___
\___/___________________________________________/

Reply by Andor●October 15, 20082008-10-15

On 15 Okt., 18:16, Sampo Niskanen <spnis...@cc.hut.fi> wrote:

> Andor <andor.bari...@gmail.com> wrote:
> > This sounds a bit contradicting. Either you are intersted in 1/f^{5/3}
> > noise, or you are worried about values in the series deviating from
> > zero for long periods.
> > You can't have both.
>
> The application in my case was simulating wind turbulence. �Wind
> turbulence freqency has a power spectrum density proportional to
> 1/(1+K*f)^(5/3), which for large f is equal to 1/f^(5/3). �The
> simulation time scale, however, is quite short, so I don't want a gust
> of wind to last the whole simulation duration (that would effectively
> change the average wind speed).
>
> By choosing the number of poles suitably one can choose how much low
> frequency components to include. �The original (empirical) formula
> doesn't go to infinity at f=0 either, so one could estimate it with a
> suitable number of poles. �In my case I wanted even less low-frequency
> components, and using 2 poles with a 20Hz sampling rate yields maximum
> wind gust lengths of approximately 3-5 seconds (the spectrum turns flat
> below 0.3Hz). �So in my application, it is desireable that the noise is
> pink at the high end of the spectrum, but flat at the low end.

Ok, I see. I think I had a similar discussion with somebody here some
time ago who wanted to use "pink noise" for audio testing purposes. He
was also interested in generating series with PSD proportional to
1/(1+ K f),
with a finite cutoff, going flat towards DC, and having the 1/f
property for audible frequencies. Perhaps it would be a suitable
nomenclature to call that kind (with finite cutoff frequency towards
DC) of series "pink noise", as compared to 1/f^alpha noise. By
definition, 1/f^alpha noise has a PSD satisfying
lim_{f->0} ( P(f) / [c |f|^{-alpha}] ) = 1
for some constant c, so only the behaviour towards DC is important (as
opposed to "pink noise", where only the behaviour towards infinity is
important).
Regards,
Andor

Reply by Sampo Niskanen●October 15, 20082008-10-15

Andor <andor.bariska@gmail.com> wrote:

> This sounds a bit contradicting. Either you are intersted in 1/f^{5/3}
> noise, or you are worried about values in the series deviating from
> zero for long periods.

> You can't have both.

The application in my case was simulating wind turbulence. Wind
turbulence freqency has a power spectrum density proportional to
1/(1+K*f)^(5/3), which for large f is equal to 1/f^(5/3). The
simulation time scale, however, is quite short, so I don't want a gust
of wind to last the whole simulation duration (that would effectively
change the average wind speed).
By choosing the number of poles suitably one can choose how much low
frequency components to include. The original (empirical) formula
doesn't go to infinity at f=0 either, so one could estimate it with a
suitable number of poles. In my case I wanted even less low-frequency
components, and using 2 poles with a 20Hz sampling rate yields maximum
wind gust lengths of approximately 3-5 seconds (the spectrum turns flat
below 0.3Hz). So in my application, it is desireable that the noise is
pink at the high end of the spectrum, but flat at the low end.
--
__________________________________________________
/____\ Sampo Niskanen <=> sampo.niskanen@iki.fi \
\ http://www.iki.fi/sampo.niskanen/ \
\ ________________________________________\___
\___/___________________________________________/

Reply by Andor●October 15, 20082008-10-15

> Note that your alpha is > 1, indicating non-stationary long
> range correlation. For that case long deviations from zero
> is the rule, not the exception.

You can consider alpha = 2 (Wiener process). In that case, the
expected time between zero crossings is infinite. 5/3 isn't much
smaller than 6/3.

Reply by Andor●October 15, 20082008-10-15

Sampo Niskanen wrote:

> I haven't compared the AR and MA methods either. �I'm no DSP expert and
> once I had a working filter that produced good enough results I went
> back to my original project. �:) �In fact, I decided to use only a
> two-pole filter in my application since I didn't want the values
> deviating from zero for long periods of time.

Hi Sampo
This sounds a bit contradicting. Either you are intersted in 1/f^{5/3}
noise, or you are worried about values in the series deviating from
zero for long periods.
You can't have both.
This deviaton for longer periods from zero is one of the very
characteristics that makes 1/f^alpha noise "statistically difficult"
to handle. This phenomenon was called the "Joseph Effect" by
Mandelbrot, because of Joseph's dream that Egypt was to experience
seven years of plentiful harvests followed by seven years of famine
(long periods of deviation from the zero mean). Particularly
climatological and hydrological time series can be modeled very well
by 1/f^alpha noise. The Nile river data has an estimated Hurst
exponent of H = 0.91, which translates to an alpha of 0.82. Note that
your alpha is > 1, indicating non-stationary long range correlation.
For that case long deviations from zero is the rule, not the
exception.
Just a thought.

> Anyway, I hope the code is useful to others as well.

I'm sure it is :-).
Regards,
Andor

Reply by Sampo Niskanen●October 15, 20082008-10-15

Thanks for all the input.

> 1) When alpha = 1, how do your results compare to the Voss-McCartney
> algorithm?

I haven't compared this method to the Voss-McCartney method. In my
application I required 1/alpha^(5/3) noise, which cannot be generated by
the Voss-McCartney algorithm (or at least I couldn't figure out how to
modify it in such a way).
I haven't compared the AR and MA methods either. I'm no DSP expert and
once I had a working filter that produced good enough results I went
back to my original project. :) In fact, I decided to use only a
two-pole filter in my application since I didn't want the values
deviating from zero for long periods of time.
Anyway, I hope the code is useful to others as well.
--
__________________________________________________
/____\ Sampo Niskanen <=> sampo.niskanen@iki.fi \
\ http://www.iki.fi/sampo.niskanen/ \
\ ________________________________________\___
\___/___________________________________________/

Reply by Jerry Avins●October 14, 20082008-10-14

Stacy wrote:

> 1) When alpha = 1, how do your results compare to the Voss-McCartney
> algorithm?
>
> 2) Is this the formula for the rolloff?
> rolloff = (alpha * -10dB/decade)?

That's how I interpret it.

> 3) How do you convert from db/decade to db/octave?
> i.e. -10dB/decade is 3.0102999 dB/octave.

Multiply by .30102999 (.3 is close enough).
Jerry
--
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
** Posted from http://www.teranews.com **

Reply by Stacy●October 14, 20082008-10-14

1) When alpha = 1, how do your results compare to the Voss-McCartney
algorithm?
2) Is this the formula for the rolloff?
rolloff = (alpha * -10dB/decade)?
3) How do you convert from db/decade to db/octave?
i.e. -10dB/decade is 3.0102999 dB/octave.
Thanks for posting your code.

Reply by SteveSmith●October 13, 20082008-10-13

Hi Sampo,
Thanks for the presenting the material, very interesting and useful.
Steve