Logo coherent WaveBurst  
Library Reference Guide
Logo
gwseries.cc
Go to the documentation of this file.
1 /*
2 # Copyright (C) 2019 Gabriele Vedovato
3 #
4 # This program is free software: you can redistribute it and/or modify
5 # it under the terms of the GNU General Public License as published by
6 # the Free Software Foundation, either version 3 of the License, or
7 # (at your option) any later version.
8 #
9 # This program is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 # GNU General Public License for more details.
13 #
14 # You should have received a copy of the GNU General Public License
15 # along with this program. If not, see <https://www.gnu.org/licenses/>.
16 */
17 
18 
19 #include "gwseries.hh"
20 
21 ClassImp(gWSeries<double>)
22 
23 
24 //______________________________________________________________________________
25 /* Begin_Html
26 <center><h2>gWSeries class</h2></center>
27 This class gWSeries is derived from the WSeries class
28 It uses the methods implemented in gwavearray and the scalogram plot
29 
30 <p>
31 <h3><a name="example">Example</a></h3>
32 <p>
33 The macro <a href="./tutorials/gwat/DrawGWSeries.C.html">DrawGWSeries.C</a> is an example which shown how to use the gWSeries class.<br>
34 The pictures below gives the macro output plots.<br>
35 <p>
36 
37 End_Html
38 
39 Begin_Macro
40 DrawGWSeries.C
41 End_Macro */
42 
43 
44 using namespace std;
45 
46 
47 //______________________________________________________________________________
48 // destructor
49 template<class DataType_t>
51 {
52  if(gw!=NULL) delete gw;
53 
54  if(wsextern) {
55  this->data = NULL;
56  this->pWavelet = NULL;
57  }
58 
59 // if(wts!=NULL) delete wts;
60 }
61 
62 //______________________________________________________________________________
63 template<class DataType_t>
65 WSeries<DataType_t>(w) {this->Init();wsextern=false;}
66 
67 //______________________________________________________________________________
68 template<class DataType_t>
70 {
71  wsextern = true;
72 
73  this->data = w->data;
74  this->Rate = w->rate();
75  this->Size = w->size();
76  this->Start = w->start();
77  this->Stop = w->start() + w->size()/w->rate();
78  this->Slice = std::slice(0,w->size(),1);
79 
80  this->pWavelet = w->pWavelet;
81  this->bpp = w->getbpp();
82  this->wRate = w->wRate;
83  this->f_low = w->getlow();
84  this->f_high = w->gethigh();
85  this->w_mode = w->w_mode;
86 
87  this->Init();
88 }
89 
90 //______________________________________________________________________________
91 template<class DataType_t>
93 
94  wavearray<DataType_t>* p = this;
95  gw = new gwavearray<DataType_t>(p);
96 
97  gRandom->SetSeed(0);
98  rnID = int(gRandom->Rndm(13)*1.e9); // random name ID
99 
100 // wts=NULL;
101 }
102 
103 //______________________________________________________________________________
104 template<class DataType_t>
106 //
107 // Draw waveform scalogram
108 //
109 // Input:
110 // options - draw options (same as TH2D)
111 //
112 // return pointer to watplot object
113 //
114 
115  double tStart,tStop;
116  if(options.Contains("FULL")) {
117  options.ReplaceAll("FULL","");
118  tStart=this->start();
119  tStop=this->start()+this->size()/this->rate();
120  } else {
121  GetTimeRange(tStart, tStop);
122  }
123 
124  watplot* wts = GetWATPLOT();
125  if(wts!=NULL) delete wts;
126 
127  char name[32];sprintf(name,"sg-%d",rnID);
128  wts = new watplot(const_cast<char*>(name));
129  WSeries<DataType_t>* ws = x!=NULL ? x : this;
130  wts->plot(ws, 2, tStart, tStop,const_cast<char*>("COLZ"));
131 
132  double fLow = 32.;
133  double fHigh = this->rate()/2.;
134  wts->hist2D->GetYaxis()->SetRangeUser(fLow,fHigh);
135 
136  gw->SetWATPLOT(wts);
137 
138  return wts;
139 }
140 
141 // instantiations
142 
143 #define CLASS_INSTANTIATION(class_) template class gWSeries< class_ >;
144 
145 //CLASS_INSTANTIATION(short)
146 //CLASS_INSTANTIATION(int)
147 //CLASS_INSTANTIATION(unsigned int)
148 //CLASS_INSTANTIATION(long)
149 //CLASS_INSTANTIATION(long long)
150 //CLASS_INSTANTIATION(float)
151 CLASS_INSTANTIATION(double)
152 
double f_low
Definition: wseries.hh:464
double f_high
Definition: wseries.hh:466
double getlow() const
Definition: wseries.hh:129
int rnID
Definition: gwseries.hh:99
double fHigh
virtual void rate(double r)
Definition: wavearray.hh:141
size_t Size
data array
Definition: wavearray.hh:323
par [0] name
TString("c")
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\ layers : "<< nLAYERS<< "\ dF(hz) : "<< dF<< "\ dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1) *itime+ifreq;double time=itime *dT;double freq=(ifreq >0) ? ifreq *dF :dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
double Rate
Definition: wavearray.hh:324
double Stop
Definition: wavearray.hh:326
STL namespace.
virtual void start(double s)
Definition: wavearray.hh:137
double gethigh() const
Definition: wseries.hh:136
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
TH2F * hist2D
Definition: watplot.hh:193
double wRate
Definition: wseries.hh:462
wavearray< double > w
Definition: Test1.C:27
virtual size_t size() const
Definition: wavearray.hh:145
gWSeries< double > gw(w)
virtual double rate() const
Definition: wavearray.hh:142
size_t w_mode
Definition: wseries.hh:458
double GetTimeRange(double &tMin, double &tMax)
Definition: gwseries.hh:80
i() int(T_cor *100))
void Init()
Definition: gwseries.cc:92
virtual ~gWSeries()
Definition: gwseries.cc:50
std::slice Slice
Definition: wavearray.hh:328
virtual double start() const
Definition: wavearray.hh:138
#define CLASS_INSTANTIATION(class_)
Definition: gwseries.cc:143
double getbpp() const
Definition: wseries.hh:117
watplot * GetWATPLOT()
Definition: gwseries.hh:86
char options[256]
watplot * DrawSG(TString options="")
Definition: gwseries.hh:62
gwavearray< DataType_t > * gw
Definition: gwseries.hh:95
double fLow
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
DataType_t * data
Definition: wavearray.hh:319
WaveDWT< DataType_t > * pWavelet
Definition: wseries.hh:456
double Start
Definition: wavearray.hh:325
bool wsextern
Definition: gwseries.hh:100
double bpp
Definition: wseries.hh:460
gWSeries()
Definition: gwseries.hh:39