Logo Coherent WaveBurst  
Reference Guide 		Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
gnetwork.cc
Go to the documentation of this file.
1 
2 #include "gnetwork.hh"
3 #include "constants.hh"
4 
5 using namespace std;
6 using namespace ROOT::Math;
7 
8 ClassImp(gnetwork)
9 
10 
11 //______________________________________________________________________________
12 /* Begin_Html
13 <center><h2>gnetwork class</h2></center>
14 This class gnetwork is derived from the network class and uses the methods of the class gskymap
15 It is used to produce the skymap plots of the detector's network antenna pattern or the statistics produced in the cwb likelihood stage.
16 
17 <p>
18 <h3><a name="example">Example</a></h3>
19 <p>
20 The macro <a href="./tutorials/gwat/DrawAntennaPattern.C.html">DrawAntennaPattern.C</a> is an example which shown how to use the gnetwork class.<br>
21 The pictures below gives the macro output plot.<br>
22 <p>
23 
24 End_Html
25 
26 Begin_Macro
27 DrawAntennaPattern.C
28 End_Macro */
29 
30 
31 //______________________________________________________________________________
32 gnetwork::gnetwork(int nifo, TString* ifo, detectorParams* detParms) : network() {
33 //
34 // constructor
35 //
36 // Input: nifo - number of detectors
37 // ifo - array of built-in ifo names
38 // Ex : TString ifo[3] = {"L1","H1","V1"}
39 // detParms - array of usewr defined detector's parameters
40 // it is used when ifo entry is empty -> ""
41 //
42 
43  for(int n=0; n<nifo; n++) {
44  if((ifo!=NULL)&&(ifo[n]!="")) {
45  pD[n] = new detector(const_cast<char*>(ifo[n].Data())); // built in detector
46  } else {
47  if(detParms!=NULL) {
48  pD[n] = new detector(detParms[n]); // user define detector
49  } else {
50  cout << "gnetwork::gnetwork : Error - user define detector with no detParams" << endl;
51  exit(1);
52  }
53  }
54  }
55 
56  for(int n=0; n<nifo; n++) this->add(pD[n]);
57  Init();
58 }
59 
60 //______________________________________________________________________________
61 gnetwork::gnetwork(network* net) : network() {
62 //
63 // constructor
64 //
65 // Input: net - pointer to network
66 //
67 
68  int nifo = net->ifoListSize();
69  for(int n=0; n<nifo; n++) pD[n] = new detector(*net->getifo(n));
70 
71  for(int n=0; n<nifo; n++) this->add(pD[n]);
72  Init();
73 }
74 
75 //______________________________________________________________________________
76 gnetwork::~gnetwork() {
77 //
78 // destructor
79 //
80 
81  ClearSites();
82  ClearSitesArms();
83  ClearSitesLabel();
84  ClearCircles();
85  for(int n=0;n<NIFO_MAX;n++) if(pD[n]!=NULL) delete pD[n];
86  for(int i=0;i<2;i++) if(grP[i]!=NULL) delete grP[i];
87  if(canvasP!=NULL) delete canvasP;
88  for(size_t n=0; n<this->ifoList.size(); n++) if(this->ifoList[n]!=NULL) delete this->ifoList[n];
89 }
90 
91 //______________________________________________________________________________
92 void
93 gnetwork::Init() {
94 //
95 // initialize internal parameters
96 //
97 
98  siteLabelType=0;
99  for(int n=0;n<NIFO_MAX;n++) pD[n]=NULL;
100  for(int i=0;i<2;i++) grP[i]=NULL;
101  canvasP=NULL;
102 }
103 
104 //______________________________________________________________________________
105 void
106 gnetwork::Draw(TString smName, network* net) {
107 //
108 // Draw skymap network info
109 //
110 // Input: smName - skymap name
111 //
112 // nSensitivity : network sensitivity
113 // nAlignment : network alignment factor
114 // nCorrelation : network correlation coefficient
115 // nLikelihood : network likelihood
116 // nNullEnergy : network null energy
117 // nPenalty : signal * noise penalty factor
118 // nCorrEnergy : reduced correlated energy
119 // nNetIndex : network index
120 // nDisbalance : energy disbalance
121 // nSkyStat : sky optimization statistic
122 // nEllipticity : waveform ellipticity
123 // nPolarisation : polarisation angle
124 // nProbability : probability skymap
125 //
126 // index : theta, phi mask index array
127 // skyMask : index array for setting sky mask
128 // skyMaskCC : index array for setting sky mask Celestial Coordinates
129 // skyHole : static sky mask describing "holes"
130 // veto : veto array for pixel selection
131 // skyProb : sky probability
132 // skyENRG : energy skymap
133 //
134 // net - pointer to network object
135 // if net!=NULL then net is used to retrive the data to be plotted
136 //
137 
138  network* NET = net!=NULL ? net : this;
139 
140  if(smName=="nSensitivity") {Draw(NET->nSensitivity,smName);return;} // network sensitivity
141  if(smName=="nAlignment") {Draw(NET->nAlignment,smName);return;} // network alignment factor
142  if(smName=="nCorrelation") {Draw(NET->nCorrelation,smName);return;} // network correlation coefficient
143  if(smName=="nLikelihood") {Draw(NET->nLikelihood,smName);return;} // network likelihood
144  if(smName=="nNullEnergy") {Draw(NET->nNullEnergy,smName);return;} // network null energy
145  if(smName=="nPenalty") {Draw(NET->nPenalty,smName);return;} // signal * noise penalty factor
146  if(smName=="nCorrEnergy") {Draw(NET->nCorrEnergy,smName);return;} // reduced correlated energy
147  if(smName=="nNetIndex") {Draw(NET->nNetIndex,smName);return;} // network index
148  if(smName=="nDisbalance") {Draw(NET->nDisbalance,smName);return;} // energy disbalance
149  if(smName=="nSkyStat") {Draw(NET->nSkyStat,smName);return;} // sky optimization statistic
150  if(smName=="nEllipticity") {Draw(NET->nEllipticity,smName);return;} // waveform ellipticity
151  if(smName=="nPolarisation") {Draw(NET->nPolarisation,smName);return;} // polarisation angle
152  if(smName=="nProbability") {Draw(NET->nProbability,smName);return;} // probability skymap
153 
154  if(smName=="index") {Draw(NET->index,smName);return;} // theta, phi mask index array
155  if(smName=="skyMask") {Draw(NET->skyMask,smName);return;} // index array for setting sky mask
156  if(smName=="skyMaskCC") {Draw(NET->skyMaskCC,smName);return;} // index array for setting sky mask Celestial Coordinates
157  if(smName=="skyHole") {Draw(NET->skyHole,smName);return;} // static sky mask describing "holes"
158  if(smName=="veto") {Draw(NET->veto,smName);return;} // veto array for pixel selection
159  if(smName=="skyProb") {Draw(NET->skyProb,smName);return;} // sky probability
160  if(smName=="skyENRG") {Draw(NET->skyENRG,smName);return;} // energy skymap
161 
162  cout << "gnetwork::Draw : Warning - skymap '" << smName.Data() << "' not valid !" << endl;
163 }
164 
165 //______________________________________________________________________________
166 double
167 gnetwork::GetSite(TString ifo, TString type) {
168 //
169 // return detector infos
170 //
171 // Input: ifo - detector name
172 // type - detector info type
173 //
174 // THETA : theta full
175 // THETA_D : theta degrees
176 // THETA_M : theta minutes
177 // THETA_S : theta seconds
178 //
179 // PHI : phi full
180 // PHI_D : phi degrees
181 // PHI_M : phi minutes
182 // PHI_S : phi seconds
183 //
184 // "" : print out theta,phi infos
185 //
186 // if type invalid return 1111
187 //
188 
189  int nIFO = this->ifoListSize();
190  if(nIFO<1) {cout << "gnetwork::GetSite nIFO must be >= 1" << endl;return -1111;}
191 
192  detector *pD[NIFO];
193  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
194 
195  int ifoId=0;
196  int nifo=0;
197  for(int n=0; n<nIFO; n++) {
198  if(ifo.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId=n;}
199  }
200  if(nifo<1) {cout << "gnetwork::GetSite ifo not valid" << endl;return -1111;}
201 
202  detectorParams dP = pD[ifoId]->getDetectorParams();
203 
204  double phi=dP.longitude;
205  double theta=dP.latitude;
206 
207  char LAT;
208  double theta_t=theta;
209  if(theta_t>0) LAT='N'; else {LAT='S';theta_t=-theta_t;}
210  int theta_d = int(theta_t);
211  int theta_m = int((theta_t-theta_d)*60);
212  float theta_s = (theta_t-theta_d-theta_m/60.)*3600.;
213 
214  char LON;
215  double phi_t=phi;
216  if(phi_t>0) LON='E'; else {LON='W';phi_t=-phi_t;}
217  int phi_d = int(phi_t);
218  int phi_m = int((phi_t-phi_d)*60);
219  float phi_s = (phi_t-phi_d-phi_m/60.)*3600.;
220 
221  type.ToUpper();
222  if(type.CompareTo("THETA")==0) return theta;
223  if(type.CompareTo("THETA_D")==0) return theta_d;
224  if(type.CompareTo("THETA_M")==0) return theta_m;
225  if(type.CompareTo("THETA_S")==0) return theta_s;
226  if(type.CompareTo("PHI")==0) return phi;
227  if(type.CompareTo("PHI_D")==0) return phi_d;
228  if(type.CompareTo("PHI_M")==0) return phi_m;
229  if(type.CompareTo("PHI_S")==0) return phi_s;
230 
231  if(type.CompareTo("")==0) pD[ifoId]->print();
232 
233  return 1111;
234 }
235 
236 //______________________________________________________________________________
237 void
238 gnetwork::DrawSites(Color_t mcolor, Size_t msize, Style_t mstyle) {
239 //
240 // Draw a marker in the ifo site position (see ROOT TMarker input parameters)
241 //
242 // Input: mcolor - color of marker
243 // msize - size of marker
244 // mstyle - style of marker
245 //
246 
247  if(gSM.goff) return;
248  if(gSM.canvas==NULL) return;
249 
250  int nIFO = this->ifoListSize();
251  if(nIFO<1) return;
252 
253  ClearSites();
254 
255  detector *pD[NIFO];
256  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
257 
258  gSM.canvas->cd();
259 
260  double r2d = 180./TMath::Pi();
261  for(int n=0;n<nIFO;n++) {
262  XYZVector Rv(pD[n]->Rv[0],pD[n]->Rv[1],pD[n]->Rv[2]);
263  TVector3 vRv(Rv.X(),Rv.Y(),Rv.Z());
264 
265  double phi=vRv.Phi()*r2d;
266  double theta=vRv.Theta()*r2d;
267 
268  theta=90-theta;
269 
270  if (gSM.projection.CompareTo("hammer")==0) {
271  gSM.ProjectHammer(phi, theta, phi, theta);
272  } else if (gSM.projection.CompareTo("sinusoidal")==0) {
273  gSM.ProjectSinusoidal(phi, theta, phi, theta);
274  } else if (gSM.projection.CompareTo("parabolic")==0) {
275  gSM.ProjectParabolic(phi, theta, phi, theta);
276  }
277  if (gSM.coordinate.CompareTo("CWB")==0) GeographicToCwb(phi,theta,phi,theta);
278  TMarker* pM = new TMarker(phi,theta,mstyle);
279  pM->SetMarkerSize(msize);
280  pM->SetMarkerColor(mcolor);
281  pM->Draw();
282  siteM.push_back(pM);
283  pM = new TMarker(phi,theta,mstyle);
284  pM->SetMarkerSize(msize/2);
285  pM->SetMarkerColor(kWhite);
286  pM->Draw();
287  siteM.push_back(pM);
288  }
289 
290  return;
291 }
292 
293 //______________________________________________________________________________
294 void
295 gnetwork::DrawSitesArms(double mlength, Color_t lcolor, Size_t lwidth, Style_t lstyle) {
296 //
297 // Draw detector's arms in the ifo site position (see ROOT TPolyLine input parameters)
298 //
299 // Input: mlength - arm's length in KM
300 // mcolor - color of line
301 // msize - size of line
302 // mstyle - style of line
303 //
304 
305 
306  if(gSM.goff) return;
307  if(gSM.canvas==NULL) return;
308 
309  int nIFO = this->ifoListSize();
310  if(nIFO<1) return;
311 
312  ClearSitesArms();
313 
314  detector *pD[NIFO];
315  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
316 
317  gSM.canvas->cd();
318 
319  double mlength_deg = 360.*mlength/(TMath::TwoPi()*watconstants::EarthEquatorialRadius());
320 
321  double x[nIFO][2][2];
322  double y[nIFO][2][2];
323  double phi,theta;
324  double r2d = 180./TMath::Pi();
325  for(int n=0;n<nIFO;n++) {
326  XYZVector Rv(pD[n]->Rv[0],pD[n]->Rv[1],pD[n]->Rv[2]);
327  XYZVector Ex(pD[n]->Ex[0],pD[n]->Ex[1],pD[n]->Ex[2]);
328  XYZVector Ey(pD[n]->Ey[0],pD[n]->Ey[1],pD[n]->Ey[2]);
329  TVector3 vRv(Rv.X(),Rv.Y(),Rv.Z());
330  TVector3 vEx(Ex.X(),Ex.Y(),Ex.Z());
331  TVector3 vEy(Ey.X(),Ey.Y(),Ey.Z());
332  vEx=mlength*vEx+vRv;
333  vEy=mlength*vEy+vRv;
334 
335  phi=vRv.Phi()*r2d;
336  theta=vRv.Theta()*r2d;
337  theta=90-theta;
338  if (gSM.projection.CompareTo("hammer")==0) {
339  gSM.ProjectHammer(phi, theta, phi, theta);
340  } else if (gSM.projection.CompareTo("sinusoidal")==0) {
341  gSM.ProjectSinusoidal(phi, theta, phi, theta);
342  } else if (gSM.projection.CompareTo("parabolic")==0) {
343  gSM.ProjectParabolic(phi, theta, phi, theta);
344  }
345  if (gSM.coordinate.CompareTo("CWB")==0) GeographicToCwb(phi,theta,phi,theta);
346  x[n][0][0]=phi;y[n][0][0]=theta;
347  x[n][1][0]=phi;y[n][1][0]=theta;
348 
349  phi=vEx.Phi()*r2d;
350  theta=vEx.Theta()*r2d;
351  theta=90-theta;
352  if (gSM.projection.CompareTo("hammer")==0) {
353  gSM.ProjectHammer(phi, theta, phi, theta);
354  } else if (gSM.projection.CompareTo("sinusoidal")==0) {
355  gSM.ProjectSinusoidal(phi, theta, phi, theta);
356  } else if (gSM.projection.CompareTo("parabolic")==0) {
357  gSM.ProjectParabolic(phi, theta, phi, theta);
358  }
359  if (gSM.coordinate.CompareTo("CWB")==0) GeographicToCwb(phi,theta,phi,theta);
360  x[n][0][1]=phi;y[n][0][1]=theta;
361  // fix draw on the edge
362  if (gSM.coordinate.CompareTo("CWB")==0) {
363  if((x[n][0][0]>360-mlength_deg)&&(x[n][0][0]+x[n][0][1]-360<mlength_deg)) x[n][0][1]=360;
364  }
365  if (gSM.coordinate.CompareTo("GEOGRAPHIC")==0) {
366  if((x[n][0][0]>180-mlength_deg)&&(x[n][0][0]+x[n][0][1]<mlength_deg)) x[n][0][1]=180;
367  }
368 
369  phi=vEy.Phi()*r2d;
370  theta=vEy.Theta()*r2d;
371  theta=90-theta;
372  if (gSM.projection.CompareTo("hammer")==0) {
373  gSM.ProjectHammer(phi, theta, phi, theta);
374  } else if (gSM.projection.CompareTo("sinusoidal")==0) {
375  gSM.ProjectSinusoidal(phi, theta, phi, theta);
376  } else if (gSM.projection.CompareTo("parabolic")==0) {
377  gSM.ProjectParabolic(phi, theta, phi, theta);
378  }
379  if (gSM.coordinate.CompareTo("CWB")==0) GeographicToCwb(phi,theta,phi,theta);
380  x[n][1][1]=phi;y[n][1][1]=theta;
381  // fix draw on the edge
382  if (gSM.coordinate.CompareTo("CWB")==0) {
383  if((x[n][1][0]>360-mlength_deg)&&(x[n][1][0]+x[n][1][1]-360<mlength_deg)) x[n][1][1]=360;
384  }
385  if (gSM.coordinate.CompareTo("GEOGRAPHIC")==0) {
386  if((x[n][1][0]>180-mlength_deg)&&(x[n][1][0]+x[n][1][1]<mlength_deg)) x[n][1][1]=180;
387  }
388  }
389 
390  // Draw X,Y arms
391  for(int n=0;n<nIFO;n++) {
392  for(int m=0;m<2;m++) {
393  TPolyLine *pL = new TPolyLine(2,x[n][m],y[n][m]);
394  pL->SetLineColor(lcolor);
395  pL->SetLineStyle(lstyle);
396  pL->SetLineWidth(lwidth);
397  pL->Draw();
398  siteA.push_back(pL);
399  }
400  }
401 
402  return;
403 }
404 
405 //______________________________________________________________________________
406 void
407 gnetwork::DrawSitesShortLabel(Color_t tcolor, Size_t tsize, Font_t tfont) {
408 //
409 // Draw detector's short label in the ifo site position (see ROOT TLatex input parameters)
410 // Ex : L,H,V
411 //
412 // Input: tcolor - color of text
413 // tsize - size of text
414 // tfont - font of text
415 //
416 
417  siteLabelType=1;
418  DrawSitesLabel(tcolor, tsize, tfont);
419 }
420 
421 //______________________________________________________________________________
422 void
423 gnetwork::DrawSitesLabel(Color_t tcolor, Size_t tsize, Font_t tfont) {
424 //
425 // Draw detector's short label in the ifo site position (see ROOT TLatex input parameters)
426 // Ex : L1,H1,V1
427 //
428 // Input: tcolor - color of text
429 // tsize - size of text
430 // tfont - font of text
431 //
432 
433  if(gSM.goff) return;
434  if(gSM.canvas==NULL) return;
435 
436  int nIFO = this->ifoListSize();
437  if(nIFO<1) return;
438 
439  ClearSitesLabel();
440 
441  detector *pD[NIFO];
442  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
443 
444  gSM.canvas->cd();
445 
446  double r2d = 180./TMath::Pi();
447  for(int n=0;n<nIFO;n++) {
448  XYZVector Rv(pD[n]->Rv[0],pD[n]->Rv[1],pD[n]->Rv[2]);
449  TVector3 vRv(Rv.X(),Rv.Y(),Rv.Z());
450 
451  double phi=vRv.Phi()*r2d;
452  double theta=vRv.Theta()*r2d;
453 
454  theta=90-theta;
455 
456  if (gSM.projection.CompareTo("hammer")==0) {
457  gSM.ProjectHammer(phi, theta, phi, theta);
458  } else if (gSM.projection.CompareTo("sinusoidal")==0) {
459  gSM.ProjectSinusoidal(phi, theta, phi, theta);
460  } else if (gSM.projection.CompareTo("parabolic")==0) {
461  gSM.ProjectParabolic(phi, theta, phi, theta);
462  }
463  if (gSM.coordinate.CompareTo("CWB")==0) GeographicToCwb(phi,theta,phi,theta);
464  TString ifoName=this->getifo(n)->Name;
465  if(ifoName.CompareTo("H2")==0) {phi+=2;theta-=10;}
466  else if(ifoName.CompareTo("J1")==0) {phi+=6;theta-=1;}
467  else if(ifoName.CompareTo("H1")==0) {phi+=6;theta-=1;}
468  else if(ifoName.CompareTo("L1")==0) {phi+=3;theta+=5;}
469  else if(ifoName.CompareTo("V1")==0) {phi+=4;theta+=2;}
470  else {phi+=2;theta+=2;}
471  if(siteLabelType==1) {
472  if(ifoName.CompareTo("A1")==0) ifoName=TString("#tilde{A}");
473  else if(ifoName.CompareTo("H2")==0) ifoName=TString("#tilde{H}");
474  else if(ifoName.CompareTo("I2")==0) ifoName=TString("#tilde{I}");
475  else if(ifoName.CompareTo("R2")==0) ifoName=TString("#tilde{R}");
476  else ifoName.Resize(1);
477  }
478  TLatex *pL = new TLatex(phi+1.0,theta, ifoName);
479  pL->SetTextFont(tfont);
480  pL->SetTextSize(tsize);
481  pL->SetTextColor(tcolor);
482  pL->Draw();
483  siteL.push_back(pL);
484  }
485 
486  return;
487 }
488 
489 //______________________________________________________________________________
490 void
491 gnetwork::ClearCircles() {
492 //
493 // Clear all circles
494 //
495 
496  for(int i=0; i<(int)circleL.size(); i++) {if(circleL[i]) delete circleL[i];}
497  circleL.clear();
498  for(int i=0; i<(int)daxisM.size(); i++) {if(daxisM[i]) delete daxisM[i];}
499  daxisM.clear();
500  return;
501 }
502 
503 //______________________________________________________________________________
504 void
505 gnetwork::ClearSites() {
506 //
507 // Clear all site markers
508 //
509 
510  for(int i=0; i<(int)siteM.size(); i++) {if(siteM[i]) delete siteM[i];}
511  siteM.clear();
512  return;
513 }
514 
515 //______________________________________________________________________________
516 void
517 gnetwork::ClearSitesArms() {
518 //
519 // Clear all site arms
520 //
521 
522  for(int i=0; i<(int)siteA.size(); i++) {if(siteA[i]) delete siteA[i];}
523  siteA.clear();
524  return;
525 }
526 
527 //______________________________________________________________________________
528 void
529 gnetwork::ClearSitesLabel() {
530 //
531 // Clear all site labels
532 //
533 
534  for(int i=0; i<(int)siteT.size(); i++) {if(siteT[i]) delete siteT[i];}
535  siteT.clear();
536  for(int i=0; i<(int)siteL.size(); i++) {if(siteL[i]) delete siteL[i];}
537  siteL.clear();
538  for(int i=0; i<(int)siteP.size(); i++) {if(siteP[i]) delete siteP[i];}
539  siteP.clear();
540  return;
541 }
542 
543 //______________________________________________________________________________
544 double
545 gnetwork::GetAntennaPattern(double phi, double theta, double psi, int polarization) {
546 //
547 // Get Network Antenna Pattern value in the DPF : Dominant Polarization Frame
548 //
549 // Input: phi - phi (degrees)
550 // theta - theta (degrees)
551 // psi - polarization angle (degrees)
552 // polarization - DPF component
553 //
554 // 0 -> |Fx| - DPF
555 // 1 -> |F+| - DPF
556 // 2 -> |Fx|/|F+| - DPF
557 // 3 -> sqrt((|F+|^2+|Fx|^2)/nIFO) - DPF
558 // 4 -> |Fx|^2 - DPF
559 // 5 -> |F+|^2 - DPF
560 // 6 -> Fx - only with 1 detector
561 // 7 -> F+ - only with 1 detector
562 // 8 -> F1x/F2x - only with 2 detectors
563 // 9 -> F1+/F2+ - only with 2 detectors
564 // 10 -> sqrt(|F1+|^2+|F1x|^2)/sqrt(|F2+|^2+|F2x|^2) - only with 2 detectors
565 // 12 -> DPF Angle
566 //
567 
568  if((polarization<0||polarization>10)&&(polarization!=12))
569  {cout << "polarization must be [0-10][12]" << endl;exit(1);}
570 
571  int nIFO = this->ifoListSize();
572  if(nIFO==0)
573  {cout << "detectors are no declared" << endl;exit(1);}
574  if((nIFO!=1)&&(polarization==6||polarization==7))
575  {cout << "with polarization [6-7] works only with one detector" << endl;exit(1);}
576  if((nIFO!=2)&&(polarization>=8&&polarization<=11))
577  {cout << "with polarization [8-11] works only two one detectors" << endl;exit(1);}
578 
579  detector *pD[NIFO];
580  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
581 
582  double z=0.;
583  std::vector<wavecomplex> F;
584  F.resize(nIFO);
585 
586  for (int n=0;n<nIFO;n++) F[n] = pD[n]->antenna(theta,phi,psi);
587  double gp=0,gx=0,gI=0;
588  for (int n=0;n<nIFO;n++) {
589  gp+=F[n].real()*F[n].real();
590  gx+=F[n].imag()*F[n].imag();
591  gI+=F[n].real()*F[n].imag();
592  }
593  double gR = (gp-gx)/2.;
594  double gr = (gp+gx)/2.;
595  double gc = sqrt(gR*gR+gI*gI);
596 
597  if (polarization==0) z = sqrt(fabs(gr-gc)<1e-8?0.:gr-gc);
598  if (polarization==1) z = sqrt(gr+gc);
599  if (polarization==2&&sqrt(gr+gc)>0) z=sqrt((gr-gc)/(gr+gc));
600  if (polarization==3) z = sqrt(2*gr/nIFO);
601  if (polarization==4) z = (fabs(gr-gc)<1e-8?0.:gr-gc);
602  if (polarization==5) z = (gr+gc);
603  if (polarization==6) z = F[0].imag();
604  if (polarization==7) z = F[0].real();
605  if (polarization==8) z = (F[1].imag()!=0) ? F[0].imag()/F[1].imag() : 0.;
606  if (polarization==9) z = (F[1].real()!=0) ? F[0].real()/F[1].real() : 0.;
607  if (polarization==10) {
608  double F0=sqrt(pow(F[0].real(),2)+pow(F[0].imag(),2));
609  double F1=sqrt(pow(F[1].real(),2)+pow(F[1].imag(),2));
610  z = (F1!=0) ? F0/F1 : 0.;
611  }
612  if (polarization==12) { // DPF Angle
613  double cos2G = gc!=0 ? gR/gc : 0; // (Fp2-Fc2)/Norm
614  double sin2G = gc!=0 ? gI/gc : 0; // 2*Fpc/Norm
615  double r2d = 180./TMath::Pi();
616  z = r2d*atan2(sin2G,cos2G)/2.;
617  }
618 
619  return z;
620 }
621 
622 //______________________________________________________________________________
623 double
624 gnetwork::GetAntennaPattern(TString ifo, double phi, double theta, double psi, bool plus) {
625 //
626 // Get Detector Antenna Pattern value
627 //
628 // Input: ifo - detector name (Ex: "L1","H1","V1")
629 // phi - phi (degrees)
630 // theta - theta (degrees)
631 // psi - polarization angle (degrees)
632 // plus - true/false -> f+/fx
633 //
634 
635  int nIFO = this->ifoListSize();
636  if(nIFO<1) {cout << "gnetwork::GetAntennaPattern nIFO must be >= 1" << endl;return -100;}
637 
638  detector *pD[NIFO];
639  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
640 
641  int ifoId=0;
642  int nifo=0;
643  for(int n=0; n<nIFO; n++) {
644  if(ifo.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId=n;}
645  }
646  if(nifo<1) {cout << "gnetwork::GetAntennaPattern ifo not valid" << endl;return -100;}
647 
648  wavecomplex F = pD[ifoId]->antenna(theta,phi,psi);
649  if(plus) return F.real(); else return F.imag();
650 
651 }
652 
653 //______________________________________________________________________________
654 void
655 gnetwork::DrawAntennaPattern(int polarization, int dpaletteId, bool btitle, int order) {
656 //
657 // Draw Network Antenna Pattern sky map in the DPF : Dominant Polarization Frame
658 //
659 // Input: polarization - DPF component
660 //
661 // 0 -> |Fx| - DPF
662 // 1 -> |F+| - DPF
663 // 2 -> |Fx|/|F+| - DPF
664 // 3 -> sqrt((|F+|^2+|Fx|^2)/nIFO) - DPF
665 // 4 -> |Fx|^2 - DPF
666 // 5 -> |F+|^2 - DPF
667 // 6 -> Fx - only with 1 detector
668 // 7 -> F+ - only with 1 detector
669 // 8 -> F1x/F2x - only with 2 detectors
670 // 9 -> F1+/F2+ - only with 2 detectors
671 // 10 -> sqrt(|F1+|^2+|F1x|^2)/sqrt(|F2+|^2+|F2x|^2) - only with 2 detectors
672 // 11 -> The same as (10) but averaged over psi - only with 2 detectors
673 // 12 -> DPF Angle
674 //
675 // dpaletteId - palette ID
676 // btitle - if true then a default title is provided
677 // order - healpix order used for Antenna Pattern
678 //
679 
680  if(gSM.goff) return;
681 
682  gSM=skymap(int(order));
683 
684  if(polarization<0||polarization>12) {
685  cout << "---------------------------------------------------------------------------------------" << endl;
686  cout << "polarization must be [0-12]" << endl;
687  cout << "---------------------------------------------------------------------------------------" << endl;
688  cout << "polarization=0 -> |Fx| DPF" << endl;
689  cout << "polarization=1 -> |F+| DPF" << endl;
690  cout << "polarization=2 -> |Fx|/|F+| DPF" << endl;
691  cout << "polarization=3 -> sqrt((|F+|^2+|Fx|^2)/nIFO) DPF" << endl;
692  cout << "polarization=4 -> |Fx|^2 DPF" << endl;
693  cout << "polarization=5 -> |F+|^2 DPF" << endl;
694  cout << "polarization=6 -> Fx only with 1 detector" << endl;
695  cout << "polarization=7 -> F+ only with 1 detector" << endl;
696  cout << "polarization=8 -> F1x/F2x only with 2 detectors" << endl;
697  cout << "polarization=9 -> F1+/F2+ only with 2 detectors" << endl;
698  cout << "polarization=10 -> sqrt(|F1+|^2+|F1x|^2)/sqrt(|F2+|^2+|F2x|^2) only with 2 detectors" << endl;
699  cout << "polarization=11 -> The same as (10) but averaged over psi only with 2 detectors" << endl;
700  cout << "polarization=12 -> DPF Angle" << endl;
701  cout << "---------------------------------------------------------------------------------------" << endl;
702  return;
703  }
704 
705  int nIFO = this->ifoListSize();
706  if(nIFO==0)
707  {cout << "detectors are no declared" << endl;return;}
708  if((nIFO!=1)&&(polarization==6||polarization==7))
709  {cout << "with polarization [6-7] works only with one detector" << endl;return;}
710  if((nIFO!=2)&&(polarization>=8&&polarization<=11))
711  {cout << "with polarization [8-11] works only two one detectors" << endl;return;}
712 
713  detector *pD[NIFO];
714  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
715 
716  TString pTitle="";
717  switch(polarization) {
718  case 0:
719  pTitle=TString("|F_{x}|");
720  break;
721  case 1:
722  pTitle=TString("|F_{+}|");
723  break;
724  case 2:
725  pTitle=TString("|F_{x}|/|F_{+}|");
726  break;
727  case 3:
728  pTitle=TString("#sqrt{(|F_{+}|^{2}+|F_{x}|^{2})/nIFO}");
729  break;
730  case 4:
731  pTitle=TString("|F_{x}|^{2}");
732  break;
733  case 5:
734  pTitle=TString("|F_{+}|^{2}");
735  break;
736  case 6:
737  pTitle=TString("F_{x}");
738  break;
739  case 7:
740  pTitle=TString("F_{+}");
741  break;
742  case 8:
743  pTitle=TString("F1_{x}/F2_{x}");
744  break;
745  case 9:
746  pTitle=TString("F1_{+}/F2_{+}");
747  break;
748  case 10:
749  pTitle=TString("#sqrt{|F1_{+}|^{2}+|F1_{x}|^{2})/sqrt(|F2_{+}|^{2}+|F2_{x}|^{2}}");
750  break;
751  case 12:
752  pTitle=TString("DPF Angle");
753  break;
754  default:
755  break;
756  }
757  TString netName="";
758  for(int n=0; n<nIFO; n++) {
759  TString ifoName=this->getifo(n)->Name;
760  if(ifoName.CompareTo("A1")==0) ifoName=TString("#tilde{A}");
761  else if(ifoName.CompareTo("H2")==0) ifoName=TString("#tilde{H}");
762  else if(ifoName.CompareTo("I2")==0) ifoName=TString("#tilde{I}");
763  else if(ifoName.CompareTo("R2")==0) ifoName=TString("#tilde{R}");
764  else ifoName.Resize(1);
765  netName+=ifoName;
766  }
767  if(btitle) gSM.h2->SetTitle("Network = "+netName+
768  " "+
769  "Antenna Pattern = "+pTitle);
770 
771  std::vector<wavecomplex> F(nIFO);
772  int size=int(180*2*gSM.resolution*360*2*gSM.resolution);
773  double* x = new double[size];
774  double* y = new double[size];
775  double* z = new double[size];
776  int cnt=0;
777  for (int i=0;i<360*2*gSM.resolution;i++) {
778  for (int j=0;j<180*2*gSM.resolution;j++) {
779  double phi = (double)(i+1)/(double)(2*gSM.resolution);
780  double theta = (double)(j+1)/(double)(2*gSM.resolution);
781  double psi=0;
782  x[cnt]=phi;
783  y[cnt]=theta;
784  if (polarization!=11) {
785 /*
786  psi=0;
787  for (int n=0;n<nIFO;n++) F[n] = pD[n]->antenna(theta,phi,psi);
788  double gp=0,gx=0,gI=0;
789  for (int n=0;n<nIFO;n++) {
790  gp+=F[n].real()*F[n].real();
791  gx+=F[n].imag()*F[n].imag();
792  gI+=F[n].real()*F[n].imag();
793  }
794  double gR = (gp-gx)/2.;
795  double gr = (gp+gx)/2.;
796  double gc = sqrt(gR*gR+gI*gI);
797 
798  if (polarization==0) z[cnt] = sqrt(fabs(gr-gc)<1e-8?0.:gr-gc);
799  if (polarization==1) z[cnt] = sqrt(gr+gc);
800  if (polarization==2&&sqrt(gr+gc)>0) z[cnt]=sqrt((gr-gc)/(gr+gc));
801  if (polarization==3) z[cnt] = sqrt(2*gr/nIFO);
802  //if (polarization==3) z[cnt] = sqrt(2*gr);
803  if (polarization==4) z[cnt] = (fabs(gr-gc)<1e-8?0.:gr-gc);
804  if (polarization==5) z[cnt] = (gr+gc);
805  if (polarization==6) z[cnt] = F[0].imag();
806  if (polarization==7) z[cnt] = F[0].real();
807  if (polarization==8) z[cnt] = (F[1].imag()!=0) ? F[0].imag()/F[1].imag() : 0.;
808  if (polarization==9) z[cnt] = (F[1].real()!=0) ? F[0].real()/F[1].real() : 0.;
809  //if (polarization==8) z[cnt] = F[0].imag()/F[1].imag()>0 ? 1. : -1.;
810  //if (polarization==9) z[cnt] = F[0].real()/F[1].real()>0 ? 1. : -1.;
811  if (polarization==10) {
812  double F0=sqrt(pow(F[0].real(),2)+pow(F[0].imag(),2));
813  double F1=sqrt(pow(F[1].real(),2)+pow(F[1].imag(),2));
814  z[cnt] = (F1!=0) ? F0/F1 : 0.;
815  }
816 */
817  z[cnt] = GetAntennaPattern(phi, theta, psi, polarization);
818  } else {
819  int counts=0;
820  z[cnt] = 0.;
821  //for (int k=0;k<360;k+=8) {
822  for (int k=0;k<1;k+=1) {
823  psi=(double)k;
824  for (int n=0;n<nIFO;n++) F[n] = pD[n]->antenna(theta,phi,psi);
825  //double F0=sqrt(pow(F[0].real(),2)+pow(F[0].imag(),2));
826  //double F1=sqrt(pow(F[1].real(),2)+pow(F[1].imag(),2));
827  double F0=F[0].real()+F[0].imag();
828  double F1=F[1].real()+F[1].imag();
829  z[cnt] += (F1!=0) ? fabs(F0/F1) : 0.;
830  }
831  z[cnt] /= (double)counts;
832  }
833  // fill skymap
834  int ind = gSM.getSkyIndex(theta,phi);
835  gSM.set(ind,z[cnt]);
836  cnt++;
837  }
838  }
839 
840  if (gSM.coordinate.CompareTo("GEOGRAPHIC")==0) {
841  for (int i=0;i<size;i++) CwbToGeographic(x[i],y[i],x[i],y[i]);
842  }
843  if (gSM.coordinate.CompareTo("CELESTIAL")==0) {
844  for (int i=0;i<size;i++) CwbToCelestial(x[i],y[i],x[i],y[i]);
845  }
846  gSM.FillData(size, x, y, z);
847  gSM.Draw(dpaletteId);
848 
849  delete [] x;
850  delete [] y;
851  delete [] z;
852 }
853 
854 //______________________________________________________________________________
855 double
856 gnetwork::GetDelay(TString ifo1, TString ifo2, TString mode) {
857 //
858 // Return the delay (sec) between two detector
859 //
860 // Input: ifo1 - detector name of the first detector (Ex: "L1")
861 // ifo2 - detector name of the second detector (Ex: "V1")
862 // mode - delay mde
863 // "" : (max-min)/2 - default
864 // min : minum delay value
865 // max : maximum delay value
866 // MAX : maximum absolute delay value
867 //
868 
869  int nIFO = this->ifoListSize();
870  if(nIFO<2) {cout << "nIFO must be >= 2" << endl;return -1.;}
871 
872  if(ifo1.CompareTo(ifo2)==0) return 0.0;
873 
874  detector *pD[NIFO];
875  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
876 
877  int ifoId1=0;
878  int ifoId2=0;
879  int nifo=0;
880  for(int n=0; n<nIFO; n++) {
881  if(ifo1.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId1=n;}
882  if(ifo2.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId2=n;}
883  }
884  if(nifo<2) {cout << "gnetwork::GetDelay ifo not valid" << endl;return -1.;}
885 
886  skymap sm0, sTau;
887  double maxTau = -1.;
888  double minTau = 1.;
889  double tmax, tmin;
890 
891  tmax = tmin = 0.;
892  for(int n=0; n<nifo; n++) {
893  sTau = pD[n]->tau;
894  tmax = sTau.max();
895  tmin = sTau.min();
896  if(tmax > maxTau) maxTau = tmax;
897  if(tmin < minTau) minTau = tmin;
898  }
899  if(strstr(mode,"min")) return minTau;
900  if(strstr(mode,"max")) return maxTau;
901  if(strstr(mode,"MAX")) return fabs(maxTau)>fabs(minTau) ? fabs(maxTau) : fabs(minTau);
902  return (maxTau-minTau)/2.;
903 }
904 
905 //______________________________________________________________________________
906 double
907 gnetwork::GetDelay(TString ifo1, TString ifo2, double phi, double theta) {
908 //
909 // Return the delay (sec) between two detector from the direction phi,theta
910 //
911 // Input: ifo1 - detector name of the first detector (Ex: "L1")
912 // ifo2 - detector name of the second detector (Ex: "V1")
913 // phi - phi direction (degrees)
914 // theta - theta direction (degrees)
915 //
916 
917  int nIFO = this->ifoListSize();
918  if(nIFO<2) {cout << "gnetwork::GetDelay nIFO must be >= 2" << endl;return -1.;}
919 
920  if(ifo1.CompareTo(ifo2)==0) return 0.0;
921 
922  detector *pD[NIFO];
923  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
924 
925  int ifoId1=0;
926  int ifoId2=0;
927  int nifo=0;
928  for(int n=0; n<nIFO; n++) {
929  if(ifo1.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId1=n;}
930  if(ifo2.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId2=n;}
931  }
932  if(nifo<2) {cout << "gnetwork::GetDelay ifo not valid" << endl;return -1.;}
933 
934  return pD[ifoId1]->getTau(theta,phi)-pD[ifoId2]->getTau(theta,phi);
935 }
936 
937 //______________________________________________________________________________
938 void
939 gnetwork::DrawDelay(TString ifo1, TString ifo2, double phi, double theta, double frequency) {
940 
941  if(gSM.goff) return;
942 
943  int nIFO = this->ifoListSize();
944  if(nIFO<2) {cout << "gnetwork::DrawDelay nIFO must be >= 2" << endl;return;}
945 
946  if(ifo1.CompareTo(ifo2)==0) {cout << "gnetwork::DrawDelay ifo1 must != ifo2" << endl;return;}
947 
948  double ifo1d2_pos_delay=0.;
949  // if ifo1d2_pos_delay!=0 draw delay respect to ifo1d2_pos_delay
950  if(phi!=1000) ifo1d2_pos_delay=GetDelay(ifo1, ifo2, phi, theta);
951 
952  //double ifo1d2_max_delay=0.;
953  //if(frequency>0) ifo1d2_max_delay=GetDelay(ifo1, ifo2, "max");
954 
955  detector *pD[NIFO];
956  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
957 
958  int ifoId1=0;
959  int ifoId2=0;
960  int nifo=0;
961  for(int n=0; n<nIFO; n++) {
962  if(ifo1.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId1=n;}
963  if(ifo2.CompareTo(this->getifo(n)->Name)==0) {nifo++;ifoId2=n;}
964  }
965  if(nifo<2) {cout << "gnetwork::DrawDelay ifo not valid" << endl;return;}
966 
967  int size=int(180*2*gSM.resolution*360*2*gSM.resolution);
968  double* x = new double[size];
969  double* y = new double[size];
970  double* z = new double[size];
971  int cnt=0;
972  for (int i=0;i<360*2*gSM.resolution;i++) {
973  for (int j=0;j<180*2*gSM.resolution;j++) {
974  double phi = (double)(i+1)/(double)(2*gSM.resolution);
975  double theta = (double)(j+1)/(double)(2*gSM.resolution);
976 
977  double delay=pD[ifoId1]->getTau(theta,phi)-pD[ifoId2]->getTau(theta,phi);
978 
979  delay-=ifo1d2_pos_delay;
980  if(frequency>0) {
981  double fdelay=fabs(fmod(delay,0.5/frequency));
982  int idelay = int((fabs(delay)-fdelay)/(0.5/frequency));
983  delay = idelay%2==0 ? fdelay : (0.5/frequency)-fdelay;
984  }
985 
986  x[cnt]=phi;
987  y[cnt]=theta;
988  z[cnt]=delay;
989  cnt++;
990  }
991  }
992 
993  if (gSM.coordinate.CompareTo("GEOGRAPHIC")==0) {
994  for (int i=0;i<size;i++) CwbToGeographic(x[i],y[i],x[i],y[i]);
995  }
996  if (gSM.coordinate.CompareTo("CELESTIAL")==0) {
997  for (int i=0;i<size;i++) CwbToCelestial(x[i],y[i],x[i],y[i]);
998  }
999  gSM.FillData(size, x, y, z);
1000  gSM.Draw(0);
1001 
1002  delete [] x;
1003  delete [] y;
1004  delete [] z;
1005 }
1006 
1007 //______________________________________________________________________________
1008 int
1009 gnetwork::Delay2Coordinates(double t1, double t2, double t3, double& ph1, double& th1, double& ph2,double& th2) {
1010 //
1011 // Return the direction & mirror-direction uding the input arrival times of the 3 detectors
1012 // NOTE: Works only with 3 detectors
1013 // With 3 detector there are 2 possible solutions
1014 //
1015 // Input: t1,t2,t3 - arrival times of the 3 detectors (sec)
1016 //
1017 // Output: ph1,th1 - phi,theta coordinates of the first direction
1018 // ph2,th2 - phi,theta coordinates of the second direction
1019 //
1020 // Return 0/1 - success/error
1021 //
1022 
1023  int nIFO = this->ifoListSize();
1024  if(nIFO!=3) {cout << "gnetwork::Delay2Source nIFO must be = 3" << endl;return 1;}
1025 
1026  TString ifoName[3];
1027  for(int n=0; n<nIFO; n++) ifoName[n]=TString(this->getifo(n)->Name);
1028  //for(int n=0; n<nIFO; n++) cout << n << " " << ifoName[n].Data() << " " << ifoId[n] << endl;
1029  for(int n=0; n<nIFO; n++) for(int m=n+1; m<nIFO; m++)
1030  if(ifoName[n].CompareTo(ifoName[m])==0)
1031  {cout << "gnetwork::Delay2Source ifos must differents" << endl;return 2;}
1032 
1033  detector *pD[NIFO];
1034  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
1035 
1036  if(!gROOT->GetClass("XYZVector")) gSystem->Load("libMathCore");
1037 
1038  XYZVector vD[3];
1039  for(int n=0; n<nIFO; n++) vD[n].SetXYZ(pD[n]->Rv[0],pD[n]->Rv[1],pD[n]->Rv[2]);
1040  for(int n=0; n<nIFO; n++) vD[n]/=speedlight;
1041 
1042  double a=sqrt((vD[1]-vD[0]).Mag2());
1043  double b=sqrt((vD[2]-vD[1]).Mag2());
1044  double c=sqrt((vD[2]-vD[0]).Mag2());
1045 
1046  double s=(a+b+c)/2;
1047  double b2=(2/a)*sqrt(s*(s-a)*(s-b)*(s-c));
1048  double b1=sqrt(c*c-b2*b2);
1049 
1050  double delta=pow(b1*(t1-t2)-a*(t1-t3),2)+pow(b2*(t1-t2),2);
1051 
1052  double stheta = sqrt(delta)/(a*b2);
1053  double sphi = -b2*(t1-t2)/sqrt(delta);
1054 
1055  if (stheta>1) stheta=1;
1056  if (stheta<-1) stheta=-1;
1057 
1058  double theta = asin(stheta);
1059  double phi = asin(sphi);
1060 
1061  double ctheta = sqrt(pow(a*b2,2)-delta)/(a*b2);
1062  double cphi = (a*(t1-t3)-b1*(t1-t2))/sqrt(delta);
1063 
1064  if (TMath::IsNaN(ctheta)) ctheta=-1;
1065  if (ctheta>1) ctheta=-1;
1066  if (ctheta<-1) ctheta=-1;
1067 
1068  double thetacn = acos(-ctheta);
1069 
1070  XYZVector D1 = vD[1]-vD[0];
1071  XYZVector D3 = D1.Cross(vD[2]-vD[0]);
1072  XYZVector D2 = D1.Cross(D3);
1073 
1074  XYZVector eD1 = D1.Unit();
1075  XYZVector eD2 = D2.Unit();
1076  XYZVector eD3 = D3.Unit();
1077 
1078  Rotation3D R(eD1,eD2,eD3);
1079 
1080  //XYZVector IS;
1081  //XYZVector S;
1082  //S.SetR(1);
1083  TVector3 IS;
1084  TVector3 S(1,0,0);
1085 
1086  if (cphi<0) {
1087 
1088  S.SetTheta(theta);
1089  S.SetPhi((TMath::PiOver2()+phi));
1090  IS = R*S;
1091  th1=IS.Theta()*180./TMath::Pi();
1092  ph1=IS.Phi()*180./TMath::Pi();
1093  if (ph1<0) ph1+=360;
1094  if (ph1>360) ph1-=360;
1095 
1096  S.SetTheta(thetacn);
1097  S.SetPhi((TMath::PiOver2()+phi));
1098  IS = R*S;
1099  th2=IS.Theta()*180./TMath::Pi();
1100  ph2=IS.Phi()*180./TMath::Pi();
1101  if (ph2<0) ph2+=360;
1102  if (ph2>360) ph2-=360;
1103 
1104  } else {
1105 
1106  S.SetTheta(theta);
1107  S.SetPhi(-(TMath::PiOver2()+phi));
1108  IS = R*S;
1109  th1=IS.Theta()*180./TMath::Pi();
1110  ph1=IS.Phi()*180./TMath::Pi();
1111  if (ph1<0) ph1+=360;
1112  if (ph1>360) ph1-=360;
1113 
1114  S.SetTheta(thetacn);
1115  S.SetPhi(-(TMath::PiOver2()+phi));
1116  IS = R*S;
1117  th2=IS.Theta()*180./TMath::Pi();
1118  ph2=IS.Phi()*180./TMath::Pi();
1119  if (ph2<0) ph2+=360;
1120  if (ph2>360) ph2-=360;
1121  }
1122 
1123  return 0;
1124 }
1125 
1126 //______________________________________________________________________________
1127 void
1128 gnetwork::MakeNetworkDetectorIndex(double gamma, int loop, double snr) {
1129 //
1130 // Make Network Detector Index (1G analysis) and store it to the internal gSM sky map
1131 //
1132 // Input: gamma - 1G gamma value
1133 // loop - increase the statistic
1134 // snr - add noise to signal -> snr
1135 //
1136 
1137  if(gamma>1) {cout << "gamma must be <=1" << endl;exit(1);}
1138  if(loop<1) {cout << "gamma must be >=1" << endl;exit(1);}
1139 
1140  int nIFO = this->ifoListSize();
1141  if(nIFO==0) {cout << "detectors are no declared" << endl;exit(1);}
1142 
1143  detector *pD[NIFO];
1144  for(int n=0; n<nIFO; n++) pD[n] = this->getifo(n);
1145 
1146  std::vector<wavecomplex> F(nIFO);
1147  double* X = new double[nIFO];
1148  double* u = new double[nIFO];
1149  double* v = new double[nIFO];
1150 
1151  int L = gSM.size();
1152 
1153  double h2max=0.;
1154  gRandom->SetSeed(0);
1155  for (int l=0;l<L;l++) {
1156  double ndi=0.;
1157  for (int i=0;i<loop;i++) {
1158 
1159  double phi=gSM.getPhi(l);
1160  double theta=gSM.getTheta(l);
1161  double psi=gRandom->Uniform(0,180);
1162 
1163  double hp=gRandom->Uniform(-1,1);
1164  double hx=gRandom->Uniform(-1,1);
1165 
1166  double E=0.;
1167  for (int n=0;n<nIFO;n++) {
1168  F[n] = pD[n]->antenna(theta,phi,psi);
1169  X[n] = hp*F[n].real()+hx*F[n].imag();
1170  if(snr>=0.0) E+=X[n]*X[n];
1171  }
1172  // add noise to signal -> snr
1173  E=sqrt(E);
1174  if(snr>=0.0) for (int n=0;n<nIFO;n++) X[n] = X[n]*snr/E + gRandom->Gaus(0.,1.);
1175 
1176  double gp=0,gx=0,gI=0;
1177  for (int n=0;n<nIFO;n++) {
1178  gp+=F[n].real()*F[n].real();
1179  gx+=F[n].imag()*F[n].imag();
1180  gI+=F[n].real()*F[n].imag();
1181  }
1182  gp+=1.e-12;
1183  gx+=1.e-12;
1184 
1185  double Xp=0;
1186  double Xx=0;
1187  for (int k=0;k<nIFO;k++) {
1188  Xp+=X[k]*F[k].real();
1189  Xx+=X[k]*F[k].imag();
1190  }
1191 
1192  double uc = (Xp*gx - Xx*gI); // u cos of rotation to PCF
1193  double us = (Xx*gp - Xp*gI); // u sin of rotation to PCF
1194  double vc = (gp*uc + gI*us); // (u*f+)/|u|^2 - 'cos' for v
1195  double vs = (gx*us + gI*uc); // (u*fx)/|u|^2 - 'sin' for v
1196  double um=0;
1197  double vm=0;
1198  for (int k=0;k<nIFO;k++) {
1199  u[k]=0;
1200  v[k]=0;
1201  u[k]=F[k].real()*uc+F[k].imag()*us;
1202  um+=u[k]*u[k];
1203  v[k]=-F[k].imag()*vc+F[k].real()*vs;
1204  vm+=v[k]*v[k];
1205  }
1206  vm+=1.e-24;
1207  if(gamma>0) {
1208  double GAMMA = 1.-gamma*gamma; // network regulator
1209  ndi+=this->netx(u,um,v,vm,GAMMA);
1210  } else ndi+=this->ndi(u,um,nIFO);
1211 /*
1212  } else {
1213  double xndi=this->ndi(u,um,nIFO);
1214  if(xndi>=-gamma) ndi+=1;
1215  }
1216 */
1217  }
1218  double xndi=ndi/loop;
1219  if(h2max<xndi) h2max=xndi;
1220  gSM.set(l,xndi);
1221  }
1222 
1223  //if(gamma<0) h2->GetZaxis()->SetRangeUser(1,nIFO);
1224  if(gamma<0) gSM.h2->GetZaxis()->SetRangeUser(1,h2max);
1225 
1226  delete [] X;
1227  delete [] u;
1228  delete [] v;
1229 }
1230 
1231 //______________________________________________________________________________
1232 void
1233 gnetwork::DrawNetworkDetectorIndex(double gamma, int loop, double snr, int dpaletteId, bool btitle) {
1234 //
1235 // Draw Network Detector Index sky map (1G analysis)
1236 //
1237 // Input: gamma - 1G gamma value
1238 // loop - increase the statistic
1239 // snr - add noise to signal -> snr
1240 //
1241 // dpaletteId - palette ID
1242 // btitle - if true then a default title is provided
1243 //
1244 
1245  if(gSM.size()==0) {
1246  cout << "gnetwork::DrawNetworkDetectorIndex : Error - network skymap not initialized" << endl;
1247  exit(1);
1248  }
1249 
1250  MakeNetworkDetectorIndex(gamma, loop, snr);
1251 
1252  int nIFO = this->ifoListSize();
1253 
1254  TString netName="";
1255  for(int n=0; n<nIFO; n++) {
1256  TString ifoName=this->getifo(n)->Name;
1257  if(ifoName.CompareTo("A1")==0) ifoName=TString("#tilde{A}");
1258  else if(ifoName.CompareTo("H2")==0) ifoName=TString("#tilde{H}");
1259  else if(ifoName.CompareTo("I2")==0) ifoName=TString("#tilde{I}");
1260  else if(ifoName.CompareTo("R2")==0) ifoName=TString("#tilde{R}");
1261  else ifoName.Resize(1);
1262  netName+=ifoName;
1263  }
1264  if(btitle) gSM.h2->SetTitle("Network = "+netName+
1265  " "+
1266  "Network Detector Index");
1267 
1268  gSM.FillData();
1269  gSM.Draw(dpaletteId);
1270 }
1271 
1272 //______________________________________________________________________________
1273 void
1274 gnetwork::DrawCircles(double phi, double theta, Color_t lcolor, Width_t lwidth, Style_t lstyle, bool labels) {
1275 //
1276 // Draw iso delay circles for each couple of detectors along the direction phi,theta
1277 //
1278 // Input: phi,theta - sky direction coordinates (degrees)
1279 // lcolor - line color
1280 // lwidth - line width
1281 // lstyle - line style
1282 //
1283 
1284  DrawCircles(phi, theta, 0., lcolor, lwidth, lstyle, labels);
1285 }
1286 
1287 //______________________________________________________________________________
1288 void
1289 gnetwork::DrawCircles(double phi, double theta, double gps, Color_t lcolor, Width_t lwidth, Style_t lstyle, bool labels) {
1290 //
1291 // Draw iso delay circles for each couple of detectors along the direction phi,theta at time gps
1292 //
1293 // Input: phi,theta - sky direction coordinates (degrees)
1294 // gps - gps time (sec)
1295 // lcolor - line color
1296 // lwidth - line width
1297 // lstyle - line style
1298 //
1299 
1300  if(gSM.goff) return;
1301  if(gSM.canvas==NULL) return;
1302 
1303  int n,m,l;
1304 
1305  int nIFO = this->ifoListSize();
1306  if(nIFO<2) return;
1307 
1308  detector *pD[NIFO];
1309  for(n=0; n<nIFO; n++) pD[n] = this->getifo(n);
1310 
1311  gSM.canvas->cd();
1312 
1313  double phi_off = gps>0 ? gSM.phi2RA(phi, gps) : phi;
1314  phi_off=fmod(phi_off-phi+360,360);
1315 
1316  double ph1,th1,ph2,th2;
1317  double r2d = 180./TMath::Pi();
1318  double d2r = TMath::Pi()/180.;
1319 //Color_t colors[6] = {kBlue,kRed,kYellow,kGreen,kBlack,kWhite};
1320  Style_t marker[3] = {20,21,22}; // circle,square,triangle
1321  int ncircle=0;
1322  for(n=0;n<nIFO;n++) {
1323  for(m=n+1;m<nIFO;m++) {
1324 
1325  XYZVector D1(pD[n]->Rv[0],pD[n]->Rv[1],pD[n]->Rv[2]);
1326  XYZVector D2(pD[m]->Rv[0],pD[m]->Rv[1],pD[m]->Rv[2]);
1327 
1328  D1=D1/speedlight;
1329  D2=D2/speedlight;
1330 
1331  XYZVector D12 = D1-D2;
1332  if(D12.Mag2()==0) continue; // the 2 detectors are located in the same site
1333  TVector3 vD12(D12.X(),D12.Y(),D12.Z());
1334 
1335  TVector3 vSD(1,0,0);
1336  // coordinates are transformed in the CWB frame
1337  if (gSM.coordinate.CompareTo("GEOGRAPHIC")==0) {
1338  GeographicToCwb(phi,theta,ph1,th1);
1339  vSD.SetTheta(th1*d2r);
1340  vSD.SetPhi(ph1*d2r);
1341  } else if (gSM.coordinate.CompareTo("CELESTIAL")==0) {
1342  CelestialToCwb(phi,theta,ph1,th1);
1343  ph1=fmod(ph1-phi_off+360,360);
1344  vSD.SetTheta(th1*d2r);
1345  vSD.SetPhi(ph1*d2r);
1346  } else {
1347  vSD.SetTheta(theta*d2r);
1348  vSD.SetPhi(phi*d2r);
1349  }
1350 
1351  TVector3 rvSD;
1352  double as;
1353  int L=4*360;
1354  wavearray<float> th(L);
1355  wavearray<float> ph(L);
1356  for (l=0;l<L;l++) {
1357  as = l*d2r/4.;
1358  TRotation vR;
1359  vR.Rotate(as,vD12);
1360  rvSD = vR*vSD;
1361  th[l] = rvSD.Theta()*r2d;
1362  ph[l] = rvSD.Phi()*r2d+phi_off;
1363  ph[l]=fmod(ph[l]+360,360);
1364  }
1365 
1366  wavearray<double> x(L);
1367  wavearray<double> y(L);
1368  int P=0;
1369  for (l=0;l<L;l++) {
1370  // the main circle is shown in sub-circles to avoid to draw fake lines from ph=0 to 360
1371  // if((fabs(ph[l]-ph[l==0?L-1:l-1])<180)&&(l<L-1)) {
1372  // TVector3 is in radians [phi=-Pi,+Pi] [theta=0,Pi]
1373  ph1=ph[l];
1374  ph2=ph[l==0?L-1:l-1];
1375  if (gSM.coordinate=="GEOGRAPHIC") {
1376  CwbToGeographic(ph1,th1,ph1,th1);
1377  CwbToGeographic(ph2,th2,ph2,th2);
1378  }
1379  if((fabs(ph1-ph2)<180)&&(l<L-1)) {
1380  x[P]=ph[l]; y[P]=th[l];
1381 
1382  double X=x[P];
1383  if (gSM.projection=="hammer") {
1384  if(X>180) {x[P]+=360; y[P]-=90;}
1385  else {y[P]-=90;}
1386  gSM.ProjectHammer(x[P], y[P], x[P], y[P]);
1387  y[P]=y[P]+90;
1388  } else if (gSM.projection=="sinusoidal") {
1389  if(X>180) {x[P]-=360; y[P]+=90;}
1390  else {y[P]-=90;}
1391  gSM.ProjectSinusoidal(x[P], y[P], x[P], y[P]);
1392  if(X>180) {x[P]=-x[P]; y[P]-=90;}
1393  else {y[P]+=90;}
1394  } else if (gSM.projection=="parabolic") {
1395  if(X>180) {x[P]-=360; y[P]-=90;}
1396  else {y[P]-=90;}
1397  gSM.ProjectParabolic(x[P], y[P], x[P], y[P]);
1398  y[P]+=90;
1399  } else {
1400 // x[P]=x[P]-180; y[P]=y[P]-90;
1401  }
1402  if (gSM.coordinate=="GEOGRAPHIC") CwbToGeographic(x[P],y[P],x[P],y[P]);
1403  if (gSM.coordinate=="CELESTIAL") CwbToCelestial(x[P],y[P],x[P],y[P]);
1404  P++;
1405  } else {
1406  TPolyLine *pL = new TPolyLine(P,x.data,y.data);
1407  //pL->SetLineColor(colors[ncircle%6]);
1408  pL->SetLineColor(lcolor);
1409  pL->SetLineStyle(lstyle);
1410  pL->SetLineWidth(lwidth);
1411  pL->Draw();
1412  circleL.push_back(pL);
1413  P=0;
1414  }
1415  }
1416  // draw detector axis
1417  double vph,vth;
1418  if(vSD.Dot(vD12)>0) {
1419  vph=D12.Phi()*r2d;
1420  vth=D12.Theta()*r2d;
1421  } else {
1422  vph=(-D12).Phi()*r2d;
1423  vth=(-D12).Theta()*r2d;
1424  }
1425  vph=fmod(vph+phi_off+360,360);
1426 
1427  double X=vph;
1428  if (gSM.projection=="hammer") {
1429  if(X>180) {vph+=360; vth-=90;}
1430  if(X<180) {vth-=90;}
1431  gSM.ProjectHammer(vph, vth, vph, vth);
1432  vth=vth+90;
1433  } else if (gSM.projection=="sinusoidal") {
1434  if(X>180) {vph-=360; vth+=90;}
1435  else {vth-=90;}
1436  gSM.ProjectSinusoidal(vph, vth, vph, vth);
1437  if(X>180) {vph=-vph; vth-=90;}
1438  else {vth+=90;}
1439  } else if (gSM.projection=="parabolic") {
1440  if(X>180) {vph-=360; vth-=90;}
1441  else {vth-=90;}
1442  gSM.ProjectParabolic(vph, vth, vph, vth);
1443  vth+=90;
1444  } else {
1445 // vph=vph-180; vth=vth-90;
1446  }
1447  if (gSM.coordinate=="GEOGRAPHIC") CwbToGeographic(vph,vth,vph,vth);
1448  if (gSM.coordinate=="CELESTIAL") CwbToCelestial(vph,vth,vph,vth);
1449  //if (coordinate.CompareTo("CWB")==0) {vph+=180;vth+=90;}
1450  TMarker* pM = new TMarker(vph,vth,marker[ncircle%3]);
1451  pM->SetMarkerSize(1.0);
1452  pM->SetMarkerColor(kBlack);
1453  if(labels) pM->Draw();
1454  daxisM.push_back(pM);
1455 
1456  // draw ifos pair labels on the detector axis
1457  char ifoPair[16];sprintf(ifoPair,"%c%c",pD[n]->Name[0],pD[m]->Name[0]);
1458  TLatex *pL = new TLatex(vph+4.0,vth,ifoPair);
1459  pL->SetTextFont(32);
1460  pL->SetTextSize(0.045);
1461  pL->SetTextColor(kBlack);
1462  if(labels) pL->Draw();
1463  siteP.push_back(pL);
1464 
1465  ncircle++;
1466  }
1467  }
1468 
1469  return;
1470 }
1471 
1472 //______________________________________________________________________________
1473 void gnetwork::DumpObject(char* file)
1474 {
1475 //
1476 // dump gnetwork object into root file
1477 //
1478 // Input: file - root file name
1479 //
1480 
1481  TFile* rfile = new TFile(file, "RECREATE");
1482  this->Write("gnetwork"); // write this object
1483  rfile->Close();
1484 }
1485 
1486 //______________________________________________________________________________
1487 void gnetwork::LoadObject(char* file)
1488 {
1489 //
1490 // load gnetwork object from root file
1491 //
1492 // Input: file - root file name
1493 //
1494 
1495  TFile* rfile = new TFile(file);
1496  gnetwork* gnet = (gnetwork*)rfile->Get("gnetwork;1");
1497  if(gnet==NULL) {
1498  cout << "gnetwork::LoadObject : Error - input file don't contains object gnetwork" << endl;
1499  exit(1);
1500  }
1501  *this=*gnet;
1502 
1503  siteM.clear();
1504  siteA.clear();
1505  siteT.clear();
1506  siteL.clear();
1507  siteP.clear();
1508  circleL.clear();
1509  daxisM.clear();
1510 
1511  gSM.FillData();
1512  rfile->Close();
1513 }
1514 
1515 //______________________________________________________________________________
1516 TCanvas* gnetwork::DrawPolargram(int ptype, network* net)
1517 {
1518 //
1519 // Plot polargrams of the event pixels in the network plane
1520 //
1521 // Input: ptype - plot type [0/1]
1522 // 0 : projection on network plane of 00/90 amplitudes
1523 // 1 : standard response of 00/90 amplitudes
1524 // net - pointer to network object
1525 // if net!=NULL then net is used to retrive the data pixels to be plotted
1526 //
1527 // Output: return canvas object
1528 //
1529 
1530  if(ptype!=0 && ptype!=1) {
1531  cout << "gnetwork::DrawPolargram : Error - wrong input parameter, must be [0/1]" << endl;
1532  exit(1);
1533  }
1534 
1535  network* NET = net!=NULL ? net : this;
1536 
1537  // get size
1538  int size = ptype==0 ? NET->p00_POL[0].size() : NET->r00_POL[0].size();
1539  if(size==0) {
1540  cout << "gnetwork::DrawPolargram : Warning - event pixels size=0" << endl;
1541  return NULL;
1542  }
1543 
1544  // find max radius
1545  double rmax=0;
1546  if(ptype==0) {
1547  for(int i=0;i<size;i++) if(NET->p00_POL[0][i]>rmax) rmax=NET->p00_POL[0][i];
1548  for(int i=0;i<size;i++) if(NET->p90_POL[0][i]>rmax) rmax=NET->p90_POL[0][i];
1549  }
1550  if(ptype==1) {
1551  for(int i=0;i<size;i++) if(NET->r00_POL[0][i]>rmax) rmax=NET->r00_POL[0][i];
1552  for(int i=0;i<size;i++) if(NET->r90_POL[0][i]>rmax) rmax=NET->r90_POL[0][i];
1553  }
1554  rmax=TMath::Nint(1.1*rmax+0.5);
1555 
1556  for(int n=0;n<2;n++) if(grP[n]!=NULL) delete grP[n];
1557  if(canvasP!=NULL) delete canvasP;
1558  canvasP = new TCanvas("canvasP","polargram",800,800);
1559  gStyle->SetLineColor(kBlack);
1560  double* rpol;
1561  double* apol;
1562  for(int n=0;n<2;n++) {
1563  // projection on network plane 00 amplitudes
1564  if(ptype==0 && n==0) {rpol = NET->p00_POL[0].data; apol = NET->p00_POL[1].data;}
1565  // projection on network plane 90 amplitudes
1566  if(ptype==0 && n==1) {rpol = NET->p90_POL[0].data; apol = NET->p90_POL[1].data;}
1567  // standard response 00 amplitudes
1568  if(ptype==1 && n==0) {rpol = NET->r00_POL[0].data; apol = NET->r00_POL[1].data;}
1569  // standard response 90 amplitudes
1570  if(ptype==1 && n==1) {rpol = NET->r90_POL[0].data; apol = NET->r90_POL[1].data;}
1571  grP[n] = new TGraphPolar(size,apol,rpol);
1572 
1573  grP[n]->SetMarkerStyle(20);
1574  grP[n]->SetMarkerSize(0.8);
1575  grP[n]->SetTitle(" ");
1576 
1577  if(n==0) {
1578  grP[n]->SetMarkerColor(kBlue);
1579  grP[n]->SetLineColor(kBlue);
1580  grP[n]->Draw("P");
1581  }
1582 
1583  if(n==1) {
1584  grP[n]->SetMarkerColor(kRed);
1585  grP[n]->SetLineColor(kRed);
1586  grP[n]->Draw("PSAME");
1587  }
1588 
1589  canvasP->Update();
1590  grP[n]->GetPolargram()->SetTextColor(8);
1591  grP[n]->GetPolargram()->SetRangePolar(-TMath::Pi(),TMath::Pi());
1592  grP[n]->SetMinRadial(0);
1593  grP[n]->SetMaxRadial(rmax);
1594  grP[n]->GetPolargram()->SetNdivPolar(612);
1595  grP[n]->GetPolargram()->SetNdivRadial(504);
1596  grP[n]->GetPolargram()->SetToDegree();
1597  grP[n]->GetPolargram()->SetPolarLabelSize(0.03);
1598  grP[n]->GetPolargram()->SetRadialLabelSize(0.03);
1599  }
1600 
1601  canvasP->Update();
1602 
1603  return canvasP;
1604 }
1605 
network::ifoName
std::vector< char * > ifoName
Definition: network.hh:591
F1
#define F1
Definition: Regression_H1.C:13
gnetwork::DrawPolargram
TCanvas * DrawPolargram(int ptype, network *net=NULL)
Definition: gnetwork.cc:1516
network::veto
wavearray< short > veto
Definition: network.hh:626
network::getifo
detector * getifo(size_t n)
param: detector index
Definition: network.hh:418
gnetwork::circleL
std::vector< TPolyLine * > circleL
Definition: gnetwork.hh:99
detector::getDetectorParams
detectorParams getDetectorParams()
Definition: detector.cc:201
wavearray::size
virtual size_t size() const
Definition: wavearray.hh:127
NIFO
#define NIFO
Definition: wat.hh:56
wavecomplex::imag
double imag() const
Definition: wavecomplex.hh:52
network::skyHole
wavearray< double > skyHole
Definition: network.hh:625
Init
void Init()
Definition: ChirpMass.C:280
network::add
size_t add(detector *)
param: detector structure return number of detectors in the network
Definition: network.cc:2528
gnetwork::DrawSitesArms
void DrawSitesArms(double mlength=600000., Color_t lcolor=kBlack, Size_t lwidth=1.0, Style_t lstyle=1)
Definition: gnetwork.cc:295
network::p90_POL
wavearray< double > p90_POL[2]
buffer for projection on network plane 00 ampl
Definition: network.hh:642
gnetwork::DrawAntennaPattern
void DrawAntennaPattern(int polarization=-1, int dpaletteId=0, bool btitle=true, int order=6)
Definition: gnetwork.cc:655
gnetwork::DrawDelay
void DrawDelay(TString ifo1, TString ifo2, double phi=1000., double theta=1000., double frequency=0.)
Definition: gnetwork.cc:939
gnetwork::canvasP
TCanvas * canvasP
Definition: gnetwork.hh:104
CwbToCelestial
void CwbToCelestial(double ilongitude, double ilatitude, double &olongitude, double &olatitude, double gps=0)
Definition: skycoord.hh:410
gnetwork::DrawCircles
void DrawCircles(double phi, double theta, double gps, Color_t lcolor=kBlack, Width_t lwidth=1, Style_t lstyle=1, bool labels=false)
Definition: gnetwork.cc:1289
ptype
TString ptype[nPLOT]
Definition: cwb_mkfad.C:94
a
wavearray< double > a(hp.size())
v
WSeries< float > v[nIFO]
Definition: cwb_net.C:62
n
int n
Definition: cwb_net.C:10
gnetwork::GetSite
double GetSite(TString ifo, TString type="")
Definition: gnetwork.cc:167
network::nSkyStat
skymap nSkyStat
Definition: network.hh:610
z
wavearray< double > z
Definition: Test10.C:32
const_cast< char * >
cout<< "skymap size : "<< L<< endl;for(int l=0;l< L;l++) sm.set(l, l);sm > const_cast< char * >("skymap.dat")
gnetwork::pD
detector * pD[NIFO_MAX]
Definition: gnetwork.hh:102
TString
TString("c")
Definition: cwb_report_skymap.C:111
gnetwork::~gnetwork
virtual ~gnetwork()
Definition: gnetwork.cc:76
gskymap::set
void set(size_t i, double a)
Definition: gskymap.hh:110
gnetwork::Delay2Coordinates
int Delay2Coordinates(double t1, double t2, double t3, double &ph1, double &th1, double &ph2, double &th2)
Definition: gnetwork.cc:1009
network::nPenalty
skymap nPenalty
Definition: network.hh:606
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
frequency
double frequency
Definition: DrawClusterWithSkyplot.C:20
theta
float theta
Definition: cwb_condor_create_ced.C:51
gskymap::coordinate
TString coordinate
Definition: gskymap.hh:211
ph
TH2F * ph
Definition: ReadHistFromJobFile.C:17
watconstants::EarthEquatorialRadius
double EarthEquatorialRadius()
Definition: constants.hh:168
P
TRandom3 P
Definition: compare_bkg.C:296
skymap::min
double min()
Definition: skymap.cc:441
wavecomplex::real
double real() const
Definition: wavecomplex.hh:51
gnetwork::DrawSitesLabel
void DrawSitesLabel(Color_t tcolor=kBlack, Size_t tsize=0.045, Font_t tfont=32)
Definition: gnetwork.cc:423
network::p00_POL
wavearray< double > p00_POL[2]
buffers for cluster MRA energy
Definition: network.hh:641
std
STL namespace.
skymap::getTheta
double getTheta(size_t i)
Definition: skymap.hh:206
network::nPolarisation
skymap nPolarisation
Definition: network.hh:612
gnetwork::DumpObject
void DumpObject(char *file)
Definition: gnetwork.cc:1473
polarization
int polarization
Definition: cwb_draw_antpat.C:23
size
Long_t size
Definition: cwb_condor_check.C:48
hp
wavearray< double > hp
Definition: DrawInspiral.C:43
gskymap::canvas
TCanvas * canvas
Definition: gskymap.hh:194
gnetwork::siteP
std::vector< TLatex * > siteP
Definition: gnetwork.hh:98
m
int m
Definition: cwb_net.C:10
skymap::max
double max()
Definition: skymap.cc:420
j
int j
Definition: cwb_net.C:10
i
i drho i
Definition: cwb_epparameters.C:85
detectorParams::longitude
double longitude
Definition: detector.hh:34
detector::tau
skymap tau
Definition: detector.hh:328
gskymap::Draw
void Draw(int dpaletteId=0, Option_t *option="colfz")
Definition: gskymap.cc:442
network::ifoList
std::vector< detector * > ifoList
Definition: network.hh:590
gskymap::h2
TH2D * h2
`
Definition: gskymap.hh:195
gnetwork::siteM
std::vector< TMarker * > siteM
Definition: gnetwork.hh:94
network::skyENRG
wavearray< double > skyENRG
Definition: network.hh:628
network::nCorrEnergy
skymap nCorrEnergy
Definition: network.hh:607
ifo
char ifo[NIFO_MAX][8]
Definition: cwb1G_parameters.C:17
net
network ** net
NOISE_MDC_SIMULATION.
Definition: CWB_Plugin_HEN_BKG_Config.C:4
network::ifoListSize
size_t ifoListSize()
Definition: network.hh:413
network::netx
static int netx(double *, double, double *, double, double)
Definition: network.hh:812
network::nLikelihood
skymap nLikelihood
Definition: network.hh:604
mode
size_t mode
Definition: cwb1G_parameters.C:134
nIFO
#define nIFO
Definition: DrawSensitivitiesS5.C:11
skymap::getSkyIndex
size_t getSkyIndex(double th, double ph)
param: theta param: phi
Definition: skymap.cc:702
gnetwork::DrawNetworkDetectorIndex
void DrawNetworkDetectorIndex(double gamma, int loop, double snr=-1.0, int dpaletteId=DUMMY_PALETTE_ID, bool btitle=true)
Definition: gnetwork.cc:1233
gskymap::FillData
void FillData(int size, double *phi, double *theta, double *binc)
Definition: gskymap.cc:376
gnetwork::DrawSitesShortLabel
void DrawSitesShortLabel(Color_t tcolor=kBlack, Size_t tsize=0.052, Font_t tfont=32)
Definition: gnetwork.cc:407
Write
tlive_fix Write()
phi
float phi
Definition: cwb_condor_create_ced.C:51
detector::antenna
wavecomplex antenna(double, double, double=0.)
param: source theta,phi, polarization angle psi in degrees
Definition: detector.cc:473
psi
float psi
Definition: cwb_condor_create_ced.C:51
gr
TGraph * gr
Definition: cwb_draw_sensitivity.C:4
hx
wavearray< double > hx
Definition: DrawInspiral.C:44
D1
detector D1
Definition: TestReadWriteDetectorObject.C:9
int
i() int(T_cor *100))
NET
network NET
Definition: cwb_dump_inj.C:12
Pi
double Pi
Definition: DrawPhaseShift.C:12
gx
gwavearray< double > * gx
Definition: CompareWDMvsTime.C:18
NIFO_MAX
const int NIFO_MAX
Definition: wat.hh:4
gnetwork::DrawSites
void DrawSites(Color_t mcolor=kBlack, Size_t msize=2.0, Style_t mstyle=20)
Definition: gnetwork.cc:238
gnetwork::gnetwork
gnetwork()
Definition: gnetwork.hh:21
gnetwork::siteLabelType
unsigned int siteLabelType
Definition: gnetwork.hh:92
speedlight
#define speedlight
Definition: watfun.hh:15
k
int k
Definition: cwb_report_prod_2.C:150
gnetwork::ndi
static double ndi(double *, double, int)
Definition: gnetwork.hh:110
add
r add(tfmap,"hchannel")
gnetwork::LoadObject
void LoadObject(char *file)
Definition: gnetwork.cc:1487
network::index
wavearray< int > index
Definition: network.hh:622
F
double F
Definition: DrawPhaseShift.C:14
skymap
Definition: skymap.hh:45
detectorParams::latitude
double latitude
Definition: detector.hh:33
vR
std::vector< double > vR
Definition: cbc_plots_sim4_ifar.C:1479
e
double e
Definition: CreateListEBBH.C:35
network::skyProb
wavearray< double > skyProb
Definition: network.hh:627
skymap::phi2RA
double phi2RA(double ph, double gps)
Definition: skymap.hh:194
ind
int * ind
Definition: cwb_fix_live_slag_missed.C:89
gnetwork::gSM
gskymap gSM
Definition: gnetwork.hh:90
gnetwork::GetDelay
double GetDelay(TString ifo1, TString ifo2, double phi, double theta)
Definition: gnetwork.cc:907
s
s s
Definition: cwb_net.C:137
gskymap::resolution
double resolution
Definition: gskymap.hh:206
nifo
int nifo
Definition: cwb_dump_events.C:32
skymap::getPhi
double getPhi(size_t i)
Definition: skymap.hh:146
gnetwork::siteT
std::vector< TText * > siteT
Definition: gnetwork.hh:96
CelestialToCwb
void CelestialToCwb(double ilongitude, double ilatitude, double &olongitude, double &olatitude, double gps=0)
Definition: skycoord.hh:419
gps
double gps
Definition: cwb_report_prod_2.C:194
network::nAlignment
skymap nAlignment
Definition: network.hh:602
network::print
void print()
Definition: network.cc:8150
gnetwork::MakeNetworkDetectorIndex
void MakeNetworkDetectorIndex(double gamma, int loop, double snr=-1.0)
Definition: gnetwork.cc:1128
gskymap::ProjectParabolic
void ProjectParabolic(Double_t l, Double_t b, Double_t &Al, Double_t &Ab)
Definition: gskymap.cc:875
network::nDisbalance
skymap nDisbalance
Definition: network.hh:609
detector::Name
char Name[16]
Definition: detector.hh:309
fabs
double fabs(const Complex &x)
Definition: numpy.cc:37
cwb_online.file
string file
Definition: cwb_online.py:385
gskymap::ProjectHammer
void ProjectHammer(Double_t l, Double_t b, Double_t &Al, Double_t &Ab)
Definition: gskymap.cc:818
network::gamma
double gamma
Definition: network.hh:574
CwbToGeographic
void CwbToGeographic(double ilongitude, double ilatitude, double &olongitude, double &olatitude)
Definition: skycoord.hh:396
cnt
int cnt
Definition: cwb_condor_cleanup.C:19
network::skyMaskCC
wavearray< double > skyMaskCC
Definition: network.hh:624
S
Meyer< double > S(1024, 2)
l
int l
Definition: cbc_plots.C:434
network::nNetIndex
skymap nNetIndex
Definition: network.hh:608
sprintf
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
network::nSensitivity
skymap nSensitivity
list of wdm tranformations
Definition: network.hh:601
detector::getTau
double getTau(double, double)
param: source theta,phi angles in degrees
Definition: detector.cc:681
gnetwork::ClearSites
void ClearSites()
Definition: gnetwork.cc:505
gnetwork::siteA
std::vector< TPolyLine * > siteA
Definition: gnetwork.hh:95
gnetwork::ClearSitesLabel
void ClearSitesLabel()
Definition: gnetwork.cc:529
gnetwork::GetAntennaPattern
double GetAntennaPattern(double phi, double theta, double psi=0., int polarization=1)
Definition: gnetwork.cc:545
btitle
bool btitle
Definition: DrawGnetwork2.C:16
gskymap::ProjectSinusoidal
void ProjectSinusoidal(Double_t l, Double_t b, Double_t &Al, Double_t &Ab)
Definition: gskymap.cc:857
wavearray::data
DataType_t * data
Definition: wavearray.hh:301
gnetwork::Draw
void Draw(char *smName, network *net=NULL)
Definition: gnetwork.hh:59
gnetwork::ClearSitesArms
void ClearSitesArms()
Definition: gnetwork.cc:517
t1
TH1 * t1
Definition: cbc_plots_sim4_ifar.C:1896
SetXYZ
xyz SetXYZ(xgc, ygc, zgc)
gnetwork::ClearCircles
void ClearCircles()
Definition: gnetwork.cc:491
gnetwork::grP
TGraphPolar * grP[2]
canvas used for the polargrams
Definition: gnetwork.hh:105
network::nNullEnergy
skymap nNullEnergy
Definition: network.hh:605
network::skyMask
wavearray< short > skyMask
Definition: network.hh:623
GeographicToCwb
void GeographicToCwb(double ilongitude, double ilatitude, double &olongitude, double &olatitude)
Definition: skycoord.hh:403
network::delay
void delay(double theta, double phi)
Definition: network.cc:7898
snr
snr * snr
Definition: ComputeSNR.C:71
network::nCorrelation
skymap nCorrelation
Definition: network.hh:603
gskymap::goff
bool goff
Definition: gskymap.hh:207
detParms
detectorParams detParms[4]
Definition: DrawInspiralUserDet.C:19
skymap::size
size_t size()
Definition: skymap.hh:118
gnetwork::siteL
std::vector< TLatex * > siteL
Definition: gnetwork.hh:97
y
wavearray< double > y
Definition: Test10.C:31
gnetwork::Init
void Init()
Definition: gnetwork.cc:93
network::delta
double delta
Definition: network.hh:573
gskymap::projection
TString projection
Definition: gskymap.hh:212
gnetwork::daxisM
std::vector< TMarker * > daxisM
Definition: gnetwork.hh:100
network::r90_POL
wavearray< double > r90_POL[2]
buffer for standard response 00 ampl
Definition: network.hh:644
network::nProbability
skymap nProbability
Definition: network.hh:613
pD
detector ** pD
Definition: cwb_dump_events.C:34
exit
exit(0)
GAMMA
#define GAMMA(TYPE)
Definition: xroot.hh:5
network::r00_POL
wavearray< double > r00_POL[2]
buffer for projection on network plane 90 ampl
Definition: network.hh:643
network::nEllipticity
skymap nEllipticity
Definition: network.hh:611