Logo Coherent WaveBurst  
Reference Guide
Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Toolbox.hh
Go to the documentation of this file.
1 /**********************************************************
2  * Package: Toolbox Class Library
3  * File name: Toolbox.hh
4  * Author: Gabriele Vedovato (vedovato@lnl.infn.it)
5  **********************************************************/
6 
7 
8 #ifndef TOOLBOX_HH
9 #define TOOLBOX_HH
10 
11 #include "TROOT.h"
12 #include "TMath.h"
13 #include "TSystem.h"
14 #include "TString.h"
15 #include "TTree.h"
16 #include "TFile.h"
17 #include "TRandom.h"
18 #include "TRandom3.h"
19 #include "TObjArray.h"
20 #include "TObjString.h"
21 #include "TSystemDirectory.h"
22 #include "TChain.h"
23 #include "TApplication.h"
24 #include "TH1I.h"
25 #include "TF1.h"
26 #include "TH1D.h"
27 #include "TStyle.h"
28 #include "TLegend.h"
29 #include "TCanvas.h"
30 #include "TGraphSmooth.h"
31 
32 #include <string>
33 #include <iostream>
34 #include <fstream>
35 #include <stdlib.h>
36 #include <math.h>
37 #include <ctype.h>
38 #include <vector>
39 #include <string.h>
40 
41 #include "wavecomplex.hh"
42 #include "wavearray.hh"
43 #include "network.hh"
44 #include "Meyer.hh"
45 #include "netevent.hh"
46 
47 #include "FrameL.h"
48 
49 #include "History.hh"
50 #include "Toolfun.hh"
51 
52 #define LST_TREE_NAME "frl"
53 
54 enum CWB_CAT {
55  CWB_CAT0 = 0,
56  CWB_CAT1 = 1,
57  CWB_CAT2 = 2,
58  CWB_CAT3 = 3,
59  CWB_HVETO = 10,
60  CWB_PEM = 11,
61  CWB_EXC = 12,
62  CWB_USER = 13
63 };
64 
65 struct dqfile {
66  char ifo[32];
67  char file[1024];
69  double shift;
70  bool invert;
71  bool c4;
72 };
73 
74 struct frfile {
75  int start;
76  int stop;
77  int length;
78  vector<TString> file;
79 };
80 
81 struct slag {
82  int jobId; // job id : sequential progressive number
83  vector<int> slagId; // slag id vector : [0]=jobId - [1]=1/0 1=header slag - [2,..,nIFO+1] ifo slag
84  vector<int> segId; // seg id vector : [0,..,nIFO-1] ifo segment number
85 };
86 
87 struct mdcshift {
88  double startMDC;
89  double stopMDC;
90  double offset;
91 };
92 
93 struct ifoparms {
95  double latitude;
96  double longitude;
97  double elevation;
98  double AltX; // elevation of the x arm
99  double AzX; // azimut of the x arm (angle-deg from nord)
100  double AltY; // elevation of the y arm
101  double AzY; // azimut of the y arm (angle-deg from nord)
102 };
103 
104 static vector<double> DEFAULT_DOUBLE_VECTOR;
105 
106 using namespace std;
107 
108 namespace CWB {
109 
110 class Toolbox {
111 
112 public:
113 
114  /* ************************ */
115  /* * segment methods * */
116  /* ************************ */
117 
118  static vector<waveSegment>
119  readSegments(TString ifile);
120  static vector<waveSegment>
121  unionSegments(vector<waveSegment>& ilist);
122  static vector<waveSegment>
123  sortSegments(vector<waveSegment>& ilist);
124  static vector<waveSegment>
125  invertSegments(vector<waveSegment>& ilist);
126  static vector<waveSegment>
127  mergeSegLists(vector<waveSegment>& ilist1, vector<waveSegment>& ilist2);
128  static waveSegment
129  getMaxSeg(vector<waveSegment> list);
130  static double
131  getTimeSegList(vector<waveSegment> list);
132  static int
133  dumpSegList(vector<waveSegment> list, TString fName, bool c4=false);
134  static vector<waveSegment>
135  getSlagJobList(vector<waveSegment> ilist, int seglen=600);
136 
137  /* ************************ */
138  /* * dq methods * */
139  /* ************************ */
140 
141  static vector<waveSegment>
142  readSegList(dqfile DQF);
143  static vector<waveSegment>
144  readSegList(int nDQF, dqfile* DQF, CWB_CAT dqcat);
145  static void
146  setSlagShifts(slag SLAG, vector<TString> ifos, double segLen, int nDQF, dqfile* DQF);
147 
148  /* ************************ */
149  /* * job methods * */
150  /* ************************ */
151 
152  static vector<waveSegment>
153  getJobList(vector<waveSegment> ilist, double segLen=600.,
154  double segMLS=300., double segEdge=8.);
155  static vector<waveSegment>
156  getJobList(vector<waveSegment> cat1List, vector<waveSegment> cat2List,
157  double segLen=600., double segMLS=300., double segTHR=30., double segEdge=8.);
158  static void
159  dumpJobList(vector<waveSegment> ilist, TString fName, double segLen=600.,
160  double segMLS=300., double segEdge=8.);
161  static vector<waveSegment>
162  getSegList(int jobId, int nIFO, double segLen, double segMLS, double segEdge,
163  vector<waveSegment> dqList);
164  static double
165  getLiveTime(vector<waveSegment>& jobList, vector<waveSegment>& dqList);
166  static double
167  getLiveTime(vector<waveSegment>& jobList, vector<waveSegment>& dqList,
168  vector<double> shiftList);
169 
170  /* ************************ */
171  /* * slag methods * */
172  /* ************************ */
173 
174  static vector<slag>
175  getSlagList(size_t nIFO, size_t slagSize, int slagSegs,
176  int slagOff, size_t slagMin, size_t slagMax,
177  size_t* slagSite, char* slagFile=NULL);
178  static void
179  dumpSlagList(vector<slag> slagList, TString slagFile, bool slagOnly=false);
180  static slag
181  getSlag(vector<slag> slagList, int jobid);
182  static vector<slag>
183  getSlagList(vector<slag> slagList, vector<TString> ifos, double segLen,
184  double segMin, double segEdge, int nDQF, dqfile* iDQF, CWB_CAT dqcat);
185  static void
186  getSlagList(vector<slag>& slagList, int slagSize,
187  int slagRank, int nifo, vector<int>& id);
188  static vector<waveSegment>
189  getSegList(slag SLAG, vector<waveSegment> jobList, double segLen,
190  double segMLS, double segEdge, vector<waveSegment> dqList);
191 
192  /* ************************ */
193  /* * MDC methods * */
194  /* ************************ */
195 
196  static double
197  getMDCShift(mdcshift mshift, double time);
198  static int
199  shiftBurstMDCLog(std::vector<std::string>& mdcList,
200  vector<TString> ifos, double mdc_shift);
201  static TString
202  SetMDCLog(TString log, int pos, TString val);
203  static TString
204  SetMDCLog(TString log, int pos, double val);
205  static TString
206  GetMDCLog(TString log, int pos);
207  static int
208  GetMDCLogSize(TString log);
209 
210  /* ************************ */
211  /* * condor methods * */
212  /* ************************ */
213 
214  static int
215  createDagFile(vector<waveSegment> jobList, TString condor_dir, TString label,
216  int jobmin=0, int jobmax=0);
217  static int
218  createDagFile(vector<slag> slagList, TString condor_dir, TString label,
219  int jobmin=0, int jobmax=0);
220  static int
221  createDagFile(vector<int> jobList, TString condor_dir, TString label,
222  vector<TString> jobFiles, TString stage="CWB_STAGE_FULL", int jobmin=0, int jobmax=0);
223  static int
224  createDagFile(vector<waveSegment> jobList, TString condor_dir, TString label,
225  vector<TString> jobFiles, TString stage="CWB_STAGE_FULL", int jobmin=0, int jobmax=0);
226  static int
227  createDagFile(vector<slag> slagList, TString condor_dir, TString label,
228  vector<TString> jobFiles, TString stage="CWB_STAGE_FULL", int jobmin=0, int jobmax=0);
229  static int
230  createSubFile(TString label, TString condor_dir, TString out_dir, TString err_dir,
231  TString log_dir, TString ext="", TString condor_tag="");
232  static vector<int>
233  getCondorJobList(TString condor_dir, TString label);
234  static vector<float>
235  getJobBenchmark(TString ifName, int stageID, TString bench);
236  static vector<float>
237  getJobBenchmark(TString ifName, int stageID, int resID, int factorID, TString bench);
238  static TString
239  DAG2LSF(char* dagFile, char* data_label, char* nodedir, char* data_dir,
240  char* condor_dir, char* log_dir, char* output_dir, char* work_dir);
241 
242  /* ************************ */
243  /* * merge methods * */
244  /* ************************ */
245 
246  static vector<TString>
247  mergeCWBTrees(TString dir_name, bool simulation, TString odir,
248  TString label, bool brms=false, bool bvar=false, bool bpsm=false);
249  static void
250  mergeCWBTrees(vector<TString> fileList, bool simulation, TString odir,
251  TString label, bool brms=false, bool bvar=false, bool bpsm=false);
252  static vector<TString>
253  mergeCWBTrees(int nthreads, TString dir_name, bool simulation, TString odir,
254  TString label, bool brms=false, bool bvar=false, bool bpsm=false);
255  static void
256  mergeCWBTrees(int nthreads, vector<TString> fileList, bool simulation, TString odir,
257  TString label, bool brms=false, bool bvar=false, bool bpsm=false);
258  static void
259  mergeTrees(vector<TString> fileList, TString treeName, TString odir, TString ofName, bool bhistory);
260  static vector<TString>
261  readFileList(TString ifName);
262  static void
263  dumpFileList(vector<TString> fileList, TString ofName);
264  static vector<int>
265  getMergeJobList(TString merge_dir, TString label, int version);
266  static vector<int>
267  getMergeJobList(TString ifname, vector<TString>& jobFileList);
268  static vector<int>
269  getMergeJobList(TString ifname);
270 
271  /* ************************ */
272  /* * history methods * */
273  /* ************************ */
274 
275  static char*
276  readFile(TString ifName);
277  static char*
278  getEnvCWB();
279 
280  /* ************************ */
281  /* * spectra methods * */
282  /* ************************ */
283 
284  static void
285  makeSpectrum(wavearray<double>& psd, wavearray<double> x, double chuncklen=8,
286  double scratchlen=0, bool oneside=true);
287  static void
288  makeSpectrum(TString ofname, wavearray<double> x, double chuncklen=8,
289  double scratchlen=0, bool oneside=true);
290  static void
291  getSimNoise(wavearray<double>& u, TString fName, int seed, int run);
292  static void
293  convertSampleRate(wavearray<double> iw, wavearray<double> ow);
294  static void
295  convertSampleRate(wavearray<double> ix, wavearray<double> iy,
297  static wavearray<double>
298  GetDetectorPSD(TString fName, double fWidth=8192., double dFreq=1.);
299 
300  /* ************************ */
301  /* * system methods * */
302  /* ************************ */
303 
304  static vector<TString>
305  getFileListFromDir(TString dir_name, TString endString="",
306  TString beginString="", TString containString="",bool fast=false);
307  static std::map<int,TString>
308  getJobFileMapFromDir(TString dir_name, TString endString="",
309  TString beginString="", TString containString="",bool fast=false);
310  static bool
311  isFileExisting(TString fName);
312  static bool
313  checkFile(TString fName, bool question=false, TString message="");
314  static void
315  mkDir(TString dir, bool question=false, bool remove=true);
316  static bool
317  rmDir(TString dir, bool question=true);
318  static void
319  PrintProcInfo(TString str="");
320  static TString
321  getTemporaryFileName(TString label="temporary", TString extension="tmp");
322  static int
323  mksTemp(char *fTemplate) {return mkstemp(fTemplate);} // wrapper to system funcion (used in CINT)
324  static vector<TString>
325  sortStrings(vector<TString> ilist);
326 
327  /* ******************************** */
328  /* * post processing methods * */
329  /* ******************************** */
330 
331  static int
332  setVeto(TString ifName, TString idir, TString odir, vector<TString> ifos,
333  int nVDQF, dqfile* VDQF, int nDQF, dqfile* DQF, double segLen,
334  double segMLS, double segEdge);
335  static int
336  setCuts(TString ifName, TString idir, TString odir,
337  TString trname, TString cuts, TString olabel);
338  static int
339  setIFAR(TString ifName, TString idir, TString odir, TString trname,
340  TString sels, TString farFile, int irho, TString olabel, bool inclusive=true);
341  static int
342  setChunk(TString ifName, TString idir, TString odir, TString trname, int chunk);
343  static void
344  setUniqueEvents(TString ifwave, TString ofwave, int nIFO, int pp_irho);
345  static double
346  getLiveTime(int nIFO, TChain& liv, wavearray<double>& Trun,
348  wavearray<double>& Tdlag, int lag_number=-1, int slag_number=-1, int dummy=0);
349  static double
350  getZeroLiveTime(int nIFO, TChain& liv, int dummy=0);
351  static void
352  doPoissonPlot(int nIFO, wavearray<double>* Wlag,
354  static bool
355  isLeafInTree(TTree* itree, TString leaf);
356  static int
357  setMultiplicity(TString ifName, TString idir, TString odir, int nIFO, double dTime);
358  static wavecomplex
359  getRate(double rho, double Tgap, int nIFO, TChain& wav, wavearray<int>& Wsel,
361  static int GetStepFunction(TString fName, vector<double>& x, vector<double>& y,
362  vector<double>& ex = DEFAULT_DOUBLE_VECTOR, vector<double>& ey = DEFAULT_DOUBLE_VECTOR);
363  static double GetStepFunction(TString option, TString fName, double V=0,
364  vector<double>& x = DEFAULT_DOUBLE_VECTOR, vector<double>& y = DEFAULT_DOUBLE_VECTOR,
365  vector<double>& ex = DEFAULT_DOUBLE_VECTOR, vector<double>& ey = DEFAULT_DOUBLE_VECTOR);
366 
367  /* ************************ */
368  /* * utility methods * */
369  /* ************************ */
370 
371  static void
372  resampleToPowerOfTwo(wavearray<double>& w);
373  static int
374  getJobId(TString file, TString fext="root");
375  static TString
376  getParameter(TString options, TString param="");
377  static TString
378  getFileName(FILE* fp);
379  static TString
380  getFileName(char* symlink);
381  static void
383  static TString
384  WriteFrameFile(wavearray<double> x, TString chName, TString frName,
386  static bool
387  question(TString question);
388 
389  /* ************************************************** */
390  /* * Hilbert & Wigner-Ville Transforms methods * */
391  /* ************************************************** */
392 
393  static wavearray<double>
394  getHilbertTransform(wavearray<double> x);
395  static wavearray<double>
396  getHilbertEnvelope(wavearray<double> x);
397  static wavearray<double>
398  getHilbertIFrequency(wavearray<double> x);
399  static wavearray<double>
400  getWignerVilleTransform(wavearray<double> x);
401  static void
402  unWrapPhase(wavearray<float>& p);
403  static void
404  getSineFittingParams(double a, double b, double c, double rate,
405  double& amplitude, double& omega, double& phase);
406 
407 private:
408 
409  static void
410  blackmanharris (double* window, int n);
411 
412  ClassDef(Toolbox,1)
413 };
414 
415 } // end namespace
416 
417 #endif
double rho
double segMLS
Definition: test_config1.C:47
int slagSize
Definition: test_config1.C:65
char val[20]
double startMDC
Definition: Toolbox.hh:88
TString ofName
size_t * slagSite
Definition: test_config1.C:69
int slagOff
Definition: test_config1.C:68
int stop
Definition: Toolbox.hh:76
mdcshift mdc_shift
Definition: test_config1.C:93
TB createDagFile(rslagList, full_condor_dir, data_label, jobFiles, cwb_stage_name)
Definition: ced.hh:24
CWB run(runID)
wavearray< double > a(hp.size())
wavearray< double > Trun(500000)
cout<< "slagList size : "<< slagList.size()<< endl;cout<< endl<< "Start segments selection from dq cat1 list ..."<< endl<< endl;rslagList=TB.getSlagList(slagList, ifos, segLen, segMLS, segEdge, nDQF, DQF, CWB_CAT1);cout<< "Number of selected jobs after cat1 : "<< rslagList.size()<< endl;cout<< endl<< "Start segments selection from dq cat2 list ..."<< endl<< endl;rslagList=TB.getSlagList(rslagList, ifos, segLen, segTHR, segEdge, nDQF, DQF, CWB_CAT2);cout<< "Number of selected jobs after cat2 : "<< rslagList.size()<< endl;vector< TString > jobFiles
int n
Definition: cwb_net.C:10
double longitude
Definition: Toolbox.hh:96
dqfile VDQF[100]
bool invert
Definition: Toolbox.hh:70
TString("c")
char ifo[32]
Definition: Toolbox.hh:66
int start
Definition: Toolbox.hh:75
double stopMDC
Definition: Toolbox.hh:89
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<< "\t layers : "<< nLAYERS<< "\t dF(hz) : "<< dF<< "\t 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 phase
bool c4
Definition: Toolbox.hh:71
CWB_CAT dqcat
char odir[1024]
double amplitude
wavearray< double > psd(33)
int slagMax
Definition: test_config1.C:67
static int mksTemp(char *fTemplate)
Definition: Toolbox.hh:323
vector< int > segId
Definition: Toolbox.hh:84
CWB_CAT
Definition: Toolbox.hh:54
STL namespace.
TB createSubFile(data_label, full_condor_dir, full_condor_out_dir, full_condor_err_dir, condor_log, extention, condor_tag)
static vector< double > DEFAULT_DOUBLE_VECTOR
Definition: Toolbox.hh:104
char ofile[512]
double segEdge
Definition: test_config1.C:49
int nVDQF
wavearray< int > Wsel(ntrg)
TList * list
TB mergeTrees(fileList, TREE_NAME, OUTPUT_MERGE_DIR, OFILE_NAME, false)
int jobId
Definition: Toolbox.hh:82
double fWidth
Definition: DrawPSD.C:16
double Tgap
Definition: test_config1.C:23
TString symlink
nDQF
Definition: cwb_eced.C:92
double segTHR
Definition: test_config1.C:48
wavearray< double > w
Definition: Test1.C:27
#define nIFO
CWB_CAT cat
Definition: Toolbox.hh:68
char data_label[512]
Definition: test_config1.C:160
TString chName[NIFO_MAX]
char str[1024]
char file[1024]
Definition: Toolbox.hh:67
double time[6]
Definition: cbc_plots.C:435
wavearray< double > Rlag
int jobmax
int length
Definition: Toolbox.hh:77
vector< int > jobList
vector< slag > slagList
TString label
Definition: MergeTrees.C:21
bool log
Definition: WaveMDC.C:41
dqfile DQF[12]
Definition: test_config1.C:171
TString frDir[NIFO_MAX]
TB dumpSegList(dq1List, dq1ListFile, false)
wavearray< double > Tlag
segLen
Definition: cwb_eced.C:7
double shift
Definition: Toolbox.hh:69
char merge_dir[512]
Definition: test_config1.C:147
wavearray< double > Tdlag
char data_dir[512]
Definition: test_config1.C:152
double offset
Definition: Toolbox.hh:90
vector< TString > getFileListFromDir(TString dir_name, TString endString="", TString beginString="")
Definition: cwb_csh2sh.C:148
vector< waveSegment > cat1List
char log_dir[512]
Definition: test_config1.C:151
TFile * ifile
double elevation
Definition: Toolbox.hh:97
double AltY
Definition: Toolbox.hh:100
float irho
char options[256]
TString frName[NIFO_MAX]
int nifo
TB mkDir(eced_dir, check, remove)
double AltX
Definition: Toolbox.hh:98
bool inclusive
Definition: cwb_setifar.C:139
cout<< "Starting reading output directory ..."<< endl;vector< TString > fileList
int jobmin
vector< int > slagId
Definition: Toolbox.hh:83
double AzY
Definition: Toolbox.hh:101
TString dir_name[NDIR]
Definition: cwb_mkdir.C:52
string file
Definition: cwb_online.py:385
double latitude
Definition: Toolbox.hh:95
Definition: Toolbox.hh:81
cout<< "total cat1 livetime : "<< int(cat1_time)<< " sec "<< cat1_time/3600.<< " h "<< cat1_time/86400.<< " day"<< endl;cout<< endl;vector< waveSegment > cat2List
Definition: cwb_dump_job.C:22
double omega
Definition: testWDM_4.C:12
TTree * itree
wavearray< double > Wlag[NIFO_MAX+1]
slagFile
Definition: cwb_tune_slag.C:7
int jobId
char nodedir[1024]
Definition: test_config1.C:187
char condor_dir[512]
Definition: test_config1.C:148
double dFreq
Definition: DrawPSD.C:17
vector< TString > file
Definition: Toolbox.hh:78
char work_dir[512]
Definition: test_config1.C:143
int slagMin
Definition: test_config1.C:66
string version
Definition: cWB_conf.py:136
simulation
Definition: cwb_eced.C:9
char fName[256]
char output_dir[512]
Definition: test_config1.C:146
void getUniqueFileList(TString ifile, TString ofile)
Definition: FrDisplay.cc:1404
wavearray< double > y
Definition: Test10.C:31
TString frLabel[NIFO_MAX]
TString name
Definition: Toolbox.hh:94
i drho pp_irho
TB checkFile(gSystem->Getenv("CWB_ROOTLOGON_FILE"))
TString ifos[60]
double AzX
Definition: Toolbox.hh:99
wavearray< double > Wslag[NIFO_MAX+1]
TString treeName
Definition: MergeTrees.C:15