Logo coherent WaveBurst  
Library Reference Guide 		Logo
CWB_Plugin_QLWveto.C
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 #define XIFO 4
20 
21 #pragma GCC system_header
22 
23 #include "cwb.hh"
24 #include "config.hh"
25 #include "network.hh"
26 #include "wavearray.hh"
27 #include "TString.h"
28 #include "TObjArray.h"
29 #include "TObjString.h"
30 #include "TRandom.h"
31 #include "TComplex.h"
32 #include "TMath.h"
33 #include "mdc.hh"
34 #include "watplot.hh"
35 #include "gwavearray.hh"
36 #include <vector>
37 
38 //#define PLOT_LIKELIHOOD
39 //#define PLOT_WHITENED_WAVEFORMS
40 
41 #define NTHR 1
42 #define ATHR 7.58859
43 
44 float GetQveto(wavearray<double>* wf);
45 void GetLveto(netcluster* pwc, int cid, int nifo, float* Lveto);
46 void GetWveto(netcluster* pwc, int cid, int nifo, float* Wveto);
47 void PlotWaveform(TString ifo, wavearray<double>* wfREC,
48  CWB::config* cfg, bool fft=false, bool strain=false);
49 void ClearWaveforms(detector* ifo);
50 
51 
52 void
53 CWB_Plugin(TFile* jfile, CWB::config* cfg, network* NET, WSeries<double>* x, TString ifo, int type) {
54 //!MISCELLANEA
55 // Extract whitened reconstructed waveforms, and compute the Qveto, Lveto & Wveto parameters
56 
57  cout << endl;
58  cout << "-----> CWB_Plugin_QLWveto.C" << endl;
59  cout << "ifo " << ifo.Data() << endl;
60  cout << "type " << type << endl;
61  cout << endl;
62 
63  float Qveto[2*NIFO_MAX]; // Qveto
64  float Lveto[3]; // Lveto
65  float Wveto[2]; // Wveto
66 
67  if(type==CWB_PLUGIN_CONFIG) {
68  cfg->outPlugin=true; // disable built-in output root file
69  }
70 
71  if(type==CWB_PLUGIN_ILIKELIHOOD) {
72  NET->wfsave=true; // enable waveform save
73 
74  // search output root file in the system list
75  TFile* froot = NULL;
76  TList *files = (TList*)gROOT->GetListOfFiles();
77  TString outDump="";
78  netevent* EVT;
79  int nIFO = NET->ifoListSize(); // number of detectors
80  if (files) {
81  TIter next(files);
82  TSystemFile *file;
83  TString fname;
84  bool check=false;
85  while ((file=(TSystemFile*)next())) {
86  fname = file->GetName();
87  // set output root file as the current file
88  if(fname.Contains("wave_")) {
89  froot=(TFile*)file;froot->cd();
90  outDump=fname;
91  outDump.ReplaceAll(".root.tmp",".txt");
92  //cout << "output file name : " << fname << endl;
93  }
94  }
95  if(!froot) {
96  cout << "CWB_Plugin_QLWveto.C : Error - output root file not found" << endl;
97  gSystem->Exit(1);
98  }
99  } else {
100  cout << "CWB_Plugin_QLWveto.C : Error - output root file not found" << endl;
101  gSystem->Exit(1);
102  }
103 
104  TTree* net_tree = (TTree *) froot->Get("waveburst");
105  if(net_tree==NULL) {
106  EVT = new netevent(nIFO);
107  net_tree = EVT->setTree();
108  net_tree->Branch("Qveto",Qveto,TString::Format("Qveto[%i]/F",2*cfg->nIFO));
109  net_tree->Branch("Lveto",Lveto,TString::Format("Lveto[%i]/F",3));
110  net_tree->Branch("Wveto",Wveto,TString::Format("Wveto[%i]/F",2));
111  }
112  }
113 
114  if(type==CWB_PLUGIN_OLIKELIHOOD) {
115 
116  if(TString(cfg->analysis)!="2G") {
117  cout << "CWB_Plugin_QLWveto.C -> "
118  << "CWB_PLUGIN_OLIKELIHOOD implemented only for 2G" << endl;
119  gSystem->Exit(1);
120  }
121 
122  // import ifactor
123  int gIFACTOR=-1; IMPORT(int,gIFACTOR)
124  cout << "-----> CWB_Plugin_QLWveto.C -> "
125  << " gIFACTOR : " << gIFACTOR << endl;
126 
127  // import slagShift
128  float* gSLAGSHIFT=NULL; IMPORT(float*,gSLAGSHIFT)
129 
130  int nIFO = NET->ifoListSize(); // number of detectors
131  int K = NET->nLag; // number of time lag
132  netevent* EVT;
133  wavearray<double> id;
134  //double factor = cfg->simulation==3||cfg->simulation==4 ? -gIFACTOR : cfg->factors[gIFACTOR];
135  double factor = cfg->factors[gIFACTOR];
136  int rate = 0; // select all resolutions
137 
138  // search output root file in the system list
139  TFile* froot = NULL;
140  TList *files = (TList*)gROOT->GetListOfFiles();
141  TString outDump="";
142  if (files) {
143  TIter next(files);
144  TSystemFile *file;
145  TString fname;
146  bool check=false;
147  while ((file=(TSystemFile*)next())) {
148  fname = file->GetName();
149  // set output root file as the current file
150  if(fname.Contains("wave_")) {
151  froot=(TFile*)file;froot->cd();
152  outDump=fname;
153  outDump.ReplaceAll(".root.tmp",".txt");
154  //cout << "output file name : " << fname << endl;
155  }
156  }
157  if(!froot) {
158  cout << "CWB_Plugin_QLWveto.C : Error - output root file not found" << endl;
159  gSystem->Exit(1);
160  }
161  } else {
162  cout << "CWB_Plugin_QLWveto.C : Error - output root file not found" << endl;
163  gSystem->Exit(1);
164  }
165 
166  TTree* net_tree = (TTree *) froot->Get("waveburst");
167  if(net_tree!=NULL) {
168  EVT = new netevent(net_tree,nIFO);
169  net_tree->SetBranchAddress("Qveto",Qveto);
170  net_tree->SetBranchAddress("Lveto",Lveto);
171  net_tree->SetBranchAddress("Wveto",Wveto);
172  } else {
173  EVT = new netevent(nIFO);
174  net_tree = EVT->setTree();
175  net_tree->Branch("Qveto",Qveto,TString::Format("Qveto[%i]/F",2*cfg->nIFO));
176  net_tree->Branch("Lveto",Lveto,TString::Format("Lveto[%i]/F",3));
177  net_tree->Branch("Wveto",Wveto,TString::Format("Wveto[%i]/F",2));
178  }
179  EVT->setSLags(gSLAGSHIFT); // set slags into netevent
180 
181  for(int k=0; k<K; k++) { // loop over the lags
182 
183  id = NET->getwc(k)->get(const_cast<char*>("ID"), 0, 'L', rate);
184 
185  for(int j=0; j<(int)id.size(); j++) { // loop over cluster index
186 
187  int ID = size_t(id.data[j]+0.5);
188 
189  if(NET->getwc(k)->sCuts[ID-1]!=-1) continue; // skip rejected/processed clusters
190 
191  double ofactor=0;
192  if(cfg->simulation==4) ofactor=-factor;
193  else if(cfg->simulation==3) ofactor=-gIFACTOR;
194  else ofactor=factor;
195 
196  EVT->output2G(NULL,NET,ID,k,ofactor); // get reconstructed parameters
197 
198  wavearray<double>** pwfREC = new wavearray<double>*[nIFO];
199  detector* pd = NET->getifo(0);
200  int idSize = pd->RWFID.size();
201 
202  int wfIndex=-1;
203  for (int mm=0; mm<idSize; mm++) if (pd->RWFID[mm]==ID) wfIndex=mm;
204  if(wfIndex==-1) continue;
205 
206  netcluster* pwc = NET->getwc(k);
207  cout << endl << "----------------------------------------------------------------" << endl;
208  GetLveto(pwc, ID, nIFO, Lveto);
209  cout << "Lveto : " << "fmean : " << Lveto[0] << " frms : " << Lveto[1]
210  << " Energy Ratio : " << Lveto[2] << endl << endl;
211  GetWveto(pwc, ID, nIFO, Wveto);
212  cout << "Wveto : " << " Slope : " << Wveto[0] << " Correlation : " << Wveto[1] << endl << endl;
213 
214  // extract whitened reconstructed waveforms
215  for(int n=0; n<nIFO; n++) {
216 
217  pd = NET->getifo(n);
218 
219  pwfREC[n] = pd->RWFP[wfIndex];
220  wavearray<double>* wfREC = pwfREC[n]; // array of reconstructed waveforms
221 
222 #ifdef PLOT_WHITENED_WAVEFORMS
223  //PlotWaveform(NET->ifoName[n], wfREC, cfg, false, false);
224  PlotWaveform(NET->ifoName[n], wfREC, cfg, true, false);
225 #endif
226  // reconstructed whitened waveform
227  NET->getMRAwave(ID,k,'S',0,true);
228  Qveto[n] = GetQveto(&(pd->waveForm));
229  // whitened waveform
230  NET->getMRAwave(ID,k,'W',0,true);
231  Qveto[n+nIFO] = GetQveto(&(pd->waveBand));
232 
233  //Qveto[n] = GetQveto(wfREC);
234  cout << "Qveto : " << pd->Name << " Qveto[R] = " << Qveto[n]
235  << " Qveto[W] = " << Qveto[n+nIFO] << endl;
236 
237  if(!cfg->simulation) ClearWaveforms(pd); // release waveform memory
238  }
239  cout << "----------------------------------------------------------------" << endl;
240  delete [] pwfREC;
241 
242  std::vector<int> sCuts = NET->getwc(k)->sCuts; // save cCuts
243  // set sCuts=1 to the events which must be not copied with cps to pwc
244  for(int i=0; i<(int)sCuts.size(); i++) if(i!=ID-1) NET->getwc(k)->sCuts[i]=1;
245 
246  // ID can not be used to get the event, to get event use ID=1 (ex: for watplot)
247  NET->getwc(k)->sCuts = sCuts; // restore cCuts
248 
249  if(cfg->dump) EVT->dopen(outDump.Data(),const_cast<char*>("a"),false);
250  EVT->output2G(net_tree,NET,ID,k,ofactor); // get reconstructed parameters
251  if(cfg->dump) {
252  // add Qveto to dump file
253  fprintf(EVT->fP,"Qveto: ");
254  for(int i=0; i<2*nIFO; i++) fprintf(EVT->fP,"%f ",Qveto[i]);
255  fprintf(EVT->fP,"\n");
256  // add Lveto to dump file
257  fprintf(EVT->fP,"Lveto: ");
258  for(int i=0; i<3; i++) fprintf(EVT->fP,"%f ",Lveto[i]);
259  fprintf(EVT->fP,"\n");
260  // add Wveto to dump file
261  fprintf(EVT->fP,"Wveto: ");
262  for(int i=0; i<2; i++) fprintf(EVT->fP,"%f ",Wveto[i]);
263  fprintf(EVT->fP,"\n");
264  }
265  if(cfg->dump) EVT->dclose();
266  if(!cfg->cedDump) NET->getwc(k)->sCuts[ID-1]=1; // mark as processed
267  }
268  }
269 
270  jfile->cd();
271  if(EVT) delete EVT;
272  }
273  return;
274 }
275 
276 float
277 GetQveto(wavearray<double>* wf) {
278 
279  wavearray<double> x = *wf;;
280 
281  // resample data by a factor 4
282  int xsize=x.size();
283  x.FFTW(1);
284  x.resize(4*x.size());
285  x.rate(4*x.rate());
286  for(int j=xsize;j<x.size();j++) x[j]=0;
287  x.FFTW(-1);
288 
289  // extract max/min values and save the absolute values in the array 'a'
290  wavearray<double> a(x.size());
291  int size=0;
292  double dt = 1./x.rate();
293  double prev=x[0];
294  double xmax=0;
295  for (int i=1;i<x.size();i++) {
296  if(fabs(x[i])>xmax) xmax=fabs(x[i]);
297  if(prev*x[i]<0) {
298  a[size]=xmax;
299  size++;
300  xmax=0;
301  }
302  prev=x[i];
303  }
304 
305  // find max value/index ans save on amax/imax
306  int imax=-1;
307  double amax=0;
308  for (int i=1;i<size;i++) {
309  if(a[i]>amax) {amax=a[i];imax=i;}
310  }
311 
312 /*
313  cout << endl;
314  cout << "a[imax-2] " << a[imax-2] << endl;
315  cout << "a[imax-1] " << a[imax-1] << endl;
316  cout << "a[imax] " << a[imax] << endl;
317  cout << "a[imax+1] " << a[imax+1] << endl;
318  cout << "a[imax+2] " << a[imax+2] << endl;
319  cout << endl;
320 */
321 
322  // compute Qveto
323  double ein=0; // energy of max values inside NTHR
324  double eout=0; // energy of max values outside NTHR
325  for (int i=0;i<size;i++) {
326  if(abs(imax-i)<=NTHR) {
327  ein+=a[i]*a[i];
328  //cout << i << " ein " << a[i] << " " << amax << endl;
329  } else {
330  if(a[i]>amax/ATHR) eout+=a[i]*a[i];
331  //if(a[i]>amax/ATHR) cout << i << " eout " << a[i] << " " << amax << endl;
332  }
333  }
334  float Qveto = ein>0 ? eout/ein : 0.;
335  //cout << "Qveto : " << Qveto << " ein : " << ein << " eout : " << eout << endl;
336 
337  return Qveto;
338 }
339 
340 void
341 GetLveto(netcluster* pwc, int cid, int nifo, float* Lveto) {
342 //
343 // input
344 // pwc : pointer to netcluster object
345 // cid : cluster id
346 // nifo : number of detectors
347 // output
348 // Lveto[0] : line frequency
349 // Lveto[1] : line bandwitdh
350 // Lveto[2] : line enery ratio (line_energy / total_energy)
351 //
352 
353  Lveto[0] = Lveto[1] = Lveto[2] = 0;
354 
355  std::vector<int>* vint = &(pwc->cList[cid-1]); // pixel list
356  int V = vint->size(); // cluster size
357  if(!V) return;
358 
359  // ------------------------------------------------------------------
360  // Find max pixel parameters
361  // ------------------------------------------------------------------
362 
363  double likeMax=0; // maximum pixel's energy
364  double likeTot=0; // total cluster energy
365  double freqMax; // frequency of the pixel with max energy
366  double dfreqMax; // df of the pixel with max energy
367  for(int n=0; n<V; n++) {
368  netpixel* pix = pwc->getPixel(cid,n);
369  if(pix->layers%2==0) {
370  cout << "CWB_Plugin_QLWveto.C - Error : is enabled only for WDM (2G)" << endl;
371  exit(1);
372  }
373  if(!pix->core) continue; // select only the principal components pixels
374 
375  double likePix=0;
376  for(int m=0; m<nifo; m++) {
377  likePix += pow(pix->getdata('S',m),2); // like whitened reconstructed signal 00
378  likePix += pow(pix->getdata('P',m),2); // like whitened reconstructed signal 90
379  }
380 
381  double freq = pix->frequency*pix->rate/2.;
382  double df = pix->rate/2.;
383 
384  likeTot+=likePix;
385  if(likePix>likeMax) {likeMax=likePix;freqMax=freq;dfreqMax=df;}
386  }
387  //cout << "likeMax : " << likeMax << " likeTot : " << likeTot
388  // << " freqMax : " << freqMax << " dfreqMax : " << dfreqMax << endl;
389 
390  // ------------------------------------------------------------------
391  // Compute Lveto parameters
392  // ------------------------------------------------------------------
393 
394  double fmean=0; // line mean frequency
395  double frms=0; // line bandwidth
396  double likeLin=0; // line energy
397  for(int n=0; n<V; n++) {
398  netpixel* pix = pwc->getPixel(cid,n);
399  if(!pix->core) continue; // select only the principal components pixels
400 
401  double likePix=0;
402  for(int m=0; m<nifo; m++) {
403  likePix += pow(pix->getdata('S',m),2); // like whitened reconstructed signal 00
404  likePix += pow(pix->getdata('P',m),2); // like whitened reconstructed signal 90
405  }
406 
407  // the estimation is done for all pixels
408  // with freq in the range [freqMax-dfreqMax, freqMax+dfreqMax]
409  double freq = pix->frequency*pix->rate/2.;
410  if(fabs(freq-freqMax)<=dfreqMax) {
411  likeLin += likePix;
412  fmean += likePix*freq;
413  frms += likePix*freq*freq;
414  }
415  }
416 
417  fmean = fmean/likeLin;
418  frms = frms/likeLin-fmean*fmean;
419  frms = frms>0 ? sqrt(frms) : 0.;
420 
421  if(frms<dfreqMax/2.) frms=dfreqMax/2.;
422 
423  // ------------------------------------------------------------------
424  // Save Lveto parameters
425  // ------------------------------------------------------------------
426 
427  Lveto[0] = fmean; // line mean frequency
428  Lveto[1] = frms; // line bandwidth
429  Lveto[2] = likeTot>0. ? likeLin/likeTot : 0.; // energy ratio energy inside_line/total
430 
431  // ------------------------------------------------------------------
432  // plot time-frequency energy
433  // ------------------------------------------------------------------
434 
435 #if defined PLOT_LIKELIHOOD
436  watplot WTS(const_cast<char*>("wts"));
437  WTS.plot(pwc, cid, nifo, 'L', 0, const_cast<char*>("COLZ"));
438  WTS.canvas->Print("l_tfmap_scalogram.png");
439 #endif
440 
441  return;
442 }
443 
444 void
445 GetWveto(netcluster* pwc, int cid, int nifo, float* Wveto) {
446 //
447 // input
448 // pwc : pointer to netcluster object
449 // cid : cluster id
450 // nifo : number of detectors
451 // output
452 // Wveto[0] : whistle slope
453 // Wveto[1] : whistle correlation
454 //
455 
456  Wveto[0] = Wveto[1] = 0;
457 
458  std::vector<int>* vint = &(pwc->cList[cid-1]); // pixel list
459  int V = vint->size(); // cluster size
460  if(!V) return;
461 
462  std::vector<double> x;
463  std::vector<double> y;
464  std::vector<double> w;
465 
466  // extract pixels
467  for(int j=0; j<V; j++) {
468  netpixel* pix = pwc->getPixel(cid,j);
469  if(pix->layers%2==0) {
470  cout << "CWB_Plugin_QLWveto.C - Error : is enabled only for WDM (2G)" << endl;
471  exit(1);
472  }
473  if(pix->likelihood<1. || pix->frequency==0) continue;
474  //if(!pix->core) continue; // select only the principal components pixels
475 
476  double time = int(double(pix->time)/pix->layers)/pix->rate; // time in seconds from the start
477  double freq = pix->frequency*pix->rate/2.;
478 
479  x.push_back(time);
480  y.push_back(freq);
481  w.push_back(pix->likelihood);
482  }
483  int size = x.size();
484  if(size<5) return;
485 
486  double xcm, ycm, qxx, qyy, qxy, ew;
487  xcm = ycm = qxx = qyy = qxy = ew = 0;
488 
489  for(int i=0; i<size; ++i) {
490  xcm += x[i]*w[i];
491  ycm += y[i]*w[i];
492  ew += w[i];
493  }
494  xcm /= ew;
495  ycm /= ew;
496 
497  for(int i=0; i<size; ++i) {
498  qxx += (x[i] - xcm)*(x[i] - xcm)*w[i];
499  qyy += (y[i] - ycm)*(y[i] - ycm)*w[i];
500  qxy += (x[i] - xcm)*(y[i] - ycm)*w[i];
501  }
502 
503  double beta = qxy/qxx; // slope
504  double alpha = ycm-beta*xcm; // intercept
505  double corr = qxy/sqrt(qxx*qyy); // correlation
506  double duration = sqrt(qxx/ew); // duration
507  double bandwidth = sqrt(qyy/ew); // bandwidth
508 
509  //printf("alpha : %lf , beta : %lf , corr : %lf, dur : %lf , bw : %lf \n",
510  // alpha, beta, corr, duration, bandwidth);
511 
512  Wveto[0] = beta; // slope
513  Wveto[1] = fabs(corr); // correlation
514 
515  return;
516 }
517 
518 void
519 PlotWaveform(TString ifo, wavearray<double>* wfREC,
520  CWB::config* cfg, bool fft, bool strain) {
521 
522  watplot PTS(const_cast<char*>("ptswrc"),200,20,800,500);
523 
524  //cout << "Print " << fname << endl;
525  double tmin = wfREC->start();
526  wfREC->start(wfREC->start()-tmin);
527  if(fft) {
528  PTS.plot(wfREC, const_cast<char*>("ALP"), 1, 0, 0, true, cfg->fLow, cfg->fHigh);
529  } else {
530  PTS.plot(wfREC, const_cast<char*>("ALP"), 1, 0, 0);
531  }
532  PTS.graph[0]->SetLineWidth(1);
533  wfREC->start(wfREC->start()+tmin);
534 
535  char label[64]="";
536  if(fft) sprintf(label,"%s","fft");
537  else sprintf(label,"%s","time");
538  if(strain) sprintf(label,"%s_%s",label,"strain");
539  else sprintf(label,"%s_%s",label,"white");
540 
541  char fname[1024];
542  //sprintf(fname, "%s_wf_%s_rec_gps_%d.root",ifo.Data(),label,int(tmin));
543  sprintf(fname, "%s_wf_%s_rec_gps_%d.png",ifo.Data(),label,int(tmin));
544  PTS.canvas->Print(fname);
545  cout << "write : " << fname << endl;
546  //PTS.canvas->Write(REPLACE(fname,dirCED,gtype));
547 }
548 
549 void
550 ClearWaveforms(detector* ifo) {
551 
552  int n;
553 
554  n = ifo->IWFP.size();
555  for (int i=0;i<n;i++) {
556  wavearray<double>* wf = (wavearray<double>*)ifo->IWFP[i];
557  delete wf;
558  }
559  ifo->IWFP.clear();
560  ifo->IWFID.clear();
561 
562  n = ifo->RWFP.size();
563  for (int i=0;i<n;i++) {
564  wavearray<double>* wf = (wavearray<double>*)ifo->RWFP[i];
565  delete wf;
566  }
567  ifo->RWFP.clear();
568  ifo->RWFID.clear();
569 }
570 
network::ifoName
std::vector< char * > ifoName
Definition: network.hh:609
network::getifo
detector * getifo(size_t n)
param: detector index
Definition: network.hh:436
cfg
CWB::config * cfg
Definition: TestCWB_Plugin.C:23
CWB::config::analysis
char analysis[8]
Definition: config.hh:117
outDump
TString outDump
Definition: CWB_Plugin_xWRC.C:115
CWB_Plugin
void CWB_Plugin(TFile *jfile, CWB::config *cfg, network *NET, WSeries< double > *x, TString ifo, int type)
COHERENCE.
Definition: CWB_Plugin_QLWveto.C:53
network::nLag
size_t nLag
Definition: network.hh:573
duration
double duration
Definition: DrawClusterWithSkyplot.C:21
netevent::setSLags
void setSLags(float *slag)
Definition: netevent.cc:426
wavearray::rate
virtual void rate(double r)
Definition: wavearray.hh:141
factor
float factor
Definition: cwb_report_cbc.C:777
a
wavearray< double > a(hp.size())
netpixel::frequency
size_t frequency
Definition: netpixel.hh:111
netpixel::likelihood
float likelihood
Definition: netpixel.hh:114
n
int n
Definition: cwb_net.C:28
GetQveto
float GetQveto(wavearray< double > *wf)
Definition: CWB_Plugin_QLWveto.C:277
beta
double beta
Definition: CreateListEBBH.C:20
netevent::setTree
TTree * setTree()
Definition: netevent.cc:434
TString
TString("c")
Definition: cwb_report_skymap.C:129
ID
int ID
Definition: TestMDC.C:70
x
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
Definition: BaseFunctionWDM.C:51
netcluster::get
wavearray< double > get(char *name, size_t index=0, char atype='R', int type=1, bool=true)
param: string with parameter name param: index in the amplitude array, which define detector param: c...
Definition: netcluster.cc:2207
detector::RWFP
std::vector< wavearray< double > * > RWFP
Definition: detector.hh:382
CWB::config::fLow
double fLow
Definition: config.hh:140
pix
netpixel pix(nifo)
pwc
netcluster * pwc
Definition: cwb_job_obj.C:38
watplot::graph
std::vector< TGraph * > graph
Definition: watplot.hh:194
CWB::config::cedDump
bool cedDump
Definition: config.hh:297
PlotWaveform
void PlotWaveform(TString ifo, wavearray< double > *wfREC, CWB::config *cfg, bool fft=false, bool strain=false)
Definition: CWB_Plugin_QLWveto.C:519
wf
waveform wf
Definition: CWB_Plugin_HEN_SIM_Config.C:88
size
Long_t size
Definition: cwb_condor_check.C:66
detector::waveBand
WSeries< double > waveBand
Definition: detector.hh:356
netpixel::layers
size_t layers
Definition: netpixel.hh:112
m
int m
Definition: cwb_net.C:28
netcluster::cList
std::vector< vector_int > cList
Definition: netcluster.hh:397
wavearray::start
virtual void start(double s)
Definition: wavearray.hh:137
j
int j
Definition: cwb_net.C:28
network::getMRAwave
bool getMRAwave(size_t ID, size_t lag, char atype='S', int mode=0, bool tof=false)
Definition: network.cc:3666
i
i drho i
Definition: cwb_epparameters.C:88
watplot::plot
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
CWB::config::outPlugin
bool outPlugin
Definition: config.hh:369
netpixel::core
bool core
Definition: netpixel.hh:120
ifo
char ifo[NIFO_MAX][8]
Definition: cwb1G_parameters.C:17
netevent::dopen
void dopen(const char *fname, char *mode, bool header=true)
Definition: netevent.hh:312
network::ifoListSize
size_t ifoListSize()
Definition: network.hh:431
fprintf
fprintf(stdout,"start=%f duration=%f rate=%f\, x.start(), x.size()/x.rate(), x.rate())
w
wavearray< double > w
Definition: Test1.C:27
nIFO
#define nIFO
Definition: DrawSensitivitiesS5.C:11
wavearray::size
virtual size_t size() const
Definition: wavearray.hh:145
watplot::canvas
TCanvas * canvas
Definition: watplot.hh:192
GetWveto
void GetWveto(netcluster *pwc, int cid, int nifo, float *Wveto)
Definition: CWB_Plugin_QLWveto.C:445
freq
wavearray< double > freq
Definition: Regression_H1.C:79
jfile
jfile
Definition: cwb_job_obj.C:43
netevent::output2G
void output2G(TTree *, network *, size_t, int, double)
Definition: netevent.cc:711
net_tree
TTree * net_tree
Definition: CWB_Plugin_xWRC.C:114
CWB::config::simulation
int simulation
Definition: config.hh:199
IMPORT
#define IMPORT(TYPE, VAR)
Definition: cwb.hh:69
WTS
watplot * WTS
Definition: ChirpMass.C:119
int
i() int(T_cor *100))
NET
network NET
Definition: cwb_dump_inj.C:30
label
TString label
Definition: MergeTrees.C:21
NIFO_MAX
const int NIFO_MAX
Definition: wat.hh:22
next
TIter next(twave->GetListOfBranches())
fname
char fname[1024]
Definition: cwb_report_prod_1.C:128
ATHR
#define ATHR
Definition: CWB_Plugin_QLWveto.C:42
detector::RWFID
std::vector< int > RWFID
Definition: detector.hh:381
k
int k
Definition: cwb_report_prod_2.C:168
detector::IWFP
std::vector< wavearray< double > * > IWFP
Definition: detector.hh:379
netpixel::time
size_t time
Definition: netpixel.hh:110
detector::IWFID
std::vector< int > IWFID
Definition: detector.hh:378
froot
TFile * froot
Definition: TestReadWriteDetectorObject.C:11
netcluster::getPixel
netpixel * getPixel(size_t n, size_t i)
Definition: netcluster.hh:413
dt
double dt
Definition: DrawPhaseShift.C:13
wavearray::FFTW
virtual void FFTW(int=1)
Definition: wavearray.cc:896
nifo
int nifo
Definition: cwb_dump_events.C:50
netevent::fP
FILE * fP
injected reconstructed xcor waveform
Definition: netevent.hh:141
network::wfsave
bool wfsave
Definition: network.hh:600
EVT
netevent EVT(itree, nifo)
netcluster::sCuts
std::vector< int > sCuts
Definition: netcluster.hh:392
detector::waveForm
WSeries< double > waveForm
Definition: detector.hh:355
netevent::dclose
void dclose()
Definition: netevent.hh:321
CWB::config::fHigh
double fHigh
Definition: config.hh:141
detector::Name
char Name[16]
Definition: detector.hh:327
fabs
double fabs(const Complex &x)
Definition: numpy.cc:55
gIFACTOR
int gIFACTOR
Definition: cwb_mkhtml_cbc.C:138
df
double df
Definition: HowToAccessToWSeries.C:36
network::getwc
netcluster * getwc(size_t n)
param: delay index
Definition: network.hh:439
ClearWaveforms
void ClearWaveforms(detector *ifo)
Definition: CWB_Plugin_QLWveto.C:550
sprintf
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
CWB::config::nIFO
int nIFO
Definition: config.hh:120
K
int K
Definition: CWB_Plugin_Sim4_BRST_Config.C:87
id
Long_t id
Definition: cwb_condor_check.C:66
CWB::config::factors
double factors[FACTORS_MAX]
Definition: config.hh:202
netpixel::rate
float rate
Definition: netpixel.hh:113
netpixel::getdata
double getdata(char type='R', size_t n=0)
Definition: netpixel.hh:74
GetLveto
void GetLveto(netcluster *pwc, int cid, int nifo, float *Lveto)
Definition: CWB_Plugin_QLWveto.C:341
NTHR
#define NTHR
Definition: CWB_Plugin_QLWveto.C:41
wavearray::resize
virtual void resize(unsigned int)
Definition: wavearray.cc:463
check
int check
Definition: CWB_configPlugin_Test.C:5
y
wavearray< double > y
Definition: Test10.C:31
exit
exit(0)
CWB::config::dump
bool dump
Definition: config.hh:295