Logo coherent WaveBurst  
Library Reference Guide
Logo
livetime.cc
Go to the documentation of this file.
1 /*
2 # Copyright (C) 2019 Gabriele Vedovato, Sergey Klimenko
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 #include "livetime.hh"
20 #include "TH2.h"
21 #include "TStyle.h"
22 #include "TCanvas.h"
23 
24 ClassImp(livetime) // used by THtml doc
25 
26 livetime& livetime::operator=(livetime& a)
27 {
28  this->fChain= NULL;
29  this->fCurrent= a.fCurrent;
30  this->run= a.run;
31  this->gps= a.gps;
32  this->live= a.live;
33  this->lag[NIFO_MAX+1]= a.lag[NIFO_MAX+1];
34  this->slag[NIFO_MAX+1]= a.slag[NIFO_MAX+1];
35  for(size_t n=0; n<(NIFO_MAX+1); n++) {
36  this->lag[n]= a.lag[n];
37  this->slag[n]= a.slag[n];
38  }
39  for(size_t n=0; n<NIFO_MAX; n++) {
40  this->start[n]= a.start[n];
41  this->stop[n]= a.stop[n];
42  }
43  return *this;
44 }
45 
46 // Set branch addresses
47 void livetime::Init(TTree *tree)
48 {
49  if (tree == 0) return;
50  fChain = tree;
51  fCurrent = -1;
52  fChain->SetMakeClass(1);
53 
54  fChain->SetBranchAddress("run",&(this->run));
55  fChain->SetBranchAddress("gps",&(this->gps));
56  fChain->SetBranchAddress("live",&(this->live));
57  fChain->SetBranchAddress("lag",(this->lag));
58  fChain->SetBranchAddress("slag",(this->slag));
59  fChain->SetBranchAddress("start",(this->start));
60  fChain->SetBranchAddress("stop",(this->stop));
61 
62  Notify();
63 }
64 
66 {
67  if (!lag) lag= (Float_t*)malloc((NIFO_MAX+1)*sizeof(Int_t));
68  else lag= (Float_t*)realloc(lag,(NIFO_MAX+1)*sizeof(Int_t));
69  if (!slag) slag= (Float_t*)malloc((NIFO_MAX+1)*sizeof(Int_t));
70  else slag= (Float_t*)realloc(slag,(NIFO_MAX+1)*sizeof(Int_t));
71  if (!start) start= (Double_t*)malloc((NIFO_MAX)*sizeof(Double_t));
72  else start= (Double_t*)realloc(start,(NIFO_MAX)*sizeof(Double_t));
73  if (!stop) stop= (Double_t*)malloc((NIFO_MAX)*sizeof(Double_t));
74  else stop= (Double_t*)realloc(start,(NIFO_MAX)*sizeof(Double_t));
75  return;
76 }
77 
78 //++++++++++++++++++++++++++++++++++++++++++++++
79 // set noise livetime tree
80 //++++++++++++++++++++++++++++++++++++++++++++++
82 {
83  TTree* waveTree = new TTree("liveTime","liveTime");
84 
85  //==================================
86  // Define livetime tree
87  //==================================
88 
89  char clag[16];sprintf(clag,"lag[%d]/F",NIFO_MAX+1);
90  char cslag[16];sprintf(cslag,"slag[%d]/F",NIFO_MAX+1);
91  char cstart[16];sprintf(cstart,"start[%d]/D",NIFO_MAX);
92  char cstop[16];sprintf(cstop,"stop[%d]/D",NIFO_MAX);
93 
94  waveTree->Branch("run", &(this->run), "run/I");
95  waveTree->Branch("gps", &(this->gps), "gps/D");
96  waveTree->Branch("live", &(this->live), "live/D");
97  waveTree->Branch("lag", (this->lag), clag);
98  waveTree->Branch("slag", (this->slag), cslag);
99  waveTree->Branch("start", (this->start), cstart);
100  waveTree->Branch("stop", (this->stop), cstop);
101 
102  return waveTree;
103 }
104 
105 
106 //++++++++++++++++++++++++++++++++++++++++++++++++++++
107 // dump livetime into tree
108 //++++++++++++++++++++++++++++++++++++++++++++++++++++
109 void livetime::output(TTree* waveTree, network* net, float* slag, vector<waveSegment> detSegs)
110 {
111  size_t n,m;
112  size_t M = net->ifoListSize();
113  if(!M) return;
114 
115  // add detectors to tree user info
116 // if(waveTree!=NULL) for(m=0;m<M;m++) waveTree->GetUserInfo()->Add(net->getifo(m));
117 
118 //Fill tree
119  for(m=0; m<(NIFO_MAX+1); m++) {
120  this->lag[m] = -1;
121  this->slag[m] = -1;
122  }
123  for(m=0; m<NIFO_MAX; m++) {
124  this->start[m] = 0;
125  this->stop[m] = 0;
126  }
127  for(m=0; m<detSegs.size(); m++) {
128  this->start[m] = detSegs[m].start;
129  this->stop[m] = detSegs[m].stop;
130  }
131  this->run = net->nRun;
132  this->gps = net->getifo(0)->getTFmap()->start();
133 
134  if(slag!=NULL) for(m=0; m<M+1; m++) this->slag[m] = slag[m];
135 
136  for(n=0; n<net->nLag; n++){ // loop on time lags
137 
138  this->live = net->getliveTime(n);
139 
140  for(m=0; m<NIFO_MAX; m++){ // loop on detectors
141  if(m<M) {
142  this->lag[m] = net->getifo(m)->lagShift.data[n];
143  }
144  }
145  this->lag[M] = net->getwc(n)->shift;
146  waveTree->Fill();
147  }
148 }
149 
150 
151 
153 {
154  // Called when loading a new file.
155  // Get branch pointers.
156  b_run = fChain->GetBranch("run");
157  b_gps = fChain->GetBranch("gps");
158  b_live = fChain->GetBranch("live");
159  b_lag = fChain->GetBranch("lag");
160  b_slag = fChain->GetBranch("slag");
161  b_start= fChain->GetBranch("start");
162  b_stop = fChain->GetBranch("stop");
163  return kTRUE;
164 }
165 
167 {
168  if (!fChain) return 0;
169  return fChain->GetEntry(entry);
170 };
171 
173 {
174 // Print contents of entry.
175 // If entry is not specified, print current entry
176  if (!fChain) return;
177  fChain->Show(entry);
178 }
179 
180 
181 /*
182 Int_t livetime::LoadTree(Int_t entry)
183 {
184 // Set the environment to read one entry
185  if (!fChain) return -5;
186  Int_t centry = fChain->LoadTree(entry);
187  if (centry < 0) return centry;
188  if (fChain->IsA() != TChain::Class()) return centry;
189  TChain *chain = (TChain*)fChain;
190  if (chain->GetTreeNumber() != fCurrent) {
191  fCurrent = chain->GetTreeNumber();
192  Notify();
193  }
194  return centry;
195 }
196 
197 Int_t livetime::Cut(Int_t entry)
198 {
199 // This function may be called from Loop.
200 // returns 1 if entry is accepted.
201 // returns -1 otherwise.
202  return 1;
203 }
204 
205 void livetime::Loop()
206 {
207 // In a ROOT session, you can do:
208 // Root > .L livetime.C
209 // Root > livetime t
210 // Root > t.GetEntry(12); // Fill t data members with entry number 12
211 // Root > t.Show(); // Show rmss of entry 12
212 // Root > t.Show(16); // Read and show values of entry 16
213 // Root > t.Loop(); // Loop on all entries
214 //
215 
216 // This is the loop skeleton where:
217 // jentry is the global entry number in the chain
218 // ientry is the entry number in the current Tree
219 // Note that the argument to GetEntry must be:
220 // jentry for TChain::GetEntry
221 // ientry for TTree::GetEntry and TBranch::GetEntry
222 //
223 // To read only selected branches, Insert statements like:
224 // METHOD1:
225 // fChain->SetBranchStatus("*",0); // disable all branches
226 // fChain->SetBranchStatus("branchname",1); // activate branchname
227 // METHOD2: replace line
228 // fChain->GetEntry(jentry); //read all branches
229 //by b_branchname->GetEntry(ientry); //read only this branch
230  if (fChain == 0) return;
231 
232  Int_t nentries = Int_t(fChain->GetEntriesFast());
233 
234  Int_t nbytes = 0, nb = 0;
235  for (Int_t jentry=0; jentry<nentries;jentry++) {
236  Int_t ientry = LoadTree(jentry);
237  if (ientry < 0) break;
238  nb = fChain->GetEntry(jentry); nbytes += nb;
239  // if (Cut(ientry) < 0) continue;
240  }
241 }
242 */
TTree * tree
Definition: TimeSortTree.C:20
detector * getifo(size_t n)
param: detector index
Definition: network.hh:436
void Init(TTree *)
Definition: livetime.cc:47
size_t nLag
Definition: network.hh:573
double M
Definition: DrawEBHH.C:13
TString live
CWB run(runID)
wavearray< double > a(hp.size())
int n
Definition: cwb_net.C:28
TBranch * b_live
Definition: livetime.hh:60
Double_t * start
time slag [sec]
Definition: livetime.hh:53
size_t nRun
Definition: network.hh:572
double getliveTime(size_t n)
Definition: network.hh:428
Double_t gps
Definition: livetime.hh:49
TBranch * b_start
Definition: livetime.hh:63
int m
Definition: cwb_net.C:28
virtual void start(double s)
Definition: wavearray.hh:137
network ** net
NOISE_MDC_SIMULATION.
size_t ifoListSize()
Definition: network.hh:431
Int_t fCurrent
pointer to the analyzed TTree or TChain
Definition: livetime.hh:43
Double_t live
Definition: livetime.hh:50
double shift
Definition: netcluster.hh:382
Int_t run
current Tree number in a TChain
Definition: livetime.hh:48
TTree * fChain
Definition: livetime.hh:42
TBranch * b_gps
Definition: livetime.hh:59
const int NIFO_MAX
Definition: wat.hh:22
Float_t * lag
Definition: livetime.hh:51
TBranch * b_lag
Definition: livetime.hh:61
void allocate()
Definition: livetime.cc:65
void output(TTree *waveTree, network *net, float *slag=NULL, vector< waveSegment > detSegs=DEFAULT_WAVESEGMENT)
Definition: livetime.cc:109
wavearray< double > lagShift
Definition: detector.hh:369
double * entry
Definition: cwb_setcuts.C:224
Int_t GetEntry(Int_t)
Definition: livetime.cc:166
vector< waveSegment > detSegs
Definition: cwb_net.C:193
Bool_t Notify()
Definition: livetime.cc:152
double gps
Double_t * stop
array of gps start time segment for each detector
Definition: livetime.hh:54
WSeries< double > * getTFmap()
param: no parameters
Definition: detector.hh:179
Definition: Toolbox.hh:99
netcluster * getwc(size_t n)
param: delay index
Definition: network.hh:439
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
DataType_t * data
Definition: wavearray.hh:319
TTree * setTree()
Definition: livetime.cc:81
void Show(Int_t entry=-1)
Definition: livetime.cc:172
TBranch * b_run
array of gps stop time segment for each detector
Definition: livetime.hh:58
TBranch * b_slag
Definition: livetime.hh:62
TBranch * b_stop
Definition: livetime.hh:64