Fast 2d brown/pink noise generation?

Discussion in 'Python' started by Mathias, Mar 4, 2005.

  1. Mathias

    Mathias Guest

    Dear NG,

    I'm looking for a fast way to produce 2d-noise images with 1/f or 1/f^2
    spectrum. I currently generate the noise via inverse FFT, but since I
    need lots of data (~10^7 for a monte carlo simulation) it needs to be
    really fast. Does someone know a faster way than my approach?

    - Dimensionality is between 20x20 and 100x100
    - The spectrum doesn't need to be exactly pink/brown, an approximation
    is fine.
    - Implementation in either matlab or scientific python (LAPACK anyway)


    Thanks a lot for hints, literature, algorithms!
    Mathias
    Mathias, Mar 4, 2005
    #1
    1. Advertising

  2. Mathias

    Robert Kern Guest

    Mathias wrote:
    > Dear NG,
    >
    > I'm looking for a fast way to produce 2d-noise images with 1/f or 1/f^2
    > spectrum. I currently generate the noise via inverse FFT, but since I
    > need lots of data (~10^7 for a monte carlo simulation) it needs to be
    > really fast. Does someone know a faster way than my approach?
    >
    > - Dimensionality is between 20x20 and 100x100
    > - The spectrum doesn't need to be exactly pink/brown, an approximation
    > is fine.
    > - Implementation in either matlab or scientific python (LAPACK anyway)


    This is a 1D version that I have using scipy. It's naive, so I'm sure
    that it is slower. However, I believe the general technique can be
    implemented on a larger scale.

    The basic idea is to sum up a bunch of white time series with different
    time steps. The first level is white noise at every time step. The
    second level changes at every second time step. The third changes at
    every fourth, etc.

    I think you can replicate this by generating a few white noise arrays of
    the appropriate sizes, judiciously using repeat(), and summing them
    together. I got this scheme from an article I found by googling for pink
    noise algorithms, I believe.

    from scipy import *

    class PinkGenerator(object):
    updateTable = [0,1,0,2,0,1,0,3,0,1,0,2,0,1,0,4]
    updateTable.extend(updateTable[:-1])
    del updateTable[-1]

    def __init__(self, rng=stats.norm):
    self.key = 0
    self.rng = rng
    self.whiteValues = self.rng.rvs(size=5)

    def getNextValue(self):
    self.key += 1
    self.key = self.key % len(self.updateTable)
    self.whiteValues[self.updateTable[self.key]] = self.rng.rvs()[0]
    return (sum(self.whiteValues) + self.rng.rvs()[0])/6

    def getManyValues(self, size):
    data = zeros((size,), Float)
    for i in range(size):
    data = self.getNextValue()
    return data

    def sampleData(self, size=1024):
    data = self.getManyValues(size)
    p = power(absolute(fftshift(fft(data))), 2)/size
    f = fftshift(fftfreq(size))
    return data, f, p

    --
    Robert Kern


    "In the fields of hell where the grass grows high
    Are the graves of dreams allowed to die."
    -- Richard Harter
    Robert Kern, Mar 5, 2005
    #2
    1. Advertising

Want to reply to this thread or ask your own question?

It takes just 2 minutes to sign up (and it's free!). Just click the sign up button to choose a username and then you can ask your own questions on the forum.
Similar Threads
  1. Peter Strøiman
    Replies:
    1
    Views:
    2,070
    Peter Strøiman
    Aug 23, 2005
  2. Mathias
    Replies:
    0
    Views:
    368
    Mathias
    Mar 4, 2005
  3. Adam Hartshorne
    Replies:
    1
    Views:
    645
    Victor Bazarov
    Jul 9, 2005
  4. braskor
    Replies:
    0
    Views:
    418
    braskor
    Sep 24, 2008
  5. Jim Weirich
    Replies:
    0
    Views:
    141
    Jim Weirich
    Aug 21, 2003
Loading...

Share This Page