Logo Coherent WaveBurst  
Reference Guide 		Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
netevent.cc
Go to the documentation of this file.
1 // netevent class to process and store cWB triggers in root file
2 // S.Klimenko, University of Florida, Gainesville, FL
3 // G.Vedovato, University of Padova, Padova, Italy
4 
5 #include "netevent.hh"
6 #include "wseries.hh"
7 
8 #include "TH2.h"
9 #include "TStyle.h"
10 #include "TCanvas.h"
11 #include "TSystem.h"
12 #include "TMath.h"
13 #include "TMarker.h"
14 #include "TVector3.h"
15 #include "TRotation.h"
16 #include "TPolyLine.h"
17 #include "Math/Rotation3D.h"
18 #include "Math/Vector3Dfwd.h"
19 
20 #include "Meyer.hh"
21 #include <string.h>
22 
23 #define WAVE_TREE_NAME "waveburst"
24 
25 ClassImp(netevent) // used by THtml doc
26 
27 using namespace ROOT::Math;
28 
29 TTree* netevent::Init(TString fName, int n)
30 {
31  iFile = TFile::Open(fName);
32  if((iFile==NULL) || (iFile!=NULL && !iFile->IsOpen())) {
33  cout << "netevent::Init : Error opening root file " << fName.Data() << endl;
34  exit(1);
35  }
36 
37  TTree* tree = (TTree *) iFile->Get(WAVE_TREE_NAME);
38  if(tree) {
39  ndim = tree->GetUserInfo()->GetSize(); // get number of detectors
40  if(ndim>0) {
41  if(n>0 && ndim!=n) {
42  cout << "netevent::Init : number of detectors declared in the constructor (" << n
43  << ") are not equals to the one ("<<ndim<<") declared in the root file : "
44  << fName.Data() << endl;
45  exit(1);
46  }
47  } else ndim=n;
48  } else {
49  cout << "netevent::Init : object tree " << WAVE_TREE_NAME
50  << " not present in the root file " << fName.Data() << endl;
51  exit(1);
52  }
53 
54 
55  if(ndim==0) {
56  cout << "netevent::Init : detector number is not declared in the constructor or"
57  << " not present in the root file " << fName.Data() << endl;
58  exit(1);
59  }
60 
61 
62  return tree;
63 }
64 
65 // Set branch addresses
66 void netevent::Init(TTree *tree)
67 {
68  if (tree == 0) return;
69  fChain = tree;
70  fCurrent = -1;
71  fChain->SetMakeClass(1);
72 
73  fChain->SetBranchAddress("ndim",&ndim);
74  fChain->SetBranchAddress("run",&run);
75  fChain->SetBranchAddress("nevent",&nevent);
76  fChain->SetBranchAddress("eventID",eventID);
77  fChain->SetBranchAddress("type",type);
78  fChain->SetBranchAddress("name",&name);
79  fChain->SetBranchAddress("rate",rate);
80 
81  fChain->SetBranchAddress("volume",volume);
82  fChain->SetBranchAddress("size",size);
83  fChain->SetBranchAddress("usize",&usize);
84 
85  fChain->SetBranchAddress("gap",gap);
86  fChain->SetBranchAddress("lag",lag);
87  fChain->SetBranchAddress("slag",slag);
88  fChain->SetBranchAddress("strain",strain);
89  fChain->SetBranchAddress("phi",phi);
90  fChain->SetBranchAddress("theta",theta);
91  fChain->SetBranchAddress("psi",psi);
92  fChain->SetBranchAddress("iota",iota);
93  fChain->SetBranchAddress("bp",bp);
94  fChain->SetBranchAddress("bx",bx);
95 
96 
97  fChain->SetBranchAddress("time",time);
98  fChain->SetBranchAddress("gps",gps);
99  fChain->SetBranchAddress("right",right);
100  fChain->SetBranchAddress("left",left);
101  fChain->SetBranchAddress("duration",duration);
102  fChain->SetBranchAddress("start",start);
103  fChain->SetBranchAddress("stop",stop);
104 
105  fChain->SetBranchAddress("frequency",frequency);
106  fChain->SetBranchAddress("low",low);
107  fChain->SetBranchAddress("high",high);
108  fChain->SetBranchAddress("bandwidth",bandwidth);
109 
110  fChain->SetBranchAddress("hrss",hrss);
111  fChain->SetBranchAddress("noise",noise);
112  fChain->SetBranchAddress("erA",erA);
113  if(fChain->GetBranch("Psm")!=NULL) fChain->SetBranchAddress("Psm",&Psm);
114  fChain->SetBranchAddress("null",null);
115  fChain->SetBranchAddress("nill",nill);
116  fChain->SetBranchAddress("netcc",netcc);
117  fChain->SetBranchAddress("neted",neted);
118  fChain->SetBranchAddress("rho",rho);
119 
120  fChain->SetBranchAddress("gnet",&gnet);
121  fChain->SetBranchAddress("anet",&anet);
122  fChain->SetBranchAddress("inet",&inet);
123  fChain->SetBranchAddress("ecor",&ecor);
124  fChain->SetBranchAddress("norm",&norm);
125  fChain->SetBranchAddress("ECOR",&ECOR);
126  fChain->SetBranchAddress("penalty",&penalty);
127  fChain->SetBranchAddress("likelihood",&likelihood);
128 
129  fChain->SetBranchAddress("factor",&factor);
130  fChain->SetBranchAddress("range",range);
131  fChain->SetBranchAddress("chirp",chirp);
132  fChain->SetBranchAddress("eBBH",eBBH);
133  fChain->SetBranchAddress("Deff",Deff);
134  fChain->SetBranchAddress("mass",mass);
135  fChain->SetBranchAddress("spin",spin);
136 
137  fChain->SetBranchAddress("snr",snr);
138  fChain->SetBranchAddress("xSNR",xSNR);
139  fChain->SetBranchAddress("sSNR",sSNR);
140  fChain->SetBranchAddress("iSNR",iSNR);
141  fChain->SetBranchAddress("oSNR",oSNR);
142  fChain->SetBranchAddress("ioSNR",ioSNR);
143 
144  Notify();
145 }
146 
147 // allocate memory
148 void netevent::allocate()
149 {
150  if(ndim==0) {
151  cout << "netevent::allocate : Error - number of detectors must be > 0" << endl;
152  exit(1);
153  }
154 
155  if (!eventID) eventID= (Int_t*)malloc(2*sizeof(Int_t));
156  else eventID= (Int_t*)realloc(eventID,2*sizeof(Int_t));
157  if (!type) type= (Int_t*)malloc(2*sizeof(Int_t));
158  else type= (Int_t*)realloc(type,2*sizeof(Int_t));
159  if (!name) name= new string();
160  else {delete name;name= new string();}
161  if (!rate) rate= (Int_t*)malloc(ndim*sizeof(Int_t));
162  else rate= (Int_t*)realloc(rate,ndim*sizeof(Int_t));
163 
164  if (!volume) volume= (Int_t*)malloc(ndim*sizeof(Int_t));
165  else volume= (Int_t*)realloc(volume,ndim*sizeof(Int_t));
166  if (!size) size= (Int_t*)malloc(ndim*sizeof(Int_t));
167  else size= (Int_t*)realloc(size,ndim*sizeof(Int_t));
168 
169  if (!gap) gap= (Float_t*)malloc(ndim*sizeof(Float_t));
170  else gap= (Float_t*)realloc(gap,ndim*sizeof(Float_t));
171  if (!lag) lag= (Float_t*)malloc((ndim+1)*sizeof(Float_t));
172  else lag= (Float_t*)realloc(lag,(ndim+1)*sizeof(Float_t));
173  if (!slag) slag= (Float_t*)malloc((ndim+1)*sizeof(Float_t));
174  else slag= (Float_t*)realloc(slag,(ndim+1)*sizeof(Float_t));
175  if (!gps) gps= (Double_t*)malloc((ndim)*sizeof(Double_t));
176  else gps= (Double_t*)realloc(gps,(ndim)*sizeof(Double_t));
177  if (!strain) strain= (Double_t*)malloc(2*sizeof(Double_t));
178  else strain= (Double_t*)realloc(strain,2*sizeof(Double_t));
179  if (!phi) phi= (Float_t*)malloc(4*sizeof(Float_t));
180  else phi= (Float_t*)realloc(phi,4*sizeof(Float_t));
181  if (!theta) theta= (Float_t*)malloc(4*sizeof(Float_t));
182  else theta= (Float_t*)realloc(theta,4*sizeof(Float_t));
183  if (!psi) psi= (Float_t*)malloc(2*sizeof(Float_t));
184  else psi= (Float_t*)realloc(psi,2*sizeof(Float_t));
185  if (!iota) iota= (Float_t*)malloc(2*sizeof(Float_t));
186  else iota= (Float_t*)realloc(iota,2*sizeof(Float_t));
187  if (!bp) bp= (Float_t*)malloc(ndim*2*sizeof(Float_t));
188  else bp= (Float_t*)realloc(bp,ndim*2*sizeof(Float_t));
189  if (!bx) bx= (Float_t*)malloc(ndim*2*sizeof(Float_t));
190  else bx= (Float_t*)realloc(bx,ndim*2*sizeof(Float_t));
191 
192  if (!time) time= (Double_t*)malloc(ndim*2*sizeof(Double_t));
193  else time= (Double_t*)realloc(time,ndim*2*sizeof(Double_t));
194  if (!right) right= (Float_t*)malloc(ndim*sizeof(Float_t));
195  else right= (Float_t*)realloc(right,ndim*sizeof(Float_t));
196  if (!left) left= (Float_t*)malloc(ndim*sizeof(Float_t));
197  else left= (Float_t*)realloc(left,ndim*sizeof(Float_t));
198  if (!duration) duration= (Float_t*)malloc(ndim*sizeof(Float_t));
199  else duration= (Float_t*)realloc(duration,ndim*sizeof(Float_t));
200  if (!start) start= (Double_t*)malloc(ndim*sizeof(Double_t));
201  else start= (Double_t*)realloc(start,ndim*sizeof(Double_t));
202  if (!stop) stop= (Double_t*)malloc(ndim*sizeof(Double_t));
203  else stop= (Double_t*)realloc(stop,ndim*sizeof(Double_t));
204 
205  if (!frequency) frequency=(Float_t*)malloc(ndim*sizeof(Float_t));
206  else frequency=(Float_t*)realloc(frequency,ndim*sizeof(Float_t));
207  if (!low) low= (Float_t*)malloc(ndim*sizeof(Float_t));
208  else low= (Float_t*)realloc(low,ndim*sizeof(Float_t));
209  if (!high) high= (Float_t*)malloc(ndim*sizeof(Float_t));
210  else high= (Float_t*)realloc(high,ndim*sizeof(Float_t));
211  if (!bandwidth) bandwidth=(Float_t*)malloc(ndim*sizeof(Float_t));
212  else bandwidth=(Float_t*)realloc(bandwidth,ndim*sizeof(Float_t));
213  if (!hrss) hrss= (Double_t*)malloc(ndim*2*sizeof(Double_t));
214  else hrss= (Double_t*)realloc(hrss,ndim*2*sizeof(Double_t));
215  if (!noise) noise= (Double_t*)malloc(ndim*sizeof(Double_t));
216  else noise= (Double_t*)realloc(noise,ndim*sizeof(Double_t));
217  if (!null) null= (Float_t*)malloc(ndim*sizeof(Float_t));
218  else null= (Float_t*)realloc(null,ndim*sizeof(Float_t));
219  if (!nill) nill= (Float_t*)malloc(ndim*sizeof(Float_t));
220  else nill= (Float_t*)realloc(nill,ndim*sizeof(Float_t));
221  if (!netcc) netcc= (Float_t*)malloc(4*sizeof(Float_t));
222  else netcc= (Float_t*)realloc(netcc,4*sizeof(Float_t));
223  if (!neted) neted= (Float_t*)malloc(5*sizeof(Float_t));
224  else neted= (Float_t*)realloc(neted,5*sizeof(Float_t));
225  if (!rho) rho= (Float_t*)malloc(2*sizeof(Float_t));
226  else rho= (Float_t*)realloc(rho,2*sizeof(Float_t));
227  if (!erA) erA= (Float_t*)malloc(11*sizeof(Float_t));
228  else erA= (Float_t*)realloc(erA,11*sizeof(Float_t));
229  if (!range) range= (Float_t*)malloc(2*sizeof(Float_t));
230  else range= (Float_t*)realloc(range,2*sizeof(Float_t));
231  if (!chirp) chirp= (Float_t*)malloc(6*sizeof(Float_t));
232  else chirp= (Float_t*)realloc(chirp,6*sizeof(Float_t));
233  if (!eBBH) eBBH= (Float_t*)malloc(4*sizeof(Float_t));
234  else eBBH= (Float_t*)realloc(eBBH,4*sizeof(Float_t));
235  if (!Deff) Deff= (Float_t*)malloc(ndim*sizeof(Float_t));
236  else Deff= (Float_t*)realloc(Deff,ndim*sizeof(Float_t));
237  if (!mass) mass= (Float_t*)malloc(2*sizeof(Float_t));
238  else mass= (Float_t*)realloc(mass,ndim*sizeof(Float_t));
239  if (!spin) spin= (Float_t*)malloc(6*sizeof(Float_t));
240  else spin= (Float_t*)realloc(spin,6*sizeof(Float_t));
241 
242  if (!snr) snr= (Float_t*)malloc(ndim*sizeof(Float_t));
243  else snr= (Float_t*)realloc(snr,ndim*sizeof(Float_t));
244  if (!xSNR) xSNR= (Float_t*)malloc(ndim*sizeof(Float_t));
245  else xSNR= (Float_t*)realloc(xSNR,ndim*sizeof(Float_t));
246  if (!sSNR) sSNR= (Float_t*)malloc(ndim*sizeof(Float_t));
247  else sSNR= (Float_t*)realloc(sSNR,ndim*sizeof(Float_t));
248  if (!iSNR) iSNR= (Float_t*)malloc(ndim*sizeof(Float_t));
249  else iSNR= (Float_t*)realloc(iSNR,ndim*sizeof(Float_t));
250  if (!oSNR) oSNR= (Float_t*)malloc(ndim*sizeof(Float_t));
251  else oSNR= (Float_t*)realloc(oSNR,ndim*sizeof(Float_t));
252  if (!ioSNR) ioSNR= (Float_t*)malloc(ndim*sizeof(Float_t));
253  else ioSNR= (Float_t*)realloc(ioSNR,ndim*sizeof(Float_t));
254 
255  if (!Psm) Psm= WAT::USE_HEALPIX() ? new skymap(int(0)) : new skymap(0.);
256  else {delete Psm; Psm= WAT::USE_HEALPIX() ? new skymap(int(0)) : new skymap(0.);}
257 
258  return;
259 }
260 
261 // init array
262 void netevent::init()
263 {
264  for(int i=0;i<2;i++) eventID[i]=0;
265  for(int i=0;i<2;i++) type[i]=0;
266  for(int i=0;i<2;i++) strain[i]=0;
267  for(int i=0;i<4;i++) phi[i]=0;
268  for(int i=0;i<2;i++) psi[i]=0;
269  for(int i=0;i<4;i++) theta[i]=0;
270  for(int i=0;i<2;i++) psi[i]=0;
271  for(int i=0;i<2;i++) iota[i]=0;
272  for(int i=0;i<4;i++) netcc[i]=0;
273  for(int i=0;i<5;i++) neted[i]=0;
274  for(int i=0;i<2;i++) rho[i]=0;
275  for(int i=0;i<11;i++) erA[i]=0;
276  for(int i=0;i<2;i++) mass[i]=0;
277  for(int i=0;i<2;i++) range[i]=0;
278  for(int i=0;i<6;i++) chirp[i]=0;
279  for(int i=0;i<4;i++) eBBH[i]=0;
280  for(int i=0;i<6;i++) spin[i]=0;
281 
282  if(Psm) *Psm=0;
283 
284  for(int i=0;i<ndim+1;i++) lag[i]=0;
285  for(int i=0;i<ndim+1;i++) slag[i]=0;
286 
287  for(int i=0;i<2*ndim;i++) bp[i]=0;
288  for(int i=0;i<2*ndim;i++) bx[i]=0;
289  for(int i=0;i<2*ndim;i++) hrss[i]=0;
290  for(int i=0;i<2*ndim;i++) time[i]=0;
291 
292  for(int i=0;i<ndim;i++) {
293  rate[i]=0; volume[i]=0; size[i]=0; gap[i]=0;
294  gps[i]=0; snr[i]=0; right[i]=0; left[i]=0;
295  duration[i]=0; start[i]=0; stop[i]=0; frequency[i]=0; low[i]=0;
296  high[i]=0; bandwidth[i]=0; noise[i]=0; null[i]=0; nill[i]=0;
297  sSNR[i]=0; Deff[i]=0; iSNR[i]=0; oSNR[i]=0; ioSNR[i]=0; xSNR[i]=0;
298  }
299 
300  return;
301 }
302 
303 
304 Bool_t netevent::Notify()
305 {
306  // Called when loading a new file.
307  // Get branch pointers.
308  b_ndim = fChain->GetBranch("ndim");
309  b_run = fChain->GetBranch("run");
310  b_nevent = fChain->GetBranch("nevent");
311  b_eventID = fChain->GetBranch("eventID");
312  b_type = fChain->GetBranch("type");
313  b_name = fChain->GetBranch("name");
314  b_rate = fChain->GetBranch("rate");
315 
316  b_volume = fChain->GetBranch("volume");
317  b_size = fChain->GetBranch("size");
318  b_usize = fChain->GetBranch("usize");
319 
320  b_gap = fChain->GetBranch("gap");
321  b_lag = fChain->GetBranch("lag");
322  b_slag = fChain->GetBranch("slag");
323  b_strain = fChain->GetBranch("strain");
324  b_phi = fChain->GetBranch("phi");
325  b_theta = fChain->GetBranch("theta");
326  b_psi = fChain->GetBranch("psi");
327  b_iota = fChain->GetBranch("iota");
328  b_bp= fChain->GetBranch("bx");
329  b_bx = fChain->GetBranch("bp");
330 
331  b_time = fChain->GetBranch("time");
332  b_gps = fChain->GetBranch("gps");
333  b_right = fChain->GetBranch("right");
334  b_left = fChain->GetBranch("left");
335  b_start = fChain->GetBranch("start");
336  b_stop = fChain->GetBranch("stop");
337  b_duration = fChain->GetBranch("duration");
338 
339  b_frequency = fChain->GetBranch("frequency");
340  b_low = fChain->GetBranch("low");
341  b_high = fChain->GetBranch("high");
342  b_bandwidth = fChain->GetBranch("bandwidth");
343 
344  b_hrss = fChain->GetBranch("hrss");
345  b_noise = fChain->GetBranch("noise");
346  b_erA = fChain->GetBranch("erA");
347  b_Psm = fChain->GetBranch("Psm");
348  b_null = fChain->GetBranch("null");
349  b_nill = fChain->GetBranch("nill");
350  b_netcc = fChain->GetBranch("netcc");
351  b_neted = fChain->GetBranch("neted");
352  b_rho = fChain->GetBranch("rho");
353 
354  b_gnet = fChain->GetBranch("gnet");
355  b_anet = fChain->GetBranch("anet");
356  b_inet = fChain->GetBranch("inet");
357  b_ecor = fChain->GetBranch("ecor");
358  b_norm = fChain->GetBranch("norm");
359  b_ECOR = fChain->GetBranch("ECOR");
360  b_penalty = fChain->GetBranch("penalty");
361  b_likelihood = fChain->GetBranch("likelihood");
362 
363  b_factor = fChain->GetBranch("factor");
364  b_range = fChain->GetBranch("range");
365  b_chirp = fChain->GetBranch("chirp");
366  b_eBBH = fChain->GetBranch("eBBH");
367  b_Deff = fChain->GetBranch("Deff");
368  b_mass = fChain->GetBranch("mass");
369  b_spin = fChain->GetBranch("spin");
370 
371  b_snr = fChain->GetBranch("snr");
372  b_xSNR = fChain->GetBranch("xSNR");
373  b_sSNR = fChain->GetBranch("sSNR");
374  b_iSNR = fChain->GetBranch("iSNR");
375  b_oSNR = fChain->GetBranch("oSNR");
376  b_ioSNR = fChain->GetBranch("ioSNR");
377 
378  return kTRUE;
379 }
380 
381 Int_t netevent::GetEntries()
382 {
383  if (!fChain) return 0;
384  return fChain->GetEntries();
385 };
386 
387 Int_t netevent::GetEntry(Int_t entry)
388 {
389  if (!fChain) return 0;
390  return fChain->GetEntry(entry);
391 };
392 
393 void netevent::Show(Int_t entry)
394 {
395 // Print contents of entry.
396 // If entry is not specified, print current entry
397  if (!fChain) return;
398  fChain->Show(entry);
399 }
400 
401 //++++++++++++++++++++++++++++++++++++++++++++++
402 // set super lags
403 //++++++++++++++++++++++++++++++++++++++++++++++
404 void netevent::setSLags(float* slag)
405 {
406  for(int n=0;n<=ndim;n++) this->slag[n]=slag[n];
407 }
408 
409 //++++++++++++++++++++++++++++++++++++++++++++++
410 // set single event tree
411 //++++++++++++++++++++++++++++++++++++++++++++++
412 TTree* netevent::setTree()
413 {
414  TTree* waveTree = new TTree(WAVE_TREE_NAME,WAVE_TREE_NAME);
415 
416  int i;
417 
418  char crate[16];
419 
420  char cvolume[16];
421  char csize[16];
422 
423  char cgap[16];
424  char clag[16];
425  char cslag[16];
426  char cgps[16];
427  char cbp[16];
428  char cbx[16];
429 
430  char ctime[16];
431  char cright[16];
432  char cleft[16];
433  char cduration[16];
434  char cstart[16];
435  char cstop[16];
436 
437  char cfrequency[16];
438  char clow[16];
439  char chigh[16];
440  char cbandwidth[16];
441  char chrss[16];
442  char cnoise[16];
443  char cnull[16];
444  char cnill[16];
445  char crange[16];
446  char cchirp[16];
447  char ceBBH[16];
448  char cDeff[16];
449  char cmass[16];
450  char cspin[16];
451  char csnr[16];
452  char csSNR[16];
453  char cxSNR[16];
454  char ciSNR[16];
455  char coSNR[16];
456  char cioSNR[16];
457 
458 
459  sprintf(crate, "rate[%1d]/I",ndim);
460 
461  sprintf(cvolume, "volume[%1d]/I",ndim);
462  sprintf(csize, "size[%1d]/I",ndim);
463 
464  sprintf(cgap, "gap[%1d]/F",ndim);
465  sprintf(clag, "lag[%1d]/F",ndim+1);
466  sprintf(cslag, "slag[%1d]/F",ndim+1);
467  sprintf(cgps, "gps[%1d]/D",ndim);
468 
469  if(ndim<5) sprintf(cbp, "bp[%1d]/F",ndim*2);
470  else sprintf(cbp, "bp[%2d]/F",ndim*2);
471  if(ndim<5) sprintf(cbx, "bx[%1d]/F",ndim*2);
472  else sprintf(cbx, "bx[%2d]/F",ndim*2);
473 
474  if(ndim<5) sprintf(ctime, "time[%1d]/D",ndim*2);
475  else sprintf(ctime, "time[%2d]/D",ndim*2);
476 
477  sprintf(cright, "right[%1d]/F",ndim);
478  sprintf(cleft, "left[%1d]/F",ndim);
479  sprintf(cduration, "duration[%1d]/F",ndim);
480  sprintf(cstart, "start[%1d]/D",ndim);
481  sprintf(cstop, "stop[%1d]/D",ndim);
482 
483  sprintf(cfrequency, "frequency[%1d]/F",ndim);
484  sprintf(clow, "low[%1d]/F",ndim);
485  sprintf(chigh, "high[%1d]/F",ndim);
486  sprintf(cbandwidth, "bandwidth[%1d]/F",ndim);
487 
488  if(ndim<5) sprintf(chrss, "hrss[%1d]/D",ndim*2);
489  else sprintf(chrss, "hrss[%2d]/D",ndim*2);
490 
491  sprintf(cnoise, "noise[%1d]/D",ndim);
492  sprintf(cnull, "null[%1d]/F",ndim);
493  sprintf(cnill, "nill[%1d]/F",ndim);
494  sprintf(cDeff, "Deff[%1d]/F",ndim);
495  sprintf(crange, "range[2]/F");
496  sprintf(ceBBH, "eBBH[4]/F");
497  sprintf(cchirp, "chirp[6]/F");
498  sprintf(cmass, "mass[2]/F");
499  sprintf(cspin, "spin[6]/F");
500  sprintf(csnr, "snr[%1d]/F",ndim);
501  sprintf(csSNR, "sSNR[%1d]/F",ndim);
502  sprintf(cxSNR, "xSNR[%1d]/F",ndim);
503  sprintf(ciSNR, "iSNR[%1d]/F",ndim);
504  sprintf(coSNR, "oSNR[%1d]/F",ndim);
505  sprintf(cioSNR, "ioSNR[%1d]/F",ndim);
506 
507  //==================================
508  // Define trigger tree
509  //==================================
510 
511  waveTree->Branch("ndim", &ndim, "ndim/I");
512  waveTree->Branch("run", &run, "run/I");
513  waveTree->Branch("nevent", &nevent, "nevent/I");
514  waveTree->Branch("eventID", eventID, "eventID[2]/I");
515  waveTree->Branch("type", type, "type[2]/I");
516  waveTree->Branch("name", name);
517  waveTree->Branch("rate", rate, crate);
518 
519  waveTree->Branch("usize", &usize, "usize/I");
520  waveTree->Branch("volume", volume, cvolume);
521  waveTree->Branch("size", size, csize);
522 
523  waveTree->Branch("gap", gap, cgap);
524  waveTree->Branch("lag", lag, clag);
525  waveTree->Branch("slag", slag, cslag);
526  waveTree->Branch("strain", strain, "strain[2]/D");
527  waveTree->Branch("phi", phi, "phi[4]/F");
528  waveTree->Branch("theta", theta, "theta[4]/F");
529  waveTree->Branch("psi", psi, "psi[2]/F");
530  waveTree->Branch("iota", iota, "iota[2]/F");
531  waveTree->Branch("bp", bp, cbp);
532  waveTree->Branch("bx", bx, cbx);
533 
534  waveTree->Branch("time", time, ctime);
535  waveTree->Branch("gps", gps, cgps);
536  waveTree->Branch("left", left, cleft);
537  waveTree->Branch("right", right, cright);
538  waveTree->Branch("start", start, cstart);
539  waveTree->Branch("stop", stop, cstop);
540  waveTree->Branch("duration", duration, cduration);
541 
542  waveTree->Branch("frequency", frequency, cfrequency);
543  waveTree->Branch("low", low, clow);
544  waveTree->Branch("high", high, chigh);
545  waveTree->Branch("bandwidth", bandwidth, cbandwidth);
546 
547  waveTree->Branch("hrss", hrss, chrss);
548  waveTree->Branch("noise", noise, cnoise);
549  //waveTree->Branch("ndm", ndm, cndm);
550  waveTree->Branch("erA", erA, "erA[11]/F");
551 
552  waveTree->Branch("Psm","skymap",&Psm, 32000,0);
553 
554  waveTree->Branch("null", null, cnull);
555  waveTree->Branch("nill", nill, cnill);
556  waveTree->Branch("netcc", netcc, "netcc[4]/F");
557  waveTree->Branch("neted", neted, "neted[5]/F");
558  waveTree->Branch("rho", rho, "rho[2]/F");
559 
560  waveTree->Branch("gnet", &gnet, "gnet/F");
561  waveTree->Branch("anet", &anet, "anet/F");
562  waveTree->Branch("inet", &inet, "inet/F");
563  waveTree->Branch("ecor", &ecor, "ecor/F");
564  waveTree->Branch("norm", &norm, "norm/F");
565  waveTree->Branch("ECOR", &ECOR, "ECOR/F");
566  waveTree->Branch("penalty", &penalty, "penalty/F");
567  waveTree->Branch("likelihood", &likelihood, "likelihood/F");
568 
569  waveTree->Branch("factor", &factor, "factor/F");
570  waveTree->Branch("chirp", chirp, cchirp);
571  waveTree->Branch("range", range, crange);
572  waveTree->Branch("eBBH", eBBH, ceBBH);
573  waveTree->Branch("Deff", Deff, cDeff);
574  waveTree->Branch("mass", mass, cmass);
575  waveTree->Branch("spin", spin, cspin);
576 
577  waveTree->Branch("snr", snr, csnr);
578  waveTree->Branch("xSNR", xSNR, cxSNR);
579  waveTree->Branch("sSNR", sSNR, csSNR);
580  waveTree->Branch("iSNR", iSNR, ciSNR);
581  waveTree->Branch("oSNR", oSNR, coSNR);
582  waveTree->Branch("ioSNR", ioSNR, cioSNR);
583 
584  return waveTree;
585 }
586 
587 netevent& netevent::operator=(const netevent& a)
588 {
589  int i,j;
590 
591  if(ndim != a.ndim) { ndim=a.ndim; this->allocate(); }
592 
593  ifoList= a.ifoList;
594 
595  run= a.run;
596  nevent= a.nevent;
597  eventID[0]= a.eventID[0];
598  eventID[1]= a.eventID[1];
599  type[0]= a.type[0];
600  type[1]= a.type[1];
601  *name= *a.name;
602  strain[0]= a.strain[0];
603  strain[1]= a.strain[1];
604  usize= a.usize;
605  psi[0]= a.psi[0];
606  psi[1]= a.psi[1];
607  iota[0]= a.iota[0];
608  iota[1]= a.iota[1];
609  mass[0]= a.mass[0];
610  mass[1]= a.mass[1];
611  gnet= a.gnet;
612  anet= a.anet;
613  inet= a.inet;
614  ecor= a.ecor;
615  norm= a.norm;
616  ECOR= a.ECOR;
617  penalty= a.penalty;
618  likelihood= a.likelihood;
619  factor= a.factor;
620  Psave= a.Psave;
621  *Psm = *(a.Psm);
622 
623  for(i=0; i<6; i++) spin[i] = a.spin[i];
624  for(i=0; i<4; i++) netcc[i] = a.netcc[i];
625  for(i=0; i<5; i++) neted[i] = a.neted[i];
626  for(i=0; i<2; i++) rho[i] = a.rho[i];
627  for(i=0; i<4; i++) theta[i] = a.theta[i];
628  for(i=0; i<4; i++) phi[i] = a.phi[i];
629  for(i=0; i<11; i++) erA[i] = a.erA[i];
630  for(i=0; i<2; i++) range[i] = a.range[i];
631  for(i=0; i<4; i++) eBBH[i] = a.eBBH[i];
632  for(i=0; i<6; i++) chirp[i] = a.chirp[i];
633 
634  lag[ndim]= a.lag[ndim];
635  slag[ndim]= a.slag[ndim];
636  for(i=0; i<ndim; i++){
637 
638  gps[i]= a.gps[i];
639 
640  rate[i]= a.rate[i];
641 
642  gap[i]= a.gap[i];
643  lag[i]= a.lag[i];
644  slag[i]= a.slag[i];
645  volume[i]= a.volume[i];
646  size[i]= a.size[i];
647  bp[i]= a.bp[i];
648  bp[i+ndim]= a.bp[i+ndim];
649  bx[i]= a.bx[i];
650  bx[i+ndim]= a.bx[i+ndim];
651 
652  time[i]= a.time[i];
653  time[i+ndim]= a.time[i+ndim];
654  right[i]= a.right[i];
655  left[i]= a.left[i];
656  duration[i]= a.duration[i];
657  start[i]= a.start[i];
658  stop[i]= a.stop[i];
659 
660  frequency[i]= a.frequency[i];
661  low[i]= a.low[i];
662  high[i]= a.high[i];
663  bandwidth[i]= a.bandwidth[i];
664 
665  hrss[i]= a.hrss[i];
666  hrss[i+ndim]= a.hrss[i+ndim];
667  noise[i]= a.noise[i];
668  null[i]= a.null[i];
669  nill[i]= a.nill[i];
670  Deff[i]= a.Deff[i];
671 
672  snr[i]= a.snr[i];
673  sSNR[i]= a.sSNR[i];
674  xSNR[i]= a.xSNR[i];
675  iSNR[i]= a.iSNR[i];
676  oSNR[i]= a.oSNR[i];
677  ioSNR[i]= a.ioSNR[i];
678 
679  }
680  return *this;
681 }
682 
683 
684 //++++++++++++++++++++++++++++++++++++++++++++++++++++
685 // output network events for 2G analysis
686 //++++++++++++++++++++++++++++++++++++++++++++++++++++
687 void netevent::output2G(TTree* waveTree, network* net, size_t ID, int LAG, double factor)
688 {
689 // output event ID in lag LAG to root file
690  if(!net || !net->nLag) exit(1);
691 
692  int i,j,k,ind,kid;
693  int n,m,K,M,injID;
694  detector* pd;
695  netcluster* p;
696  int nIFO = (int)net->ifoListSize(); // number of detectors
697  wavecomplex Aa, gC;
698  double Pi = 3.14159265358979312;
699  this->factor = fabs(factor); // used to tag events in root file
700  factor = factor<=0 ? 1 : fabs(factor); // used to rescale injected events
701  double inRate = net->getifo(0)->rate; // get original rate
702  bool pat0 = net->pattern==0 ? true : false; // pattern flag
703 
704  if(!nIFO) return;
705 
706 // injections
707 
708  injection INJ(nIFO);
709  size_t mdcID, id;
710  double T, injTime, TAU, pcc;
711  double x, ndmLike, Etot;
712  skymap* psm;
713  netcluster* pwc;
714  std::vector<float>* vP;
715  std::vector<int>* vI;
716 
717  bool ellips = net->tYPe=='i' || net->tYPe=='I' ||
718  net->tYPe=='g' || net->tYPe=='G' ||
719  net->tYPe=='s' || net->tYPe=='S' ||
720  net->tYPe=='r' || net->tYPe=='R';
721 
722  bool burst = net->tYPe=='b' || net->tYPe=='B';
723 
724  // arrays for cluster parameters
725 
726  wavearray<double> clusterID_net;
727  wavearray<double> vol0_net;
728  wavearray<double> vol1_net;
729  wavearray<double> size_net;
730  wavearray<double> start_net;
731  wavearray<double> stop_net;
732  wavearray<double> low_net;
733  wavearray<double> high_net;
734  wavearray<double> noise_net;
735  wavearray<double> NOISE_net;
736  wavearray<double> cFreq_net;
737  wavearray<double> rate_net;
738  wavearray<double> duration_net;
739  wavearray<double> bandwidth_net;
740 
741  // why do we need all this?
742 
743  this->ifoList = net->ifoList; // add detectors to tree user info
744  if(waveTree!=NULL) {
745  int dsize = waveTree->GetUserInfo()->GetSize();
746  if(dsize!=0 && dsize!=nIFO) {
747  cout<<"netevent::output2G(): wrong user detector list in header tree"<<endl; exit(1);
748  }
749  if(dsize==0) {
750  for(int n=0;n<nIFO;n++) {
751  // must be done a copy because detector object
752  // is destroyed when waveTree is destroyed
753  detector* pD = new detector(*net->getifo(n));
754  waveTree->GetUserInfo()->Add(pD);
755  }
756  }
757  }
758 
759  this->run = net->nRun;
760  this->nevent = 0;
761 
762  pwc = net->getwc(LAG); // pointer to netcluster
763  if(!pwc->size()) return;
764 
765  clusterID_net = pwc->get((char*)"ID",0,'S',0);
766  K = clusterID_net.size();
767  kid = -1;
768  for(k=0; k<K; k++) { // find cluster
769  id = size_t(clusterID_net.data[k]+0.1);
770  if(ID&&(ID==id)) {kid=k; break;}
771  }
772  if(kid<0) return;
773 
774 // read cluster parameters
775 
776  start_net.resize(0);
777  stop_net.resize(0);
778  noise_net.resize(0);
779  NOISE_net.resize(0);
780 
781  rate_net = pwc->get((char*)"rate",0,'R',0);
782  vol0_net = pwc->get((char*)"volume",0,'R',0); // stored in volume[0]
783  vol1_net = pwc->get((char*)"VOLUME",0,'R',0); // stored in volume[1]
784  size_net = pwc->get((char*)"size",0,'R',0); // stored in size[0]
785  low_net = pwc->get((char*)"low",0,'R',0);
786  high_net = pwc->get((char*)"high",0,'R',0);
787  cFreq_net = pwc->get((char*)"freq",0,'L',0,false);
788  duration_net = pwc->get((char*)"duration",0,'L',0,false);
789  bandwidth_net = pwc->get((char*)"bandwidth",0,'L',0,false);
790 
791  for(i=1; i<=int(nIFO); i++) { // loop on detectors
792  start_net.append(pwc->get((char*)"start",i,'L',0));
793  stop_net.append(pwc->get((char*)"stop",i,'L',0));
794  noise_net.append(pwc->get((char*)"noise",i,'S',0));
795  NOISE_net.append(pwc->get((char*)"NOISE",i,'S',0));
796  }
797 
798  psm = &(net->getifo(0)->tau);
799  vI = &(net->wc_List[LAG].p_Ind[ID-1]);
800  ind = (*vI)[0]; // reconstructed sky index
801 
802  for(i=0; i<int(nIFO); i++) // store gps time of data intervals
803  this->gps[i] = pwc->start+(slag[i]-slag[0]);
804 
805  clusterdata* pcd = &(pwc->cData[ID-1]);
806 
807  this->ndim = nIFO;
808  psm->gps = pcd->cTime + this->gps[0];
809  this->ecor = pcd->netecor;
810  this->nevent += 1;
811  this->eventID[0] = ID;
812  this->eventID[1] = 0;
813  //this->iota[0] = pcd->iota; // ellipticity : 2*cos(iota)/(1+cos(iota)^2)
814  this->iota[0] = pcd->skyChi2; // sky chi2 for all resolutions (temporary)
815  //this->psi[0] = pcd->psi;
816  this->psi[0] = pcd->Gnoise; // estimated Gaussian noise (temporary)
817  this->phi[0] = psm->getPhi(ind); // reconstructed phi
818  this->phi[2] = psm->getRA(ind);
819  this->phi[3] = pcd->phi; // detection phi
820  this->theta[0] = psm->getTheta(ind); // reconstructed theta
821  this->theta[2] = psm->getDEC(ind);
822  this->theta[3] = pcd->theta; // detection theta
823  this->gnet = pcd->gNET;
824  this->anet = pcd->aNET;
825  this->inet = pcd->iNET; // network index
826  this->norm = pcd->norm; // packet norm
827  this->likelihood = pcd->likenet; // total likelihood in sky loop
828  this->volume[0] = int(vol0_net.data[kid]+0.5); // event volume
829  this->volume[1] = int(vol1_net.data[kid]+0.5); // selected pixels volume
830  this->size[0] = int(size_net.data[kid]+0.5); // event size
831  this->size[1] = pcd->skySize; // signal size
832  this->chirp[1] = pcd->mchirp; // reconstructed chirp mass
833  this->chirp[2] = pcd->mchirperr; // reconstructed chirp mass error
834  this->chirp[3] = pcd->chirpEllip; // ellipticity parameter
835  this->chirp[4] = pcd->chirpPfrac; // pixel fraction
836  this->chirp[5] = pcd->chirpEfrac; // energy fraction
837  //this->chirp[6] = pcd->chi2chirp; // chirp chi2/DoF
838  //this->chirp[7] = pcd->tmrgr; // t merger parameter
839  //this->chirp[8] = pcd->tmrgrerr; // t merger error
840 
841  this->range[0] = 0.; // reconstructed distance to source
842 
843  TAU = psm->get(this->theta[0],this->phi[0]);
844 
845  M = 0; gC = 0.;
846  this->strain[0] = 0.;
847  this->penalty = 0.;
848 
849  for(i=0; i<5; i++) { this->neted[i] = 0.; } // zero neted array
850 
851  this->lag[nIFO] = pwc->shift;
852 
853  net->getMRAwave(ID,LAG,'s',net->optim?1:0); // get signal strain
854 
855  for(i=0; i<int(nIFO); i++) { // loop on detectors
856  pd = net->getifo(i);
857  Aa = pd->antenna(this->theta[0],this->phi[0],this->psi[0]); // antenna pattern
858  m = i*K+kid;
859 
860  this->type[0] = 1;
861  this->rate[i] = net->optim ? int(rate_net.data[kid]+0.1) : 0;
862 
863  this->gap[i] = 0.;
864  this->lag[i] = pd->lagShift.data[LAG];
865 
866  this->snr[i] = pd->enrg;
867  this->nill[i] = pd->xSNR-pd->sSNR;
868  this->null[i] = pd->null; // null per detector
869  this->xSNR[i] = pd->xSNR;
870  this->sSNR[i] = pd->sSNR;
871  this->time[i] = pcd->cTime + this->gps[i];
872 
873  if(i) { // take delays into account
874  psm = &(net->getifo(i)->tau);
875  this->time[i] += psm->get(this->theta[0],this->phi[0])-TAU;
876  }
877 
878  this->left[i] = start_net.data[m];
879  this->right[i] = pwc->stop-pwc->start-stop_net.data[m];
880  this->duration[i] = stop_net.data[m] - start_net.data[m];
881  this->start[i] = start_net.data[m] + this->gps[i];
882  this->stop[i] = stop_net.data[m] + this->gps[i];
883 
884  // take lag shift into account
885  double xstart = this->gps[i]+net->Edge; // start data
886  double xstop = this->gps[i]+pwc->stop-pwc->start-net->Edge; // end data
887  this->time[i] += lag[i];
888  if(this->time[i]>xstop) this->time[i] = xstart+(this->time[i]-xstop); // circular buffer
889 
890  this->frequency[i] = cFreq_net.data[k];
891  this->low[i] = low_net.data[k];
892  this->high[i] = high_net.data[k];
893  this->bandwidth[i] = high_net.data[k] - low_net.data[k];
894 
895  this->hrss[i] = sqrt(pd->get_SS()/inRate);
896  this->noise[i] = pow(10.,noise_net.data[m])/sqrt(inRate);
897  this->bp[i]= Aa.real();
898  this->bx[i]= Aa.imag();
899  this->strain[0] += this->hrss[i]*this->hrss[i];
900 
901  Aa /= pow(10.,NOISE_net.data[m]);
902  gC += Aa*Aa;
903 
904  psm->gps = pcd->cTime+this->gps[0];
905  }
906 
907  // Fix start,stop,duration,left,right when simulation!=0.
908  // Due to circular buffer the events detected on the edges
909  // of the segment could have a wrong start,stop values
910  // The start,stop of det with lag=0 are used to fix the wrong values of the other detectors
911  int mdet=0; for(i=0; i<int(nIFO); i++) if(this->lag[i]==0) mdet=i;
912  for(i=0; i<int(nIFO); i++) {
913  if(this->duration[i]!=this->duration[mdet]) {
914  double xstart = this->gps[i]+net->Edge; // start data
915  double xstop = this->gps[i]+pwc->stop-pwc->start-net->Edge; // end data
916 
917  this->start[i] = this->start[mdet]+this->lag[i]+(slag[i]-slag[mdet]);
918  if(this->start[i]>xstop) this->start[i] = xstart+(this->start[i]-xstop); // circular buffer
919 
920  this->stop[i] = this->stop[mdet]+this->lag[i]+(slag[i]-slag[mdet]);
921  if(this->stop[i]>xstop) this->stop[i] = xstart+(this->stop[i]-xstop); // circular buffer
922 
923  this->duration[i] = this->duration[mdet];
924  this->left[i] = this->start[i]-pwc->start;
925  this->right[i] = pwc->stop-this->stop[i];
926  }
927  }
928  this->duration[0] = duration_net.data[kid];
929  this->bandwidth[0] = bandwidth_net.data[kid];
930  this->frequency[0] = pcd->cFreq;
931 
932  ind = pwc->sArea[ID-1].size();
933  for(i=0; i<11; i++)
934  this->erA[i] = i<ind ? pwc->sArea[ID-1][i] : 0.;
935 
936  this->ECOR = pcd->normcor; // normalized coherent energy
937  this->netcc[0] = pcd->netcc; // MRA or SRA cc statistic
938  this->netcc[1] = pcd->skycc; // all-resolution cc statistic
939  this->netcc[2] = pcd->subnet; // MRA or SRA sub-network statistic
940  this->netcc[3] = pcd->SUBNET; // all-resolution sub-network statistic
941 
942  this->neted[0] = pcd->netED; // network ED
943  this->neted[1] = pcd->netnull; // total null energy with Gaussian bias correction
944  this->neted[2] = pcd->energy; // total event energy
945  this->neted[3] = pcd->likesky; // total likelihood at all resolutions
946  this->neted[4] = pcd->enrgsky; // total energy at all resolutions
947  //this->neted[5] = pcd->Gnoise; // estimated contribution of Gaussian noise
948  //this->neted[6] = pcd->skyChi2; // sky chi2 for all resolutions
949 
950  double chrho = this->chirp[3]*sqrt(this->chirp[5]); // reduction factor for chirp events
951  this->rho[0] = pcd->netRHO; // reduced coherent SNR per detector
952  this->rho[1] = pat0 ?pcd->netrho:pcd->netRHO*chrho; // reduced coherent SNR per detector for chirp events
953 
954  this->strain[0] = sqrt(this->strain[0]);
955  if(!ellips) this->psi[0] = gC.arg()*180/Pi;
956  this->penalty = pcd->netnull/nIFO;
957  this->penalty /= pat0 ? this->size[0] : pcd->nDoF; // cluster chi2/nDoF
958 
959 
960 // set injections if there are any
961 
962  M = net->mdc__IDSize();
963  if(!LAG) { // only for zero lag
964  injTime = 1.e12;
965  injID = -1;
966  for(m=0; m<M; m++) {
967  mdcID = net->getmdc__ID(m);
968  T = fabs(this->time[0] - net->getmdcTime(mdcID));
969  if(T<injTime && INJ.fill_in(net,mdcID)) {
970  injTime = T;
971  injID = mdcID;
972  }
973  // printf("%d %12.3f %12.3f\n",mdcID,net->getmdcTime(mdcID),T);
974  }
975 
976  if(INJ.fill_in(net,injID)) { // set injections
977  this->range[1] = INJ.distance/factor;
978  this->chirp[0] = INJ.mchirp;
979  this->eBBH[1] = INJ.e0;
980  this->eBBH[2] = INJ.rp0;
981  this->eBBH[3] = INJ.redshift;
982  this->strain[1] = INJ.strain*factor;
983  this->type[1] = INJ.type;
984  *this->name = net->getmdcType(this->type[1]-1);
985  this->theta[1] = INJ.theta[0];
986  this->phi[1] = INJ.phi[0];
987  this->psi[1] = INJ.psi[0];
988  this->iota[1] = 2.*INJ.iota[0]/(1+INJ.iota[0]*INJ.iota[0]); // ellipticity
989  this->mass[0] = INJ.mass[0];
990  this->mass[1] = INJ.mass[1];
991 
992  for(i=0; i<6; i++) this->spin[i] = INJ.spin[i];
993 
994 // printf("injection type: %d %12.3f\n",INJ.type,INJ.time[0]);
995 
996  wavearray<double>** pwfINJ = new wavearray<double>*[nIFO];
997  wavearray<double>** pwfREC = new wavearray<double>*[nIFO];
998  pd = net->getifo(0);
999  int idSize = pd->RWFID.size();
1000  int wfIndex=-1;
1001  for (int mm=0; mm<idSize; mm++) if (pd->RWFID[mm]==ID) wfIndex=mm;
1002 
1003  for(j=0; j<nIFO; j++) {
1004  pd = net->getifo(j);
1005  Aa = pd->antenna(this->theta[1],this->phi[1],this->psi[1]); // inj antenna pattern
1006  this->hrss[j+nIFO] = INJ.hrss[j]*factor;
1007  this->bp[j+nIFO] = Aa.real();
1008  this->bx[j+nIFO] = Aa.imag();
1009  this->time[j+nIFO] = INJ.time[j];
1010  this->Deff[j] = INJ.Deff[j]/factor;
1011 
1012  pwfINJ[j] = INJ.pwf[j];
1013  if (pwfINJ[j]==NULL) {
1014  cout << "Error : Injected waveform not saved !!! : detector "
1015  << net->ifoName[j] << endl;
1016  continue;
1017  }
1018  if (wfIndex<0) {
1019  cout << "Error : Reconstructed waveform not saved !!! : ID -> "
1020  << ID << " : detector " << net->ifoName[j] << endl;
1021  continue;
1022  }
1023 
1024  if (wfIndex>=0) pwfREC[j] = pd->RWFP[wfIndex];
1025  double R = pd->TFmap.rate();
1026  //double rFactor = log(2.)/log(pd->rate/R); // rescale waveform
1027  double rFactor = 1.;
1028  rFactor *= factor;
1029  wavearray<double>* wfINJ = pwfINJ[j];
1030  *wfINJ*=rFactor;
1031  wavearray<double>* wfREC = pwfREC[j];
1032 
1033  double bINJ = wfINJ->start();
1034  double eINJ = wfINJ->start()+wfINJ->size()/R;
1035  double bREC = wfREC->start();
1036  double eREC = wfREC->start()+wfREC->size()/R;
1037  //cout.precision(14);
1038  //cout << "bINJ : " << bINJ << " eINJ : " << eINJ << endl;
1039  //cout << "bREC : " << bREC << " eREC : " << eREC << endl;
1040 
1041  int oINJ = bINJ>bREC ? 0 : int((bREC-bINJ)*R+0.5);
1042  int oREC = bINJ<bREC ? 0 : int((bINJ-bREC)*R+0.5);
1043  //cout << "oINJ : " << oINJ << " oREC : " << oREC << endl;
1044 
1045  double startXCOR = bINJ>bREC ? bINJ : bREC;
1046  double endXCOR = eINJ<eREC ? eINJ : eREC;
1047  int sizeXCOR = int((endXCOR-startXCOR)*R+0.5);
1048  //cout << "startXCOR : " << startXCOR << " endXCOR : " << endXCOR << " sizeXCOR :" << sizeXCOR << endl;
1049 
1050  if (sizeXCOR<=0) {*wfINJ*=1./rFactor; continue;}
1051 
1052  // the enINJ, enREC, xcorINJ_REC are computed in the INJ range
1053 
1054  double enINJ=0;
1055  for (int i=0;i<wfINJ->size();i++) enINJ+=wfINJ->data[i]*wfINJ->data[i];
1056  //for (int i=0;i<sizeXCOR;i++) enINJ+=wfINJ->data[i+oINJ]*wfINJ->data[i+oINJ];
1057 
1058  double enREC=0;
1059  for (int i=0;i<wfREC->size();i++) enREC+=wfREC->data[i]*wfREC->data[i];
1060  //for (int i=0;i<sizeXCOR;i++) enREC+=wfREC->data[i+oREC]*wfREC->data[i+oREC];
1061 
1062  double xcorINJ_REC=0;
1063  for (int i=0;i<sizeXCOR;i++) xcorINJ_REC+=wfINJ->data[i+oINJ]*wfREC->data[i+oREC];
1064 
1065  WSeries<double> wfREC_SUB_INJ;
1066  wfREC_SUB_INJ.resize(sizeXCOR);
1067  for (int i=0;i<sizeXCOR;i++) wfREC_SUB_INJ.data[i]=wfREC->data[i+oREC]-wfINJ->data[i+oINJ];
1068  wfREC_SUB_INJ.start(startXCOR);
1069  wfREC_SUB_INJ.rate(wfREC->rate());
1070 
1071  this->iSNR[j] = enINJ;
1072  this->oSNR[j] = enREC;
1073  this->ioSNR[j] = xcorINJ_REC;
1074 
1075  //double erINJ_REC = enINJ+enREC-2*xcorINJ_REC;
1076  //cout << "enINJ : " << enINJ << " enREC : " << enREC << " xcorINJ_REC : " << xcorINJ_REC << endl;
1077  //cout << "erINJ_REC/enINJ : " << erINJ_REC/enINJ << endl;
1078 
1079  *wfINJ*=1./rFactor;
1080  }
1081  delete [] pwfINJ;
1082  delete [] pwfREC;
1083  }
1084  else { // no injections
1085  this->range[1] = 0.;
1086  this->chirp[0] = 0.;
1087  this->eBBH[1] = 0.;
1088  this->eBBH[2] = 0.;
1089  this->eBBH[3] = 0.;
1090  this->strain[1] = 0.;
1091  this->type[1] = 0;
1092  this->theta[1] = 0.;
1093  this->phi[1] = 0.;
1094  this->psi[1] = 0.;
1095  this->iota[1] = 0.;
1096  this->mass[0] = 0.;
1097  this->mass[1] = 0.;
1098 
1099  for(i=0; i<6; i++) this->spin[i] = 0.;
1100 
1101  for(j=0; j<nIFO; j++) {
1102  this->hrss[j+nIFO] = 0.;
1103  this->bp[j+nIFO] = 0.;
1104  this->bx[j+nIFO] = 0.;
1105  this->time[j+nIFO] = 0.;
1106  this->Deff[j] = 0.;
1107  this->iSNR[j] = 0.;
1108  this->oSNR[j] = 0.;
1109  this->ioSNR[j] = 0.;
1110  }
1111  }
1112  }
1113 
1114  if(this->fP!=NULL) {
1115  fprintf(fP,"\n# trigger %d in lag %d for \n",int(ID),int(LAG));
1116  this->Dump("2G");
1117  vP = &(net->wc_List[LAG].p_Map[ID-1]);
1118  vI = &(net->wc_List[LAG].p_Ind[ID-1]);
1119  x = cos(psm->theta_1*PI/180.)-cos(psm->theta_2*PI/180.);
1120  x*= (psm->phi_2-psm->phi_1)*180/PI/psm->size();
1121  fprintf(fP,"sky_res: %f\n",sqrt(fabs(x)));
1122  fprintf(fP,"map_lenght: %d\n",int(vP->size()));
1123  fprintf(fP,"#skyID theta DEC step phi R.A step probability cumulative\n");
1124  x = 0;
1125  for(j=0; j<int(vP->size()); j++) {
1126  i = (*vI)[j];
1127  if(net->mdc__IDSize()) { // simulation mode
1128  x = (j==int(vP->size())-1) ? 0 : x+(*vP)[j]; // last value is the inj sky index (x=0)
1129  } else x+=(*vP)[j];
1130  fprintf(fP,"%6d %5.1f %5.1f %6.2f %5.1f %5.1f %6.2f %e %e\n",
1131  int(i),psm->getTheta(i),psm->getDEC(i),psm->getThetaStep(i),
1132  psm->getPhi(i),psm->getRA(i),psm->getPhiStep(i),(*vP)[j],x);
1133  }
1134  }
1135 
1136  // save the probability skymap to tree
1137  if(waveTree!=NULL && net->wc_List[LAG].p_Map.size()) {
1138 
1139  vP = &(net->wc_List[LAG].p_Map[ID-1]);
1140  vI = &(net->wc_List[LAG].p_Ind[ID-1]);
1141  if(this->Psave) {
1142 
1143  *Psm = *psm;
1144  *Psm = 0.;
1145 
1146  int k;
1147  double th,ph;
1148  for(j=0; j<int(vP->size()); j++) {
1149  i = (*vI)[j];
1150  th = Psm->getTheta(i);
1151  ph = Psm->getPhi(i);
1152  k=Psm->getSkyIndex(th, ph);
1153  Psm->set(k,(*vP)[j]);
1154  }
1155  }
1156  }
1157 
1158  if(waveTree!=NULL) waveTree->Fill();
1159 
1160 }
1161 
1162 
1163 //++++++++++++++++++++++++++++++++++++++++++++++++++++
1164 // output network events
1165 //++++++++++++++++++++++++++++++++++++++++++++++++++++
1166 void netevent::output(TTree* waveTree, network* net, double factor, size_t iID, int LAG)
1167 {
1168  if(!net || !net->nLag) exit(1);
1169 
1170  int i,j,k,ind;
1171  int n,m,K,M,injID;
1172  detector* pd;
1173  netcluster* p;
1174  int N = (int)net->ifoListSize(); // number of detectors
1175  int bLag = LAG<0 ? 0 : LAG;
1176  int eLag = LAG<0 ? net->nLag : LAG+1;
1177  wavecomplex Aa, gC;
1178  double Pi = 3.14159265358979312;
1179  this->factor = fabs(factor); // used to tag events in root file
1180  factor = factor<=0 ? 1 : fabs(factor); // used to rescale injected events
1181  double inRate = net->getifo(0)->rate; // get original rate
1182 
1183  if(!N) return;
1184  double* ndm = (double*)malloc(N*N*sizeof(double));
1185  int iTYPE = net->MRA ? 0 : 1;
1186 
1187 // injections
1188 
1189  injection INJ(N);
1190  size_t mdcID, ID;
1191  double T, injTime, TAU, pcc;
1192  double x, tot_null, ndmLike, Etot;
1193  skymap* psm;
1194  std::vector<float>* vP;
1195  std::vector<int>* vI;
1196 
1197  bool ellips = net->tYPe=='i' || net->tYPe=='I' ||
1198  net->tYPe=='g' || net->tYPe=='G' ||
1199  net->tYPe=='s' || net->tYPe=='S' ||
1200  net->tYPe=='r' || net->tYPe=='R';
1201 
1202  bool burst = net->tYPe=='b' || net->tYPe=='B';
1203 
1204  wavearray<double> skSNR(N);
1205  wavearray<double> xkSNR(N);
1206 
1207 // arrays for cluster parameters
1208 
1209  wavearray<double> clusterID_net;
1210  wavearray<double> volume_net;
1211  wavearray<double> size_net;
1212  wavearray<double> start_net;
1213  wavearray<double> stop_net;
1214  wavearray<double> time_net;
1215  wavearray<double> frequency_net;
1216  wavearray<double> low_net;
1217  wavearray<double> high_net;
1218  wavearray<double> LH_net;
1219  wavearray<double> null_net;
1220  wavearray<double> rSNR_net;
1221  wavearray<double> gSNR_net;
1222  wavearray<double> gSNR_NET;
1223  wavearray<double> rate_net;
1224  wavearray<double> hrss_net;
1225  wavearray<double> hrss_NET;
1226  wavearray<double> noise_net;
1227  wavearray<double> NOISE_net;
1228  wavearray<double> psi_net;
1229  wavearray<double> ell_net;
1230  wavearray<double> phi_net;
1231  wavearray<double> theta_net;
1232  wavearray<double> energy_net;
1233  wavearray<double> energy_NET;
1234 
1235  this->ifoList = net->ifoList;
1236  // add detectors to tree user info
1237  if(waveTree!=NULL) {
1238  for(int n=0;n<N;n++) {
1239  // must be done a copy because detector object
1240  // is destroyed when waveTree is destroyed
1241  detector* pD = new detector(*net->getifo(n));
1242  waveTree->GetUserInfo()->Add(pD);
1243  }
1244  }
1245 
1246  this->run = net->nRun;
1247  this->nevent = 0;
1248 
1249  for(n=bLag; n<eLag; n++){ // loop on time lags
1250 
1251  p = net->getwc(n); // pointer to netcluster
1252  if(!p->size()) continue;
1253 
1254  clusterID_net = p->get((char*)"ID",0,'S',iTYPE);
1255  K = clusterID_net.size();
1256 
1257  if(!K) continue;
1258 
1259 // read cluster parameters
1260 
1261  time_net.resize(0);
1262  start_net.resize(0);
1263  stop_net.resize(0);
1264  frequency_net.resize(0);
1265  energy_net.resize(0);
1266  energy_NET.resize(0);
1267  rSNR_net.resize(0);
1268  gSNR_net.resize(0);
1269  gSNR_NET.resize(0);
1270  hrss_net.resize(0);
1271  hrss_NET.resize(0);
1272  noise_net.resize(0);
1273  NOISE_net.resize(0);
1274  null_net.resize(0);
1275 
1276  LH_net = p->get((char*)"likelihood",0,'R',iTYPE);
1277  rate_net = p->get((char*)"rate",0,'R',iTYPE);
1278  ell_net = p->get((char*)"ellipticity",0,'R',iTYPE);
1279  psi_net = p->get((char*)"psi",0,'R',iTYPE);
1280  phi_net = p->get((char*)"phi",0,'R',iTYPE);
1281  theta_net = p->get((char*)"theta",0,'R',iTYPE);
1282  size_net = p->get((char*)"size",0,'R',iTYPE);
1283  volume_net = p->get((char*)"volume",0,'R',iTYPE);
1284  low_net = p->get((char*)"low",0,'R',iTYPE);
1285  high_net = p->get((char*)"high",0,'R',iTYPE);
1286 
1287  for(i=1; i<int(N+1); i++) // loop on detectors
1288  { // read cluster parameters
1289  time_net.append(p->get((char*)"time",i,'L',0));
1290  start_net.append(p->get((char*)"start",i,'L',iTYPE));
1291  stop_net.append(p->get((char*)"stop",i,'L',iTYPE));
1292  frequency_net.append(p->get((char*)"FREQUENCY",i,'L',0));
1293  rSNR_net.append(p->get((char*)"SNR",i,'R',iTYPE));
1294  gSNR_net.append(p->get((char*)"SNR",i,'S',iTYPE));
1295  gSNR_NET.append(p->get((char*)"SNR",i,'P',iTYPE));
1296  hrss_net.append(p->get((char*)"hrss",i,'W',iTYPE));
1297  hrss_NET.append(p->get((char*)"hrss",i,'U',iTYPE));
1298  noise_net.append(p->get((char*)"noise",i,'S',iTYPE));
1299  NOISE_net.append(p->get((char*)"NOISE",i,'S',iTYPE));
1300  energy_net.append(p->get((char*)"energy",i,'S',iTYPE));
1301  energy_NET.append(p->get((char*)"energy",i,'P',iTYPE));
1302  null_net.append(p->get((char*)"null",i,'W',iTYPE));
1303  }
1304 
1305  if(ellips) {
1306  energy_net += energy_NET;
1307  energy_net *= net->MRA ? 1.0 : 0.5;
1308  }
1309 
1310  for(k=0; k<K; k++) // loop on events
1311  {
1312  ID = size_t(clusterID_net.data[k]+0.1);
1313  if((iID)&&(ID!=iID)) continue;
1314  if(ellips) {
1315  net->SETNDM(ID,n,true,net->MRA?0:1);
1316  }
1317  else if(burst)
1318  net->setndm(ID,n,true,1);
1319 // else
1320 // cout<<"netevent::output(): incorrect search option"<<endl; //exit(1);}
1321 
1322  psm = &(net->getifo(0)->tau);
1323  vI = &(net->wc_List[n].p_Ind[ID-1]);
1324  ind = (*vI)[0]; // reconstructed sky index
1325 
1326  this->ndim = N;
1327  for(i=0; i<N; i++) // loop on detectors
1328  this->gps[i] = net->getifo(i)->getTFmap()->start()+(slag[i]-slag[0]);
1329  psm->gps = time_net.data[k]+this->gps[0];
1330  this->ecor = net->eCOR;
1331  this->nevent += 1;
1332  this->eventID[0] = ID;
1333  this->eventID[1] = n;
1334  this->psi[0] = psi_net.data[k];
1335  this->phi[0] = psm->getPhi(ind); // reconstructed phi
1336  this->phi[2] = psm->getRA(ind);
1337  this->phi[3] = phi_net.data[k]; // detection phi
1338  this->theta[0] = psm->getTheta(ind); // reconstructed theta
1339  this->theta[2] = psm->getDEC(ind);
1340  this->theta[3] = theta_net.data[k]; // detection theta
1341  this->gnet = net->gNET;
1342  this->anet = net->aNET;
1343  this->inet = net->iNET;
1344  this->iota[0] = net->norm;
1345  this->norm = net->norm;
1346  this->likelihood = 0.;
1347 
1348  TAU = psm->get(this->theta[0],this->phi[0]);
1349 
1350  M = 0; gC = 0.;
1351  this->strain[0] = 0.;
1352  this->penalty = 0.;
1353  tot_null = 0.;
1354 
1355  Etot = 0.;
1356  for(i=0; i<N; i++) Etot += energy_net.data[i*K+k];
1357 
1358  for(i=0; i<5; i++) { this->neted[i] = 0.; }
1359 
1360  this->lag[N] = p->shift;
1361 
1362  for(i=0; i<N; i++) // loop on detectors
1363  {
1364  pd = net->getifo(i);
1365  Aa = pd->antenna(this->theta[0],this->phi[0]); // antenna pattern
1366  m = i*K+k;
1367 
1368  if(i) { // take delays into account
1369  psm = &(net->getifo(i)->tau);
1370  time_net.data[m] += psm->get(this->theta[0],this->phi[0])-TAU;
1371  }
1372 
1373  this->type[0] = 1;
1374  this->rate[i] = int(rate_net.data[k]+0.1);
1375 
1376  this->gap[i] = 0.;
1377  this->lag[i] = pd->lagShift.data[n];
1378 
1379  this->volume[i] = int(volume_net.data[k]+0.5);
1380  this->size[i] = int(size_net.data[k]+0.5);
1381 
1382  this->snr[i] = energy_net.data[m];
1383  this->nill[i] = pd->xSNR-pd->sSNR;
1384  this->null[i] = pd->null;
1385  this->xSNR[i] = pd->xSNR;
1386  this->sSNR[i] = pd->sSNR;
1387  this->neted[0] += fabs(pd->xSNR-pd->sSNR);
1388  this->neted[1] += fabs(pd->ED[1]);
1389  this->neted[2] += fabs(pd->ED[2]);
1390  this->neted[3] += fabs(pd->ED[3]);
1391  this->neted[4] += fabs(pd->ED[4]);
1392  this->likelihood += snr[i] - pd->null;
1393 
1394  this->time[i] = time_net.data[m] + this->gps[i];
1395  this->left[i] = start_net.data[m];
1396  this->right[i] = p->stop-p->start-stop_net.data[m];
1397  this->duration[i] = stop_net.data[m] - start_net.data[m];
1398  this->start[i] = start_net.data[m] + this->gps[i];
1399  this->stop[i] = stop_net.data[m] + this->gps[i];
1400 
1401  this->frequency[i] = frequency_net.data[m];
1402  this->low[i] = low_net.data[k];
1403  this->high[i] = high_net.data[k];
1404  this->bandwidth[i] = high_net.data[k] - low_net.data[k];
1405 
1406  if(net->MRA) {
1407  this->hrss[i] = pow(pow(10.,hrss_net.data[m]),2)+pow(pow(10.,hrss_NET.data[m]),2);
1408  this->hrss[i] = sqrt(this->hrss[i]/inRate);
1409  } else {
1410  this->hrss[i] = pow(10.,hrss_net.data[m])/sqrt(inRate);
1411  }
1412  this->noise[i] = pow(10.,noise_net.data[m])/sqrt(inRate);
1413  this->bp[i]= Aa.real();
1414  this->bx[i]= Aa.imag();
1415  this->strain[0]+= this->hrss[i]*this->hrss[i];
1416 
1417  tot_null += this->null[i];
1418  Aa /= pow(10.,NOISE_net.data[m]);
1419  gC += Aa*Aa;
1420 
1421  skSNR.data[i] = pd->sSNR;
1422  xkSNR.data[i] = pd->xSNR;
1423 
1424  x = 1. - pd->sSNR/(this->snr[i]+net->precision*Etot);
1425  if(x<this->penalty) this->penalty = x;
1426 
1427  psm->gps = time_net.data[m]+this->gps[0];
1428 
1429  }
1430 
1431  ndmLike = 0.;
1432  this->ECOR = 0.;
1433  this->netcc[2] = 0.;
1434  this->netcc[3] = 0.;
1435 
1436  for(i=0; i<N; i++) { // loop on detectors
1437  for(j=i; j<int(N); j++) {
1438  ndm[M] = net->getNDM(i,j);
1439 
1440  if(i!=j) {
1441  ndm[M] *= 2.;
1442  pcc = 2*sqrt(net->getNDM(i,i)*net->getNDM(j,j));
1443  pcc = pcc>0. ? ndm[M]/pcc : 0.;
1444  this->ECOR += ndm[M]*fabs(pcc); // effective correlated energy
1445  x = skSNR.data[i] + skSNR.data[j];
1446  this->netcc[2] += x*pcc;
1447  x = xkSNR.data[i] + xkSNR.data[j];
1448  this->netcc[3] += x*pcc;
1449  }
1450 
1451  ndmLike += ndm[M++];
1452  }
1453  }
1454 
1455  if(ndmLike>0.) { // temporary 2G condition
1456  ind = p->sArea[ID-1].size();
1457  for(i=0; i<11; i++)
1458  this->erA[i] = i<ind ? p->sArea[ID-1][i] : 0.;
1459 
1460  x = this->ecor;
1461  this->netcc[0] = x/(tot_null+x); // network (ecor) correlation
1462 
1463  x = this->ECOR;
1464  this->netcc[1] = x/(tot_null+x); // network (ECOR) correlation
1465 
1466  this->netcc[2] /= (N-1)*ndmLike; // network Pearson's correlation
1467  this->netcc[3] /= (N-1)*ndmLike; // network Pearson's correlation
1468  }
1469 
1470  this->rho[0] = sqrt(ecor*this->netcc[0]/N);
1471  this->rho[1] = sqrt(ECOR*this->netcc[0]/N);
1472 
1473  this->strain[0] = sqrt(this->strain[0]);
1474  this->penalty = sqrt(1./(1-this->penalty));
1475  if(!ellips) this->psi[0] = gC.arg()*180/Pi;
1476 
1477 // set injections if there are any
1478 
1479  M = net->mdc__IDSize();
1480  if(!n) { // only for zero lag
1481  injTime = 1.e12;
1482  injID = -1;
1483  for(m=0; m<M; m++) {
1484  mdcID = net->getmdc__ID(m);
1485  T = fabs(this->time[0] - net->getmdcTime(mdcID));
1486  if(T<injTime && INJ.fill_in(net,mdcID)) {
1487  injTime = T;
1488  injID = mdcID;
1489  }
1490 // printf("%d %12.3f %12.3f\n",mdcID,net->getmdcTime(mdcID),T);
1491  }
1492 
1493  if(INJ.fill_in(net,injID)) { // set injections
1494  this->range[1] = INJ.distance/factor;
1495  this->chirp[0] = INJ.mchirp;
1496  this->eBBH[1] = INJ.e0;
1497  this->eBBH[2] = INJ.rp0;
1498  this->eBBH[3] = INJ.redshift;
1499  this->strain[1] = INJ.strain*factor;
1500  this->type[1] = INJ.type;
1501  *this->name = net->getmdcType(this->type[1]-1);
1502  this->theta[1] = INJ.theta[0];
1503  this->phi[1] = INJ.phi[0];
1504  this->psi[1] = INJ.psi[0];
1505  this->iota[1] = INJ.iota[0];
1506  this->mass[0] = INJ.mass[0];
1507  this->mass[1] = INJ.mass[1];
1508 
1509  for(i=0; i<6; i++) this->spin[i] = INJ.spin[i];
1510 
1511 // printf("injection type: %d %12.3f\n",INJ.type,INJ.time[0]);
1512 
1513  wavearray<double>** pwfINJ = new wavearray<double>*[N];
1514  wavearray<double>** pwfREC = new wavearray<double>*[N];
1515  pd = net->getifo(0);
1516  int idSize = pd->RWFID.size();
1517  int wfIndex=-1;
1518  for (int mm=0; mm<idSize; mm++) if (pd->RWFID[mm]==ID) wfIndex=mm;
1519 
1520  for(j=0; j<int(N); j++) {
1521  pd = net->getifo(j);
1522  Aa = pd->antenna(this->theta[1],this->phi[1]); // inj antenna pattern
1523  this->hrss[j+N] = INJ.hrss[j]*factor;
1524  this->bp[j+N] = Aa.real();
1525  this->bx[j+N] = Aa.imag();
1526  this->time[j+N] = INJ.time[j];
1527  this->Deff[j] = INJ.Deff[j]/factor;
1528 
1529  pwfINJ[j] = INJ.pwf[j];
1530  if (pwfINJ[j]==NULL) {
1531  cout << "Error : Injected waveform not saved !!! : detector "
1532  << net->ifoName[j] << endl;
1533  continue;
1534  }
1535  if (wfIndex<0) {
1536  cout << "Error : Reconstructed waveform not saved !!! : ID -> "
1537  << ID << " : detector " << net->ifoName[j] << endl;
1538  continue;
1539  }
1540 
1541  if (wfIndex>=0) pwfREC[j] = pd->RWFP[wfIndex];
1542  double R = pd->TFmap.rate();
1543  //double rFactor = log(2.)/log(pd->rate/R); // rescale waveform
1544  double rFactor = 1.;
1545  rFactor *= factor;
1546  wavearray<double>* wfINJ = pwfINJ[j];
1547  *wfINJ*=rFactor;
1548  wavearray<double>* wfREC = pwfREC[j];
1549 
1550  double bINJ = wfINJ->start();
1551  double eINJ = wfINJ->start()+wfINJ->size()/R;
1552  double bREC = wfREC->start();
1553  double eREC = wfREC->start()+wfREC->size()/R;
1554  //cout.precision(14);
1555  //cout << "bINJ : " << bINJ << " eINJ : " << eINJ << endl;
1556  //cout << "bREC : " << bREC << " eREC : " << eREC << endl;
1557 
1558  int oINJ = bINJ>bREC ? 0 : int((bREC-bINJ)*R+0.5);
1559  int oREC = bINJ<bREC ? 0 : int((bINJ-bREC)*R+0.5);
1560  //cout << "oINJ : " << oINJ << " oREC : " << oREC << endl;
1561 
1562  double startXCOR = bINJ>bREC ? bINJ : bREC;
1563  double endXCOR = eINJ<eREC ? eINJ : eREC;
1564  int sizeXCOR = int((endXCOR-startXCOR)*R+0.5);
1565  //cout << "startXCOR : " << startXCOR << " endXCOR : " << endXCOR << " sizeXCOR :" << sizeXCOR << endl;
1566 
1567  if (sizeXCOR<=0) continue;
1568 
1569  // the enINJ, enREC, xcorINJ_REC are computed in the INJ range
1570 
1571  double enINJ=0;
1572  for (int i=0;i<wfINJ->size();i++) enINJ+=wfINJ->data[i]*wfINJ->data[i];
1573  //for (int i=0;i<sizeXCOR;i++) enINJ+=wfINJ->data[i+oINJ]*wfINJ->data[i+oINJ];
1574 
1575  double enREC=0;
1576  //for (int i=0;i<wfREC->size();i++) enREC+=wfREC->data[i]*wfREC->data[i];
1577  for (int i=0;i<sizeXCOR;i++) enREC+=wfREC->data[i+oREC]*wfREC->data[i+oREC];
1578 
1579  double xcorINJ_REC=0;
1580  for (int i=0;i<sizeXCOR;i++) xcorINJ_REC+=wfINJ->data[i+oINJ]*wfREC->data[i+oREC];
1581 
1582  WSeries<double> wfREC_SUB_INJ;
1583  wfREC_SUB_INJ.resize(sizeXCOR);
1584  for (int i=0;i<sizeXCOR;i++) wfREC_SUB_INJ.data[i]=wfREC->data[i+oREC]-wfINJ->data[i+oINJ];
1585  wfREC_SUB_INJ.start(startXCOR);
1586  wfREC_SUB_INJ.rate(wfREC->rate());
1587  this->iSNR[j] = enINJ;
1588  this->oSNR[j] = enREC;
1589  this->ioSNR[j] = xcorINJ_REC;
1590 
1591  //double erINJ_REC = enINJ+enREC-2*xcorINJ_REC;
1592  //cout << "enINJ : " << enINJ << " enREC : " << enREC << " xcorINJ_REC : " << xcorINJ_REC << endl;
1593  //cout << "erINJ_REC/enINJ : " << erINJ_REC/enINJ << endl;
1594 
1595  *wfINJ*=1./rFactor;
1596  }
1597  delete [] pwfINJ;
1598  delete [] pwfREC;
1599  }
1600  else { // no injections
1601  this->range[1] = 0.;
1602  this->chirp[0] = 0.;
1603  this->eBBH[1] = 0.;
1604  this->eBBH[2] = 0.;
1605  this->eBBH[3] = 0.;
1606  this->strain[1] = 0.;
1607  this->type[1] = 0;
1608  this->theta[1] = 0.;
1609  this->phi[1] = 0.;
1610  this->psi[1] = 0.;
1611  this->iota[1] = 0.;
1612  this->mass[0] = 0.;
1613  this->mass[1] = 0.;
1614 
1615  for(i=0; i<6; i++) this->spin[i] = 0.;
1616 
1617  for(j=0; j<int(N); j++) {
1618  this->hrss[j+N] = 0.;
1619  this->bp[j+N] = 0.;
1620  this->bx[j+N] = 0.;
1621  this->time[j+N] = 0.;
1622  this->Deff[j] = 0.;
1623  this->iSNR[j] = 0.;
1624  this->oSNR[j] = 0.;
1625  this->ioSNR[j] = 0.;
1626  }
1627  }
1628  }
1629 
1630  if(this->fP!=NULL) {
1631  fprintf(fP,"\n# trigger %d in lag %d for \n",int(ID),int(n));
1632  this->Dump("1G");
1633  vP = &(net->wc_List[n].p_Map[ID-1]);
1634  vI = &(net->wc_List[n].p_Ind[ID-1]);
1635  x = cos(psm->theta_1*PI/180.)-cos(psm->theta_2*PI/180.);
1636  x*= (psm->phi_2-psm->phi_1)*180/PI/psm->size();
1637  fprintf(fP,"sky_res: %f\n",sqrt(fabs(x)));
1638  fprintf(fP,"map_lenght: %d\n",int(vP->size()));
1639  fprintf(fP,"#skyID theta DEC step phi R.A step probability cumulative\n");
1640  x = 0.;
1641  for(j=0; j<int(vP->size()); j++) {
1642  i = (*vI)[j];
1643  if(net->mdc__IDSize()) { // simulation mode
1644  x = (j==int(vP->size())-1) ? 0 : x+(*vP)[j]; // last value is the inj sky index (x=0)
1645  } else x+=(*vP)[j];
1646  fprintf(fP,"%6d %5.1f %5.1f %6.2f %5.1f %5.1f %6.2f %e %e\n",
1647  int(i),psm->getTheta(i),psm->getDEC(i),psm->getThetaStep(i),
1648  psm->getPhi(i),psm->getRA(i),psm->getPhiStep(i),(*vP)[j],x);
1649  }
1650  }
1651 
1652  // save the probability skymap to tree
1653  if(waveTree!=NULL && net->wc_List[n].p_Map.size()) {
1654 
1655  vP = &(net->wc_List[n].p_Map[ID-1]);
1656  vI = &(net->wc_List[n].p_Ind[ID-1]);
1657  if(this->Psave) {
1658 
1659  *Psm = *psm;
1660  *Psm = 0.;
1661 
1662  int k;
1663  double th,ra;
1664  for(j=0; j<int(vP->size()); j++) {
1665  i = (*vI)[j];
1666  th = Psm->getTheta(i);
1667  ra = Psm->getRA(i);
1668  k=Psm->getSkyIndex(th, ra);
1669  Psm->set(k,(*vP)[j]);
1670  }
1671  }
1672  }
1673 
1674  if(waveTree!=NULL) waveTree->Fill();
1675 
1676  }
1677  }
1678  if(ndm) free(ndm);
1679 }
1680 
1681 
tree
TTree * tree
Definition: TimeSortTree.C:20
network::ifoName
std::vector< char * > ifoName
Definition: network.hh:591
netevent::anet
Float_t anet
Definition: netevent.hh:99
network::getifo
detector * getifo(size_t n)
param: detector index
Definition: network.hh:418
netevent::b_Psm
TBranch * b_Psm
Definition: netevent.hh:167
netevent::b_high
TBranch * b_high
Definition: netevent.hh:161
netevent::mass
Float_t * mass
effective range for each detector
Definition: netevent.hh:112
wavearray::append
size_t append(const wavearray< DataType_t > &)
Definition: wavearray.cc:775
clusterdata::SUBNET
float SUBNET
Definition: netcluster.hh:48
netcluster::stop
double stop
Definition: netcluster.hh:362
WSeries::resize
virtual void resize(unsigned int)
Definition: wseries.cc:883
netevent::b_ECOR
TBranch * b_ECOR
Definition: netevent.hh:179
wavearray::size
virtual size_t size() const
Definition: wavearray.hh:127
netevent::eBBH
Float_t * eBBH
chirp array: 0-injmass,1-recmass,2-merr,3-tmrgr,4-terr,5-chi2
Definition: netevent.hh:110
injection::time
Double_t * time
beam pattern coefficients for hx
Definition: injection.hh:70
netevent::b_hrss
TBranch * b_hrss
Definition: netevent.hh:164
netevent::ndim
Int_t ndim
current Tree number in a TChain
Definition: netevent.hh:47
injection::rp0
Float_t rp0
Definition: injection.hh:57
wavecomplex::imag
double imag() const
Definition: wavecomplex.hh:52
network::nLag
size_t nLag
Definition: network.hh:555
netevent::rho
Float_t * rho
biased null statistics
Definition: netevent.hh:93
netevent::high
Float_t * high
min frequency
Definition: netevent.hh:85
injection::e0
Float_t e0
Definition: injection.hh:58
Init
void Init()
Definition: ChirpMass.C:280
netevent::range
Float_t * range
Definition: netevent.hh:108
netevent::b_volume
TBranch * b_volume
Definition: netevent.hh:136
network::wc_List
std::vector< netcluster > wc_List
Definition: network.hh:592
netevent::ecor
Float_t ecor
Definition: netevent.hh:101
netevent::setSLags
void setSLags(float *slag)
Definition: netevent.cc:404
netevent::b_range
TBranch * b_range
Definition: netevent.hh:184
netevent::Notify
Bool_t Notify()
Definition: netevent.cc:304
netevent::start
Double_t * start
cluster duration = stopW-startW
Definition: netevent.hh:80
network::getmdcTime
std::vector< double > * getmdcTime()
Definition: network.hh:404
fprintf
fprintf(stdout,"start=%f duration=%f rate=%f\n", x.start(), x.size()/x.rate(), x.rate())
clusterdata::likenet
float likenet
Definition: netcluster.hh:37
clusterdata::gNET
float gNET
Definition: netcluster.hh:59
wavearray::rate
virtual void rate(double r)
Definition: wavearray.hh:123
netevent::duration
Float_t * duration
max cluster time relative to segment start
Definition: netevent.hh:79
netevent::low
Float_t * low
average center_of_snr frequency
Definition: netevent.hh:84
netevent::b_theta
TBranch * b_theta
Definition: netevent.hh:145
detector::sSNR
double sSNR
Definition: detector.hh:320
a
wavearray< double > a(hp.size())
netevent::b_type
TBranch * b_type
Definition: netevent.hh:132
clusterdata::energy
float energy
Definition: netcluster.hh:35
n
int n
Definition: cwb_net.C:10
netevent::setTree
TTree * setTree()
Definition: netevent.cc:412
TString
TString("c")
Definition: cwb_report_skymap.C:111
clusterdata::netED
float netED
Definition: netcluster.hh:41
network::nRun
size_t nRun
Definition: network.hh:554
netevent::b_ndim
TBranch * b_ndim
Definition: netevent.hh:128
clusterdata::chirpEfrac
float chirpEfrac
Definition: netcluster.hh:69
netevent::inet
Float_t inet
Definition: netevent.hh:100
ID
int ID
Definition: TestMDC.C:70
netevent::size
Int_t * size
cluster volume
Definition: netevent.hh:61
skymap::gps
double gps
Definition: skymap.hh:307
netevent::right
Float_t * right
segment start GPS time
Definition: netevent.hh:77
clusterdata::mchirp
float mchirp
Definition: netcluster.hh:66
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<< "\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
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:2188
detector::RWFP
std::vector< wavearray< double > * > RWFP
Definition: detector.hh:364
netevent::b_neted
TBranch * b_neted
Definition: netevent.hh:171
netevent::GetEntry
Int_t GetEntry(Int_t)
Definition: netevent.cc:387
netevent::ECOR
Float_t ECOR
Definition: netevent.hh:103
clusterdata::chirpPfrac
float chirpPfrac
Definition: netcluster.hh:70
netevent::b_duration
TBranch * b_duration
Definition: netevent.hh:155
netevent::gap
Float_t * gap
cluster union size
Definition: netevent.hh:64
clusterdata::normcor
float normcor
Definition: netcluster.hh:39
netcluster::sArea
std::vector< vector_float > sArea
Definition: netcluster.hh:383
pwc
netcluster * pwc
Definition: cwb_job_obj.C:20
netevent::b_rho
TBranch * b_rho
Definition: netevent.hh:172
ph
TH2F * ph
Definition: ReadHistFromJobFile.C:17
network::getmdcType
string getmdcType(size_t n)
Definition: network.hh:402
detector::get_SS
double get_SS()
Definition: detector.hh:291
netevent::left
Float_t * left
min cluster time relative to segment start
Definition: netevent.hh:78
netevent::fCurrent
Int_t fCurrent
pointer to the analyzed TTree or TChain
Definition: netevent.hh:46
wavecomplex::real
double real() const
Definition: wavecomplex.hh:51
netevent::b_norm
TBranch * b_norm
Definition: netevent.hh:178
skymap::getTheta
double getTheta(size_t i)
Definition: skymap.hh:206
netevent::ioSNR
Float_t * ioSNR
reconstructed snr waveform
Definition: netevent.hh:120
injection::iota
Float_t * iota
source psi angle
Definition: injection.hh:64
netevent::fChain
TTree * fChain
root input file cointainig the analyzed TTree
Definition: netevent.hh:45
skymap::getThetaStep
double getThetaStep(size_t i)
Definition: skymap.hh:224
skymap::theta_1
double theta_1
Definition: skymap.hh:303
netevent::b_ioSNR
TBranch * b_ioSNR
Definition: netevent.hh:196
netevent::Init
TTree * Init(TString fName, int n)
Definition: netevent.cc:29
netevent::b_snr
TBranch * b_snr
Definition: netevent.hh:191
netevent::type
Int_t * type
event ID
Definition: netevent.hh:56
M
#define M
Definition: UniqSLagsList.C:3
m
int m
Definition: cwb_net.C:10
skymap::getPhiStep
double getPhiStep(size_t i)
Definition: skymap.hh:164
wavearray::start
virtual void start(double s)
Definition: wavearray.hh:119
j
int j
Definition: cwb_net.C:10
detector::rate
double rate
Definition: detector.hh:323
network::getMRAwave
bool getMRAwave(size_t ID, size_t lag, char atype='S', int mode=0, bool tof=false)
Definition: network.cc:3635
injection::strain
Double_t strain
Definition: injection.hh:62
i
i drho i
Definition: cwb_epparameters.C:85
detector::tau
skymap tau
Definition: detector.hh:328
network::ifoList
std::vector< detector * > ifoList
Definition: network.hh:590
N
#define N
netevent::b_ecor
TBranch * b_ecor
Definition: netevent.hh:177
netevent::b_inet
TBranch * b_inet
Definition: netevent.hh:176
injection::theta
Float_t * theta
source phi angle
Definition: injection.hh:66
netevent::b_iota
TBranch * b_iota
Definition: netevent.hh:147
network::Edge
double Edge
Definition: network.hh:560
net
network ** net
NOISE_MDC_SIMULATION.
Definition: CWB_Plugin_HEN_BKG_Config.C:4
skymap::theta_2
double theta_2
Definition: skymap.hh:304
network::ifoListSize
size_t ifoListSize()
Definition: network.hh:413
netevent::run
Int_t run
max size used by allocate() for the probability maps
Definition: netevent.hh:53
PI
#define PI
Definition: watfun.hh:14
clusterdata::likesky
float likesky
Definition: netcluster.hh:43
network::aNET
double aNET
Definition: network.hh:562
network::optim
bool optim
Definition: network.hh:572
nIFO
#define nIFO
Definition: DrawSensitivitiesS5.C:11
skymap::getSkyIndex
size_t getSkyIndex(double th, double ph)
param: theta param: phi
Definition: skymap.cc:702
netevent::b_psi
TBranch * b_psi
Definition: netevent.hh:146
netevent::b_phi
TBranch * b_phi
Definition: netevent.hh:144
netevent::frequency
Float_t * frequency
GPS stop time of the cluster.
Definition: netevent.hh:83
netevent::b_rate
TBranch * b_rate
Definition: netevent.hh:134
network::precision
double precision
Definition: network.hh:575
netevent::b_netcc
TBranch * b_netcc
Definition: netevent.hh:170
factor
double factor
Definition: ReadCStrainFromJobFile.C:18
ra
double ra
Definition: ConvertGWGC.C:46
injection::redshift
Float_t redshift
Definition: injection.hh:59
netevent::factor
Float_t factor
Definition: netevent.hh:107
detector::antenna
wavecomplex antenna(double, double, double=0.)
param: source theta,phi, polarization angle psi in degrees
Definition: detector.cc:473
netevent::output2G
void output2G(TTree *, network *, size_t, int, double)
Definition: netevent.cc:687
injection::psi
Float_t * psi
Definition: injection.hh:63
netevent::nevent
Int_t nevent
run ID
Definition: netevent.hh:54
netevent::b_low
TBranch * b_low
Definition: netevent.hh:160
netcluster::shift
double shift
Definition: netcluster.hh:364
netevent::Psave
Int_t Psave
number of detectors
Definition: netevent.hh:48
netevent::b_left
TBranch * b_left
Definition: netevent.hh:154
injection::hrss
Double_t * hrss
estimated bandwidth
Definition: injection.hh:76
WAVE_TREE_NAME
#define WAVE_TREE_NAME
Definition: netevent.cc:23
netevent::volume
Int_t * volume
1/rate - wavelet time resolution
Definition: netevent.hh:60
netevent::b_spin
TBranch * b_spin
Definition: netevent.hh:189
netevent::b_penalty
TBranch * b_penalty
Definition: netevent.hh:180
detector::ED
double ED[5]
Definition: detector.hh:316
netevent::Deff
Float_t * Deff
eBBH array
Definition: netevent.hh:111
netevent::gps
Double_t * gps
average center_of_gravity time
Definition: netevent.hh:76
skymap::getRA
double getRA(size_t i)
Definition: skymap.hh:197
netevent::b_eBBH
TBranch * b_eBBH
Definition: netevent.hh:186
netevent::b_eventID
TBranch * b_eventID
Definition: netevent.hh:131
netevent::b_nevent
TBranch * b_nevent
Definition: netevent.hh:130
int
i() int(T_cor *100))
netevent::b_erA
TBranch * b_erA
Definition: netevent.hh:166
Pi
double Pi
Definition: DrawPhaseShift.C:12
injection::mchirp
Float_t mchirp
Definition: injection.hh:55
netevent::b_sSNR
TBranch * b_sSNR
Definition: netevent.hh:193
netevent::penalty
Float_t penalty
Definition: netevent.hh:104
netevent::stop
Double_t * stop
GPS start time of the cluster.
Definition: netevent.hh:81
clusterdata::enrgsky
float enrgsky
Definition: netcluster.hh:36
network::eCOR
double eCOR
Definition: network.hh:564
netevent::b_null
TBranch * b_null
Definition: netevent.hh:168
netevent::gnet
Float_t gnet
network energy disbalance: 0 - total, 1 - 00-phase, 2 - 90-phase
Definition: netevent.hh:98
netcluster::start
double start
Definition: netcluster.hh:361
network::gNET
double gNET
Definition: network.hh:561
netevent::Dump
void Dump(TString analysis="2G")
Definition: netevent.hh:306
netevent::b_time
TBranch * b_time
Definition: netevent.hh:151
skymap::getDEC
double getDEC(size_t i)
Definition: skymap.hh:233
netevent::lag
Float_t * lag
time between consecutive events
Definition: netevent.hh:65
detector::RWFID
std::vector< int > RWFID
Definition: detector.hh:363
clusterdata::netcc
float netcc
Definition: netcluster.hh:45
injection::mass
Float_t * mass
detector specific effective distance
Definition: injection.hh:79
netevent::xSNR
Float_t * xSNR
energy/noise_variance
Definition: netevent.hh:116
netevent::b_lag
TBranch * b_lag
Definition: netevent.hh:141
network::mdc__IDSize
size_t mdc__IDSize()
Definition: network.hh:396
netevent::b_name
TBranch * b_name
Definition: netevent.hh:133
netevent::ifoList
std::vector< detector * > ifoList
dump file
Definition: netevent.hh:124
netevent::Show
void Show(Int_t entry=-1)
Definition: netevent.cc:393
k
int k
Definition: cwb_report_prod_2.C:150
network::pattern
int pattern
Definition: network.hh:583
netevent::b_likelihood
TBranch * b_likelihood
Definition: netevent.hh:181
netevent::strain
Double_t * strain
time slag [sec]
Definition: netevent.hh:67
netevent::netcc
Float_t * netcc
effective correlated SNR
Definition: netevent.hh:94
netevent::slag
Float_t * slag
time lag [sec]
Definition: netevent.hh:66
netevent::b_gnet
TBranch * b_gnet
Definition: netevent.hh:174
netevent::b_anet
TBranch * b_anet
Definition: netevent.hh:175
netevent::output
void output(TTree *=NULL, network *=NULL, double=0., size_t=0, int=-1)
Definition: netevent.cc:1166
netevent::b_bandwidth
TBranch * b_bandwidth
Definition: netevent.hh:162
netevent::b_usize
TBranch * b_usize
Definition: netevent.hh:138
skymap
Definition: skymap.hh:45
detector::lagShift
wavearray< double > lagShift
Definition: detector.hh:351
entry
double * entry
Definition: cwb_setcuts.C:206
netevent::b_bp
TBranch * b_bp
Definition: netevent.hh:148
injection::phi
Float_t * phi
source iota angle
Definition: injection.hh:65
netevent::b_slag
TBranch * b_slag
Definition: netevent.hh:142
netevent::iota
Float_t * iota
[0]-reconstructed psi or phase of gc, [1]-injected psi angle
Definition: netevent.hh:71
netcluster::size
size_t size()
Definition: netcluster.hh:129
wavecomplex::arg
double arg() const
Definition: wavecomplex.hh:53
WAT::USE_HEALPIX
bool USE_HEALPIX()
Definition: wat.hh:13
netevent::b_strain
TBranch * b_strain
Definition: netevent.hh:143
netevent::b_frequency
TBranch * b_frequency
Definition: netevent.hh:159
netevent::psi
Float_t * psi
[0]-reconstructed, [1]-injected theta angle, [2]-DEC
Definition: netevent.hh:70
netevent::hrss
Double_t * hrss
high-low
Definition: netevent.hh:87
injection::type
Int_t type
Definition: injection.hh:50
netevent::b_gap
TBranch * b_gap
Definition: netevent.hh:140
ind
int * ind
Definition: cwb_fix_live_slag_missed.C:89
netevent::b_start
TBranch * b_start
Definition: netevent.hh:156
netevent::b_size
TBranch * b_size
Definition: netevent.hh:137
clusterdata::norm
float norm
Definition: netcluster.hh:62
clusterdata::aNET
float aNET
Definition: netcluster.hh:60
netevent::b_gps
TBranch * b_gps
Definition: netevent.hh:152
clusterdata::mchirperr
float mchirperr
Definition: netcluster.hh:67
netevent::b_oSNR
TBranch * b_oSNR
Definition: netevent.hh:195
netevent::fP
FILE * fP
injected reconstructed xcor waveform
Definition: netevent.hh:122
skymap::getPhi
double getPhi(size_t i)
Definition: skymap.hh:146
detector::TFmap
WSeries< double > TFmap
Definition: detector.hh:336
netevent::phi
Float_t * phi
sqrt(h+*h+ + hx*hx)
Definition: netevent.hh:68
netevent::norm
Float_t norm
Definition: netevent.hh:102
netevent::b_Deff
TBranch * b_Deff
Definition: netevent.hh:187
clusterdata::Gnoise
float Gnoise
Definition: netcluster.hh:42
clusterdata::iNET
float iNET
Definition: netcluster.hh:61
netevent::b_noise
TBranch * b_noise
Definition: netevent.hh:165
netevent::theta
Float_t * theta
[0]-reconstructed, [1]-injected phi angle, [2]-RA
Definition: netevent.hh:69
netevent::time
Double_t * time
beam pattern coefficients for hx
Definition: netevent.hh:75
T
double T
Definition: testWDM_4.C:11
netevent::b_right
TBranch * b_right
Definition: netevent.hh:153
netevent::rate
Int_t * rate
event name: "" - prod, mdc_name - sim
Definition: netevent.hh:58
clusterdata::netRHO
float netRHO
Definition: netcluster.hh:50
netevent::b_iSNR
TBranch * b_iSNR
Definition: netevent.hh:194
detector::getTFmap
WSeries< double > * getTFmap()
param: no parameters
Definition: detector.hh:161
injection::distance
Float_t distance
Definition: injection.hh:54
detector::xSNR
double xSNR
Definition: detector.hh:321
network::getNDM
double getNDM(size_t i, size_t j)
param: first detector param: second detector
Definition: network.hh:122
fabs
double fabs(const Complex &x)
Definition: numpy.cc:37
netevent::nill
Float_t * nill
un-biased null statistics
Definition: netevent.hh:92
netevent::spin
Float_t * spin
mass[2], binary mass parameters
Definition: netevent.hh:113
netevent::b_mass
TBranch * b_mass
Definition: netevent.hh:188
netevent::sSNR
Float_t * sSNR
data-signal correlation Xk*Sk
Definition: netevent.hh:117
network::norm
double norm
Definition: network.hh:565
netevent::erA
Float_t * erA
noise rms
Definition: netevent.hh:89
netevent::b_stop
TBranch * b_stop
Definition: netevent.hh:157
netevent::null
Float_t * null
probability cc skymap
Definition: netevent.hh:91
netevent::likelihood
Float_t likelihood
Definition: netevent.hh:105
slag
Definition: Toolbox.hh:81
network::getwc
netcluster * getwc(size_t n)
param: delay index
Definition: network.hh:421
netevent::bx
Float_t * bx
beam pattern coefficients for hp
Definition: netevent.hh:73
clusterdata::subnet
float subnet
Definition: netcluster.hh:47
sprintf
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
clusterdata::cFreq
float cFreq
Definition: netcluster.hh:58
netevent::b_bx
TBranch * b_bx
Definition: netevent.hh:149
clusterdata::skycc
float skycc
Definition: netcluster.hh:44
netcluster::cData
std::vector< clusterdata > cData
Definition: netcluster.hh:373
netevent::operator=
virtual netevent & operator=(const netevent &)
Definition: netevent.cc:587
clusterdata::cTime
float cTime
Definition: netcluster.hh:57
clusterdata::netnull
float netnull
Definition: netcluster.hh:40
clusterdata::skyChi2
float skyChi2
Definition: netcluster.hh:46
netevent::noise
Double_t * noise
hrss
Definition: netevent.hh:88
skymap::set
void set(size_t i, double a)
param: sky index param: value to set
Definition: skymap.hh:104
skymap::phi_2
double phi_2
Definition: skymap.hh:306
netevent::name
string * name
event type: [0] - prod, [1]-sim
Definition: netevent.hh:57
wavearray::data
DataType_t * data
Definition: wavearray.hh:301
netevent::snr
Float_t * snr
spin[6], binary spin parameters
Definition: netevent.hh:115
detector::null
double null
Definition: detector.hh:318
netevent::b_factor
TBranch * b_factor
Definition: netevent.hh:183
K
int K
Definition: CWB_Plugin_Sim4_BRST_Config.C:87
id
Long_t id
Definition: cwb_condor_check.C:48
clusterdata::nDoF
float nDoF
Definition: netcluster.hh:63
skymap::phi_1
double phi_1
Definition: skymap.hh:305
skymap::get
double get(size_t i)
param: sky index
Definition: skymap.cc:681
ctime
double ctime
Definition: TestCreateFrList.C:32
netevent::GetEntries
Int_t GetEntries()
Definition: netevent.cc:381
injection::fill_in
Bool_t fill_in(network *, int, bool=true)
Definition: injection.cc:325
netevent::neted
Float_t * neted
network correlation coefficients: 0-net,1-pc,2-cc,3-net2
Definition: netevent.hh:95
skymap::size
size_t size()
Definition: skymap.hh:118
netevent::b_xSNR
TBranch * b_xSNR
Definition: netevent.hh:192
network::getmdc__ID
size_t getmdc__ID(size_t n)
Definition: network.hh:408
fName
char fName[256]
Definition: CreateListEBBH.C:25
for
for(int i=0;i< 101;++i) Cos2[2][i]=0
Definition: CWB_Plugin_MDC_OTF_Config_BRST.C:36
injection::pwf
wavearray< double > ** pwf
[x1,y1,z1,x2,y2,z2] components of spin vector
Definition: injection.hh:82
wavearray::resize
virtual void resize(unsigned int)
Definition: wavearray.cc:445
network::MRA
bool MRA
Definition: network.hh:581
network::tYPe
char tYPe
Definition: network.hh:570
network::setndm
bool setndm(size_t, size_t, bool=true, int=1)
param: cluster ID param: lag index param: statistic identificator param: resolution idenificator retu...
Definition: network.cc:6494
injection::spin
Float_t * spin
[m1,m2], binary mass parameters
Definition: injection.hh:80
clusterdata::chirpEllip
float chirpEllip
Definition: netcluster.hh:71
network::SETNDM
bool SETNDM(size_t, size_t, bool=true, int=1)
Definition: network.cc:6793
clusterdata::netecor
float netecor
Definition: netcluster.hh:38
netevent::chirp
Float_t * chirp
range to source: [0/1]-rec/inj
Definition: netevent.hh:109
clusterdata::theta
float theta
Definition: netcluster.hh:52
TAU
double TAU[nWF]
Definition: wmdc_config_rde.C:28
pD
detector ** pD
Definition: cwb_dump_events.C:34
netevent::b_nill
TBranch * b_nill
Definition: netevent.hh:169
netevent::b_chirp
TBranch * b_chirp
Definition: netevent.hh:185
network::iNET
double iNET
Definition: network.hh:563
netevent::b_run
TBranch * b_run
Definition: netevent.hh:129
netevent::Psm
skymap * Psm
error angle
Definition: netevent.hh:90
injection::Deff
Float_t * Deff
injected snr in the detectors
Definition: injection.hh:78
clusterdata::netrho
float netrho
Definition: netcluster.hh:51
exit
exit(0)
detector::enrg
double enrg
Definition: detector.hh:319
netevent::bp
Float_t * bp
[0]-reconstructed iota angle, [1]-injected iota angle
Definition: netevent.hh:72
netevent::iSNR
Float_t * iSNR
energy of reconstructed responses Sk*Sk
Definition: netevent.hh:118
netevent::eventID
Int_t * eventID
event count
Definition: netevent.hh:55
netevent::bandwidth
Float_t * bandwidth
max frequency
Definition: netevent.hh:86
netevent::usize
Int_t usize
cluster size (black pixels only)
Definition: netevent.hh:62
netevent::oSNR
Float_t * oSNR
injected snr waveform
Definition: netevent.hh:119