Logo coherent WaveBurst  
Library Reference Guide
Logo
testWDM_1.C
Go to the documentation of this file.
1 { // apply transformation on white noise and plot it
2 
3  int Rate = 1024; // sampling rate (Hz)
4  int Duration = 100; // duration (seconds)
5  int N = Rate*Duration; // number of samples
6  wavearray<double> ts(N); //time series container
7  ts.rate(Rate);
8 
9  // time series is filled with white noise data:
10  TRandom3 rnd(0);
11  for(int i=0; i<N; i++) ts[i] = rnd.Gaus();
12 
13  // produce the TF map:
14  WDM<double> wdm(32, 64, 4, 8); // define a WDM transform (32 bands)
15  WSeries<double> tfmap; // TF map container
16  tfmap.Forward(ts, wdm); // apply the WDM to the time series
17 
18  // another example: add a sine wave @ 100Hz:
19  double T = 1./Rate;
20  double omega = TMath::TwoPi()*100;
21  for(int i=0; i<N; i++) ts[i] += sin(omega*i*T);
22  tfmap.Forward(ts, wdm);
23  watplot wmap;
24 
25  // plot the TF map as energy average of 0 and 90 degree phases (quadratures)
26  wmap.plot(tfmap);
27 
28  return wmap.canvas; // used by THtml doc
29 }
virtual void rate(double r)
Definition: wavearray.hh:141
WDM< double > wdm(nLAYERS, nLAYERS, 6, 10)
i drho i
void plot(wavearray< double > &, char *=NULL, int=1, double=0., double=0., bool=false, float=0., float=0., bool=false, float=0., bool=false)
Definition: watplot.cc:150
TCanvas * canvas
Definition: watplot.hh:192
int Duration
Definition: testWDM_1.C:4
WSeries< double > tfmap
watplot * wmap
wavearray< double > ts(N)
double T
Definition: testWDM_4.C:11
TRandom3 rnd(0)
int N
Definition: testWDM_1.C:5
void Forward(int n=-1)
param: wavelet - n is number of steps (-1 means full decomposition)
Definition: wseries.cc:246
double omega
Definition: testWDM_4.C:12