Logo coherent WaveBurst  
Library Reference Guide
Logo
gnetwork.hh
Go to the documentation of this file.
1 /*
2 # Copyright (C) 2019 Gabriele Vedovato
3 #
4 # This program is free software: you can redistribute it and/or modify
5 # it under the terms of the GNU General Public License as published by
6 # the Free Software Foundation, either version 3 of the License, or
7 # (at your option) any later version.
8 #
9 # This program is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 # GNU General Public License for more details.
13 #
14 # You should have received a copy of the GNU General Public License
15 # along with this program. If not, see <https://www.gnu.org/licenses/>.
16 */
17 
18 
19 /**********************************************************
20  * Package: graphic wat Class Library
21  * File name: gnetwork.hh
22  * Author: Gabriele Vedovato (vedovato@lnl.infn.it)
23  **********************************************************/
24 
25 
26 #ifndef GNETWORK_HH
27 #define GNETWORK_HH
28 
29 #include "gskymap.hh"
30 #include "network.hh"
31 #include "TGraphPolar.h"
32 
33 class gnetwork : public network
34 {
35 
36 public:
37 
38  // Constructors
39  gnetwork() : network() {Init();}
42  virtual ~gnetwork();
43 
44  gskymap* GetGskymap() {return &gSM;}
45  void SetGskymap(gskymap& gSM) {this->gSM = gSM;}
46 
47  void setSkyMaps(double sms,double t1,double t2,double p1,double p2) {
48  network::setSkyMaps(sms,t1,t2,p1,p2);gSM=gskymap(sms,t1,t2,p1,p2);}
49  void setSkyMaps(int healpix_order) {
50  network::setSkyMaps(int(healpix_order));gSM=gskymap(int(healpix_order));}
51 
52  void DrawAntennaPattern(int polarization=-1, int dpaletteId=0, bool btitle=true, int order=6);
53  void DrawSites(Color_t mcolor=kBlack, Size_t msize=2.0, Style_t mstyle=20);
54  void DrawSitesArms(double mlength=600000., Color_t lcolor=kBlack, Size_t lwidth=1.0, Style_t lstyle=1);
55  void DrawSitesShortLabel(Color_t tcolor=kBlack, Size_t tsize=0.052, Font_t tfont=32);
56  void DrawSitesLabel(Color_t tcolor=kBlack, Size_t tsize=0.045, Font_t tfont=32);
57  void DrawDelay(TString ifo1, TString ifo2,
58  double phi=1000., double theta=1000., double frequency=0.);
59  double GetDelay(TString ifo1, TString ifo2, double phi, double theta);
60  double GetDelay(TString ifo1, TString ifo2, TString mode="");
61  double GetAntennaPattern(double phi, double theta, double psi=0., int polarization=1);
62  double GetAntennaPattern(TString ifo, double phi, double theta, double psi=0., bool plus=true);
63 
64  int Delay2Coordinates(double t1, double t2, double t3, double& ph1, double& th1, double& ph2,double& th2);
65 
66  double GetSite(TString ifo, TString type="");
67 
68  void MakeNetworkDetectorIndex(double gamma, int loop, double snr=-1.0);
69  void DrawNetworkDetectorIndex(double gamma, int loop, double snr=-1.0,
70  int dpaletteId=DUMMY_PALETTE_ID, bool btitle=true);
71 
72  void DrawCircles(double phi, double theta, double gps,
73  Color_t lcolor=kBlack, Width_t lwidth=1, Style_t lstyle=1, bool labels=false);
74  void DrawCircles(double phi, double theta,
75  Color_t lcolor=kBlack, Width_t lwidth=1, Style_t lstyle=1, bool labels=false);
76 
77  void Draw(char* smName, network* net=NULL) {Draw(TString(smName),net);}
78  void Draw(TString smName, network* net=NULL);
79 
80  void SetScale(TString ifo, double scale);
81  double GetScale(TString ifo);
82 
83  void ClearSites();
84  void ClearSitesArms();
85  void ClearSitesLabel();
86  void ClearCircles();
87 
88  //: dump gnetwork into root file
89  void DumpObject(char* file);
90  //: load gnetwork object from root file
91  void LoadObject(char* file);
92 
93  static inline double ndi(double*, double, int);
94 
95  TCanvas* DrawPolargram(int ptype, network* net=NULL);
96 
97 private:
98 
99  void Init();
100  void Draw(skymap sm, TString title="") {gSM=sm;gSM.SetTitle(TString("network ")+title);gSM.Draw();}
102  skymap sm=nSensitivity;for(int l=0;l<(int)sm.size();l++) sm.set(l,double(w[l]));
103  Draw(sm,title);return;}
105  skymap sm=nSensitivity;for(int l=0;l<(int)sm.size();l++) sm.set(l,double(w[l]));
106  Draw(sm,title);return;}
108  skymap sm=nSensitivity;for(int l=0;l<(int)sm.size();l++) sm.set(l,double(w[l]));
109  Draw(sm,title);return;}
110 
112 
113  double scale[NIFO_MAX]; //! default=1, used to scale the antenna pattern acccording the PSD contribution of each detector (at a fixed frequency) to the network
114 
115  unsigned int siteLabelType;
116 
117  std::vector<TMarker*> siteM; //!
118  std::vector<TPolyLine*> siteA; //!
119  std::vector<TText*> siteT; //!
120  std::vector<TLatex*> siteL; //!
121  std::vector<TLatex*> siteP; //!
122  std::vector<TPolyLine*> circleL; //!
123  std::vector<TMarker*> daxisM; //!
124 
126 
127  TCanvas* canvasP; //! canvas used for the polargrams
128  TGraphPolar* grP[2]; //! polargraphs used for the polargrams
129 
130  ClassDef(gnetwork,3)
131 };
132 
133 inline double gnetwork::ndi(double* u, double um, int nIFO) {
134  double d = 0.;
135  for(int i=0;i<nIFO;i++) d+= (u[i]*u[i]*u[i]*u[i]);
136  d = um > 0. ? 1./(d/(um*um)) : 0.;
137  return d;
138 }
139 
140 #endif
TCanvas * DrawPolargram(int ptype, network *net=NULL)
Definition: gnetwork.cc:1535
std::vector< TPolyLine * > circleL
Definition: gnetwork.hh:122
double GetScale(TString ifo)
Definition: gnetwork.cc:1654
TH1 * t1
void DrawSitesArms(double mlength=600000., Color_t lcolor=kBlack, Size_t lwidth=1.0, Style_t lstyle=1)
Definition: gnetwork.cc:314
void DrawAntennaPattern(int polarization=-1, int dpaletteId=0, bool btitle=true, int order=6)
Definition: gnetwork.cc:674
void DrawDelay(TString ifo1, TString ifo2, double phi=1000., double theta=1000., double frequency=0.)
Definition: gnetwork.cc:958
TCanvas * canvasP
Definition: gnetwork.hh:127
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:1308
TString ptype[nPLOT]
Definition: cwb_mkfad.C:112
void SetGskymap(gskymap &gSM)
Definition: gnetwork.hh:45
double GetSite(TString ifo, TString type="")
Definition: gnetwork.cc:186
detector * pD[NIFO_MAX]
Definition: gnetwork.hh:125
TString("c")
virtual ~gnetwork()
Definition: gnetwork.cc:94
int Delay2Coordinates(double t1, double t2, double t3, double &ph1, double &th1, double &ph2, double &th2)
Definition: gnetwork.cc:1028
void Draw(wavearray< int > w, TString title="")
Definition: gnetwork.hh:101
double frequency
float theta
void DrawSitesLabel(Color_t tcolor=kBlack, Size_t tsize=0.045, Font_t tfont=32)
Definition: gnetwork.cc:442
void DumpObject(char *file)
Definition: gnetwork.cc:1492
void Draw(int dpaletteId=1, Option_t *option="colfz")
Definition: gskymap.cc:460
int polarization
std::vector< TLatex * > siteP
Definition: gnetwork.hh:121
i drho i
std::vector< TMarker * > siteM
Definition: gnetwork.hh:117
char ifo[NIFO_MAX][8]
network ** net
NOISE_MDC_SIMULATION.
#define DUMMY_PALETTE_ID
Definition: gskymap.hh:65
size_t mode
watplot p2(const_cast< char *>("TFMap2"))
wavearray< double > w
Definition: Test1.C:27
#define nIFO
void DrawNetworkDetectorIndex(double gamma, int loop, double snr=-1.0, int dpaletteId=DUMMY_PALETTE_ID, bool btitle=true)
Definition: gnetwork.cc:1252
void Draw(wavearray< double > w, TString title="")
Definition: gnetwork.hh:107
void DrawSitesShortLabel(Color_t tcolor=kBlack, Size_t tsize=0.052, Font_t tfont=32)
Definition: gnetwork.cc:426
float phi
float psi
void setSkyMaps(double sms, double t1, double t2, double p1, double p2)
Definition: gnetwork.hh:47
i() int(T_cor *100))
void setSkyMaps(int healpix_order)
Definition: gnetwork.hh:49
watplot p1(const_cast< char *>("TFMap1"))
const int NIFO_MAX
Definition: wat.hh:22
void DrawSites(Color_t mcolor=kBlack, Size_t msize=2.0, Style_t mstyle=20)
Definition: gnetwork.cc:257
int l
gnetwork()
Definition: gnetwork.hh:39
unsigned int siteLabelType
default=1, used to scale the antenna pattern acccording the PSD contribution of each detector (at a f...
Definition: gnetwork.hh:115
static double ndi(double *, double, int)
Definition: gnetwork.hh:133
void LoadObject(char *file)
Definition: gnetwork.cc:1506
Definition: skymap.hh:63
void SetScale(TString ifo, double scale)
Definition: gnetwork.cc:1627
double scale[NIFO_MAX]
Definition: gnetwork.hh:113
void SetTitle(TString title)
Definition: gskymap.hh:152
gskymap gSM
Definition: gnetwork.hh:111
double GetDelay(TString ifo1, TString ifo2, double phi, double theta)
Definition: gnetwork.cc:926
gskymap * GetGskymap()
Definition: gnetwork.hh:44
int nifo
std::vector< TText * > siteT
Definition: gnetwork.hh:119
char title[256]
Definition: SSeriesExample.C:1
void Draw(skymap sm, TString title="")
Definition: gnetwork.hh:100
double gps
void MakeNetworkDetectorIndex(double gamma, int loop, double snr=-1.0)
Definition: gnetwork.cc:1147
double gamma
Definition: network.hh:592
skymap nSensitivity
list of wdm tranformations
Definition: network.hh:619
void ClearSites()
Definition: gnetwork.cc:524
std::vector< TPolyLine * > siteA
Definition: gnetwork.hh:118
void ClearSitesLabel()
Definition: gnetwork.cc:548
double GetAntennaPattern(double phi, double theta, double psi=0., int polarization=1)
Definition: gnetwork.cc:564
void set(size_t i, double a)
param: sky index param: value to set
Definition: skymap.hh:122
bool btitle
Definition: DrawGnetwork2.C:16
void Draw(char *smName, network *net=NULL)
Definition: gnetwork.hh:77
void ClearSitesArms()
Definition: gnetwork.cc:536
void ClearCircles()
Definition: gnetwork.cc:510
TGraphPolar * grP[2]
canvas used for the polargrams
Definition: gnetwork.hh:128
void Draw(wavearray< short > w, TString title="")
Definition: gnetwork.hh:104
detectorParams detParms[4]
size_t size()
Definition: skymap.hh:136
std::vector< TLatex * > siteL
Definition: gnetwork.hh:120
void Init()
Definition: gnetwork.cc:111
std::vector< TMarker * > daxisM
Definition: gnetwork.hh:123
void setSkyMaps(double, double=0., double=180., double=0., double=360.)
param: sky map granularity step, degrees param: theta begin, degrees param: theta end...
Definition: network.cc:2687