Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NBRun.hh
Go to the documentation of this file.
1 #ifndef __NBRUN__
2 #define __NBRUN__
3 
5 #include "TFile.h"
6 #include "TVectorD.h"
7 #include "SUPERMODULE.hh"
8 #include "PMT_choice.hh"
9 
10 using namespace std ;
11 
12 using namespace JDETECTOR ;
13 
14 /**
15  * \author rgruiz
16  */
17 
18 
19 /**
20  * Class dedicated to the nanobeacon analyses, where the Modules in the detector are not regarded as single entities. Instead, they are related
21  * to the rest of the DOMs in a DU through the emission of the nanobeacon pulses. A DOM can be emitter (source) of light that is detected by other DOMS, and
22  * it can also be a receiver (target) of the light emitted by the other DOMS. The SUPERRUN class is meant to own the different SUPERMODULEs in a DU, to coherently set their mutual rerefences and to analyze the signal produced by the different nanobeacons in the different PMTs of the DU.
23  */
24 class NBRun {
25 
26  vector<int> topPMTs , bottomPMTs , srcPMTs , tgtPMTs;
27 
29 
31 
33 
34  string filename ;
35 
36  int run_number ;
37 
38  double voltage ;
39 
41 
43 
44 public :
45 
46 
47  /**
48  * Default constructor.
49  */
50  NBRun (){} ;
51 
52 
53 
54  /**
55  * Constructor.
56  * It initializes the Structure of SuperModules and leaves it ready for analyzing it.
57  *
58  * \param filename_ The name of a .root file containing the time vs ToT hit distributions for different PMTs in the DU
59  * \param detector_ A detector file with the DU that is going to be analyzed.
60  * \param string_ String number of the string to be analyzed (as it appears in the detector file)
61  * \param pmt_top Integer number representing the user's choice
62  * \param pmt_bottom Integer number representing the user's choice
63  * \param max_distance_ Integer number representing the maximum number of neighbors to for each module.
64  */
65  NBRun (string filename_ , JDetector detector_ , int string_ , int pmt_top , int pmt_bottom , int max_distance_):
66  nSuperModules (getNumberOfFloors(detector_))
67  {
68 
69  filename = filename_ ;
70 
71  string_number = string_ ;
72 
73  max_distance = max_distance_ ;
74 
75  setDetector (detector_) ;
76 
77  initializeSuperModules () ;
78 
79  setPMTs (pmt_top , pmt_bottom) ;
80 
81  readFile() ;
82 
83  }
84 
85 
86 
87  /**
88  * Destructor.
89  */
90  ~NBRun(){
91 
92  for (auto & sm : SuperMods){
93 
94  delete(sm) ;
95 
96  }
97 
98  };
99 
100 
101 
102  /**
103  * Loops over all the SUPERMODULES to compute the mean ToT of the hits from the different nanobeacon pulses recorded by the different PMTs.
104  */
106 
107  for (auto & sm : SuperMods){
108 
109  sm->compute_mean_tot_ref() ;
110 
111  sm->compute_mean_tot_tgt() ;
112 
113  }
114 
115  }
116 
117 
118 
119  /**
120  * Loops over all the SUPERMODULES. For each of them, it searches for good sources and good targets .
121  */
123 
124  for (auto & sm : SuperMods){
125 
126  sm->select_good_srcs() ;
127 
128  sm->select_good_tgts() ;
129 
130  }
131 
132  }
133 
134 
135  /**
136  * Loops over all the SUPERMODULES. Analyzes the nanobeacon pulses detected by the different PMTs. After they have been analyzed, the method find_good_couples is called.
137  *
138  * \param option Level of analysis of the nanobeacon pulses. 0 for a fast analysis 1 for a fit to a model
139  */
140  void analyze(int option){
141 
142  for (auto & sm : SuperMods){
143 
144  switch(option){
145 
146  case 0 :
147 
148  sm->analyzeAll() ;
149 
150  break ;
151 
152  case 1 :
153 
154  sm->fitAll() ;
155 
156  break ;
157 
158  }
159 
160  }
161 
162  find_good_couples();
163 
164  }
165 
166 
167  /**
168  * Reads a .root file containing the nanobeacon pulses observed by the different PMTs in the detector.
169  */
170  void readFile (){
171 
172  TFile* file = TFile::Open(filename.c_str()) ;
173 
174  readBasicInfo (file) ;
175 
176  setRefs (file) ;
177 
178  setSrcs (file) ;
179 
180  }
181 
182 
183 
184  /**
185  * Searches in the .root file for histograms corresponding to pulses produced in one target supermodule, by all the possible sources.
186  *
187  * \param sm supermodule
188  * \param file root file
189  */
190  void findSrcs(SuperModule* sm , TFile* file){
191 
192  int thisfloor = sm->getFloor () ;
193 
194  int thisstring = sm->getString () ;
195 
196  vector<SuperModule*> possible_srcs = sm->get_possible_srcs () ;
197 
198  for (auto & src : possible_srcs){
199 
200  int thatfloor = src->getFloor () ;
201 
202  int thatstring = src->getString () ;
203 
204  vector<SuperPMT*> spms ;
205 
206  char name[250] ;
207 
208  for (auto & pm : bottomPMTs){
209 
210  sprintf(name , "S%dF%d/S%dF%dPMT%d_beaconS%dF%d" , thisstring , thisfloor, thisstring , thisfloor , pm , thatstring , thatfloor);
211 
212  if (file -> Get (name) ){
213 
214  TH2D* h = (TH2D*)file -> Get(name) ;
215 
216  h->SetDirectory(0) ;
217 
218  JPMT pmt = sm->getPMT ( pm ) ;
219 
220  SuperPMT* p = new SuperPMT ( pmt , h , pm) ;
221 
222  spms.push_back (p) ;
223 
224  }
225 
226  }
227 
228  if(spms.size()>0){
229 
230  pair < SuperModule* , vector<SuperPMT*> > src_pair (src , spms) ;
231 
232  pair < SuperModule* , vector<SuperPMT*> > tgt_pair (sm , spms) ;
233 
234  sm->add_src (src_pair) ;
235 
236  src->add_tgt (tgt_pair) ;
237 
238  }
239 
240  }
241 
242  }
243 
244 
245 
246  /**
247  * Loops over all the SuoperModules in the DU. For each supermodule, it calls the findSrcs method.
248  *
249  * \param file root file
250  */
251  void setSrcs (TFile* file){
252 
253  for(auto & sm : SuperMods){
254 
255  findSrcs (sm , file) ;
256 
257  }
258 
259  }
260 
261 
262 
263  /**
264  * Loops over all the SuoperModules in the DU. For each supermodule, it calls the findRefs method.
265  *
266  * \param file root file
267  */
268  void setRefs (TFile* file){
269 
270  for (auto & sm : SuperMods){
271 
272  findRefs (sm , file) ;
273 
274  }
275 
276  }
277 
278 
279 
280  /**
281  * Searches in the .root file for histograms corresponding to pulses produced in one supermodule, by its own nanobeacon.
282  *
283  * \param sm supermodule
284  * \param file root file
285  */
286  void findRefs (SuperModule* sm , TFile* file ){
287 
288  int floor = sm->getFloor () ;
289 
290  int string = sm->getString () ;
291 
292  char name[250] ;
293 
294  for (auto & pm : topPMTs){
295 
296  sprintf(name , "S%dF%d/S%dF%dPMT%d_beaconS%dF%d" , string , floor, string , floor , pm , string , floor);
297 
298  if (file -> Get (name) ){
299 
300  TH2D* h = (TH2D*)file -> Get(name) ;
301 
302  h->SetDirectory(0) ;
303 
304  JPMT pmt = sm->getPMT ( pm ) ;
305 
306  SuperPMT* p = new SuperPMT (pmt , h , pm) ;
307 
308  sm->add_ref_pmt(p) ;
309 
310  }
311 
312  }
313 
314  }
315 
316 
317 
318  /**
319  * Get the SuperModules in the DU.
320  *
321  * \return Vector of SuperModules
322  */
324 
325  return SuperMods ;
326 
327  }
328 
329 
330 
331  /**
332  * Select the PMTs in the upper and lower hemispheres of a DOM.
333  *
334  * \param top_option user choice
335  * \param bottom_option user choice
336  */
337  void setPMTs (int top_option , int bottom_option){
338 
339  topPMTs = setTopPMTs (top_option) ;
340 
341  bottomPMTs = setBottomPMTs (bottom_option) ;
342 
343  }
344 
345 
346  /**
347  * Sets the references between the different SuperModules in the DU
348  */
350 
351  SuperMods.resize(nSuperModules) ;
352 
353  for (int i = 0 ; i < (int)SuperMods.size() ; i++){
354 
355  JModule m = detector.getModule(JModuleLocation(string_number , i+1)) ;
356 
357  SuperModule* s_m = new SuperModule (m) ;
358 
359  for (int j = i+1 ; (j < i + 1 + max_distance) && (j<(int)nSuperModules) ; j++){
360 
361  s_m->add_possible_tgt(SuperMods[j]) ;
362 
363  }
364 
365  for (int j = max(i - max_distance , 0) ; j<i ; j++){
366 
367  s_m->add_possible_src(SuperMods[j]) ;
368 
369  }
370 
371  SuperMods[i] = s_m ;
372 
373  }
374 
375  }
376 
377 
378 
379  /**
380  * Sets the detector
381  *
382  * \param detector_ a JDetector
383  */
384  void setDetector(JDetector detector_){
385 
386  detector = detector_ ;
387 
388  }
389 
390 
391  /**
392  * Get the detector
393  *
394  * \return a JDetector
395  */
397 
398  return detector ;
399 
400  }
401 
402 
403 
404  /**
405  * Reads the basic info from the .root file such as the run number and the nanobeacon voltage
406  *
407  * \param file root file
408  */
409  void readBasicInfo(TFile* file){
410 
411  filename = file->GetName() ;
412 
413  TVectorD* run_info = (TVectorD*)file->Get("Run_Info") ;
414 
415  run_number = (*run_info)(0) ;
416 
417  voltage = (*run_info)(1) ;
418 
419  }
420 
421 
422  /**
423  * Get run number
424  *
425  * \return run number
426  */
428 
429  return run_number ;
430 
431  }
432 
433 
434  /**
435  * Get nanobeacon voltage
436  *
437  * \return nanobeacon voltage
438  */
439  double getVoltage(){
440 
441  return voltage ;
442 
443  }
444 
445 
446 };
447 
448 #endif
JDetector getDetector()
Get the detector.
Definition: NBRun.hh:396
int getString()
Returns the string of this supermodule.
Definition: SUPERMODULE.hh:651
vector< SuperModule * > SuperMods
Definition: NBRun.hh:42
Class dedicated to the nanobeacon analyses, where the Modules in the detector are not regarded as sin...
Definition: SUPERMODULE.hh:26
Data structure for a composite optical module.
Definition: JModule.hh:47
vector< SuperModule * > getSuperModules()
Get the SuperModules in the DU.
Definition: NBRun.hh:323
~NBRun()
Destructor.
Definition: NBRun.hh:90
void readFile()
Reads a .root file containing the nanobeacon pulses observed by the different PMTs in the detector...
Definition: NBRun.hh:170
vector< int > setTopPMTs(int option)
Select the PMTs in the upper hemisphere of a DOM.
Definition: PMT_choice.hh:48
int max_distance
Definition: NBRun.hh:30
Detector data structure.
Definition: JDetector.hh:77
int string_number
Definition: NBRun.hh:32
void setSrcs(TFile *file)
Loops over all the SuoperModules in the DU.
Definition: NBRun.hh:251
string filename
Definition: NBRun.hh:34
void add_possible_src(SuperModule *super_m)
Adds a SUPERMODULE to the list of possible sources of this supermodule.
Definition: SUPERMODULE.hh:729
vector< int > setBottomPMTs(int option)
Select the PMTs in the lower hemisphere of a DOM.
Definition: PMT_choice.hh:100
void setDetector(JDetector detector_)
Sets the detector.
Definition: NBRun.hh:384
int getFloor()
Returns the floor of this supermodule.
Definition: SUPERMODULE.hh:638
void setRefs(TFile *file)
Loops over all the SuoperModules in the DU.
Definition: NBRun.hh:268
void computeMeanToTs()
Loops over all the SUPERMODULES to compute the mean ToT of the hits from the different nanobeacon pul...
Definition: NBRun.hh:105
JDetector detector
Definition: NBRun.hh:40
void add_possible_tgt(SuperModule *super_m)
Adds a SUPERMODULE to the list of possible targets of this supermodule.
Definition: SUPERMODULE.hh:742
Class containing a JPMT and a NBPulse object.
Definition: SUPERPMT.hh:22
void findRefs(SuperModule *sm, TFile *file)
Searches in the .root file for histograms corresponding to pulses produced in one supermodule...
Definition: NBRun.hh:286
void add_ref_pmt(SuperPMT *superPM)
Adds a reference SUPERPMT to this supermodule.
Definition: SUPERMODULE.hh:716
int nSuperModules
Definition: NBRun.hh:28
Data structure for PMT geometry and calibration.
Definition: JPMT.hh:52
int getNumberOfFloors(const JDetector &detector)
Get number of floors.
void setPMTs(int top_option, int bottom_option)
Select the PMTs in the upper and lower hemispheres of a DOM.
Definition: NBRun.hh:337
void findSrcs(SuperModule *sm, TFile *file)
Searches in the .root file for histograms corresponding to pulses produced in one target supermodule...
Definition: NBRun.hh:190
double getVoltage()
Get nanobeacon voltage.
Definition: NBRun.hh:439
void readBasicInfo(TFile *file)
Reads the basic info from the .root file such as the run number and the nanobeacon voltage...
Definition: NBRun.hh:409
int getRunNumber()
Get run number.
Definition: NBRun.hh:427
int run_number
Definition: NBRun.hh:36
void analyze(int option)
Loops over all the SUPERMODULES.
Definition: NBRun.hh:140
void find_good_couples()
Loops over all the SUPERMODULES.
Definition: NBRun.hh:122
vector< int > topPMTs
Definition: NBRun.hh:26
JPMT getPMT(int i)
Returns a JPMT from this supermodule.
Definition: SUPERMODULE.hh:625
double voltage
Definition: NBRun.hh:38
Class dedicated to the nanobeacon analyses, where the Modules in the detector are not regarded as sin...
Definition: NBRun.hh:24
void add_src(pair< SuperModule *, vector< SuperPMT * > > p)
Adds a source.
Definition: SUPERMODULE.hh:599
NBRun()
Default constructor.
Definition: NBRun.hh:50
Logical location of module.
vector< SuperModule * > get_possible_srcs()
Returns vector with pointers to all the supermodules below this one in the DU.
Definition: SUPERMODULE.hh:690
NBRun(string filename_, JDetector detector_, int string_, int pmt_top, int pmt_bottom, int max_distance_)
Constructor.
Definition: NBRun.hh:65
void initializeSuperModules()
Sets the references between the different SuperModules in the DU.
Definition: NBRun.hh:349