On Mar 23, 2:46 pm, Piergiorgio Sartor <piergiorgio.sartor.this.should.not.be.u...@nexgo.REMOVETHIS.de> wrote:> maxlittle2...@googlemail.com wrote: > > After having spent a frustrating hour searching around to find some > > simple code that would generate 1/f noise, or more generally 1/f^alpha > > noise, I decided to write it myself. So here goes: > > >http://www.eng.ox.ac.uk/samp/powernoise_soft.html > > > I hope you enjoy it. Comments and feedback welcome. > > Nice tool, thanks! > > I do, anyway, have a question. > > How did you test the function in order to be sure > it is generating the noise shape you wanted? > > Thanks again, > > bye, > > -- > > piergiorgioHi piergiorgio Thanks for the comments. Good question about testing -- I did two tests: (a) plots of PSD on log-log scales for a variety of alpha values, (b) independent tests of estimated scaling exponent using detrended fluctuation analysis: http://en.wikipedia.org/wiki/Detrended_fluctuation_analysis http://www.eng.ox.ac.uk/samp/dfa_soft.html Also, it's been tested by a couple of colleagues. No problems so far. Thanks, Max
Generating 1/f (flicker/pink) noise
Started by ●March 21, 2008
Reply by ●March 23, 20082008-03-23
Reply by ●March 23, 20082008-03-23
On Mar 22, 4:03 pm, maxlittle2...@googlemail.com wrote:> ... > Ah, well I did find those resources, but no simple "one pager" of > Matlab code for generating perfect 1/f^alpha noise. > ... > MaxThere can be no page that is 'simple' and 'perfect' for all applications. One characteristic of your implementation choice is that the power spectrum does not just have the desired expectation value in any frequency region, but it equals the expectation value. While this is a possible power spectum value, it is not a likely one, and is unlikely to repeat. This can be useful in certain testing situations. It might be inappropriate in others. An alternative is to use independent Gaussian samples for real and imaginary frequency components and weight them with the desired frequency weighting. An example implementation of this, modified from your implementation, follows below. This could provide more realistic statistical behavior in some applications. Dale B. Dalrymple http://dbdimages.com function x = powernoisen(alpha, N) % Generate samples of power law noise. The power spectrum % of the signal scales as f^(-alpha), and the phases are uniformly % distributed in the range -pi to +pi. % 20080323: % The power spectrum is now Chi-squared distributed. % % Useage: % x = powernoisen(alpha, N) % % Inputs: % alpha - power law scaling exponent % N - number of samples to generate % % Output: % x - N x 1 vector of power law samples % % Original: % (cc) Max Little, 2008. This software is licensed under the % Attribution-Share Alike 2.5 Generic Creative Commons license: % http://creativecommons.org/licenses/by-sa/2.5/ % If you use this work, please cite: % Little MA et al. (2007), "Exploiting nonlinear recurrence and fractal % scaling properties for voice disorder detection", Biomed Eng Online, 6:23 % % As of 20080323 markup % If you use this work, consider saying hi on comp.dsp % Dale B. Dalrymple N2 = floor(N/2)-1; f = (2:(N2+1))'; A2 = 1./(f.^(alpha/2)); %p2 = (rand(N2,1)-0.5)*2*pi; % 20080323 p2 = randn(N2,1) + i * randn(N2,1); % 20080323 %d2 = A2.*exp(i*p2); % 20080323 d2 = A2.*p2; % 20080323 d = [1; d2; 1/((N2+2)^alpha); flipud(conj(d2))]; x = real(ifft(d));
Reply by ●March 24, 20082008-03-24
On 23 Mar, 18:33, dbd <d...@ieee.org> wrote:> On Mar 22, 4:03 pm, maxlittle2...@googlemail.com wrote: > > > ... > > Ah, well I did find those resources, but no simple "one pager" of > > Matlab code for generating perfect 1/f^alpha noise. > > ... > > Max > > There can be no page that is 'simple' and 'perfect' for all > applications. > > One characteristic of your implementation choice is that the power > spectrum does not just have the desired expectation value in any > frequency region, but it equals the expectation value. While this is a > possible power spectum value, it is not a likely one, and is unlikely > to repeat. This can be useful in certain testing situations. It might > be inappropriate in others. An alternative is to use independent > Gaussian samples for real and imaginary frequency components and > weight them with the desired frequency weighting. An example > implementation of this, modified from your implementation, follows > below. This could provide more realistic statistical behavior in some > applications. > > Dale B. Dalrymplehttp://dbdimages.comDale Thanks, these are excellent points. I've added an option to the routine so that a choice can be made between deterministic and stochastic power. I've also added an option that normalises the signal amplitude so that it lies in the range -1 to +1. This may be useful in certain circumstances. Grab the code from: http://www.eng.ox.ac.uk/samp/powernoise_soft.html Max
Reply by ●March 25, 20082008-03-25
On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2007@googlemail.com wrote:>Dear All > >After having spent a frustrating hour searching around to find some >simple code that would generate 1/f noise, or more generally 1/f^alpha >noise, I decided to write it myself. So here goes: > >http://www.eng.ox.ac.uk/samp/powernoise_soft.html > >I hope you enjoy it. Comments and feedback welcome. > >All the best, >MaxYou might like to point out on your webpage that this method of generating pink noise has has O(nlog(n)) time complexity and O(n) storage requirements for generating n samples. This limits it to relatively small sample sizes (which, nevertheless, may be large enough for a particular need). Other methods (e.g. additive or subtractive synthesis, mentioned elsewhere in this thread) have O(1) time and storage requirements, and are suitable for streaming continuous data. Regards, Allan
Reply by ●March 25, 20082008-03-25
On 25 Mrz., 14:49, Allan Herriman <allanherri...@hotmail.com> wrote:> On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2...@googlemail.com > wrote: > > >Dear All > > >After having spent a frustrating hour searching around to find some > >simple code that would generate 1/f noise, or more generally 1/f^alpha > >noise, I decided to write it myself. So here goes: > > >http://www.eng.ox.ac.uk/samp/powernoise_soft.html > > >I hope you enjoy it. Comments and feedback welcome. > > >All the best, > >Max > > You might like to point out on your webpage that this method of > generating pink noise has has O(nlog(n)) time complexity and O(n) > storage requirements for generating n samples. �This limits it to > relatively small sample sizes (which, nevertheless, may be large > enough for a particular need). > > Other methods (e.g. additive or subtractive synthesis, mentioned > elsewhere in this thread) have O(1) time and storage requirements, and > are suitable for streaming continuous data.I find that hard to believe. The storage and the processing time is independent of the number of samples to be computed? Regards, Andor
Reply by ●March 25, 20082008-03-25
On Tue, 25 Mar 2008 07:23:21 -0700 (PDT), Andor <andor.bariska@gmail.com> wrote:>On 25 Mrz., 14:49, Allan Herriman <allanherri...@hotmail.com> wrote: >> On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2...@googlemail.com >> wrote: >> >> >Dear All >> >> >After having spent a frustrating hour searching around to find some >> >simple code that would generate 1/f noise, or more generally 1/f^alpha >> >noise, I decided to write it myself. So here goes: >> >> >http://www.eng.ox.ac.uk/samp/powernoise_soft.html >> >> >I hope you enjoy it. Comments and feedback welcome. >> >> >All the best, >> >Max >> >> You might like to point out on your webpage that this method of >> generating pink noise has has O(nlog(n)) time complexity and O(n) >> storage requirements for generating n samples. �This limits it to >> relatively small sample sizes (which, nevertheless, may be large >> enough for a particular need). >> >> Other methods (e.g. additive or subtractive synthesis, mentioned >> elsewhere in this thread) have O(1) time and storage requirements, and >> are suitable for streaming continuous data. > >I find that hard to believe. The storage and the processing time is >independent of the number of samples to be computed?Short term brain malfunction, sorry. I meant to say: O(n) time, and O(1) storage (if streaming) or O(n) storage (if saving the data somewhere). Regards, Allan
Reply by ●March 25, 20082008-03-25
On Mar 24, 4:51 pm, maxlittle2...@googlemail.com wrote:> On 23 Mar, 18:33, dbd <d...@ieee.org> wrote: > > > > > On Mar 22, 4:03 pm, maxlittle2...@googlemail.com wrote: > > > > ... > > > Ah, well I did find those resources, but no simple "one pager" of > > > Matlab code for generating perfect 1/f^alpha noise. > > > ... > > > Max > > > There can be no page that is 'simple' and 'perfect' for all > > applications. > > > One characteristic of your implementation choice is that the power > > spectrum does not just have the desired expectation value in any > > frequency region, but it equals the expectation value. While this is a > > possible power spectum value, it is not a likely one, and is unlikely > > to repeat. This can be useful in certain testing situations. It might > > be inappropriate in others. An alternative is to use independent > > Gaussian samples for real and imaginary frequency components and > > weight them with the desired frequency weighting. An example > > implementation of this, modified from your implementation, follows > > below. This could provide more realistic statistical behavior in some > > applications. > > > Dale B. Dalrymplehttp://dbdimages.com > > Dale > > Thanks, these are excellent points. I've added an option to the > routine so that a choice can be made between deterministic and > stochastic power. I've also added an option that normalises the > signal amplitude so that it lies in the range -1 to +1. This may be > useful in certain circumstances. Grab the code from: > > http://www.eng.ox.ac.uk/samp/powernoise_soft.html > > MaxMax If you really think a normalization and a peak limiting capability should be included inside the function, I suggest you do so with a gain and a clipping capability instead of putting the mean and scale of the output under the control of the extremal behavior of the noise generated. There are test purposes where the original function or my addition could be called multiple times to generate longer sequences of test data. Even these special cases no longer work when gain and offset change from block to block. Dale B. Dalrymple
Reply by ●March 25, 20082008-03-25
On 25 Mrz., 15:43, Allan Herriman <allanherri...@hotmail.com> wrote:> On Tue, 25 Mar 2008 07:23:21 -0700 (PDT), Andor > > > > > > <andor.bari...@gmail.com> wrote: > >On 25 Mrz., 14:49, Allan Herriman <allanherri...@hotmail.com> wrote: > >> On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2...@googlemail.com > >> wrote: > > >> >Dear All > > >> >After having spent a frustrating hour searching around to find some > >> >simple code that would generate 1/f noise, or more generally 1/f^alpha > >> >noise, I decided to write it myself. So here goes: > > >> >http://www.eng.ox.ac.uk/samp/powernoise_soft.html > > >> >I hope you enjoy it. Comments and feedback welcome. > > >> >All the best, > >> >Max > > >> You might like to point out on your webpage that this method of > >> generating pink noise has has O(nlog(n)) time complexity and O(n) > >> storage requirements for generating n samples. �This limits it to > >> relatively small sample sizes (which, nevertheless, may be large > >> enough for a particular need). > > >> Other methods (e.g. additive or subtractive synthesis, mentioned > >> elsewhere in this thread) have O(1) time and storage requirements, and > >> are suitable for streaming continuous data. > > >I find that hard to believe. The storage and the processing time is > >independent of the number of samples to be computed? > > Short term brain malfunction, sorry. � > > I meant to say: O(n) time, and O(1) storage (if streaming) or O(n) > storage (if saving the data somewhere).Ok, that makes more sense. Inspired by this thread (and other circumstances), I have read the Kasdin reference that I posted above somewhere. There it is shown that the LTI filtering techniques to construct 1/f noise either 1. have a finite number of coefficients (either transversal or recursive) and consequently, the 1/f-law breaks down for frequencies below some cutoff frequency or 2. use filters with growing number of filter coefficients (either transversal or recursive) as the simulation proceeds, resulting in sequences with asymptotically unbiased PSDs (true 1/f spectra). This seems natural with my understanding of the autocorrelation function of 1/f-noises, which only decays hyperbolically. Using filtered sequences with a constant number of coefficients, the autocorrelation function is either truncated (FIR filters) or decays exponentially (stable rational transfer function filter) or does not decay at all (unstable rational transfer function filters, for example to construct 1/f^2 noise). Are any of the methods that only require O(n) computation time asymptotically unbiased? Regards, Andor
Reply by ●March 25, 20082008-03-25
On Tue, 25 Mar 2008 08:03:25 -0700 (PDT), Andor <andor.bariska@gmail.com> wrote:>On 25 Mrz., 15:43, Allan Herriman <allanherri...@hotmail.com> wrote: >> On Tue, 25 Mar 2008 07:23:21 -0700 (PDT), Andor >> >> >> >> >> >> <andor.bari...@gmail.com> wrote: >> >On 25 Mrz., 14:49, Allan Herriman <allanherri...@hotmail.com> wrote: >> >> On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2...@googlemail.com >> >> wrote: >> >> >> >Dear All >> >> >> >After having spent a frustrating hour searching around to find some >> >> >simple code that would generate 1/f noise, or more generally 1/f^alpha >> >> >noise, I decided to write it myself. So here goes: >> >> >> >http://www.eng.ox.ac.uk/samp/powernoise_soft.html >> >> >> >I hope you enjoy it. Comments and feedback welcome. >> >> >> >All the best, >> >> >Max >> >> >> You might like to point out on your webpage that this method of >> >> generating pink noise has has O(nlog(n)) time complexity and O(n) >> >> storage requirements for generating n samples. �This limits it to >> >> relatively small sample sizes (which, nevertheless, may be large >> >> enough for a particular need). >> >> >> Other methods (e.g. additive or subtractive synthesis, mentioned >> >> elsewhere in this thread) have O(1) time and storage requirements, and >> >> are suitable for streaming continuous data. >> >> >I find that hard to believe. The storage and the processing time is >> >independent of the number of samples to be computed? >> >> Short term brain malfunction, sorry. � >> >> I meant to say: O(n) time, and O(1) storage (if streaming) or O(n) >> storage (if saving the data somewhere). > >Ok, that makes more sense. Inspired by this thread (and other >circumstances), I have read the Kasdin reference that I posted above >somewhere. There it is shown that the LTI filtering techniques to >construct 1/f noise either > >1. have a finite number of coefficients (either transversal or >recursive) and consequently, the 1/f-law breaks down for frequencies >below some cutoff frequency > >or > >2. use filters with growing number of filter coefficients (either >transversal or recursive) as the simulation proceeds, resulting in >sequences with asymptotically unbiased PSDs (true 1/f spectra). > >This seems natural with my understanding of the autocorrelation >function of 1/f-noises, which only decays hyperbolically. Using >filtered sequences with a constant number of coefficients, the >autocorrelation function is either truncated (FIR filters) or decays >exponentially (stable rational transfer function filter) or does not >decay at all (unstable rational transfer function filters, for example >to construct 1/f^2 noise). > >Are any of the methods that only require O(n) computation time >asymptotically unbiased?Yes, the Voss pink noise generator (additive synthesis) requires O(n) computation time (specifically: two calls to rand() and a small amount of overhead per output sample), regardless of the output bandwidth. The storage requirement is, however, directly proportional to the number of octaves generated. As an engineer (rather than a mathematician) I could claim that a modest amount of storage would get you close enough to any desired accuracy at any frequency close to (but not including) DC. Regards, Allan
Reply by ●March 25, 20082008-03-25
On Wed, 26 Mar 2008 04:29:21 +1100, Allan Herriman <allanherriman@hotmail.com> wrote:>On Tue, 25 Mar 2008 08:03:25 -0700 (PDT), Andor ><andor.bariska@gmail.com> wrote: > >>On 25 Mrz., 15:43, Allan Herriman <allanherri...@hotmail.com> wrote: >>> On Tue, 25 Mar 2008 07:23:21 -0700 (PDT), Andor >>> >>> >>> >>> >>> >>> <andor.bari...@gmail.com> wrote: >>> >On 25 Mrz., 14:49, Allan Herriman <allanherri...@hotmail.com> wrote: >>> >> On Fri, 21 Mar 2008 10:38:39 -0700 (PDT), maxlittle2...@googlemail.com >>> >> wrote: >>> >>> >> >Dear All >>> >>> >> >After having spent a frustrating hour searching around to find some >>> >> >simple code that would generate 1/f noise, or more generally 1/f^alpha >>> >> >noise, I decided to write it myself. So here goes: >>> >>> >> >http://www.eng.ox.ac.uk/samp/powernoise_soft.html >>> >>> >> >I hope you enjoy it. Comments and feedback welcome. >>> >>> >> >All the best, >>> >> >Max >>> >>> >> You might like to point out on your webpage that this method of >>> >> generating pink noise has has O(nlog(n)) time complexity and O(n) >>> >> storage requirements for generating n samples. �This limits it to >>> >> relatively small sample sizes (which, nevertheless, may be large >>> >> enough for a particular need). >>> >>> >> Other methods (e.g. additive or subtractive synthesis, mentioned >>> >> elsewhere in this thread) have O(1) time and storage requirements, and >>> >> are suitable for streaming continuous data. >>> >>> >I find that hard to believe. The storage and the processing time is >>> >independent of the number of samples to be computed? >>> >>> Short term brain malfunction, sorry. � >>> >>> I meant to say: O(n) time, and O(1) storage (if streaming) or O(n) >>> storage (if saving the data somewhere). >> >>Ok, that makes more sense. Inspired by this thread (and other >>circumstances), I have read the Kasdin reference that I posted above >>somewhere. There it is shown that the LTI filtering techniques to >>construct 1/f noise either >> >>1. have a finite number of coefficients (either transversal or >>recursive) and consequently, the 1/f-law breaks down for frequencies >>below some cutoff frequency >> >>or >> >>2. use filters with growing number of filter coefficients (either >>transversal or recursive) as the simulation proceeds, resulting in >>sequences with asymptotically unbiased PSDs (true 1/f spectra). >> >>This seems natural with my understanding of the autocorrelation >>function of 1/f-noises, which only decays hyperbolically. Using >>filtered sequences with a constant number of coefficients, the >>autocorrelation function is either truncated (FIR filters) or decays >>exponentially (stable rational transfer function filter) or does not >>decay at all (unstable rational transfer function filters, for example >>to construct 1/f^2 noise). >> >>Are any of the methods that only require O(n) computation time >>asymptotically unbiased? > >Yes, the Voss pink noise generator (additive synthesis) requires O(n) >computation time (specifically: two calls to rand() and a small amount >of overhead per output sample), regardless of the output bandwidth.Minor clarification: for each output sample, there are two calls to rand(); a small amount of overhead and the summation of a list of saved random numbers. The length of that list is proportional to the number of octaves generated. So I guess my claim of "regardless of the output bandwidth" is false if a naive implementation is used. Since only one element in the list changes each sample, one could use the old "boxcar filter recursive FIR" trick (as used in CIC filters) to give a constant amount of processing per output sample.>The storage requirement is, however, directly proportional to the >number of octaves generated. >As an engineer (rather than a mathematician) I could claim that a >modest amount of storage would get you close enough to any desired >accuracy at any frequency close to (but not including) DC. > >Regards, >Allan
