MiniballSort
Loading...
Searching...
No Matches
EventBuilder.cc
Go to the documentation of this file.
1#include "EventBuilder.hh"
2
3MiniballEventBuilder::MiniballEventBuilder( std::shared_ptr<MiniballSettings> myset ){
4
5 // First get the settings
6 set = myset;
7
8 // No calibration file by default
9 overwrite_cal = false;
10
11 // No input file at the start by default
12 flag_input_file = false;
13
14 // Progress bar starts as false
15 _prog_ = false;
16
17 // Start at MBS event 0
18 preveventid = 0;
19
20 // ------------------------------- //
21 // Initialise variables and flags //
22 // ------------------------------- //
23 build_window = set->GetEventWindow();
24
25 n_sfp.resize( set->GetNumberOfFebexSfps() );
26 n_dgf.resize( set->GetNumberOfDgfModules() );
27 n_adc.resize( set->GetNumberOfAdcModules() );
28
29 febex_time_start.resize( set->GetNumberOfFebexSfps() );
30 febex_time_stop.resize( set->GetNumberOfFebexSfps() );
31 febex_dead_time.resize( set->GetNumberOfFebexSfps() );
32 febex_time_ch.resize( set->GetNumberOfFebexSfps() );
33 sync_time.resize( set->GetNumberOfFebexSfps() );
34 pause_time.resize( set->GetNumberOfFebexSfps() );
35 resume_time.resize( set->GetNumberOfFebexSfps() );
36 n_board.resize( set->GetNumberOfFebexSfps() );
37 n_sync.resize( set->GetNumberOfFebexSfps() );
38 n_pause.resize( set->GetNumberOfFebexSfps() );
39 n_resume.resize( set->GetNumberOfFebexSfps() );
40 flag_pause.resize( set->GetNumberOfFebexSfps() );
41 flag_resume.resize( set->GetNumberOfFebexSfps() );
42
43 n_pulser.resize( set->GetNumberOfPulsers() );
44 pulser_time.resize( set->GetNumberOfPulsers() );
45 pulser_prev.resize( set->GetNumberOfPulsers() );
46
47 for( unsigned int i = 0; i < set->GetNumberOfFebexSfps(); ++i ) {
48
49 febex_time_start[i].resize( set->GetNumberOfFebexBoards() );
50 febex_time_stop[i].resize( set->GetNumberOfFebexBoards() );
51 febex_dead_time[i].resize( set->GetNumberOfFebexBoards() );
52 febex_time_ch[i].resize( set->GetNumberOfFebexBoards() );
53 sync_time[i].resize( set->GetNumberOfFebexBoards() );
54 pause_time[i].resize( set->GetNumberOfFebexBoards() );
55 resume_time[i].resize( set->GetNumberOfFebexBoards() );
56 n_board[i].resize( set->GetNumberOfFebexBoards() );
57 n_sync[i].resize( set->GetNumberOfFebexBoards() );
58 n_pause[i].resize( set->GetNumberOfFebexBoards() );
59 n_resume[i].resize( set->GetNumberOfFebexBoards() );
60 flag_pause[i].resize( set->GetNumberOfFebexBoards() );
61 flag_resume[i].resize( set->GetNumberOfFebexBoards() );
62
63 for( unsigned int j = 0; j < set->GetNumberOfFebexBoards(); ++j )
64 febex_time_ch[i][j].resize( set->GetNumberOfFebexChannels(), 0 );
65
66 }
67
68 // Histogrammer options
69 //TH1::AddDirectory(kFALSE);
70
71 // Intialise the hist list
72 histlist = new TList();
73
74}
75
77
78 // Call for every new file
79 // Reset counters etc.
80
81 time_prev = 0;
82 time_min = 0;
83 time_max = 0;
84 time_first = 0;
85 ebis_prev = 0;
86 t1_prev = 0;
87 sc_prev = 0;
88
89 n_dgf_data = 0;
90 n_adc_data = 0;
91 n_febex_data = 0;
92 n_info_data = 0;
93
94 n_ebis = 0;
95 n_rilis = 0;
96 n_t1 = 0;
97 n_sc = 0;
98
99 n_miniball = 0;
100 n_cd = 0;
101 n_pad = 0;
102 n_bd = 0;
103 n_spede = 0;
104 n_ic = 0;
105
106 gamma_ctr = 0;
107 gamma_ab_ctr = 0;
108 cd_ctr = 0;
109 bd_ctr = 0;
110 spede_ctr = 0;
111 ic_ctr = 0;
112
113 repeat_ctr = 0;
114
115
116 for( unsigned int i = 0; i < set->GetNumberOfPulsers(); ++i ) {
117
118 n_pulser[i] = 0;
119 pulser_time[i] = 0;
120 pulser_prev[i] = 0;
121
122 }
123
124 for( unsigned int i = 0; i < set->GetNumberOfFebexSfps(); ++i ) {
125
126 n_sfp[i] = 0;
127
128 for( unsigned int j = 0; j < set->GetNumberOfFebexBoards(); ++j ) {
129
130 febex_time_start[i][j] = 0;
131 febex_time_stop[i][j] = 0;
132 febex_dead_time[i][j] = 0;
133 sync_time[i][j] = 0;
134 pause_time[i][j] = 0;
135 resume_time[i][j] = 0;
136 n_board[i][j] = 0;
137 n_sync[i][j] = 0;
138 n_pause[i][j] = 0;
139 n_resume[i][j] = 0;
140 flag_pause[i][j] = false;
141 flag_resume[i][j] = false;
142
143 }
144
145 }
146
147 for( unsigned int i = 0; i < set->GetNumberOfAdcModules(); ++i )
148 n_adc[i] = 0;
149
150 for( unsigned int i = 0; i < set->GetNumberOfDgfModules(); ++i )
151 n_dgf[i] = 0;
152
153}
154
155void MiniballEventBuilder::SetInputFile( std::string input_file_name ) {
156
157 // Open next Root input file.
158 input_file = new TFile( input_file_name.data(), "read" );
159 if( input_file->IsZombie() ) {
160
161 std::cout << "Cannot open " << input_file_name << std::endl;
162 return;
163
164 }
165
166 flag_input_file = true;
167
168 // Set the input tree
169 SetInputTree( (TTree*)input_file->Get("mb_sort") );
170 SetMBSInfoTree( (TTree*)input_file->Get("mbsinfo") );
171 StartFile();
172
173 return;
174
175}
176
177void MiniballEventBuilder::SetInputTree( TTree *user_tree ){
178
179 // Find the tree and set branch addresses
180 input_tree = user_tree;
181 in_data = nullptr;
182 input_tree->SetBranchAddress( "data", &in_data );
183
184 return;
185
186}
187
189
190 // Find the tree and set branch addresses
191 mbsinfo_tree = user_tree;
192 mbs_info = nullptr;
193 mbsinfo_tree->SetBranchAddress( "mbsinfo", &mbs_info );
194
195 return;
196
197}
198
199void MiniballEventBuilder::SetOutput( std::string output_file_name, bool cWrite ) {
200
201 // These are the branches we need
202 write_evts = std::make_unique<MiniballEvts>();
203 gamma_evt = std::make_shared<GammaRayEvt>();
204 gamma_ab_evt = std::make_shared<GammaRayAddbackEvt>();
205 particle_evt = std::make_shared<ParticleEvt>();
206 spede_evt = std::make_shared<SpedeEvt>();
207 bd_evt = std::make_shared<BeamDumpEvt>();
208 ic_evt = std::make_shared<IonChamberEvt>();
209
210 // ------------------------------------------------------------------------ //
211 // Create output file and create events tree
212 // ------------------------------------------------------------------------ //
213 output_file = new TFile( output_file_name.data(), "recreate" );
214 output_tree = new TTree( "evt_tree", "evt_tree" );
215 output_tree->Branch( "MiniballEvts", "MiniballEvts", write_evts.get(), 256000, 0 );
216 output_tree->SetAutoFlush();
217
218 // Create log file.
219 std::string log_file_name = output_file_name.substr( 0, output_file_name.find_last_of(".") );
220 log_file_name += ".log";
221 log_file.open( log_file_name.data(), std::ios::app );
222
223 // Hisograms in separate function
225
226 // Write once at the start
227 if( cWrite ) output_file->Write();
228
229}
230
232
234
235 flag_close_event = false;
236 event_open = false;
237
238 hit_ctr = 0;
239
240 std::vector<float>().swap(mb_en_list);
241 std::vector<unsigned long long>().swap(mb_ts_list);
242 std::vector<unsigned char>().swap(mb_clu_list);
243 std::vector<unsigned char>().swap(mb_cry_list);
244 std::vector<unsigned char>().swap(mb_seg_list);
245
246 std::vector<float>().swap(cd_en_list);
247 std::vector<unsigned long long>().swap(cd_ts_list);
248 std::vector<unsigned char>().swap(cd_det_list);
249 std::vector<unsigned char>().swap(cd_sec_list);
250 std::vector<unsigned char>().swap(cd_side_list);
251 std::vector<unsigned char>().swap(cd_strip_list);
252
253 std::vector<float>().swap(pad_en_list);
254 std::vector<unsigned long long>().swap(pad_ts_list);
255 std::vector<unsigned char>().swap(pad_det_list);
256 std::vector<unsigned char>().swap(pad_sec_list);
257
258 std::vector<float>().swap(bd_en_list);
259 std::vector<unsigned long long>().swap(bd_ts_list);
260 std::vector<unsigned char>().swap(bd_det_list);
261
262 std::vector<float>().swap(spede_en_list);
263 std::vector<unsigned long long>().swap(spede_ts_list);
264 std::vector<unsigned char>().swap(spede_seg_list);
265
266 std::vector<float>().swap(ic_en_list);
267 std::vector<unsigned long long>().swap(ic_ts_list);
268 std::vector<unsigned char>().swap(ic_id_list);
269
270 write_evts->ClearEvt();
271
272 return;
273
274}
275
276
278
279 std::string hname, htitle;
280 std::string dirname;
281
282 // ----------------- //
283 // Timing histograms //
284 // ----------------- //
285 dirname = "timing";
286 if( !output_file->GetDirectory( dirname.data() ) )
287 output_file->mkdir( dirname.data() );
288 output_file->cd( dirname.data() );
289
290 tdiff = new TH1F( "tdiff", "Time difference to first trigger;#Delta t [ns]", 1e3, -10, 1e5 );
291 tdiff_clean = new TH1F( "tdiff_clean", "Time difference to first trigger without noise;#Delta t [ns]", 1e3, -10, 1e5 );
292 histlist->Add(tdiff);
293 histlist->Add(tdiff_clean);
294
295 pulser_freq = new TProfile( "pulser_freq", "Frequency of pulser in FEBEX DAQ as a function of time;time [ns];f [Hz]", 10.8e4, 0, 10.8e12 );
296 pulser_period = new TH1F( "pulser_period", "Period of pulser in FEBEX DAQ;T [ns]", 10e3, 0, 10e9 );
297 int npulserbins = set->GetNumberOfPulsers()-1;
298 if( npulserbins <= 0 ) npulserbins = 1;
299 pulser_tdiff = new TH2F( "pulser_tdiff", "Time difference of pulser 0 to all other pulsers in FEBEX DAQ;Pulser ID;{#Delta}t [ns]",
300 npulserbins, 0.5, set->GetNumberOfPulsers()-0.5, 2001, -10005, 10005 );
301 ebis_freq = new TProfile( "ebis_freq", "Frequency of EBIS events as a function of time;time [ns];f [Hz]", 10.8e4, 0, 10.8e12 );
302 ebis_period = new TH1F( "ebis_period", "Period of EBIS events;T [ns]", 10e3, 0, 10e9 );
303 t1_freq = new TProfile( "t1_freq", "Frequency of T1 events (p+ on ISOLDE target) as a function of time;time [ns];f [Hz]", 10.8e4, 0, 10.8e12 );
304 t1_period = new TH1F( "t1_period", "Period of T1 events (p+ on ISOLDE target);T [ns]", 10e3, 0, 10e9 );
305 sc_freq = new TProfile( "sc_freq", "Frequency of SuperCycle events as a function of time;time [ns];f [Hz]", 10.8e4, 0, 10.8e12 );
306 sc_period = new TH1F( "sc_period", "Period of SuperCycle events;T [ns]", 10e3, 0, 10e9 );
307 histlist->Add(pulser_freq);
310 histlist->Add(ebis_freq);
311 histlist->Add(ebis_period);
312 histlist->Add(t1_freq);
313 histlist->Add(t1_period);
314 histlist->Add(sc_freq);
315 histlist->Add(sc_period);
316
317 // ------------------- //
318 // Miniball histograms //
319 // ------------------- //
320 dirname = "miniball";
321 if( !output_file->GetDirectory( dirname.data() ) )
322 output_file->mkdir( dirname.data() );
323 output_file->cd( dirname.data() );
324
325 mb_td_core_seg = new TH1F( "mb_td_core_seg", "Time difference between core and segment in same crystal;#Delta t [ns]", 499, -2495, 2495 );
326 mb_td_core_core = new TH1F( "mb_td_core_core", "Time difference between two cores in same cluster;#Delta t [ns]", 499, -2495, 2495 );
329
330 mb_en_core_seg.resize( set->GetNumberOfMiniballClusters() );
331 mb_en_core_seg_ebis_on.resize( set->GetNumberOfMiniballClusters() );
332
333 for( unsigned int i = 0; i < set->GetNumberOfMiniballClusters(); ++i ) {
334
335 dirname = "miniball/cluster_" + std::to_string(i);
336 if( !output_file->GetDirectory( dirname.data() ) )
337 output_file->mkdir( dirname.data() );
338 output_file->cd( dirname.data() );
339
340 mb_en_core_seg[i].resize( set->GetNumberOfMiniballCrystals() );
341 mb_en_core_seg_ebis_on[i].resize( set->GetNumberOfMiniballCrystals() );
342
343 for( unsigned int j = 0; j < set->GetNumberOfMiniballCrystals(); ++j ) {
344
345 hname = "mb_en_core_seg_" + std::to_string(i) + "_";
346 hname += std::to_string(j);
347 htitle = "Gamma-ray spectrum from cluster " + std::to_string(i);
348 htitle += " core " + std::to_string(j) + ", gated by segment ";
349 htitle += " (multiplicity = 1 only);segment ID;Energy (keV)";
350 mb_en_core_seg[i][j] = new TH2F( hname.data(), htitle.data(), 7, -0.5, 6.5, 4096, -0.5, 4095.5 );
351 histlist->Add(mb_en_core_seg[i][j]);
352
353 hname = "mb_en_core_seg_" + std::to_string(i) + "_";
354 hname += std::to_string(j) + "_ebis_on";
355 htitle = "Gamma-ray spectrum from cluster " + std::to_string(i);
356 htitle += " core " + std::to_string(j) + ", gated by segment ";
357 htitle += " (multiplicity = 1 only) gated by EBIS time (1.5 ms);segment ID;Energy (keV)";
358 mb_en_core_seg_ebis_on[i][j] = new TH2F( hname.data(), htitle.data(), 7, -0.5, 6.5, 4096, -0.5, 4095.5 );
360
361 } // j
362
363 } // i
364
365 // ------------- //
366 // CD histograms //
367 // ------------- //
368 dirname = "cd";
369 if( !output_file->GetDirectory( dirname.data() ) )
370 output_file->mkdir( dirname.data() );
371 output_file->cd( dirname.data() );
372
373 cd_pen_id.resize( set->GetNumberOfCDDetectors() );
374 cd_nen_id.resize( set->GetNumberOfCDDetectors() );
375 cd_pn_1v1.resize( set->GetNumberOfCDDetectors() );
376 cd_pn_1v2.resize( set->GetNumberOfCDDetectors() );
377 cd_pn_2v1.resize( set->GetNumberOfCDDetectors() );
378 cd_pn_2v2.resize( set->GetNumberOfCDDetectors() );
379 cd_pn_td.resize( set->GetNumberOfCDDetectors() );
380 cd_pp_td.resize( set->GetNumberOfCDDetectors() );
381 cd_ppad_td.resize( set->GetNumberOfCDDetectors() );
382 cd_nn_td.resize( set->GetNumberOfCDDetectors() );
383 cd_pn_mult.resize( set->GetNumberOfCDDetectors() );
384 cd_ppad_mult.resize( set->GetNumberOfCDDetectors() );
385 pad_en_id.resize( set->GetNumberOfCDDetectors() );
386
387 for( unsigned int i = 0; i < set->GetNumberOfCDDetectors(); ++i ) {
388
389 cd_pen_id[i].resize( set->GetNumberOfCDSectors() );
390 cd_nen_id[i].resize( set->GetNumberOfCDSectors() );
391 cd_pn_1v1[i].resize( set->GetNumberOfCDSectors() );
392 cd_pn_1v2[i].resize( set->GetNumberOfCDSectors() );
393 cd_pn_2v1[i].resize( set->GetNumberOfCDSectors() );
394 cd_pn_2v2[i].resize( set->GetNumberOfCDSectors() );
395 cd_pn_td[i].resize( set->GetNumberOfCDSectors() );
396 cd_pp_td[i].resize( set->GetNumberOfCDSectors() );
397 cd_ppad_td[i].resize( set->GetNumberOfCDSectors() );
398 cd_nn_td[i].resize( set->GetNumberOfCDSectors() );
399 cd_pn_mult[i].resize( set->GetNumberOfCDSectors() );
400 cd_ppad_mult[i].resize( set->GetNumberOfCDSectors() );
401
402 for( unsigned int j = 0; j < set->GetNumberOfCDSectors(); ++j ) {
403
404 hname = "cd_pen_id_" + std::to_string(i) + "_" + std::to_string(j);
405 htitle = "CD p-side energy for detector " + std::to_string(i);
406 htitle += ", sector " + std::to_string(j);
407 htitle += ";Strip ID;Energy (keV);Counts per strip, per 500 keV";
408 cd_pen_id[i][j] = new TH2F( hname.data(), htitle.data(),
409 set->GetNumberOfCDPStrips(), -0.5, set->GetNumberOfCDPStrips() - 0.5,
410 4000, 0, 2000e3 );
411 histlist->Add(cd_pen_id[i][j]);
412
413 hname = "cd_nen_id_" + std::to_string(i) + "_" + std::to_string(j);
414 htitle = "CD n-side energy for detector " + std::to_string(i);
415 htitle += ", sector " + std::to_string(j);
416 htitle += ";Strip ID;Energy (keV);Counts per strip, per 500 keV";
417 cd_nen_id[i][j] = new TH2F( hname.data(), htitle.data(),
418 set->GetNumberOfCDNStrips(), -0.5, set->GetNumberOfCDNStrips() - 0.5,
419 4000, 0, 2000e3 );
420 histlist->Add(cd_nen_id[i][j]);
421
422 hname = "cd_pn_1v1_" + std::to_string(i) + "_" + std::to_string(j);
423 htitle = "CD p-side vs n-side energy, multiplicity 1v1";
424 htitle += "for detector " + std::to_string(i);
425 htitle += ", sector " + std::to_string(j);
426 htitle += ";p-side Energy (keV);n-side Energy (keV);Counts";
427 cd_pn_1v1[i][j] = new TH2F( hname.data(), htitle.data(), 4000, 0, 200e3, 400, 0, 200e3 );
428 histlist->Add(cd_pn_1v1[i][j]);
429
430 hname = "cd_pn_1v2_" + std::to_string(i) + "_" + std::to_string(j);
431 htitle = "CD p-side vs n-side energy, multiplicity 1v2";
432 htitle += "for detector " + std::to_string(i);
433 htitle += ", sector " + std::to_string(j);
434 htitle += ";p-side Energy (keV);n-side Energy (keV);Counts";
435 cd_pn_1v2[i][j] = new TH2F( hname.data(), htitle.data(), 4000, 0, 200e3, 400, 0, 200e3 );
436 histlist->Add(cd_pn_1v2[i][j]);
437
438 hname = "cd_pn_2v1_" + std::to_string(i) + "_" + std::to_string(j);
439 htitle = "CD p-side vs n-side energy, multiplicity 2v1";
440 htitle += "for detector " + std::to_string(i);
441 htitle += ", sector " + std::to_string(j);
442 htitle += ";p-side Energy (keV);n-side Energy (keV);Counts";
443 cd_pn_2v1[i][j] = new TH2F( hname.data(), htitle.data(), 4000, 0, 200e3, 400, 0, 200e3 );
444 histlist->Add(cd_pn_2v1[i][j]);
445
446 hname = "cd_pn_2v2_" + std::to_string(i) + "_" + std::to_string(j);
447 htitle = "CD p-side vs n-side energy, multiplicity 2v2";
448 htitle += "for detector " + std::to_string(i);
449 htitle += ", sector " + std::to_string(j);
450 htitle += ";p-side Energy (keV);n-side Energy (keV);Counts";
451 cd_pn_2v2[i][j] = new TH2F( hname.data(), htitle.data(), 4000, 0, 200e3, 400, 0, 200e3 );
452 histlist->Add(cd_pn_2v2[i][j]);
453
454 hname = "cd_pn_td_" + std::to_string(i) + "_" + std::to_string(j);
455 htitle = "CD p-side vs n-side time difference ";
456 htitle += "for detector " + std::to_string(i);
457 htitle += ", sector " + std::to_string(j);
458 htitle += ";time difference (ns);Counts per 10 ns";
459 cd_pn_td[i][j] = new TH1F( hname.data(), htitle.data(), 800, -4e3, 4e3 );
460 histlist->Add(cd_pn_td[i][j]);
461
462 hname = "cd_pp_td_" + std::to_string(i) + "_" + std::to_string(j);
463 htitle = "CD p-side vs p-side time difference ";
464 htitle += "for detector " + std::to_string(i);
465 htitle += ", sector " + std::to_string(j);
466 htitle += ";time difference (ns);Counts per 10 ns";
467 cd_pp_td[i][j] = new TH1F( hname.data(), htitle.data(), 800, -4e3, 4e3 );
468 histlist->Add(cd_pp_td[i][j]);
469
470 hname = "cd_nn_td_" + std::to_string(i) + "_" + std::to_string(j);
471 htitle = "CD n-side vs n-side time difference ";
472 htitle += "for detector " + std::to_string(i);
473 htitle += ", sector " + std::to_string(j);
474 htitle += ";time difference (ns);Counts per 10 ns";
475 cd_nn_td[i][j] = new TH1F( hname.data(), htitle.data(), 800, -4e3, 4e3 );
476 histlist->Add(cd_nn_td[i][j]);
477
478 hname = "cd_ppad_td_" + std::to_string(i) + "_" + std::to_string(j);
479 htitle = "CD p-side vs pad time difference ";
480 htitle += "for detector " + std::to_string(i);
481 htitle += ", sector " + std::to_string(j);
482 htitle += ";time difference (ns);Counts per 10 ns";
483 cd_ppad_td[i][j] = new TH1F( hname.data(), htitle.data(), 800, -4e3, 4e3 );
484 histlist->Add(cd_ppad_td[i][j]);
485
486 hname = "cd_pn_mult_" + std::to_string(i) + "_" + std::to_string(j);
487 htitle = "CD p-side vs n-side multiplicity ";
488 htitle += "for detector " + std::to_string(i);
489 htitle += ", sector " + std::to_string(j);
490 htitle += ";p-side multiplicity;n-side multiplicity";
491 cd_pn_mult[i][j] = new TH2F( hname.data(), htitle.data(), 10, -0.5, 9.5, 10, -0.5, 9.5 );
492 histlist->Add(cd_pn_mult[i][j]);
493
494 hname = "cd_ppad_mult_" + std::to_string(i) + "_" + std::to_string(j);
495 htitle = "CD vs Pad multiplicity ";
496 htitle += "for detector " + std::to_string(i);
497 htitle += ", sector " + std::to_string(j);
498 htitle += ";CD multiplicity;Pad multiplicity";
499 cd_ppad_mult[i][j] = new TH2F( hname.data(), htitle.data(), 10, -0.5, 9.5, 10, -0.5, 9.5 );
500 histlist->Add(cd_ppad_mult[i][j]);
501
502 } // j
503
504 hname = "pad_en_" + std::to_string(i);
505 htitle = "Pad energy for each sector in detector " + std::to_string(i);
506 htitle += ";Sector;Energy (keV);Counts per quadrant, per 25 keV";
507 pad_en_id[i] = new TH2F( hname.data(), htitle.data(),
508 set->GetNumberOfCDSectors(), -0.5, set->GetNumberOfCDSectors() - 0.5,
509 8000, 0, 200e3 );
510 histlist->Add(pad_en_id[i]);
511
512 } // i
513
514
515 // --------------------- //
516 // IonChamber histograms //
517 // --------------------- //
518 dirname = "ic";
519 if( !output_file->GetDirectory( dirname.data() ) )
520 output_file->mkdir( dirname.data() );
521 output_file->cd( dirname.data() );
522
523 ic_td = new TH1F( "ic_td", "Time difference between signals in the ion chamber;#Delta t [ns]", 499, -2495, 2495 );
524 ic_dE = new TH1F( "ic_dE", "Ionisation chamber;Energy in first layer (Gas) (arb. units);Counts", 4096, 0, 10000 );
525 ic_E = new TH1F( "ic_E", "Ionisation chamber;Energy in final layer (Si) (arb. units);Counts", 4096, 0, 10000 );
526 ic_dE_E = new TH2F( "ic_dE_E", "Ionisation chamber;Rest energy, E (arb. units);Energy Loss, dE (arb. units);Counts", 4096, 0, 10000, 4096, 0, 10000 );
527 histlist->Add(ic_td);
528 histlist->Add(ic_dE);
529 histlist->Add(ic_E);
530 histlist->Add(ic_dE_E);
531
532
533 // flag to denote that hists are ready (used for spy)
534 hists_ready = true;
535
536 return;
537
538}
539
540// Reset histograms in the DataSpy
542
543 // Loop over the hist list
544 TIter next( histlist->MakeIterator() );
545 while( TObject *obj = next() ) {
546
547 if( obj->InheritsFrom( "TH2" ) )
548 ( (TH2*)obj )->Reset("ICESM");
549 else if( obj->InheritsFrom( "TH1" ) )
550 ( (TH1*)obj )->Reset("ICESM");
551
552 }
553
554 return;
555
556}
557
558
560
561 // Temporary variables for addback
562 unsigned long long MaxTime; // time of event with maximum energy
563 unsigned char MaxSegId; // segment with maximum energy
564 unsigned char MaxCryId; // crystal with maximum energy
565 float MaxEnergy; // maximum core energy
566 float MaxSegEnergy; // maximum segment energy
567 float SegSumEnergy; // add segment energies
568 float AbSumEnergy; // add core energies for addback
569 unsigned char seg_mul; // segment multiplicity
570 unsigned char ab_mul; // addback multiplicity
571 std::vector<unsigned char> ab_index; // index of addback already used
572
573 // Loop over all the events in Miniball detectors
574 for( unsigned int i = 0; i < mb_en_list.size(); ++i ) {
575
576 // Check if it's a core event
577 if( mb_seg_list.at(i) != 0 ) continue;
578
579 // Segment veto start as false
580 bool segment_veto = false;
581
582 // Reset addback variables
583 MaxSegId = 0; // initialise as core (if no segment hit (dead), use core!)
584 MaxSegEnergy = 0.;
585 SegSumEnergy = 0.;
586 seg_mul = 0;
587
588 // Loop again to find the matching segments
589 for( unsigned int j = 0; j < mb_en_list.size(); ++j ) {
590
591 // Skip the same event
592 if( i == j ) continue;
593
594 // Skip if it's a core again, also fill time diff plot
595 if( mb_seg_list.at(j) == 0 ) {
596
597 // Fill the time difference spectrum
598 mb_td_core_core->Fill( (long long)mb_ts_list.at(i) - (long long)mb_ts_list.at(j) );
599 continue;
600
601 }
602
603 // Skip if it's not the same crystal and cluster
604 if( mb_clu_list.at(i) != mb_clu_list.at(j) ||
605 mb_cry_list.at(i) != mb_cry_list.at(j) ) continue;
606
607 // Check for a vetoed segment
608 if( set->IsMiniballSegmentVetoed( mb_clu_list.at(i), mb_cry_list.at(i), mb_seg_list.at(j) ) ) {
609
610 segment_veto = true;
611 break;
612
613 }
614
615 // Fill the time difference spectrum
616 mb_td_core_seg->Fill( (long long)mb_ts_list.at(i) - (long long)mb_ts_list.at(j) );
617
618 // Skip if we are outside of the hit window
619 if( TMath::Abs( (double)mb_ts_list.at(i) - (double)mb_ts_list.at(j) )
620 > set->GetMiniballCrystalHitWindow() ) continue;
621
622 // Increment the segment multiplicity and sum energy
623 seg_mul++;
624 SegSumEnergy += mb_en_list.at(j);
625
626 // Is this bigger than the current maximum energy?
627 if( mb_en_list.at(j) > MaxSegEnergy ){
628
629 MaxSegEnergy = mb_en_list.at(j);
630 MaxSegId = mb_seg_list.at(j);
631
632 }
633
634 } // j: matching segments
635
636
637 // If any one of the segments that triggered are being vetoed, skip this gamma ray
638 if( segment_veto ) continue;
639
640 // Fill the segment spectra with core energies
641 mb_en_core_seg[mb_clu_list.at(i)][mb_cry_list.at(i)]->Fill( MaxSegId, mb_en_list.at(i) );
642 if( mb_ts_list.at(i) - ebis_time < 1.5e6 )
643 mb_en_core_seg_ebis_on[mb_clu_list.at(i)][mb_cry_list.at(i)]->Fill( MaxSegId, mb_en_list.at(i) );
644
645
646 //if( MaxSegId == 0 && mb_en_list.at(i) > 150.0 )
647 // std::cout << std::endl << mb_en_list.at(i) << "\t" << MaxSegEnergy << std::endl;
648
649 // Build the single crystal gamma-ray event
650 gamma_ctr++;
651 gamma_evt->SetEnergy( mb_en_list.at(i) );
652 gamma_evt->SetSegmentMaxEnergy( MaxSegEnergy );
653 gamma_evt->SetSegmentSumEnergy( SegSumEnergy );
654 gamma_evt->SetSegmentMultiplicity( seg_mul );
655 gamma_evt->SetAddbackMultiplicity( 1 );
656 gamma_evt->SetCluster( mb_clu_list.at(i) );
657 gamma_evt->SetCrystal( mb_cry_list.at(i) );
658 gamma_evt->SetSegment( MaxSegId );
659 gamma_evt->SetTime( mb_ts_list.at(i) );
660 write_evts->AddEvt( gamma_evt );
661
662 } // i: core events
663
664
665 // Loop over all the gamma-ray singles for addback
666 for( unsigned int i = 0; i < write_evts->GetGammaRayMultiplicity(); ++i ) {
667
668 // Check we haven't already used this event
669 if( std::find( ab_index.begin(), ab_index.end(), i ) != ab_index.end() )
670 continue;
671
672 // Reset addback variables
673 AbSumEnergy = write_evts->GetGammaRayEvt(i)->GetEnergy();
674 MaxCryId = write_evts->GetGammaRayEvt(i)->GetCrystal();
675 MaxSegId = write_evts->GetGammaRayEvt(i)->GetSegment();
676 MaxEnergy = AbSumEnergy;
677 MaxSegEnergy = write_evts->GetGammaRayEvt(i)->GetSegmentMaxEnergy();
678 SegSumEnergy = write_evts->GetGammaRayEvt(i)->GetSegmentSumEnergy();
679 MaxTime = write_evts->GetGammaRayEvt(i)->GetTime();
680 seg_mul = write_evts->GetGammaRayEvt(i)->GetSegmentMultiplicity();
681 ab_mul = 1; // this is already the first event
682
683 // Loop to find a matching event for addback
684 for( unsigned int j = i+1; j < write_evts->GetGammaRayMultiplicity(); ++j ) {
685
686 // Make sure we are in the same cluster
687 // In the future we might consider a more intelligent
688 // algorithm, which uses the line-of-sight idea
689 if( write_evts->GetGammaRayEvt(i)->GetCluster() !=
690 write_evts->GetGammaRayEvt(j)->GetCluster() ) continue;
691
692 // Skip if we are outside of the hit window
693 if( TMath::Abs( (double)write_evts->GetGammaRayEvt(i)->GetTime() -
694 (double)write_evts->GetGammaRayEvt(j)->GetTime() )
695 > set->GetMiniballAddbackHitWindow() ) continue;
696
697 // Check we haven't already used this event
698 if( std::find( ab_index.begin(), ab_index.end(), j ) != ab_index.end() )
699 continue;
700
701 // Then we can add them back
702 ab_mul++;
703 AbSumEnergy += write_evts->GetGammaRayEvt(j)->GetEnergy();
704 SegSumEnergy += write_evts->GetGammaRayEvt(j)->GetSegmentSumEnergy();
705 seg_mul += write_evts->GetGammaRayEvt(j)->GetSegmentMultiplicity();
706 ab_index.push_back(j);
707
708 // Is this bigger than the current maximum energy?
709 if( write_evts->GetGammaRayEvt(j)->GetEnergy() > MaxEnergy ){
710
711 MaxEnergy = write_evts->GetGammaRayEvt(j)->GetEnergy();
712 MaxSegEnergy = write_evts->GetGammaRayEvt(j)->GetSegmentMaxEnergy();
713 MaxCryId = write_evts->GetGammaRayEvt(j)->GetCrystal();
714 MaxSegId = write_evts->GetGammaRayEvt(j)->GetSegment();
715 MaxTime = write_evts->GetGammaRayEvt(j)->GetTime();
716
717 }
718
719 } // j: loop for matching addback
720
721 // Build the single crystal gamma-ray event
722 gamma_ab_ctr++;
723 gamma_ab_evt->SetEnergy( AbSumEnergy );
724 gamma_ab_evt->SetSegmentMaxEnergy( MaxSegEnergy );
725 gamma_ab_evt->SetSegmentSumEnergy( SegSumEnergy );
726 gamma_ab_evt->SetSegmentMultiplicity( seg_mul );
727 gamma_ab_evt->SetAddbackMultiplicity( ab_mul );
728 gamma_ab_evt->SetCluster( write_evts->GetGammaRayEvt(i)->GetCluster() );
729 gamma_ab_evt->SetCrystal( MaxCryId );
730 gamma_ab_evt->SetSegment( MaxSegId );
731 gamma_ab_evt->SetTime( MaxTime );
732 write_evts->AddEvt( gamma_ab_evt );
733
734 } // i: gamma-ray singles
735
736 return;
737
738}
739
740
742
743 // Variables for the finder algorithm
744 std::vector<unsigned char> pindex;
745 std::vector<unsigned char> nindex;
746
747 // Loop over each detector and sector
748 for( unsigned int i = 0; i < set->GetNumberOfCDDetectors(); ++i ){
749
750 for( unsigned int j = 0; j < set->GetNumberOfCDSectors(); ++j ){
751
752 // Reset variables for a new detector element
753 pindex.clear();
754 nindex.clear();
755 std::vector<unsigned char>().swap(pindex);
756 std::vector<unsigned char>().swap(nindex);
757 int pmax_idx = -1, nmax_idx = -1;
758 float pmax_en = -999., nmax_en = -999.;
759 float pad_coinc_en = 0.0;
760 float psum_en, nsum_en;
761 unsigned long long pad_coinc_ts = 0;
762 unsigned int padmult = 0;
763
764 // Calculate p/n side multiplicities and get indicies
765 for( unsigned int k = 0; k < cd_en_list.size(); ++k ){
766
767 // Test that we have the correct detector and quadrant
768 if( i != cd_det_list.at(k) || j != cd_sec_list.at(k) )
769 continue;
770
771 // Check max energy and push back the multiplicity
772 if( cd_side_list.at(k) == 0 ) {
773
774 pindex.push_back(k);
775
776 // Check if it is max energy
777 if( cd_en_list.at(k) > pmax_en ){
778
779 pmax_en = cd_en_list.at(k);
780 pmax_idx = k;
781
782 }
783
784 } // p-side
785
786 else if( cd_side_list.at(k) == 1 ) {
787
788 nindex.push_back(k);
789
790 // Check if it is max energy
791 if( cd_en_list.at(k) > nmax_en ){
792
793 nmax_en = cd_en_list.at(k);
794 nmax_idx = k;
795
796 }
797
798 } // n-side
799
800 } // k: all CD events
801
802 // Look for pad events
803 for( unsigned int k = 0; k < pad_en_list.size(); ++k ){
804
805 // Test that we have the correct detector and quadrant
806 //if( i != pad_det_list.at(k) || j != pad_sec_list.at(k) )
807 // continue;
808
809 // The following is a hack because of the cabling of the Pad
810 // detector in September 2023 for the IS656 run
811 //if( i != pad_det_list.at(k) ) continue;
812 //if( ( j == 0 || j == 3 ) && pad_sec_list.at(k) != 0 ) continue;
813 //if( ( j == 1 || j == 2 ) && pad_sec_list.at(k) != 1 ) continue;
814
815 // Count the pad multiplicity (panic if it is >1)
816 padmult++;
817
818 // Plot time differences
819 for( unsigned int p1 = 0; p1 < pindex.size(); ++p1 ){
820
821 cd_ppad_td[i][j]->Fill( (double)cd_ts_list.at( pindex[p1] ) -
822 (double)pad_ts_list.at(k) );
823
824 }
825
827 if( pmax_idx >= 0 ) {
828
829 if( TMath::Abs( (long long)pad_ts_list.at(k) - (long long)cd_ts_list.at( pmax_idx ) )
830 < set->GetPadHitWindow() ){
831
832 pad_coinc_en = pad_en_list.at(k);
833 pad_coinc_ts = pad_ts_list.at(k);
834
835 pad_en_id[i]->Fill( j, pad_coinc_en );
836
837 }
838
839 }
840
841 // Hack to recalibrate the pads that are coupled
842 // IS595 - 10th October 2023
843 //if( i == 0 && j == 1 ) pad_coinc_en *= 9.8;
844
845 } // k: all pad events
846
847
848 // Plot multiplcities
849 if( pindex.size() || nindex.size() )
850 cd_pn_mult[i][j]->Fill( pindex.size(), nindex.size() );
851
852 // Plot time differences
853 for( unsigned int p1 = 0; p1 < pindex.size(); ++p1 ){
854
855 for( unsigned int n1 = 0; n1 < nindex.size(); ++n1 ){
856
857 cd_pn_td[i][j]->Fill( (double)cd_ts_list.at( pindex[p1] ) -
858 (double)cd_ts_list.at( nindex[n1] ) );
859
860 } // n1
861
862 for( unsigned int p2 = p1+1; p2 < pindex.size(); ++p2 ){
863
864 cd_pp_td[i][j]->Fill( (double)cd_ts_list.at( pindex[p1] ) -
865 (double)cd_ts_list.at( pindex[p2] ) );
866
867 } // p2
868
869 } // p1
870
871 for( unsigned int n1 = 0; n1 < nindex.size(); ++n1 ){
872
873 for( unsigned int n2 = n1+1; n2 < nindex.size(); ++n2 ){
874
875 cd_nn_td[i][j]->Fill( (double)cd_ts_list.at( nindex[n1] ) -
876 (double)cd_ts_list.at( nindex[n2] ) );
877
878 } // n2
879
880 } // n1
881
882
883 // ----------------------- //
884 // Particle reconstruction //
885 // ----------------------- //
886 // 1 vs 1 - easiest situation
887 if( pindex.size() == 1 && nindex.size() == 1 ) {
888
889 // Set event
890 particle_evt->SetEnergyP( cd_en_list.at( pindex[0] ) );
891 particle_evt->SetEnergyN( cd_en_list.at( nindex[0] ) );
892 particle_evt->SetTimeP( cd_ts_list.at( pindex[0] ) );
893 particle_evt->SetTimeN( cd_ts_list.at( nindex[0] ) );
894 particle_evt->SetDetector( i );
895 particle_evt->SetSector( j );
896 particle_evt->SetStripP( cd_strip_list.at( pindex[0] ) );
897 particle_evt->SetStripN( cd_strip_list.at( nindex[0] ) );
898 particle_evt->SetEnergyPad( pad_coinc_en );
899 particle_evt->SetTimePad( pad_coinc_ts );
900
901 // Fill tree
902 write_evts->AddEvt( particle_evt );
903 cd_ctr++;
904
905 // Fill histograms
906 cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex[0] ),
907 cd_en_list.at( pindex[0] ) );
908 cd_nen_id[i][j]->Fill( cd_strip_list.at( nindex[0] ),
909 cd_en_list.at( nindex[0] ) );
910 cd_pn_1v1[i][j]->Fill( cd_en_list.at( pindex[0] ),
911 cd_en_list.at( nindex[0] ) );
912 cd_ppad_mult[i][j]->Fill( 1, padmult );
913
914 } // 1 vs 1
915
916 // 1 vs 2 - n-side charge sharing?
917 else if( pindex.size() == 1 && nindex.size() == 2 ) {
918
919 // Neighbour strips
920 if( TMath::Abs( cd_strip_list.at( nindex[0] ) - cd_strip_list.at( nindex[1] ) ) == 1 ) {
921
922 // Simple sum of both energies, cross-talk not included yet
923 nsum_en = cd_en_list.at( nindex[0] );
924 nsum_en += cd_en_list.at( nindex[1] );
925
926 // Set event
927 particle_evt->SetEnergyP( cd_en_list.at( pindex[0] ) );
928 particle_evt->SetEnergyN( nsum_en );
929 particle_evt->SetTimeP( cd_ts_list.at( pindex[0] ) );
930 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
931 particle_evt->SetDetector( i );
932 particle_evt->SetSector( j );
933 particle_evt->SetStripP( cd_strip_list.at( pindex[0] ) );
934 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
935 particle_evt->SetEnergyPad( pad_coinc_en );
936 particle_evt->SetTimePad( pad_coinc_ts );
937
938 // Fill tree
939 write_evts->AddEvt( particle_evt );
940 cd_ctr++;
941
942 // Fill histograms
943 cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex[0] ),
944 cd_en_list.at( pindex[0] ) );
945 cd_nen_id[i][j]->Fill( nsum_en,
946 cd_en_list.at( nmax_idx ) );
947 cd_pn_1v2[i][j]->Fill( cd_en_list.at( pindex[0] ),
948 cd_en_list.at( nindex[0] ) );
949 cd_pn_1v2[i][j]->Fill( cd_en_list.at( pindex[0] ),
950 cd_en_list.at( nindex[1] ) );
951 cd_ppad_mult[i][i]->Fill( 1, padmult );
952
953 } // neighbour strips
954
955 // otherwise treat as 1 vs 1
956 else {
957
958 // Set event
959 particle_evt->SetEnergyP( cd_en_list.at( pindex[0] ) );
960 particle_evt->SetEnergyN( cd_en_list.at( nmax_idx ) );
961 particle_evt->SetTimeP( cd_ts_list.at( pindex[0] ) );
962 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
963 particle_evt->SetDetector( i );
964 particle_evt->SetSector( j );
965 particle_evt->SetStripP( cd_strip_list.at( pindex[0] ) );
966 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
967 particle_evt->SetEnergyPad( pad_coinc_en );
968 particle_evt->SetTimePad( pad_coinc_ts );
969
970 // Fill tree
971 write_evts->AddEvt( particle_evt );
972 cd_ctr++;
973
974 // Fill histograms
975 cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex[0] ),
976 cd_en_list.at( pindex[0] ) );
977 cd_nen_id[i][j]->Fill( cd_strip_list.at( nmax_idx ),
978 cd_en_list.at( nmax_idx ) );
979 cd_ppad_mult[i][i]->Fill( 1, padmult );
980
981
982 } // treat as 1 vs 1
983
984 } // 1 vs 2
985
986 // 2 vs 1 - p-side charge sharing?
987 else if( pindex.size() == 2 && nindex.size() == 1 ) {
988
989 // Neighbour strips
990 if( TMath::Abs( cd_strip_list.at( pindex[0] ) - cd_strip_list.at( pindex[1] ) ) == 1 ) {
991
992 // Simple sum of both energies, cross-talk not included yet
993 psum_en = cd_en_list.at( pindex[0] );
994 psum_en += cd_en_list.at( pindex[1] );
995
996 // Set event
997 particle_evt->SetEnergyP( psum_en );
998 particle_evt->SetEnergyN( cd_en_list.at( nindex[0] ) );
999 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1000 particle_evt->SetTimeN( cd_ts_list.at( nindex[0] ) );
1001 particle_evt->SetDetector( i );
1002 particle_evt->SetSector( j );
1003 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1004 particle_evt->SetStripN( cd_strip_list.at( nindex[0] ) );
1005 particle_evt->SetEnergyPad( pad_coinc_en );
1006 particle_evt->SetTimePad( pad_coinc_ts );
1007
1008 // Fill tree
1009 write_evts->AddEvt( particle_evt );
1010 cd_ctr++;
1011
1012 // Fill histograms
1013 cd_pen_id[i][j]->Fill( psum_en,
1014 cd_en_list.at( pmax_idx ) );
1015 cd_nen_id[i][j]->Fill( cd_strip_list.at( nindex[0] ),
1016 cd_en_list.at( nindex[0] ) );
1017 cd_pn_2v1[i][j]->Fill( cd_en_list.at( pindex[0] ),
1018 cd_en_list.at( nindex[0] ) );
1019 cd_pn_2v1[i][j]->Fill( cd_en_list.at( pindex[1] ),
1020 cd_en_list.at( nindex[0] ) );
1021 cd_ppad_mult[i][i]->Fill( 1, padmult );
1022
1023 } // neighbour strips
1024
1025 // otherwise treat as 1 vs 1
1026 else {
1027
1028 // Set event
1029 particle_evt->SetEnergyP( cd_en_list.at( pmax_idx ) );
1030 particle_evt->SetEnergyN( cd_en_list.at( nindex[0] ) );
1031 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1032 particle_evt->SetTimeN( cd_ts_list.at( nindex[0] ) );
1033 particle_evt->SetDetector( i );
1034 particle_evt->SetSector( j );
1035 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1036 particle_evt->SetStripN( cd_strip_list.at( nindex[0] ) );
1037 particle_evt->SetEnergyPad( pad_coinc_en );
1038 particle_evt->SetTimePad( pad_coinc_ts );
1039
1040 // Fill tree
1041 write_evts->AddEvt( particle_evt );
1042 cd_ctr++;
1043
1044 // Fill histograms
1045 cd_pen_id[i][j]->Fill( cd_strip_list.at( pmax_idx ),
1046 cd_en_list.at( pmax_idx ) );
1047 cd_nen_id[i][j]->Fill( cd_strip_list.at( nindex[0] ),
1048 cd_en_list.at( nindex[0] ) );
1049 cd_ppad_mult[i][i]->Fill( 1, padmult );
1050
1051
1052 } // treat as 1 vs 1
1053
1054 } // 2 vs 1
1055
1056 // 2 vs 2 - charge sharing on both or two particles?
1057 else if( pindex.size() == 2 && nindex.size() == 2 ) {
1058
1059 // Neighbour strips - p-side + n-side
1060 if( TMath::Abs( cd_strip_list.at( pindex[0] ) - cd_strip_list.at( pindex[1] ) ) == 1 &&
1061 TMath::Abs( cd_strip_list.at( nindex[0] ) - cd_strip_list.at( nindex[1] ) ) == 1 ) {
1062
1063 // Simple sum of both energies, cross-talk not included yet
1064 psum_en = cd_en_list.at( pindex[0] );
1065 psum_en += cd_en_list.at( pindex[1] );
1066 nsum_en = cd_en_list.at( nindex[0] );
1067 nsum_en += cd_en_list.at( nindex[1] );
1068
1069 // Set event
1070 particle_evt->SetEnergyP( psum_en );
1071 particle_evt->SetEnergyN( nsum_en );
1072 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1073 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
1074 particle_evt->SetDetector( i );
1075 particle_evt->SetSector( j );
1076 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1077 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
1078 particle_evt->SetEnergyPad( pad_coinc_en );
1079 particle_evt->SetTimePad( pad_coinc_ts );
1080
1081 // Fill tree
1082 write_evts->AddEvt( particle_evt );
1083 cd_ctr++;
1084
1085 // Fill histograms
1086 cd_pen_id[i][j]->Fill( cd_strip_list.at( pmax_idx ),
1087 psum_en );
1088 cd_nen_id[i][j]->Fill( cd_strip_list.at( nmax_idx ),
1089 nsum_en );
1090 cd_pn_2v2[i][j]->Fill( cd_en_list.at( pindex[0] ),
1091 cd_en_list.at( nindex[0] ) );
1092 cd_pn_2v2[i][j]->Fill( cd_en_list.at( pindex[0] ),
1093 cd_en_list.at( nindex[1] ) );
1094 cd_pn_2v2[i][j]->Fill( cd_en_list.at( pindex[1] ),
1095 cd_en_list.at( nindex[0] ) );
1096 cd_pn_2v2[i][j]->Fill( cd_en_list.at( pindex[1] ),
1097 cd_en_list.at( nindex[1] ) );
1098 cd_ppad_mult[i][i]->Fill( 1, padmult );
1099
1100 } // neighbour strips - p-side + n-side
1101
1102 // Neighbour strips - p-side only
1103 else if( TMath::Abs( cd_strip_list.at( pindex[0] ) - cd_strip_list.at( pindex[1] ) ) == 1 ) {
1104
1105 // Simple sum of both energies, cross-talk not included yet
1106 psum_en = cd_en_list.at( pindex.at(0) );
1107 psum_en += cd_en_list.at( pindex.at(1) );
1108
1109 // Set event
1110 particle_evt->SetEnergyP( psum_en );
1111 particle_evt->SetEnergyN( nmax_en );
1112 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1113 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
1114 particle_evt->SetDetector( i );
1115 particle_evt->SetSector( j );
1116 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1117 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
1118 particle_evt->SetEnergyPad( pad_coinc_en );
1119 particle_evt->SetTimePad( pad_coinc_ts );
1120
1121 // Fill tree
1122 write_evts->AddEvt( particle_evt );
1123 cd_ctr++;
1124
1125 // Fill histograms
1126 cd_pen_id[i][j]->Fill( cd_strip_list.at( pmax_idx ),
1127 psum_en );
1128 cd_nen_id[i][j]->Fill( cd_strip_list.at( nmax_idx ),
1129 cd_en_list.at( nmax_idx ) );
1130 cd_ppad_mult[i][i]->Fill( 1, padmult );
1131
1132 } // neighbour strips - p-side only
1133
1134 // Neighbour strips - n-side only
1135 else if( TMath::Abs( cd_strip_list.at( nindex[0] ) - cd_strip_list.at( nindex[1] ) ) == 1 ) {
1136
1137 // Simple sum of both energies, cross-talk not included yet
1138 nsum_en = cd_en_list.at( nindex.at(0) );
1139 nsum_en += cd_en_list.at( nindex.at(1) );
1140
1141 // Set event
1142 particle_evt->SetEnergyP( cd_en_list.at( pmax_idx ) );
1143 particle_evt->SetEnergyN( nsum_en );
1144 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1145 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
1146 particle_evt->SetDetector( i );
1147 particle_evt->SetSector( j );
1148 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1149 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
1150 particle_evt->SetEnergyPad( pad_coinc_en );
1151 particle_evt->SetTimePad( pad_coinc_ts );
1152
1153 // Fill tree
1154 write_evts->AddEvt( particle_evt );
1155 cd_ctr++;
1156
1157 // Fill histograms
1158 cd_pen_id[i][j]->Fill( cd_strip_list.at( pmax_idx ),
1159 cd_en_list.at( pmax_idx ) );
1160 cd_nen_id[i][j]->Fill( cd_strip_list.at( nmax_idx ),
1161 nsum_en );
1162 cd_ppad_mult[i][i]->Fill( 1, padmult );
1163
1164 } // neighbour strips - n-side only
1165
1166 // Neither of them are neighbours... 2 events?
1167 else {
1168
1169 // Event 1 is with first p-side, but which n-side?
1170 unsigned int nfriend_idx = nindex.at(0);
1171 if( TMath::Abs( cd_en_list.at( pindex.at(0) ) - cd_en_list.at( nindex.at(1) ) )
1172 < TMath::Abs( cd_en_list.at( pindex.at(0) ) - cd_en_list.at( nindex.at(0) ) ) )
1173 nfriend_idx = nindex.at(1);
1174
1175 // Set event
1176 particle_evt->SetEnergyP( cd_en_list.at( pindex.at(0) ) );
1177 particle_evt->SetEnergyN( cd_en_list.at( nfriend_idx ) );
1178 particle_evt->SetTimeP( cd_ts_list.at( pindex.at(0) ) );
1179 particle_evt->SetTimeN( cd_ts_list.at( nfriend_idx ) );
1180 particle_evt->SetDetector( i );
1181 particle_evt->SetSector( j );
1182 particle_evt->SetStripP( cd_strip_list.at( pindex.at(0) ) );
1183 particle_evt->SetStripN( cd_strip_list.at( nfriend_idx ) );
1184 particle_evt->SetEnergyPad( pad_coinc_en );
1185 particle_evt->SetTimePad( pad_coinc_ts );
1186
1187 // Fill tree for first hit
1188 write_evts->AddEvt( particle_evt );
1189 cd_ctr++;
1190
1191 // Event 2 is with first p-side, but which n-side?
1192 if( nfriend_idx == nindex.at(1) ) nfriend_idx = nindex.at(0);
1193 else nfriend_idx = nindex.at(1);
1194
1195 // Set event
1196 particle_evt->SetEnergyP( cd_en_list.at( pindex.at(1) ) );
1197 particle_evt->SetEnergyN( cd_en_list.at( nfriend_idx ) );
1198 particle_evt->SetTimeP( cd_ts_list.at( pindex.at(1) ) );
1199 particle_evt->SetTimeN( cd_ts_list.at( nfriend_idx ) );
1200 particle_evt->SetDetector( i );
1201 particle_evt->SetSector( j );
1202 particle_evt->SetStripP( cd_strip_list.at( pindex.at(1) ) );
1203 particle_evt->SetStripN( cd_strip_list.at( nfriend_idx ) );
1204 particle_evt->SetEnergyPad( pad_coinc_en );
1205 particle_evt->SetTimePad( pad_coinc_ts );
1206
1207 // Fill tree for second hit
1208 write_evts->AddEvt( particle_evt );
1209 cd_ctr++;
1210
1211
1212 // Fill histograms
1213 cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex.at(0) ),
1214 cd_en_list.at( pindex.at(0) ) );
1215 cd_nen_id[i][j]->Fill( cd_strip_list.at( nindex.at(0) ),
1216 cd_en_list.at( nindex.at(0) ) );
1217 cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex.at(1) ),
1218 cd_en_list.at( pindex.at(1) ) );
1219 cd_nen_id[i][j]->Fill( cd_strip_list.at( nindex.at(1) ),
1220 cd_en_list.at( nindex.at(1) ) );
1221 cd_ppad_mult[i][i]->Fill( 2, padmult );
1222
1223 } // neighbour strips - n-side only
1224
1225 } // 2 vs 2
1226
1227// // 1 vs 0 - p-side only, do we carry on?
1228// if( pindex.size() == 1 && nindex.size() == 0 ) {
1229//
1230// // Set event
1231// particle_evt->SetEnergyP( cd_en_list.at( pindex[0] ) );
1232// particle_evt->SetEnergyN( cd_en_list.at( 0.0 ) );
1233// particle_evt->SetTimeP( cd_ts_list.at( pindex[0] ) );
1234// particle_evt->SetTimeN( cd_ts_list.at( pindex[0] ) );
1235// particle_evt->SetDetector( i );
1236// particle_evt->SetSector( j );
1237// particle_evt->SetStripP( cd_strip_list.at( pindex[0] ) );
1238// particle_evt->SetStripN( 5.0 );
1239//
1240// // Fill tree
1241// write_evts->AddEvt( particle_evt );
1242// cd_ctr++;
1243//
1244// // Fill histograms
1245// cd_pen_id[i][j]->Fill( cd_strip_list.at( pindex[0] ),
1246// cd_en_list.at( pindex[0] ) );
1247//
1248// } // 1 vs 0
1249//
1250// // 2 vs 0 - p-side charge sharing? and no n-side
1251// else if( pindex.size() == 2 && nindex.size() == 0 ) {
1252//
1253// // Neighbour strips
1254// if( TMath::Abs( cd_strip_list.at( pindex[0] ) - cd_strip_list.at( pindex[1] ) ) == 1 ) {
1255//
1256// // Simple sum of both energies, cross-talk not included yet
1257// psum_en = cd_en_list.at( pindex[0] );
1258// psum_en += cd_en_list.at( pindex[1] );
1259//
1260// // Set event
1261// particle_evt->SetEnergyP( psum_en );
1262// particle_evt->SetEnergyN( cd_en_list.at( nindex[0] ) );
1263// particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1264// particle_evt->SetTimeN( cd_ts_list.at( nindex[0] ) );
1265// particle_evt->SetDetector( i );
1266// particle_evt->SetSector( j );
1267// particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1268// particle_evt->SetStripN( cd_strip_list.at( nindex[0] ) );
1269//
1270// // Fill tree
1271// write_evts->AddEvt( particle_evt );
1272// cd_ctr++;
1273//
1274// // Fill histograms
1275// cd_pen_id[i][j]->Fill( psum_en,
1276// cd_en_list.at( pmax_idx ) );
1277//
1278// } // neighbour strips
1279//
1280// // otherwise treat as 1 vs 0
1281// else {
1282//
1283// // Set event
1284// particle_evt->SetEnergyP( cd_en_list.at( pmax_idx ) );
1285// particle_evt->SetEnergyN( 0.0 );
1286// particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1287// particle_evt->SetTimeN( cd_ts_list.at( pmax_idx ) );
1288// particle_evt->SetDetector( i );
1289// particle_evt->SetSector( j );
1290// particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1291// particle_evt->SetStripN( 5.0 );
1292//
1293// // Fill tree
1294// write_evts->AddEvt( particle_evt );
1295// cd_ctr++;
1296//
1297// // Fill histograms
1298// cd_pen_id[i][j]->Fill( pmax_en,
1299// cd_en_list.at( pmax_idx ) );
1300//
1301// } // treat as 1 vs 0
1302//
1303// } // 2 vs 0
1304
1305 // Everything else, just take the max energy for now
1306 else if( pmax_idx >= 0 && nmax_idx >= 0 ){
1307
1308 // Set event
1309 particle_evt->SetEnergyP( pmax_en );
1310 particle_evt->SetEnergyN( nmax_en );
1311 particle_evt->SetTimeP( cd_ts_list.at( pmax_idx ) );
1312 particle_evt->SetTimeN( cd_ts_list.at( nmax_idx ) );
1313 particle_evt->SetDetector( i );
1314 particle_evt->SetSector( j );
1315 particle_evt->SetStripP( cd_strip_list.at( pmax_idx ) );
1316 particle_evt->SetStripN( cd_strip_list.at( nmax_idx ) );
1317 particle_evt->SetEnergyPad( pad_coinc_en );
1318 particle_evt->SetTimePad( pad_coinc_ts );
1319
1320 // Fill tree
1321 write_evts->AddEvt( particle_evt );
1322 cd_ctr++;
1323
1324 }
1325
1326 } // j: sector ID
1327
1328 } // i: detector ID
1329
1330 return;
1331
1332}
1333
1335
1336 // Build individual beam dump events
1337 // Loop over all the events in beam dump detectors
1338 for( unsigned int i = 0; i < bd_en_list.size(); ++i ) {
1339
1340 bd_evt->SetEnergy( bd_en_list.at(i) );
1341 bd_evt->SetTime( bd_ts_list.at(i) );
1342 bd_evt->SetDetector( bd_det_list.at(i) );
1343 write_evts->AddEvt( bd_evt );
1344 bd_ctr++;
1345
1346 }
1347
1348 return;
1349
1350}
1351
1353
1354 // Build individual Spede events
1355 // Loop over all the events in Spede detector
1356 for( unsigned int i = 0; i < spede_en_list.size(); ++i ) {
1357
1358 spede_evt->SetEnergy( spede_en_list.at(i) );
1359 spede_evt->SetTime( spede_ts_list.at(i) );
1360 spede_evt->SetSegment( spede_seg_list.at(i) );
1361 write_evts->AddEvt( spede_evt );
1362 spede_ctr++;
1363
1364 }
1365
1366 return;
1367
1368}
1369
1371
1372 // Build individual ion chamber events
1373 // Checks to prevent re-using events
1374 std::vector<unsigned int> index;
1375 std::vector<unsigned int> layer;
1376 bool flag_skip;
1377
1378 // Loop over IonChamber events
1379 for( unsigned int i = 0; i < ic_en_list.size(); ++i ) {
1380
1381 ic_evt->ClearEvt();
1382
1383 if( ic_id_list[i] == 0 ){
1384
1385 ic_evt->SetdETime( ic_ts_list[i] );
1386 ic_dE->Fill( ic_en_list[i] );
1387
1388 }
1389
1390 else if( ic_id_list[i] == set->GetNumberOfIonChamberLayers()-1 ){
1391
1392 ic_evt->SetETime( ic_ts_list[i] );
1393 ic_E->Fill( ic_en_list[i] );
1394
1395 }
1396
1397 ic_evt->AddIonChamber( ic_en_list[i], ic_id_list[i] );
1398 index.push_back( i );
1399 layer.push_back( ic_id_list[i] );
1400
1401 // Look for matching events in other layers
1402 for( unsigned int j = 0; j < ic_en_list.size(); ++j ) {
1403
1404 // Can't be coincident with itself
1405 if( i == j ) continue;
1406
1407 // Time difference plot
1408 ic_td->Fill( (double)ic_ts_list[i] - (double)ic_ts_list[j] );
1409
1410 // Check if we already used this hit
1411 flag_skip = false;
1412 for( unsigned int k = 0; k < index.size(); ++k ) {
1413 if( index[k] == j ) flag_skip = true;
1414 if( layer[k] == ic_id_list[j] ) flag_skip = true;
1415 }
1416
1417 // Found a match
1418 if( ic_id_list[j] != ic_id_list[i] && !flag_skip &&
1419 TMath::Abs( (double)ic_ts_list[i] - (double)ic_ts_list[j] ) < set->GetIonChamberHitWindow() ){
1420
1421 index.push_back( j );
1422 layer.push_back( ic_id_list[j] );
1423 ic_evt->AddIonChamber( ic_en_list[j], ic_id_list[j] );
1424
1425 if( ic_id_list[j] == 0 )
1426 ic_evt->SetdETime( ic_ts_list[j] );
1427 else if( ic_id_list[j] == set->GetNumberOfIonChamberLayers()-1 )
1428 ic_evt->SetETime( ic_ts_list[j] );
1429
1430 }
1431
1432 }
1433
1434 // Histogram the Ionchamber
1435 ic_dE_E->Fill( ic_evt->GetEnergyRest(), ic_evt->GetEnergyLoss() );
1436
1437 // Fill the tree and get ready for next ion chamber event
1438 write_evts->AddEvt( ic_evt );
1439 ic_ctr++;
1440
1441 }
1442
1443
1444
1445 return;
1446
1447}
1448
1449
1450
1452
1454
1455 // Load the full tree if possible
1456 //output_tree->SetMaxVirtualSize(1.0e9); // 1.0 GB
1457 //input_tree->SetMaxVirtualSize(2.2e9); // 2.2 GB
1458 //input_tree->LoadBaskets(2.0e9); // Load 2.0 GB of data to memory
1459
1460 if( input_tree->LoadTree(0) < 0 ){
1461
1462 std::cout << " Event Building: nothing to do" << std::endl;
1463 return 0;
1464
1465 }
1466
1467 // Get ready and go
1468 Initialise();
1469 n_entries = input_tree->GetEntries();
1470 n_mbs_entries = mbsinfo_tree->GetEntries();
1471
1472 std::cout << " Event Building: number of entries in input tree = ";
1473 std::cout << n_entries << std::endl;
1474
1475 std::cout << "\tnumber of MBS Events/triggers in input tree = ";
1476 std::cout << n_mbs_entries << std::endl;
1477
1478 // ------------------------------------------------------------------------ //
1479 // Main loop over TTree to find events
1480 // ------------------------------------------------------------------------ //
1481 for( unsigned long i = 0; i < n_entries; ++i ) {
1482
1483 // First event, yes please!
1484 if( i == 0 ){
1485
1486 input_tree->GetEntry(i);
1489
1490 // Try to get the MBS info event with the index
1491 if( mbsinfo_tree->GetEntryWithIndex( myeventid ) < 0 &&
1492 n_mbs_entries > 0 ) {
1493
1494 // Look for the matches MBS Info event if we didn't match automatically
1495 for( unsigned long j = 0; j < n_mbs_entries; ++j ){
1496
1497 mbsinfo_tree->GetEntry(j);
1498 if( mbs_info->GetEventID() == myeventid ) {
1500 break;
1501 }
1502
1503 // Panic if we failed!
1504 if( j+1 == n_mbs_entries ) {
1505 std::cerr << "Didn't find matching MBS Event IDs at start of the file: ";
1506 std::cerr << myeventid << std::endl;
1507 }
1508
1509 }
1510
1511 }
1512
1513 //std::cout << "MBS Trigger time = " << myeventtime << std::endl;
1514
1515 }
1516
1517 // Get the time of the event
1518 if( set->GetMbsEventMode() ) {
1519
1522
1523 }
1524
1525 else mytime = in_data->GetTime();
1526
1527 // check time stamp monotonically increases in time-ordered mode!
1528 if( time_prev > mytime && !set->GetMbsEventMode() ) {
1529
1530 std::cout << "Out of order event in ";
1531 std::cout << input_tree->GetName() << std::endl;
1532
1533 }
1534
1535 // check event id is increasing in MBS event ordered mode
1536 if( preveventid > myeventid && set->GetMbsEventMode() ) {
1537
1538 std::cout << "Out of order MBS event " << myeventid;
1539 std::cout << " < " << preveventid << std::endl;
1540
1541 }
1542
1543 // Get the laser status from the MBS Info tree
1544 if( mbs_info->GetPatternValue( set->GetRILISPattern() ) < 256 )
1545 mylaser = true;
1546 else mylaser = false;
1547
1548 // record time of this event
1549 time_prev = mytime;
1550
1551 // assume this is above threshold initially
1552 mythres = true;
1553
1554 // ------------------------------------------ //
1555 // Find FEBEX data
1556 // ------------------------------------------ //
1557 if( in_data->IsFebex() ) {
1558
1559 // Get the data
1561 mysfp = febex_data->GetSfp();
1562 myboard = febex_data->GetBoard();
1563 mych = febex_data->GetChannel();
1564 mypileup = febex_data->IsPileup();
1565 myclipped = febex_data->IsClipped();
1566
1567 // check for repeated timestamps in same channel
1569
1570 //std::cout << "Repeat timestamp in event " << std::dec << i << ":" << std::endl;
1571 //std::cout << "\t" << std::hex << mytime << " == ";
1572 //std::cout << std::hex << febex_time_ch[mysfp][myboard][mych];
1573 //std::cout << "\n\tSFP = " << std::dec << (int)mysfp;
1574 //std::cout << ", board = " << (int)myboard;
1575 //std::cout << ", channel = " << (int)mych << std::endl;
1576 //std::cout << std::endl;
1577
1578 repeat_ctr++;
1579
1580 } // skip repeat FEBEX event
1581
1582 // otherwise carry on
1583 else {
1584
1585 // Update calibration if necessary
1586 if( overwrite_cal ) {
1587
1588 unsigned int adc_tmp_value;
1589 if( cal->FebexType( mysfp, myboard, mych ) == "Qshort" )
1590 adc_tmp_value = febex_data->GetQshort();
1591 else if( cal->FebexType( mysfp, myboard, mych ) == "Qint" )
1592 adc_tmp_value = febex_data->GetQint();
1593 else adc_tmp_value = febex_data->GetQshort();
1594
1595 myenergy = cal->FebexEnergy( mysfp, myboard, mych, adc_tmp_value );
1596
1597 if( adc_tmp_value > cal->FebexThreshold( mysfp, myboard, mych ) )
1598 mythres = true;
1599 else mythres = false;
1600
1601 }
1602
1603 else {
1604
1605 myenergy = febex_data->GetEnergy();
1606 mythres = febex_data->IsOverThreshold();
1607
1608 }
1609
1610 // We should ignore negative energies
1611 if( myenergy < 0.0 ) mythres = false;
1612
1613 // Increment event counters
1614 n_febex_data++;
1615 n_sfp[mysfp]++;
1616 n_board[mysfp][myboard]++;
1617
1618 // Is it a gamma ray from Miniball?
1619 if( set->IsMiniball( mysfp, myboard, mych ) && mythres ) {
1620
1621 // Increment counts and open the event
1622 n_miniball++;
1623 hit_ctr++;
1624
1625 // Clipped rejection and pileup rejection
1626 if( ( !myclipped || !set->GetClippedRejection() ) &&
1627 ( !mypileup || !set->GetPileupRejection() ) ) {
1628
1629 event_open = true;
1630 mb_en_list.push_back( myenergy );
1631 mb_ts_list.push_back( mytime );
1632 mb_clu_list.push_back( set->GetMiniballCluster( mysfp, myboard, mych ) );
1633 mb_cry_list.push_back( set->GetMiniballCrystal( mysfp, myboard, mych ) );
1634 mb_seg_list.push_back( set->GetMiniballSegment( mysfp, myboard, mych ) );
1635
1636 }
1637
1638 }
1639
1640 // Is it a particle from the CD?
1641 else if( set->IsCD( mysfp, myboard, mych ) && mythres ) {
1642
1643 // Increment counts and open the event
1644 n_cd++;
1645 hit_ctr++;
1646
1647 // Clipped rejection and pileup rejection
1648 if( ( !myclipped || !set->GetClippedRejection() ) &&
1649 ( !mypileup || !set->GetPileupRejection() ) ) {
1650
1651 event_open = true;
1652 cd_en_list.push_back( myenergy );
1653 cd_ts_list.push_back( mytime );
1654 cd_det_list.push_back( set->GetCDDetector( mysfp, myboard, mych ) );
1655 cd_sec_list.push_back( set->GetCDSector( mysfp, myboard, mych ) );
1656 cd_side_list.push_back( set->GetCDSide( mysfp, myboard, mych ) );
1657 cd_strip_list.push_back( set->GetCDStrip( mysfp, myboard, mych ) );
1658
1659 }
1660
1661 }
1662
1663 // Is it a particle from the Pad?
1664 else if( set->IsPad( mysfp, myboard, mych ) && mythres ) {
1665
1666 // Increment counts and open the event
1667 n_pad++;
1668 hit_ctr++;
1669
1670 // Clipped rejection and pileup rejection
1671 if( ( !myclipped || !set->GetClippedRejection() ) &&
1672 ( !mypileup || !set->GetPileupRejection() ) ) {
1673
1674 event_open = true;
1675 pad_en_list.push_back( myenergy );
1676 pad_ts_list.push_back( mytime );
1677 pad_det_list.push_back( set->GetPadDetector( mysfp, myboard, mych ) );
1678 pad_sec_list.push_back( set->GetPadSector( mysfp, myboard, mych ) );
1679
1680 }
1681
1682 }
1683
1684 // Is it an electron from Spede?
1685 else if( set->IsSpede( mysfp, myboard, mych ) && mythres ) {
1686
1687 // Increment counts and open the event
1688 n_spede++;
1689 hit_ctr++;
1690
1691 // Clipped rejection and pileup rejection
1692 if( ( !myclipped || !set->GetClippedRejection() ) &&
1693 ( !mypileup || !set->GetPileupRejection() ) ) {
1694
1695 event_open = true;
1696 spede_en_list.push_back( myenergy );
1697 spede_ts_list.push_back( mytime );
1698 spede_seg_list.push_back( set->GetSpedeSegment( mysfp, myboard, mych ) );
1699
1700 }
1701
1702 }
1703
1704 // Is it a gamma ray from the beam dump?
1705 else if( set->IsBeamDump( mysfp, myboard, mych ) && mythres ) {
1706
1707 // Increment counts and open the event
1708 n_bd++;
1709 hit_ctr++;
1710
1711 // Clipped rejection and pileup rejection
1712 if( ( !myclipped || !set->GetClippedRejection() ) &&
1713 ( !mypileup || !set->GetPileupRejection() ) ) {
1714
1715 event_open = true;
1716 bd_en_list.push_back( myenergy );
1717 bd_ts_list.push_back( mytime );
1718 bd_det_list.push_back( set->GetBeamDumpDetector( mysfp, myboard, mych ) );
1719
1720 }
1721
1722 }
1723
1724 // Is it an IonChamber event
1725 else if( set->IsIonChamber( mysfp, myboard, mych ) && mythres ) {
1726
1727 // Increment counts and open the event
1728 n_ic++;
1729 hit_ctr++;
1730
1731 // Clipped rejection and pileup rejection
1732 if( ( !myclipped || !set->GetClippedRejection() ) &&
1733 ( !mypileup || !set->GetPileupRejection() ) ) {
1734
1735 event_open = true;
1736 ic_en_list.push_back( myenergy );
1737 ic_ts_list.push_back( mytime );
1738 ic_id_list.push_back( set->GetIonChamberLayer( mysfp, myboard, mych ) );
1739
1740 }
1741
1742 }
1743
1744 // Pulser item
1745 else if( set->IsPulser( mysfp, myboard, mych ) && mythres ) {
1746
1747 unsigned int pulserID = set->GetPulser( mysfp, myboard, mych );
1748 pulser_time[pulserID] = mytime;
1749 pulser_T = (double)pulser_time[pulserID] - (double)pulser_prev[pulserID];
1750 pulser_f = 1e9 / pulser_T;
1751 if( pulserID == 0 && pulser_prev[pulserID] != 0 ) {
1752 pulser_period->Fill( pulser_T );
1753 pulser_freq->Fill( pulser_time[pulserID], pulser_f );
1754 }
1755
1756 if( pulserID == 0 ) {
1757
1758 for( unsigned int i = 1; i < set->GetNumberOfPulsers(); i++ ) {
1759
1760 // If diff is greater than 5 ms, we have the wrong pair
1761 double tmp_tdiff = (double)pulser_time[i] - (double)pulser_time[0];
1762 if( tmp_tdiff > 1e4 ) tmp_tdiff = (double)pulser_prev[i] - (double)pulser_time[0];
1763 else if( tmp_tdiff < -1e4 ) tmp_tdiff = (double)pulser_time[i] - (double)pulser_prev[0];
1764
1765 pulser_tdiff->Fill( i, tmp_tdiff );
1766
1767 }
1768
1769 }
1770
1771 pulser_prev[pulserID] = pulser_time[pulserID];
1772 n_pulser[pulserID]++;
1773
1774 } // pulser code
1775
1776
1777 } // process febex data
1778
1779 // Update the current timestamp of this channel
1781
1782 // Is it the start event?
1783 if( febex_time_start.at( mysfp ).at( myboard ) == 0 )
1784 febex_time_start.at( mysfp ).at( myboard ) = mytime;
1785
1786 // or is it the end event (we don't know so keep updating)
1787 febex_time_stop.at( mysfp ).at( myboard ) = mytime;
1788
1789 }
1790
1791 // ------------------------------------------ //
1792 // Find DGF data
1793 // ------------------------------------------ //
1794 else if( in_data->IsDgf() ) {
1795
1796 // Get the data
1798 mydgf = dgf_data->GetModule();
1799 mych = dgf_data->GetChannel();
1800
1801 // Update calibration if necessary
1802 if( overwrite_cal ) {
1803
1804 unsigned int adc_tmp_value = dgf_data->GetQshort();
1805 myenergy = cal->DgfEnergy( mydgf, mych, adc_tmp_value );
1806
1807 if( adc_tmp_value > cal->DgfThreshold( mydgf, mych ) )
1808 mythres = true;
1809 else mythres = false;
1810
1811 }
1812
1813 else {
1814
1815 myenergy = dgf_data->GetEnergy();
1816 mythres = dgf_data->IsOverThreshold();
1817
1818 }
1819
1820 // We should ignore negative energies
1821 if( myenergy < 0.0 ) mythres = false;
1822
1823 // Increment event counters
1824 n_dgf_data++;
1825 n_dgf[mydgf]++;
1826
1827 // Is it a gamma ray from Miniball?
1828 if( set->IsMiniball( mydgf, mych ) && mythres ) {
1829
1830 // Increment counts and open the event
1831 n_miniball++;
1832 hit_ctr++;
1833 event_open = true;
1834
1835 mb_en_list.push_back( myenergy );
1836 mb_ts_list.push_back( mytime );
1837 mb_clu_list.push_back( set->GetMiniballCluster( mydgf, mych ) );
1838 mb_cry_list.push_back( set->GetMiniballCrystal( mydgf, mych ) );
1839 mb_seg_list.push_back( set->GetMiniballSegment( mydgf, mych ) );
1840
1841 }
1842
1843 // Is it a gamma ray from the beam dump?
1844 else if( set->IsBeamDump( mydgf, mych ) && mythres ) {
1845
1846 // Increment counts but do not open the event
1847 n_bd++;
1848 hit_ctr++;
1849
1850 bd_en_list.push_back( myenergy );
1851 bd_ts_list.push_back( mytime );
1852 bd_det_list.push_back( set->GetBeamDumpDetector( mydgf, mych ) );
1853
1854 }
1855
1856 }
1857
1858
1859 // ------------------------------------------ //
1860 // Find ADC data
1861 // ------------------------------------------ //
1862 if( in_data->IsAdc() ) {
1863
1864 // Get the data
1866 myadc = adc_data->GetModule();
1867 mych = adc_data->GetChannel();
1868 myclipped = adc_data->IsClipped();
1869
1870 // Update calibration if necessary
1871 if( overwrite_cal ) {
1872
1873 unsigned int adc_tmp_value = adc_data->GetQshort();
1874 myenergy = cal->AdcEnergy( myadc, mych, adc_tmp_value );
1875
1876 if( adc_tmp_value > cal->AdcThreshold( myadc, mych ) )
1877 mythres = true;
1878 else mythres = false;
1879
1880 }
1881
1882 else {
1883
1884 myenergy = adc_data->GetEnergy();
1885 mythres = adc_data->IsOverThreshold();
1886
1887 }
1888
1889 // We should ignore negative energies
1890 if( myenergy < 0.0 ) mythres = false;
1891
1892 // Increment event counters
1893 n_adc_data++;
1894 n_adc[myadc]++;
1895
1896 // Is it a particle from the CD?
1897 if( set->IsCD( myadc, mych ) && mythres ) {
1898
1899 // Increment counts and open the event
1900 n_cd++;
1901 hit_ctr++;
1902
1903 if( !myclipped || !set->GetClippedRejection() ) {
1904
1905 event_open = true;
1906 cd_en_list.push_back( myenergy );
1907 cd_ts_list.push_back( mytime );
1908 cd_det_list.push_back( set->GetCDDetector( myadc, mych ) );
1909 cd_sec_list.push_back( set->GetCDSector( myadc, mych ) );
1910 cd_side_list.push_back( set->GetCDSide( myadc, mych ) );
1911 cd_strip_list.push_back( set->GetCDStrip( myadc, mych ) );
1912
1913 }
1914
1915 }
1916
1917 // Is it a particle from the Pad?
1918 else if( set->IsPad( myadc, mych ) && mythres ) {
1919
1920 // Increment counts and open the event
1921 n_pad++;
1922 hit_ctr++;
1923
1924 if( !myclipped || !set->GetClippedRejection() ) {
1925
1926 event_open = true;
1927 pad_en_list.push_back( myenergy );
1928 pad_ts_list.push_back( mytime );
1929 pad_det_list.push_back( set->GetPadDetector( myadc, mych ) );
1930 pad_sec_list.push_back( set->GetPadSector( myadc, mych ) );
1931
1932 }
1933
1934 }
1935
1936 // Is it an IonChamber event
1937 else if( set->IsIonChamber( myadc, mych ) && mythres && !myclipped ) {
1938
1939 // Increment counts and open the event
1940 n_ic++;
1941 hit_ctr++;
1942
1943 if( !myclipped || !set->GetClippedRejection() ) {
1944
1945 event_open = true;
1946 ic_en_list.push_back( myenergy );
1947 ic_ts_list.push_back( mytime );
1948 ic_id_list.push_back( set->GetIonChamberLayer( myadc, mych ) );
1949
1950 }
1951
1952 }
1953
1954 }
1955
1956
1957
1958 // ------------------------------------------ //
1959 // Find info events, like timestamps etc
1960 // ------------------------------------------ //
1961 else if( in_data->IsInfo() ) {
1962
1963 // Increment event counter
1964 n_info_data++;
1965
1967
1968 // Update EBIS time
1969 if( info_data->GetCode() == set->GetEBISCode() &&
1970 TMath::Abs( (double)ebis_time - (double)info_data->GetTime() ) > 1e3 ) {
1971
1972 // Get the time of the EBIS
1973 ebis_time = info_data->GetTime();
1974 if( set->GetMbsEventMode() ) ebis_time += myeventtime;
1975
1976 ebis_T = (double)ebis_time - (double)ebis_prev;
1977 ebis_f = 1e9 / ebis_T;
1978 if( ebis_prev != 0 ) {
1979 ebis_period->Fill( ebis_T );
1980 ebis_freq->Fill( ebis_time, ebis_f );
1981 }
1983 n_ebis++;
1984
1985 } // EBIS code
1986
1987 // Update T1 time
1988 if( info_data->GetCode() == set->GetT1Code() &&
1989 TMath::Abs( (double)t1_time - (double)info_data->GetTime() ) > 1e3 ){
1990
1991 // Get the time of the T1
1992 t1_time = info_data->GetTime();
1993 if( set->GetMbsEventMode() ) t1_time += myeventtime;
1994
1995 t1_T = (double)t1_time - (double)t1_prev;
1996 t1_f = 1e9 / t1_T;
1997 if( t1_prev != 0 ) {
1998 t1_period->Fill( t1_T );
1999 t1_freq->Fill( t1_time, t1_f );
2000 }
2001 t1_prev = t1_time;
2002 n_t1++;
2003
2004 } // T1 code
2005
2006 // Update SuperCycle time
2007 if( info_data->GetCode() == set->GetSCCode() &&
2008 TMath::Abs( (double)sc_time - (double)info_data->GetTime() ) > 1e3 ){
2009
2010 // Get the time of the T1
2011 sc_time = info_data->GetTime();
2012 if( set->GetMbsEventMode() ) sc_time += myeventtime;
2013
2014 sc_T = (double)sc_time - (double)sc_prev;
2015 sc_f = 1e9 / sc_T;
2016 if( sc_prev != 0 ) {
2017 sc_period->Fill( sc_T );
2018 sc_freq->Fill( sc_time, sc_f );
2019 }
2020 sc_prev = sc_time;
2021 n_sc++;
2022
2023 } // SuperCycle code
2024
2025 // Update Laser time
2026 if( info_data->GetCode() == set->GetRILISCode() &&
2027 TMath::Abs( (double)laser_time - (double)info_data->GetTime() ) > 1e3 ){
2028
2029 // Get the time of the T1
2030 laser_time = info_data->GetTime();
2031 if( set->GetMbsEventMode() ) laser_time += myeventtime;
2032
2033 n_rilis++;
2034
2035 } // Laser code
2036
2037 // Update pulser time
2038 if( info_data->GetCode() >= set->GetPulserCode() &&
2039 info_data->GetCode() < set->GetPulserCode() + set->GetNumberOfPulsers() ) {
2040
2041 unsigned int pulserID = info_data->GetCode() - set->GetPulserCode();
2042 pulser_time[pulserID] = info_data->GetTime();
2043 if( set->GetMbsEventMode() ) pulser_time[pulserID] += myeventtime;
2044 pulser_T = (double)pulser_time[pulserID] - (double)pulser_prev[pulserID];
2045 pulser_f = 1e9 / pulser_T;
2046 if( pulserID == 0 && pulser_prev[pulserID] != 0 ) {
2047 pulser_period->Fill( pulser_T );
2048 pulser_freq->Fill( pulser_time[pulserID], pulser_f );
2049 }
2050
2051 if( pulserID == 0 ) {
2052
2053 for( unsigned int i = 1; i < set->GetNumberOfPulsers(); i++ ) {
2054
2055 // If diff is greater than 5 ms, we have the wrong pair
2056 double tmp_tdiff = (double)pulser_time[i] - (double)pulser_time[0];
2057 if( tmp_tdiff > 1e4 ) tmp_tdiff = (double)pulser_prev[i] - (double)pulser_time[0];
2058 else if( tmp_tdiff < -1e4 ) tmp_tdiff = (double)pulser_time[i] - (double)pulser_prev[0];
2059
2060 pulser_tdiff->Fill( i, tmp_tdiff );
2061
2062 }
2063
2064 }
2065
2066 pulser_prev[pulserID] = pulser_time[pulserID];
2067 n_pulser[pulserID]++;
2068
2069 } // pulser code
2070
2071 // Update sync time
2072 if( info_data->GetCode() == set->GetMsbSyncCode() ){
2073
2074 // Check the SFP and module number
2075 if( info_data->GetSfp() < set->GetNumberOfFebexSfps() &&
2076 info_data->GetBoard() < set->GetNumberOfFebexBoards() ) {
2077
2078 // Get the time of the sync pulse
2079 sync_time[info_data->GetSfp()][info_data->GetBoard()] = info_data->GetTime();
2080 n_sync[info_data->GetSfp()][info_data->GetBoard()]++;
2081
2082 }
2083
2084 else {
2085
2086 std::cerr << "Bad sync event in SFP " << (int)info_data->GetSfp();
2087 std::cerr << ", board " << (int)info_data->GetBoard() << std::endl;
2088
2089 }
2090
2091 } // Sync pulse code
2092
2093 // Check the pause events for each module
2094 if( info_data->GetCode() == set->GetPauseCode() ) {
2095
2096 if( info_data->GetSfp() < set->GetNumberOfFebexSfps() &&
2097 info_data->GetBoard() < set->GetNumberOfFebexBoards() ) {
2098
2099 n_pause[info_data->GetSfp()][info_data->GetBoard()]++;
2100 flag_pause[info_data->GetSfp()][info_data->GetBoard()] = true;
2101 pause_time[info_data->GetSfp()][info_data->GetBoard()] = info_data->GetTime();
2102
2103 }
2104
2105 else {
2106
2107 std::cerr << "Bad pause event in SFP " << (int)info_data->GetSfp();
2108 std::cerr << ", board " << (int)info_data->GetBoard() << std::endl;
2109
2110 }
2111
2112 } // pause code
2113
2114 // Check the resume events for each module
2115 if( info_data->GetCode() == set->GetResumeCode() ) {
2116
2117 if( info_data->GetSfp() < set->GetNumberOfFebexSfps() &&
2118 info_data->GetBoard() < set->GetNumberOfFebexBoards() ) {
2119
2120 n_resume[info_data->GetSfp()][info_data->GetBoard()]++;
2121 flag_resume[info_data->GetSfp()][info_data->GetBoard()] = true;
2122 resume_time[info_data->GetSfp()][info_data->GetBoard()] = info_data->GetTime();
2123
2124 // Work out the dead time
2125 febex_dead_time[info_data->GetSfp()][info_data->GetBoard()] += resume_time[info_data->GetSfp()][info_data->GetBoard()];
2126 febex_dead_time[info_data->GetSfp()][info_data->GetBoard()] -= pause_time[info_data->GetSfp()][info_data->GetBoard()];
2127
2128 // If we have didn't get the pause, module was stuck at start of run
2129 if( !flag_pause[info_data->GetSfp()][info_data->GetBoard()] ) {
2130
2131 std::cout << "SFP " << info_data->GetSfp();
2132 std::cout << ", board " << info_data->GetBoard();
2133 std::cout << " was blocked at start of run for ";
2134 std::cout << (double)resume_time[info_data->GetSfp()][info_data->GetBoard()]/1e9;
2135 std::cout << " seconds" << std::endl;
2136
2137 }
2138
2139 }
2140
2141 else {
2142
2143 std::cerr << "Bad resume event in SFP " << (int)info_data->GetSfp();
2144 std::cerr << ", board " << (int)info_data->GetBoard() << std::endl;
2145
2146 }
2147
2148 } // resume code
2149
2150 } // is info data
2151
2152 // Sort out the timing for the event window
2153 // but only if it isn't an info event, i.e only for real data
2154 if( !in_data->IsInfo() ) {
2155
2156 // if this is first datum included in Event
2157 if( hit_ctr == 1 && mythres ) {
2158
2159 time_min = mytime;
2160 time_max = mytime;
2162
2163 }
2164
2165 // Update min and max
2166 if( mytime > time_max ) time_max = mytime;
2167 else if( mytime < time_min ) time_min = mytime;
2168
2169 } // not info data
2170
2171 //------------------------------
2172 // check if last datum from this event and do some cleanup
2173 //------------------------------
2174
2175 if( input_tree->GetEntry(i+1) ) {
2176
2177 // Get the next MBS event ID
2180
2181 // If the next MBS event ID is the same, carry on
2182 // If not, we have to go look for the next trigger time
2183 if( myeventid != preveventid ) {
2184
2185 // Close the event
2186 flag_close_event = true;
2187
2188 // And find the next MBS event ID
2189 if( mbsinfo_tree->GetEntryWithIndex( myeventid ) < 0 &&
2190 n_mbs_entries > 0 ) {
2191
2192 std::cerr << "MBS Event " << myeventid << " not found by index, looking up manually" << std::endl;
2193
2194 // Look for the matches MBS Info event if we didn't match automatically
2195 for( unsigned long j = 0; j < n_mbs_entries; ++j ){
2196
2197 mbsinfo_tree->GetEntry(j);
2198 if( mbs_info->GetEventID() == myeventid ) {
2200 break;
2201 }
2202
2203 // Panic if we failed!
2204 if( j+1 == n_mbs_entries ) {
2205 std::cerr << "Didn't find matching MBS Event IDs at start of the file: ";
2206 std::cerr << myeventid << std::endl;
2207 }
2208 }
2209
2210 }
2211
2212 else myeventtime = mbs_info->GetTime();
2213
2214 }
2215
2216 // BELOW IS THE TIME-ORDERED METHOD!
2217
2218 // Get the time of the next event
2219 if( set->GetMbsEventMode() ) {
2220
2223
2224 }
2225
2226 else mytime = in_data->GetTime();
2227
2228 // Calculate time diff
2230
2231 // window = time_stamp_first + time_window
2232 if( time_diff > build_window )
2233 flag_close_event = true; // set flag to close this event
2234
2235 // we've gone on to the next file in the chain
2236 else if( time_diff < 0 )
2237 flag_close_event = true; // set flag to close this event
2238
2239 // Fill tdiff hist only for real data
2240 if( !in_data->IsInfo() ) {
2241
2242 tdiff->Fill( time_diff );
2243 if( !mythres )
2244 tdiff_clean->Fill( time_diff );
2245
2246 }
2247
2248 } // if next entry beyond time window: close event!
2249
2250
2251 //----------------------------
2252 // if close this event or last entry
2253 //----------------------------
2254 if( flag_close_event || (i+1) == n_entries ) {
2255
2256 // If we opened the event, then sort it out
2257 if( event_open ) {
2258
2259 //----------------------------------
2260 // Build array events, recoils, etc
2261 //----------------------------------
2262 GammaRayFinder(); // perform addback
2263 ParticleFinder(); // sort out CD n/p correlations
2264 BeamDumpFinder(); // sort out beam dump events
2265 SpedeFinder(); // sort out Spede events
2266 IonChamberFinder(); // sort out beam dump events
2267
2268 // ------------------------------------
2269 // Add timing and fill the ISSEvts tree
2270 // ------------------------------------
2271 write_evts->SetEBIS( ebis_time );
2272 write_evts->SetT1( t1_time );
2273 write_evts->SetSC( sc_time );
2274 if( TMath::Abs( (double)ebis_time - (double)laser_time ) < 1e3
2275 && laser_time > 0 ) write_evts->SetLaserStatus( true );
2276 else if( mylaser )
2277 write_evts->SetLaserStatus( true );
2278 else
2279 write_evts->SetLaserStatus( false );
2280
2281 if( write_evts->GetGammaRayMultiplicity() ||
2282 write_evts->GetGammaRayAddbackMultiplicity() ||
2283 write_evts->GetParticleMultiplicity() ||
2284 write_evts->GetSpedeMultiplicity() ||
2285 write_evts->GetIonChamberMultiplicity() ||
2286 write_evts->GetBeamDumpMultiplicity() )
2287 output_tree->Fill();
2288
2289 }
2290
2291 //--------------------------------------------------
2292 // clear values of arrays to store intermediate info
2293 //--------------------------------------------------
2294 Initialise();
2295
2296 } // if close event && hit_ctr > 0
2297
2298 // Progress bar
2299 bool update_progress = false;
2300 if( n_entries < 200 )
2301 update_progress = true;
2302 else if( i % (n_entries/100) == 0 || i+1 == n_entries )
2303 update_progress = true;
2304
2305 if( update_progress ) {
2306
2307 // Percent complete
2308 float percent = (float)(i+1)*100.0/(float)n_entries;
2309
2310 // Progress bar in GUI
2311 if( _prog_ ) {
2312
2313 prog->SetPosition( percent );
2314 gSystem->ProcessEvents();
2315
2316 }
2317
2318 // Progress bar in terminal
2319 std::cout << " " << std::setw(6) << std::setprecision(4);
2320 std::cout << percent << "% \r";
2321 std::cout.flush();
2322
2323 }
2324
2325 } // End of main loop over TTree to process raw FEBEX data entries (for n_entries)
2326
2327 //--------------------------
2328 // Clean up
2329 //--------------------------
2330
2331 std::stringstream ss_log;
2332 ss_log << "\n MiniballEventBuilder finished..." << std::endl;
2333 ss_log << " FEBEX data packets = " << n_febex_data << std::endl;
2334 for( unsigned int i = 0; i < set->GetNumberOfFebexSfps(); ++i ) {
2335 ss_log << " SFP " << i << " events = " << n_sfp[i] << std::endl;
2336 for( unsigned int j = 0; j < set->GetNumberOfFebexBoards(); ++j ) {
2337 ss_log << " Board " << j << " events = " << n_board[i][j] << std::endl;
2338 // ss_log << " pause = " << n_pause[i][j] << std::endl;
2339 // ss_log << " resume = " << n_resume[i][j] << std::endl;
2340 // ss_log << " dead time = " << (double)febex_dead_time[i][j]/1e9 << " s" << std::endl;
2341 ss_log << " sync pulses = " << n_sync[i][j] << std::endl;
2342 ss_log << " run time = " << (double)(febex_time_stop[i][j]-febex_time_start[i][j])/1e9 << " s" << std::endl;
2343 }
2344 }
2345 ss_log << " DGF data packets = " << n_dgf_data << std::endl;
2346 for( unsigned int i = 0; i < set->GetNumberOfDgfModules(); ++i ) {
2347 ss_log << " Module " << i << " events = " << n_dgf[i] << std::endl;
2348 }
2349 ss_log << " ADC data packets = " << n_adc_data << std::endl;
2350 for( unsigned int i = 0; i < set->GetNumberOfAdcModules(); ++i ) {
2351 ss_log << " Module " << i << " events = " << n_adc[i] << std::endl;
2352 }
2353 ss_log << " Info data packets = " << n_info_data << std::endl;
2354 for( unsigned int i = 0; i < set->GetNumberOfPulsers(); ++i )
2355 ss_log << " Pulser " << i << " events = " << n_pulser[i] << std::endl;
2356 ss_log << " EBIS events = " << n_ebis << std::endl;
2357 ss_log << " RILIS events = " << n_rilis << std::endl;
2358 ss_log << " T1 events = " << n_t1 << std::endl;
2359 ss_log << " SuperCycle events = " << n_sc << std::endl;
2360 ss_log << " Tree entries = " << output_tree->GetEntries() << std::endl;
2361 ss_log << " Miniball triggers = " << n_miniball << std::endl;
2362 ss_log << " Gamma singles events = " << gamma_ctr << std::endl;
2363 ss_log << " Gamma addback events = " << gamma_ab_ctr << std::endl;
2364 ss_log << " CD detector triggers = " << n_cd << std::endl;
2365 ss_log << " Pad detector triggers = " << n_pad << std::endl;
2366 ss_log << " Particle events = " << cd_ctr << std::endl;
2367 ss_log << " SPEDE triggers = " << n_spede << std::endl;
2368 ss_log << " Electron events = " << spede_ctr << std::endl;
2369 ss_log << " Beam dump triggers = " << n_bd << std::endl;
2370 ss_log << " Beam dump gamma events = " << bd_ctr << std::endl;
2371 ss_log << " IonChamber triggers = " << n_ic << std::endl;
2372 ss_log << " IonChamber ion events = " << ic_ctr << std::endl;
2373 ss_log << "\nTotal number of repeated events skipped = " << std::dec << repeat_ctr;
2374 ss_log << " / " << n_entries << " = " << (double)repeat_ctr*100.0/(double)n_entries;
2375 ss_log << "\%" << std::endl;
2376
2377 std::cout << ss_log.str();
2378 if( log_file.is_open() && flag_input_file ) log_file << ss_log.str();
2379
2380 return n_entries;
2381
2382}
unsigned int GetPatternValue(unsigned char mod, unsigned char id)
long long int GetTime()
unsigned long long int GetEventID()
unsigned long long int GetEventID() const
std::shared_ptr< InfoData > GetInfoData() const
bool IsInfo() const
std::shared_ptr< FebexData > GetFebexData() const
bool IsFebex() const
std::shared_ptr< DgfData > GetDgfData() const
std::shared_ptr< AdcData > GetAdcData() const
long long int GetTime() const
unsigned long spede_ctr
std::shared_ptr< ParticleEvt > particle_evt
std::shared_ptr< TGProgressBar > prog
std::vector< unsigned char > spede_seg_list
list of Spede segment IDs
std::vector< unsigned long long > spede_ts_list
list of Spede timestamps for ElectronFinder
std::vector< unsigned char > cd_strip_list
list of CD strip IDs
unsigned long long preveventid
previous MBS event id
unsigned long long ebis_time
TFile * input_file
Input tree.
bool mythres
above threshold?
std::vector< std::vector< TH2F * > > mb_en_core_seg
bool mylaser
laser pattern bit
std::shared_ptr< MiniballCalibration > cal
std::vector< std::vector< unsigned long > > n_resume
std::vector< std::vector< unsigned long long > > febex_dead_time
std::vector< std::vector< TH2F * > > cd_pn_1v2
std::vector< unsigned char > ic_id_list
list of IonChamber layer IDs
void SetInputFile(std::string input_file_name)
MiniballEventBuilder(std::shared_ptr< MiniballSettings > myset)
std::vector< float > mb_en_list
list of Miniball energies for GammaRayFinder
std::shared_ptr< GammaRayEvt > gamma_evt
void StartFile()
called for every file
unsigned long long myeventtime
MBS event time.
std::shared_ptr< FebexData > febex_data
unsigned long BuildEvents()
std::vector< std::vector< TH2F * > > cd_pen_id
std::vector< unsigned char > cd_sec_list
list of CD sector IDs
std::vector< std::vector< TH2F * > > cd_pn_1v1
std::vector< unsigned long long > cd_ts_list
list of CD timestamps for ParticleFinder
unsigned long long time_max
std::vector< unsigned char > cd_det_list
list of CD detector IDs
unsigned long n_spede
std::vector< unsigned char > pad_det_list
list of PAD detector IDs
std::vector< std::vector< unsigned long long > > pause_time
std::vector< unsigned char > bd_det_list
list of beam dump detector IDs
unsigned long repeat_ctr
void SetMBSInfoTree(TTree *user_tree)
unsigned long n_miniball
unsigned long long sc_time
std::vector< std::vector< TH2F * > > mb_en_core_seg_ebis_on
std::shared_ptr< AdcData > adc_data
unsigned long n_info_data
std::vector< std::vector< bool > > flag_resume
MiniballDataPackets * in_data
std::vector< std::vector< unsigned long long > > sync_time
unsigned char mydgf
DGF module number.
unsigned long gamma_ab_ctr
unsigned long long mytime
absolute timestamp
std::shared_ptr< DgfData > dgf_data
std::vector< std::vector< TH2F * > > cd_ppad_mult
unsigned char myboard
febex board number
std::vector< unsigned char > pad_sec_list
list of PAD sector IDs
std::vector< std::vector< TH2F * > > cd_pn_mult
std::vector< unsigned char > mb_cry_list
list of crystal IDs
std::vector< TH2F * > pad_en_id
bool flag_close_event
length of build window in ns
unsigned long long t1_prev
std::vector< std::vector< unsigned long long > > febex_time_start
std::vector< float > ic_en_list
list of IonChamber energies for IonChamberFinder
std::vector< std::vector< TH2F * > > cd_pn_2v1
std::vector< unsigned long long > pulser_time
std::vector< unsigned long > n_adc
std::vector< std::vector< unsigned long > > n_board
bool mypileup
pileup flag?
unsigned long n_rilis
unsigned long gamma_ctr
std::vector< float > spede_en_list
list of Spede energies for ElectronFinder
std::shared_ptr< GammaRayAddbackEvt > gamma_ab_evt
std::shared_ptr< InfoData > info_data
std::vector< unsigned long long > mb_ts_list
list of Miniball timestamps for GammaRayFinder
bool myclipped
clipped flag?
std::shared_ptr< BeamDumpEvt > bd_evt
unsigned long hit_ctr
unsigned char mysfp
sfp number
std::vector< std::vector< bool > > flag_pause
std::vector< float > bd_en_list
list of beam dump energies for BeamDumpFinder
unsigned long n_dgf_data
long myhittime
hit time with respect to event time
void Initialise()
called for every event
unsigned long long t1_time
unsigned long long time_prev
std::vector< std::vector< TH1F * > > cd_nn_td
unsigned char myadc
ADC module number.
void SetInputTree(TTree *user_tree)
std::vector< unsigned char > cd_side_list
list of CD side IDs; 0 = p, 1 = n
std::shared_ptr< MiniballSettings > set
unsigned long long ebis_prev
std::vector< std::vector< TH1F * > > cd_pp_td
std::vector< unsigned long long > pad_ts_list
list of PAD timestamps for ParticleFinder
std::vector< std::vector< TH2F * > > cd_pn_2v2
std::unique_ptr< MiniballEvts > write_evts
std::shared_ptr< IonChamberEvt > ic_evt
unsigned long long time_min
MBSInfoPackets * mbs_info
unsigned long n_adc_data
std::vector< unsigned long > n_sfp
unsigned char mych
channel number
unsigned long long myeventid
MBS event id.
unsigned long n_entries
std::vector< std::vector< unsigned long > > n_pause
std::vector< std::vector< TH2F * > > cd_nen_id
unsigned long long sc_prev
float myenergy
calibrated energy
std::shared_ptr< SpedeEvt > spede_evt
std::vector< unsigned long > n_pulser
std::vector< unsigned long long > pulser_prev
std::vector< unsigned long long > bd_ts_list
list of beam dump timestamps for BeamDumpFinder
std::vector< unsigned long long > ic_ts_list
list of IonChamber timestamps for IonChamberFinder
std::vector< unsigned char > mb_seg_list
list of segment IDs
TFile * output_file
Outputs.
std::vector< float > pad_en_list
list of PAD energies for ParticleFinder
std::ofstream log_file
Log file for recording the results of the MiniballEventBuilder.
std::vector< unsigned char > mb_clu_list
list of cluster IDs
std::vector< std::vector< unsigned long > > n_sync
std::vector< float > cd_en_list
list of CD energies for ParticleFinder
std::vector< std::vector< unsigned long long > > resume_time
void SetOutput(std::string output_file_name, bool cWrite=false)
unsigned long long time_first
unsigned long n_mbs_entries
std::vector< std::vector< TH1F * > > cd_ppad_td
std::vector< std::vector< std::vector< unsigned long long > > > febex_time_ch
unsigned long n_febex_data
std::vector< std::vector< unsigned long long > > febex_time_stop
std::vector< std::vector< TH1F * > > cd_pn_td
unsigned long long laser_time
std::vector< unsigned long > n_dgf
std::shared_ptr< MiniballSettings > myset
Definition mb_sort.cc:123