Logo coherent WaveBurst  
Library Reference Guide
Logo
Make_BRST_O1_LF.C
Go to the documentation of this file.
1 //
2 // Make BRST LF [32:1024] Hz Set for O1 (proposal)
3 // How to use it :
4 // root -b Make_BRST_O1_LF.C
5 //
6 // 1) TIME & FFT plots are generated for each waveform and saved under the ODIR/plots directory
7 // 2) A texi file html_index.texi is created under ODIR
8 // 3) texi file is converted into ODIR/html_index/index.html
9 //
10 // Author : Gabriele Vedovato
11 //
12 
13 #define ODIR "brst_o1_lf" // output plot directory
14 
16 
18 
19  gROOT->SetBatch(true);
20 
21  MDC = new CWB::mdc();
22 
23  vector<mdcpar> par;
24 
25  // -------------------------------------------
26  // GA
27  // -------------------------------------------
28 
29  par.resize(1);
30  double F_GA[4] = {0.0001,0.001,0.0025,0.004};
31  for(int n=0;n<4;n++) {
32  par[0].name="duration"; par[0].value=F_GA[n];
33  MDC->AddWaveform(MDC_GA, par);
34  }
35 
36  // -------------------------------------------
37  // SGQ3
38  // -------------------------------------------
39 
40  par.resize(2);
41  double F_Q3[3] = {70,235,849};
42  for(int n=0;n<3;n++) {
43  par[0].name="frequency"; par[0].value=F_Q3[n];
44  par[1].name="Q"; par[1].value=3.;
45  MDC->AddWaveform(MDC_SGE, par);
46  }
47 
48  // -------------------------------------------
49  // SGQ9
50  // -------------------------------------------
51 
52  par.resize(2);
53  double F_Q8d9[8] = {70,100,153,235,361,554,849,1053};
54  double Q[8] = {8.9,8.9,8.9,8.9,8.9,8.9,8.9,9};
55  for(int n=0;n<8;n++) {
56  par[0].name="frequency"; par[0].value=F_Q8d9[n];
57  par[1].name="Q"; par[1].value=Q[n];
58  if((F_Q8d9[n]==70)||(F_Q8d9[n]==100)||(F_Q8d9[n]==235)||(F_Q8d9[n]==361)||(F_Q8d9[n]==1053)) {
59  MDC->AddWaveform(MDC_SGL, par); // linear polarized for back compatibility
60  } else {
61  MDC->AddWaveform(MDC_SGE, par);
62  }
63  }
64 
65  // -------------------------------------------
66  // SGQ100
67  // -------------------------------------------
68 
69  par.resize(2);
70  double F_Q100[3] = {70,235,849};
71  for(int n=0;n<3;n++) {
72  par[0].name="frequency"; par[0].value=F_Q100[n];
73  par[1].name="Q"; par[1].value=100.;
74  MDC->AddWaveform(MDC_SGE, par);
75  }
76 
77  // -------------------------------------------
78  // WNB
79  // -------------------------------------------
80 
81  par.resize(6);
82  double F_WNB[2] = {100,250};
83  double B_WNB[2] = {100,100};
84  double D_WNB[2] = {0.1,0.1};
85  for(int n=0;n<2;n++) {
86  par[0].name="frequency"; par[0].value=F_WNB[n];
87  par[1].name="bandwidth"; par[1].value=B_WNB[n];
88  par[2].name="duration"; par[2].value=D_WNB[n];
89  for(int m=0;m<30;m++) {
90  par[3].name="pseed"; par[3].value=100000+n*100+m;
91  par[4].name="xseed"; par[4].value=100001+n*100+m;
92  par[5].name="mode"; par[5].value=0; // asymmetric
93  MDC->AddWaveform(MDC_WNB, par);
94  }
95  }
96 
97  MDC->Print();
98 
99  vector<waveform> wfList = MDC->GetWaveformList();
100 
101  // loop over the waveform list
102  for(int n=0;n<wfList.size();n++) {
103 
104  // Get the first waveform hp,hx components
105  waveform wf = MDC->GetWaveform(wfList[n].name,0);
106 
107  if(wf.status==false) {
108  cout << "Error : Waveform " << wf.name << " not exist in the MDC pool !!!" << endl << endl;
109  gSystem->Exit(1);
110  }
111 
112  cout << n << " " << wf.name << " size : " << wf.hp.size() << " rate : " << wf.hp.rate()
113  << " start : " << (int)wf.hp.start() << endl;
114 
115  wf.hp.start(0); // set start to 0 (needed by draw Method)
116  wf.hx.start(0);
117 
118  watplot* plot = NULL;
119 
120  gSystem->mkdir(TString::Format("%s",ODIR));
121  gSystem->mkdir(TString::Format("%s/plots",ODIR));
122 
123  // draw hp,hx in time domain
124  plot=MDC->Draw(wf.hp,MDC_TIME,"ALP ZOOM");
125  MDC->Draw(wf.hx,MDC_TIME,"same",kRed);
126  if(plot) plot->graph[0]->SetTitle(TString::Format("%s : h+(black), hx(red))",wf.name.Data()));
127  if(plot) plot->canvas->SaveAs(TString::Format("%s/plots/%s_TIME.png",ODIR,wf.name.Data()));
128 
129  // draw hp,hx in frequency domain
130  plot = MDC->Draw(wf.hp,MDC_FFT,"ALP NOLOGX");
131  if(wf.name=="GA0d1") MDC->Draw(wf.hx,MDC_FFT,"same NOLOGX NOLOGY",kRed);
132  else MDC->Draw(wf.hx,MDC_FFT,"same NOLOGX",kRed);
133  plot->graph[0]->GetHistogram()->GetXaxis()->SetRangeUser(32,1024);
134  if(plot) plot->graph[0]->SetTitle(TString::Format("%s : h+(black), hx(red))",wf.name.Data()));
135  if(plot) plot->canvas->SaveAs(TString::Format("%s/plots/%s_FFT.png",ODIR,wf.name.Data()));
136  }
137 
138  TString texiName=TString::Format("%s/html_index.texi",ODIR);
139  bool overwrite=CWB::Toolbox::checkFile(texiName,true);
140  if(!overwrite) gSystem->Exit(0);
141  ofstream out;
142  out.open(texiName,ios::out);
143  out << "@c include texi macros" << endl;
144  out << "@include macros.texi" << endl;
145  out << "@include mathjax.texi" << endl << endl;
146  for(int i=0;i<wfList.size();i++) {
147  out << "@center @txtfont{"<<wfList[i].name<<", red, h1}" << endl;
148  out << "@multitable @columnfractions .5 .5" << endl;
149  out << "@item @center @displayimage{../plots,"<<TString::Format("%s_TIME",wfList[i].name.Data())<<",480}"<<endl;
150  out << "@tab @center @displayimage{../plots,"<<TString::Format("%s_FFT",wfList[i].name.Data())<<",480}"<<endl;
151  out << "@end multitable" << endl;
152  }
153  out.close();
154 
155  // convert texi into html
156  TString cwb_scripts = TString(gSystem->Getenv("CWB_SCRIPTS"));
157  //TString exec_cmd = TString::Format("%s/cwb_mkhtml.csh %s wheader;rm %s",
158  // cwb_scripts.Data(),texiName.Data(),texiName.Data());
159  TString exec_cmd = TString::Format("%s/cwb_mkhtml.csh %s wheader",
160  cwb_scripts.Data(),texiName.Data());
161  int ret=gSystem->Exec(exec_cmd);
162  if(ret) {
163  cout << "Make_BRST_O1_LF.C : Error while executing cwb_mkhtml html_index.texi !!!" << endl;
164  gSystem->Exit(1);
165  }
166  cout << endl << "New html file created : " << ODIR << "/html_index/index.html" << endl << endl;
167 
168  gSystem->Exit(0);
169 }
TString cwb_scripts
Definition: mdc.hh:202
watplot * Draw(TString name, int id=0, TString polarization="hp", MDC_DRAW type=MDC_TIME, TString options="ALP", Color_t color=kBlack)
Definition: mdc.cc:2317
Definition: mdc.hh:219
bool status
Definition: mdc.hh:221
virtual void rate(double r)
Definition: wavearray.hh:141
double F_Q3[4]
Definition: mdc.hh:158
par [0] name
int n
Definition: cwb_net.C:28
Definition: mdc.hh:203
TString("c")
double B_WNB[5]
ofstream out
Definition: cwb_merge.C:214
TString mdc[4]
std::vector< TGraph * > graph
Definition: watplot.hh:194
wavearray< double > hp
Definition: mdc.hh:226
waveform wf
vector< waveform > GetWaveformList()
Definition: mdc.hh:385
int m
Definition: cwb_net.C:28
virtual void start(double s)
Definition: wavearray.hh:137
i drho i
static bool checkFile(TString fName, bool question=false, TString message="")
Definition: Toolbox.cc:4670
mdcid AddWaveform(MDC_TYPE mdc_type, vector< mdcpar > par, TString uname="")
Definition: mdc.cc:472
void Print(int level=0)
Definition: mdc.cc:2736
#define ODIR
virtual size_t size() const
Definition: wavearray.hh:145
TCanvas * canvas
Definition: watplot.hh:192
TString name
Definition: mdc.hh:222
x plot
double F_GA[4]
Definition: mdc.hh:248
double F_WNB[5]
i() int(T_cor *100))
double F_Q100[4]
double F_Q8d9[7]
double D_WNB[5]
Definition: mdc.hh:150
Definition: mdc.hh:152
vector< mdcpar > par
wavearray< double > hx
Definition: mdc.hh:227
double Q
waveform GetWaveform(int ID, int id=0)
Definition: mdc.cc:2541
bool overwrite
Definition: cwb_dump_inj.C:100
void Make_BRST_O1_LF()
CWB::mdc * MDC
Definition: mdc.hh:157