🚀 go-pugleaf

 In a previous article, rjk@sequent.uucp (Robert Kelley) writes:
>There are some obvious ways to generate 1/f time series, one way that
>comes to mind is to use an fft.  What I really want to do, is to generate
>a pseudo-random time series such that the spectral energy varies as 1/f
>(or, more generally, 1/(f**n)).
>
>Any ideas about how to do this?

Sure.  Just generate random noise (this should give you a more or less flat
frequency response when you view the FFT of it) and then run the data thru
a filter with the response you desire.  The advantage of this is that it is
done in the time domain so you can create as much data as you desire and do
it fairly fast.  You also don't need to do things in 'batches' where you would
if you tried to use the FFT in the process.  Each time you wanted to do the
conversion with an FFT, you would have to create a bunch of data points, FFT
them, Scale them, and then IFFT them.  Messy, plus lots of math which means
lots of time required.  Create an FIR filter program and give it a time domain
response that gives you the 1/f frequency domain response you desire. 
Generating the time domain impulse response of such a filter and choosing the
number of filter coefficients is left as an exercise to the interested reader.
:-)


>Robert Kelley
>rjk@sequent.com

   -dave

-- 
Dave McMahan                            mcmahan@netcom.com
					{apple,amdahl,claris}!netcom!mcmahan

rjk@sequent.uucp (Robert Kelley) writes:

> There are some obvious ways to generate 1/f time series, one way that
> comes to mind is to use an fft.  What I really want to do, is to generate
> a pseudo-random time series such that the spectral energy varies as 1/f
> (or, more generally, 1/(f**n)).

Use a random walk.  [1] x = 0.  [2] r = random integer.
[3] if r is odd, add 1 to x, else subtract 1 from x.
[4] output x. [5] goto [2].

I'm not sure how to control the exponent in 1/(f**n) tho.
Maybe increase the dimensionality of the walk?  Anybody else know?

Tom Holroyd
FAU Center for Complex Systems
tomh@bambi.ccs.fau.edu

In article <kDB7w2w163w@shark.cs.fau.edu> tomh.bbs@shark.cs.fau.edu (Tom Holroyd) writes:
>rjk@sequent.uucp (Robert Kelley) writes:

>>  What I really want to do, is to generate
>> a pseudo-random time series such that the spectral energy varies as 1/f
>> (or, more generally, 1/(f**n)).
>
>Use a random walk.  [1] x = 0.  [2] r = random integer.
>[3] if r is odd, add 1 to x, else subtract 1 from x.
>[4] output x. [5] goto [2].

I think this gives 1/f**2 rather than 1/f. Martin Gardner in Mathematical Games, Scientific American, December(?) 1978 gives an algorithm for approximating
1/f fairly easily. I'll have to go home to look it up though.

tomh.bbs@shark.cs.fau.edu (Tom Holroyd) writes:

>Use a random walk.  [1] x = 0.  [2] r = random integer.
>[3] if r is odd, add 1 to x, else subtract 1 from x.
     ^^^^^^^^^^^    
>[4] output x. [5] goto [2].

Testing if r is odd is equivalent to testing the least significant bit
of r, which in many random number generators is not very "random."
For example, I don't know what the algorithm used in ksh is, but it's
random number generator alternates between even and odd integers!  Not
a very interesting random walk.  You could test the most significant
bit instead.

>Tom Holroyd
>FAU Center for Complex Systems
>tomh@bambi.ccs.fau.edu

Steve Eddins
------------
eddins@uicbert.eecs.uic.edu 	(312) 996-5771 		FAX: (312) 413-0024
University of Illinois at Chicago, EECS Dept., M/C 154, 1120 SEO Bldg,
Box 4348, Chicago, IL  60680
-- 
Steve Eddins
------------
eddins@uicbert.eecs.uic.edu 	(312) 996-5771 		FAX: (312) 413-0024
University of Illinois at Chicago, EECS Dept., M/C 154, 1120 SEO Bldg,

The easiest way to approximate this is by taking white noise through a number
of evenly-spaced low-pass filters and adding the results, weighted according
to the frequency distribution you wish.  You can use the same noise source for
each filter, or a different one for each.  The analysis is slightly different
but both work out in the end.

The response of the filters isn't important, as long as they are identical,
with frequency scaled of course.  Simple filters will give less ripple, so a
simple RC is what you want.  That's a lossy integrator for you digital folks.
Sharper is not better in this application.

The appoximation is good because the errors cancel each other very well.  The
ripple is 0.1% even when the filters are a decade apart!  Clearly fringing
effects will be dominant here.  The filters should cover the frequency band of
interest, plus an additional one or two on each end to keep errors reasonable
at the ends.

This method works for any exponent, although there are some constraints on the
filters.  Specifically at any given frequency the output from each filter,
including weighting, should go to zero as the filter frequency goes higher or
lower.  For example for f^-3 power you need at least a second-order low-pass
filter.
--
Joe Keane, amateur mathematician