Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JPulseFinder_utils.hh
Go to the documentation of this file.
1 #ifndef __PULSE_UTILS__
2 #define __PULSE_UTILS__
3 
4 
5 // c++ standard library
6 #include <string>
7 #include <iostream>
8 #include <sstream>
9 #include <iomanip>
10 
11 // JPP
12 #include "JDAQ/JDAQ.hh"
13 #include "JDAQ/JDAQTimeslice.hh"
15 #include "JSupport/JSupport.hh"
16 #include "JDetector/JDetector.hh"
19 #include "JTrigger/JHitR0.hh"
20 #include "JTrigger/JHitL0.hh"
21 #include "JTrigger/JHitToolkit.hh"
24 
25 // root
26 #include "TH1D.h"
27 #include "TH2D.h"
28 
29 // namespaces
30 using namespace std;
31 using namespace KM3NETDAQ ; // for JDAQTimeSlice
32 using namespace JSUPPORT; // for JFileScanner
33 using namespace JLANG; // for JMultipleFileScanner
34 using namespace JDETECTOR ;
35 
36 using namespace JPP ;
37 
38 /**
39  * \author rgruiz
40  */
41 
42 
43 
44 /**
45  * Calculates an offset for a combinataion nanobeacon-pmt that allows to define the histogram range in a time window containing the corresponding pulse.
46  *
47  * \param pmt The target pmt
48  * \param NB The source module
49  * \param c_w The speed light in water
50  * \param options user options containing the stagger and pulse delay information
51  * \return The corresponding offset to apply in the definition of the hit time histogram range.
52  */
53 inline int get_offset(JPMT pmt , JModule NB , double c_w , IO options){
54 
55  double distance_pmt_nb = (pmt.getPosition() - NB.getPMT(0).getPosition()).getLength() ;
56 
57  double time_of_flight = distance_pmt_nb / c_w ;
58 
59  double delta_t0 = pmt.getT0() - NB.getPMT(0).getT0() ;
60 
61  double delta_traw = time_of_flight - delta_t0 ;
62 
63  double stagger_delay = 16 * options.stagger_16ns * (NB.getFloor() - 1) ;
64 
65  double overall_delay = 16 * options.overall_delay_16ns ;
66 
67  int total_offset = (int)round( delta_traw + stagger_delay + overall_delay);
68 
69  return total_offset ;
70 
71 }
72 
73 
74 
75 /**
76  * Allocates histograms to store the tot vs time hit distribution for all the combinations nanobeacon - pmt.
77  *
78  * \param detector The detector containing the DU to be analyzed
79  * \param options user options containing the string number and the list of used nanobeacons
80  * \param c_w The speed light in water
81  * \param top_pmts List of FPGA channels for pmts in the upper hemisphere
82  * \param bottom_pmts List of FPGA channels for pmts in the lower hemisphere.
83  * \param max_distance Maximum distance between modules
84  * \return A table of TH2D histograms defined in the correct hit time ranges.
85  */
86 inline vector < vector < TH2D* > > allocate_histograms (JDetector detector , IO options , double c_w , vector<int> top_pmts , vector<int> bottom_pmts , int max_distance){
87 
88  int nFloors = getNumberOfFloors(detector) ;
89 
90  int nHistosPerFloor = top_pmts.size() + max_distance * bottom_pmts.size() ;
91 
92  vector < vector < TH2D* > > histograms (nFloors , vector<TH2D*> (nHistosPerFloor, NULL)) ;
93 
94  char name[250] ;
95 
96  for (int tgt_floor = 1 ; tgt_floor <= nFloors ; ++tgt_floor){
97 
98  JModule tgt_module = detector.getModule( JModuleLocation(options.string_number , tgt_floor) ) ;
99 
100  int index = 0 ;
101 
102  for(auto & top_pmt : top_pmts){
103 
104  int total_offset = get_offset(tgt_module.getPMT(top_pmt) , tgt_module , c_w , options) ;
105 
106  sprintf(name , "S%iF%iPMT%i_beaconS%iF%i" , tgt_module.getString() , tgt_module.getFloor() , top_pmt , tgt_module.getString() , tgt_module.getFloor());
107 
108  histograms[tgt_floor - 1][index] = new TH2D( name , NULL , 400 , total_offset - 0.5 , 400 + total_offset - 0.5 , 256 , -0.5 , 255.5 ) ;
109 
110  histograms[tgt_floor - 1][index]->SetOption("colz") ;
111 
112  index += 1 ;
113  }
114 
115 
116  for (int src_floor = tgt_floor - max_distance ; src_floor < tgt_floor ; ++src_floor){
117 
118  if(src_floor > 0){
119 
120  JModule src_module = detector.getModule( JModuleLocation (options.string_number , src_floor) ) ;
121 
122  for (auto & bottom_pmt : bottom_pmts){
123 
124  int total_offset = get_offset(tgt_module.getPMT(bottom_pmt) , src_module , c_w , options) ;
125 
126  sprintf(name , "S%iF%iPMT%i_beaconS%iF%i" , tgt_module.getString() , tgt_module.getFloor() , bottom_pmt , src_module.getString() , src_module.getFloor());
127 
128  histograms[tgt_floor - 1][index] = new TH2D( name , NULL , 400 , total_offset - 0.5 , 400 + total_offset - 0.5 , 256 , -0.5 , 255.5 ) ;
129 
130  histograms[tgt_floor - 1][index]->SetOption("colz") ;
131 
132  index += 1 ;
133 
134  }
135 
136  }
137 
138  }
139 
140  }
141 
142  return histograms ;
143 
144 }
145 
146 
147 
148 
149 /**
150  * Fills the hit time vs ToT histograms for all the combinations nanobeacon - pmt.
151  *
152  * \param histograms The vector of histograms to be filled
153  * \param options User options containing the string number and the list files to be read
154  * \param detector The detector containing the DU to be analyzed
155  * \param top_pmts List of FPGA channels for pmts in the upper hemisphere
156  * \param bottom_pmts List of FPGA channels for pmts in the lower hemisphere.
157  * \param max_distance Maximum distance between modules
158  */
159 inline void fill_histograms (vector < vector < TH2D* > > &histograms , IO options , JDetector detector , vector<int> top_pmts , vector<int> bottom_pmts , int max_distance){
160 
161  JTriggerParameters parameters ;
162 
163  unsigned int pulse_period = 16 * options.pulse_period_16ns ;
164 
166 
167  JModuleRouter moduleRouter(detector) ;
168 
169  int slicecounter = 0 ;
170 
171  while( inputFile.hasNext() ) {
172 
173  slicecounter ++ ;
174 
175  JDAQTimesliceL0 slice = *(inputFile.next()) ;
176 
177  JTimesliceRouter timesliceRouter(parameters.numberOfBins) ;
178 
179  timesliceRouter.configure(slice);
180 
181  for (int i = 0 ; i < getNumberOfFloors(detector) ; i++){
182 
183  const JModule& module = detector.getModule(JModuleLocation(options.string_number , i+1)) ;
184 
185  if(timesliceRouter.hasSuperFrame(module.getID())){
186 
187  JDAQSuperFrame SuperFrame = timesliceRouter.getSuperFrame(module.getID()) ;
188 
189  vector<JDAQHit> buffer[getNumberOfPMTs(module)];
190 
191  for (JDAQSuperFrame::const_iterator hit = SuperFrame.begin(); hit != SuperFrame.end(); ++hit) {
192 
193  buffer[hit->getPMT()].push_back(*hit) ;
194 
195  }
196 
197  int i = 0 ;
198 
199  for (auto & top_pmt : top_pmts){
200 
201  for (int hit = 0 ; hit < (int)buffer[top_pmt].size() - 1 ; hit++){
202 
203  double tot = buffer[top_pmt][hit].getToT() ;
204 
205  double raw_time = buffer[top_pmt][hit].getT() ;
206 
207  //double t0 = module.getPMT(buffer[top_pmt][hit].getPMT()).getT0() ;
208 
209  //double calibrated_time = raw_time + t0 ;
210 
211  double mod_double = raw_time - pulse_period * int(raw_time/pulse_period) ;
212 
213  histograms[module.getFloor()-1][i]->Fill(mod_double , tot) ;
214 
215  }
216 
217  i++ ;
218 
219  }
220 
221 
222  i = 0 ;
223 
224  for (auto & bottom_pmt : bottom_pmts){
225 
226  for (int hit = 0 ; hit < (int)buffer[bottom_pmt].size() - 1 ; hit++){
227 
228  double tot = buffer[bottom_pmt][hit].getToT() ;
229 
230  double raw_time = buffer[bottom_pmt][hit].getT() ;
231 
232  //double t0 = module.getPMT(buffer[bottom_pmt][hit].getPMT()).getT0() ;
233 
234  //double calibrated_time = raw_time + t0 ;
235 
236  double mod_double = raw_time - pulse_period * int(raw_time/pulse_period) ;
237 
238  for (int j = 1 ; j <= max_distance ; j++){
239 
240  if (module.getFloor() - j <=0) continue ;
241 
242  histograms[module.getFloor() - 1][top_pmts.size() + (j-1)*bottom_pmts.size() + i]->Fill(mod_double , tot) ;
243 
244  }
245 
246  }
247 
248  i++ ;
249 
250  }
251 
252  } // end of loop over JDAQTimeSlices
253 
254  }
255 
256  }
257 
258 }
259 
260 
261 
262 
263 
264 /**
265  * Projects the hit time vs ToT histograms for all the combinations nanobeacon - pmt on the Y axis.
266  *
267  * \param ToT_vs_time The vector of histograms to be projected
268  * \return A vector with ToT distribiutions.
269  */
271 
272  vector < vector < TH1D* > > ToT (ToT_vs_time.size() , vector < TH1D* > (ToT_vs_time[0].size() , NULL)) ;
273 
274  for (int i=0 ; i<(int)ToT_vs_time.size() ; i++){
275 
276  for (int j=0 ; j<(int)ToT_vs_time[0].size() ; j++){
277 
278  if (ToT_vs_time[i][j]==NULL) continue ;
279 
280  TH1D* p_y = ToT_vs_time[i][j]->ProjectionY((string("ToT_") + ToT_vs_time[i][j]->GetName()).c_str() , 1 , ToT_vs_time[i][j]->GetNbinsX()-1) ;
281 
282  ToT[i][j] = p_y ;
283 
284  }
285 
286  }
287 
288  return ToT ;
289 
290 }
291 
292 
293 
294 /**
295  * Projects the hit time vs ToT histograms for all the combinations nanobeacon - pmt on the X axis.
296  *
297  * \param ToT_vs_time The vector of histograms to be projected
298  * \return A vector with hit time distribiutions.
299  */
301 
302  vector < vector < TH1D* > > time (ToT_vs_time.size() , vector < TH1D* > (ToT_vs_time[0].size() , NULL)) ;
303 
304  for (int i=0 ; i<(int)ToT_vs_time.size() ; i++){
305 
306  for (int j=0 ; j<(int)ToT_vs_time[0].size() ; j++){
307 
308  if (ToT_vs_time[i][j]==NULL) continue ;
309 
310  TH1D* p_x = ToT_vs_time[i][j]->ProjectionX((string("htime_") + ToT_vs_time[i][j]->GetName()).c_str() , 1 , ToT_vs_time[i][j]->GetNbinsY()-1) ;
311 
312  time[i][j] = p_x ;
313 
314  }
315 
316  }
317 
318  return time ;
319 
320 }
321 
322 
323 #endif
unsigned int pulse_period_16ns
bool hasSuperFrame(const JDAQModuleIdentifier &module) const
Check presence of module.
virtual const pointer_type & next()
Get next element.
Data structure for all trigger parameters.
int string_number
const JDAQSuperFrame & getSuperFrame(const JDAQModuleIdentifier &module) const
Get super frame.
Basic data structure for L0 hit.
Data structure for a composite optical module.
Definition: JModule.hh:47
void configure(const JDAQTimeslice &timeslice)
Configure.
int stagger_16ns
Detector data structure.
Definition: JDetector.hh:77
Structure to store the different command line arguments for JRunAnalyzer.
Definition: JRunAnalyzer.cc:40
Router for direct addressing of module data in detector data structure.
int numberOfBins
number of bins for lookup table of timeslice
int getNumberOfPMTs(const JModule &module)
Get number of PMTs.
Data structure for detector geometry and calibration.
Tools for handling different hit types.
vector< vector< TH1D * > > project_time(vector< vector< TH2D * > > ToT_vs_time)
Projects the hit time vs ToT histograms for all the combinations nanobeacon - pmt on the X axis...
int overall_delay_16ns
vector< vector< TH2D * > > allocate_histograms(JDetector detector, IO options, double c_w, vector< int > top_pmts, vector< int > bottom_pmts, int max_distance)
Allocates histograms to store the tot vs time hit distribution for all the combinations nanobeacon - ...
JMultipleFileScanner ifnames
Basic data structure for L0 hit.
int getFloor() const
Get floor number.
int getString() const
Get string number.
vector< vector< TH1D * > > project_ToT(vector< vector< TH2D * > > ToT_vs_time)
Projects the hit time vs ToT histograms for all the combinations nanobeacon - pmt on the Y axis...
Detector file.
Definition: JHead.hh:126
Hit data structure.
Definition: JDAQHit.hh:36
Router for fast addressing of hits in KM3NETDAQ::JDAQTimeslice data structure as a function of the op...
const_iterator begin() const
Definition: JDAQFrame.hh:139
int getID() const
Get identifier.
Definition: JObjectID.hh:54
Data structure for PMT geometry and calibration.
Definition: JPMT.hh:52
int getNumberOfFloors(const JDetector &detector)
Get number of floors.
const JPosition3D & getPosition() const
Get position.
Definition: JPosition3D.hh:129
const JPMT & getPMT(const int index) const
Get PMT.
Definition: JModule.hh:141
Scanning of objects from multiple files according a format that follows from the extension of each fi...
const JModule & getModule(const JModuleAddress &address) const
Get module parameters.
Definition: JDetector.hh:180
Direct access to module in detector data structure.
virtual bool hasNext()
Check availability of next element.
General purpose class for object reading from a list of file names.
ROOT TTree parameter settings.
KM3NeT DAQ constants, bit handling, etc.
void fill_histograms(vector< vector< TH2D * > > &histograms, IO options, JDetector detector, vector< int > top_pmts, vector< int > bottom_pmts, int max_distance)
Fills the hit time vs ToT histograms for all the combinations nanobeacon - pmt.
Data frame of one optical module.
int get_offset(JPMT pmt, JModule NB, double c_w, IO options)
Calculates an offset for a combinataion nanobeacon-pmt that allows to define the histogram range in a...
Timeslice data structure for L0 data.
Logical location of module.
const_iterator end() const
Definition: JDAQFrame.hh:140
double getT0() const
Get time offset.