Logo coherent WaveBurst  
Library Reference Guide
Logo
gwavearray.hh
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 /**********************************************************
20  * Package: graphic wat Class Library
21  * File name: gwavearray.hh
22  * Author: Gabriele Vedovato (vedovato@lnl.infn.it)
23  **********************************************************/
24 
25 #ifndef GWAVEARRAY_HH
26 #define GWAVEARRAY_HH
27 
28 #ifndef WAVEARRAY_HH
29 #include "wavearray.hh"
30 #endif
31 
32 #include "gwat.hh"
33 #include "watplot.hh"
34 #include "time.hh"
35 #include "STFT.hh"
36 #include "TComplex.h"
37 
38 
39 template<class DataType_t>
40 class gwavearray : public wavearray<DataType_t>
41 {
42 public:
43 
44  gwavearray() : wavearray<DataType_t>() {Init();wextern=false;} // Default constructor
45  gwavearray(const wavearray<DataType_t>&); // copy Constructor
47 
48  virtual ~gwavearray(); // Destructor
49 
50  void Draw(GWAT_DRAW type=GWAT_TIME,
51  TString options="ALP", Color_t color=kBlack);
52 
53  void PlotTime(double edge=0) {
54  SetTimeRange(this->start()+edge,this->start()+this->size()/this->rate()-edge);
55  DrawTime("CUSTOM ALP");} // *MENU*
56  void PlotFFT(double edge=0) {
57  SetTimeRange(this->start()+edge,this->start()+this->size()/this->rate()-edge);
58  DrawFFT("CUSTOM ALP");} // *MENU*
59  void PlotPSD(double edge=0) {
60  SetTimeRange(this->start()+edge,this->start()+this->size()/this->rate()-edge);
61  DrawFFT("CUSTOM ALP PSD");} // *MENU*
62  void PlotTF(double edge=0) {DrawTF();;}
63 
64  watplot* DrawTime(TString options="ALP", Color_t color=kBlack);
65  watplot* DrawFFT(TString options="ALP", Color_t color=kBlack);
67 
69  TString options="ALP", Color_t color=kBlack);
70 
71  watplot* DrawTime(wavearray<DataType_t>* x, TString options="ALP", Color_t color=kBlack);
72  watplot* DrawFFT(wavearray<DataType_t>* x, TString options="ALP", Color_t color=kBlack);
74 
75  void SetTimeRange(double tMin, double tMax) {
76  if(tMin>=tMax) return;
77  double tStart=this->start();
78  double tStop=this->start()+this->size()/this->rate();
79  this->tMin = tMin<tStart||tMin>tStop ? tStart : tMin;
80  this->tMax = tMax<tStart||tMax>tStop ? tStop : tMax;
81  }
82  double GetTimeRange(double& tMin, double& tMax, double efraction = 0.9999999);
83  double GetCentralTime();
84  void TimeShift(double tShift=0.);
85  void PhaseShift(double pShift=0.);
86 
87  CWB::STFT* GetSTFT() {return this->stft;}
88  watplot* GetWATPLOT() {return this->pts;}
89 
90  // print wseries parameters
92  virtual void Browse(TBrowser *b) {PlotTime();}
93  void SetWATPLOT(watplot* pts) {this->pts=pts;}
94 
95  void SetComment(TString comment) {this->comment=comment;}
97  void PrintComment() {cout << comment << endl;} // *MENU*
98 
99  void DumpToFile(char* fname); // *MENU*
100 
101 private:
102 
104  watplot* pts; //!
105 
106  void Init();
107 
108  bool wextern; // external wavearray
109 
110  TString comment; // store comment
111 
112  double tMin;
113  double tMax;
114  double tSave; // start time : used to manage DrawTime with option "SAME"
115 
116  ClassDef(gwavearray,3)
117 };
118 
119 #endif
void PlotTime(double edge=0)
Definition: gwavearray.hh:53
void PhaseShift(double pShift=0.)
Definition: gwavearray.cc:478
void PrintComment()
Definition: gwavearray.hh:97
void DumpToFile(char *fname)
Definition: gwavearray.cc:526
void PlotPSD(double edge=0)
Definition: gwavearray.hh:59
double tSave
Definition: gwavearray.hh:114
void print()
Definition: wavearray.cc:2226
TString comment
Definition: gwavearray.hh:110
TString("c")
void PlotTF(double edge=0)
Definition: gwavearray.hh:62
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
CWB::STFT * GetSTFT()
Definition: gwavearray.hh:87
void PlotFFT(double edge=0)
Definition: gwavearray.hh:56
double GetCentralTime()
Definition: gwavearray.cc:402
CWB::STFT * stft
Definition: gwavearray.hh:103
virtual size_t size() const
Definition: wavearray.hh:145
bool wextern
Definition: gwavearray.hh:108
virtual double rate() const
Definition: wavearray.hh:142
double tMin
Definition: gwavearray.hh:112
watplot * pts
Definition: gwavearray.hh:104
char fname[1024]
GWAT_DRAW
Definition: gwat.hh:28
void TimeShift(double tShift=0.)
Definition: gwavearray.cc:427
virtual double start() const
Definition: wavearray.hh:138
double GetTimeRange(double &tMin, double &tMax, double efraction=0.9999999)
Definition: gwavearray.cc:359
void SetComment(TString comment)
Definition: gwavearray.hh:95
void SetTimeRange(double tMin, double tMax)
Definition: gwavearray.hh:75
virtual ~gwavearray()
Definition: gwavearray.cc:48
virtual double edge() const
Definition: wavearray.hh:144
void SetWATPLOT(watplot *pts)
Definition: gwavearray.hh:93
virtual void Browse(TBrowser *b)
Definition: gwavearray.hh:92
void Init()
Definition: gwavearray.cc:78
TString GetComment()
Definition: gwavearray.hh:96
char options[256]
CWB::STFT * DrawTF(TString options="")
Definition: gwavearray.cc:291
watplot * GetWATPLOT()
Definition: gwavearray.hh:88
watplot * DrawTime(TString options="ALP", Color_t color=kBlack)
Definition: gwavearray.cc:124
double tMax
Definition: gwavearray.hh:113
void Draw(GWAT_DRAW type=GWAT_TIME, TString options="ALP", Color_t color=kBlack)
Definition: gwavearray.cc:89
watplot * DrawFFT(TString options="ALP", Color_t color=kBlack)
Definition: gwavearray.cc:203
void PrintParameters()
Definition: gwavearray.hh:91