Bug Summary

File:/home/sdobbs/work/clang/halld_recon/src/libraries/TRIGGER/DL1MCTrigger_factory.cc
Warning:line 617, column 2
Value stored to 'status' is never read

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -cc1 -triple x86_64-unknown-linux-gnu -analyze -disable-free -main-file-name DL1MCTrigger_factory.cc -analyzer-store=region -analyzer-opt-analyze-nested-blocks -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=cplusplus -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -mrelocation-model pic -pic-level 2 -fhalf-no-semantic-interposition -mframe-pointer=none -fmath-errno -fno-rounding-math -mconstructor-aliases -munwind-tables -target-cpu x86-64 -tune-cpu generic -fno-split-dwarf-inlining -debugger-tuning=gdb -resource-dir /w/halld-scifs17exp/home/sdobbs/clang/llvm-project/install/lib/clang/12.0.0 -D HAVE_CCDB -D HAVE_RCDB -D HAVE_EVIO -D HAVE_TMVA=1 -D RCDB_MYSQL=1 -D RCDB_SQLITE=1 -D SQLITE_USE_LEGACY_STRUCT=ON -I .Linux_CentOS7.7-x86_64-gcc4.8.5/libraries/TRIGGER -I libraries/TRIGGER -I . -I libraries -I libraries/include -I /w/halld-scifs17exp/home/sdobbs/clang/halld_recon/Linux_CentOS7.7-x86_64-gcc4.8.5/include -I external/xstream/include -I /usr/include/tirpc -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/root/root-6.08.06/include -I /w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/ccdb/ccdb_1.06.06/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/rcdb/rcdb_0.06.00/cpp/include -I /usr/include/mysql -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/sqlitecpp/SQLiteCpp-2.2.0^bs130/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/sqlite/sqlite-3.13.0^bs130/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/hdds/hdds-4.9.0/Linux_CentOS7.7-x86_64-gcc4.8.5/src -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/xerces-c/xerces-c-3.1.4/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/evio/evio-4.4.6/Linux-x86_64/include -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5 -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5/x86_64-redhat-linux -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5/backward -internal-isystem /usr/local/include -internal-isystem /w/halld-scifs17exp/home/sdobbs/clang/llvm-project/install/lib/clang/12.0.0/include -internal-externc-isystem /include -internal-externc-isystem /usr/include -O2 -std=c++11 -fdeprecated-macro -fdebug-compilation-dir /home/sdobbs/work/clang/halld_recon/src -ferror-limit 19 -fgnuc-version=4.2.1 -fcxx-exceptions -fexceptions -vectorize-loops -vectorize-slp -analyzer-output=html -faddrsig -o /tmp/scan-build-2021-01-21-110224-160369-1 -x c++ libraries/TRIGGER/DL1MCTrigger_factory.cc
1// Test version V0 of the L1 trigger simulation, Aug 18, 2017 (asomov)
2
3#include <iostream>
4#include <iomanip>
5#include <cmath>
6using namespace std;
7
8#include <JANA/JApplication.h>
9#include <DAQ/DCODAROCInfo.h>
10#include <DAQ/DL1Info.h>
11#include <DANA/DStatusBits.h>
12
13#include <HDDM/DEventSourceHDDM.h>
14
15using namespace jana;
16
17#include "DL1MCTrigger_factory.h"
18
19#if HAVE_RCDB1
20#include "RCDB/Connection.h"
21#include "RCDB/ConfigParser.h"
22#endif
23
24static bool print_data_message = true;
25
26//------------------
27// init
28//------------------
29jerror_t DL1MCTrigger_factory::init(void)
30{
31
32 debug = 0;
33
34 if(debug){
35 hfcal_gains = new TH1F("fcal_gains", "fcal_gains", 80, -1., 3.);
36 hfcal_gains2 = new TH2F("fcal_gains2","fcal_gains2", 71, -142., 142., 71, -142., 142.);
37 hfcal_ped = new TH1F("fcal_ped", "fcal_ped", 800, 0., 200.);
38 }
39
40 BYPASS = 0; // default is to use trigger emulation
41
42 // Default parameters for the main production trigger are taken from the
43 // spring run of 2017: 25 F + B > 45000
44
45 FCAL_ADC_PER_MEV = 3.73;
46 FCAL_CELL_THR = 65;
47 FCAL_EN_SC = 25;
48 FCAL_NSA = 10;
49 FCAL_NSB = 3;
50 FCAL_WINDOW = 10;
51
52 BCAL_ADC_PER_MEV = 34.48276; // Not corrected energy
53 BCAL_CELL_THR = 20;
54 BCAL_EN_SC = 1;
55 BCAL_NSA = 19;
56 BCAL_NSB = 3;
57 BCAL_WINDOW = 20;
58
59 FCAL_BCAL_EN = 45000;
60
61 ST_ADC_PER_MEV = 1.;
62 ST_CELL_THR = 40;
63 ST_NSA = 10;
64 ST_NSB = 3;
65 ST_WINDOW = 10;
66 ST_NHIT = 1;
67
68 BCAL_OFFSET = 2;
69
70 SIMU_BASELINE = 1;
71 SIMU_GAIN = 1;
72
73
74 simu_baseline_fcal = 1;
75 simu_baseline_bcal = 1;
76
77 simu_gain_fcal = 1;
78 simu_gain_bcal = 1;
79
80 gPARMS->SetDefaultParameter("TRIG:BYPASS", BYPASS,
81 "Bypass trigger by hard coding physics bit");
82 gPARMS->SetDefaultParameter("TRIG:FCAL_ADC_PER_MEV", FCAL_ADC_PER_MEV,
83 "FCAL energy calibration for the Trigger");
84 gPARMS->SetDefaultParameter("TRIG:FCAL_CELL_THR", FCAL_CELL_THR,
85 "FCAL energy threshold per cell");
86 gPARMS->SetDefaultParameter("TRIG:FCAL_EN_SC", FCAL_EN_SC,
87 "FCAL energy threshold");
88 gPARMS->SetDefaultParameter("TRIG:FCAL_NSA", FCAL_NSA,
89 "FCAL NSA");
90 gPARMS->SetDefaultParameter("TRIG:FCAL_NSB", FCAL_NSB,
91 "FCAL NSB");
92 gPARMS->SetDefaultParameter("TRIG:FCAL_WINDOW", FCAL_WINDOW,
93 "FCAL GTP integration window");
94
95 gPARMS->SetDefaultParameter("TRIG:BCAL_ADC_PER_MEV", BCAL_ADC_PER_MEV,
96 "BCAL energy calibration for the Trigger");
97 gPARMS->SetDefaultParameter("TRIG:BCAL_CELL_THR", BCAL_CELL_THR,
98 "BCAL energy threshold per cell");
99 gPARMS->SetDefaultParameter("TRIG:BCAL_EN_SC", BCAL_EN_SC,
100 "BCAL energy threshold");
101 gPARMS->SetDefaultParameter("TRIG:BCAL_NSA", BCAL_NSA,
102 "BCAL NSA");
103 gPARMS->SetDefaultParameter("TRIG:BCAL_NSB", BCAL_NSB,
104 "BCAL NSB");
105 gPARMS->SetDefaultParameter("TRIG:BCAL_WINDOW", BCAL_WINDOW,
106 "BCAL GTP integration window");
107
108 gPARMS->SetDefaultParameter("TRIG:ST_ADC_PER_MEV", ST_ADC_PER_MEV,
109 "ST energy calibration for the Trigger");
110 gPARMS->SetDefaultParameter("TRIG:ST_CELL_THR", ST_CELL_THR,
111 "ST energy threshold per cell");
112 gPARMS->SetDefaultParameter("TRIG:ST_NSA", ST_NSA,
113 "ST NSA");
114 gPARMS->SetDefaultParameter("TRIG:ST_NSB", ST_NSB,
115 "ST NSB");
116 gPARMS->SetDefaultParameter("TRIG:ST_WINDOW", ST_WINDOW,
117 "ST window for merging hits (GTP)");
118 gPARMS->SetDefaultParameter("TRIG:ST_NHIT", ST_NHIT,
119 "Number of hits in ST");
120
121 gPARMS->SetDefaultParameter("TRIG:FCAL_BCAL_EN", FCAL_BCAL_EN,
122 "Energy threshold for the FCAL & BCAL trigger");
123
124 gPARMS->SetDefaultParameter("TRIG:BCAL_OFFSET", BCAL_OFFSET,
125 "Timing offset between BCAL and FCAL energies at GTP (sampels)");
126
127
128 // Allows to switch off gain and baseline fluctuations
129 gPARMS->SetDefaultParameter("TRIG:SIMU_BASELINE", SIMU_BASELINE,
130 "Enable simulation of pedestal variations");
131
132 gPARMS->SetDefaultParameter("TRIG:SIMU_GAIN", SIMU_GAIN,
133 "Enable simulation of gain variations");
134
135
136 BCAL_ADC_PER_MEV_CORRECT = 22.7273;
137
138 pedestal_sigma = 1.2;
139
140 time_shift = 100;
141
142 time_min = 0;
143 time_max = (sample - 1)*max_adc_bins;
144
145 vector< vector<double > > fcal_gains_temp(DFCALGeometry::kBlocksTall,
146 vector<double>(DFCALGeometry::kBlocksWide));
147 vector< vector<double > > fcal_pedestals_temp(DFCALGeometry::kBlocksTall,
148 vector<double>(DFCALGeometry::kBlocksWide));
149
150 fcal_gains = fcal_gains_temp;
151 fcal_pedestals = fcal_pedestals_temp;
152
153 return NOERROR;
154}
155
156
157//------------------
158// brun
159//------------------
160jerror_t DL1MCTrigger_factory::brun(jana::JEventLoop *eventLoop, int32_t runnumber)
161{
162
163 int use_rcdb = 1;
164
165 int status = 0;
166
167 fcal_trig_mask.clear();
168 bcal_trig_mask.clear();
169
170 triggers_enabled.clear();
171
172 // Only print messages for one thread whenever run number change
173
174 static pthread_mutex_t print_mutex = PTHREAD_MUTEX_INITIALIZER{ { 0, 0, 0, 0, 0, 0, 0, { 0, 0 } } };
175 static set<int> runs_announced;
176 pthread_mutex_lock(&print_mutex);
177 bool print_messages = false;
178
179 if(runs_announced.find(runnumber) == runs_announced.end()){
180 print_messages = true;
181 runs_announced.insert(runnumber);
182 }
183
184 pthread_mutex_unlock(&print_mutex);
185
186 // Don't use RCDB for mc_generic: ideal MC simulation
187
188 string JANA_CALIB_CONTEXT = "";
189
190 if(getenv("JANA_CALIB_CONTEXT") != NULL__null ){
191 JANA_CALIB_CONTEXT = getenv("JANA_CALIB_CONTEXT");
192 if(JANA_CALIB_CONTEXT.find("mc_generic") != string::npos){
193 use_rcdb = 0;
194 // Don't simulate baseline fluctuations for mc_generic
195 simu_baseline_fcal = 0;
196 simu_baseline_bcal = 0;
197 // Don't simulate gain fluctuations for mc_generic
198 simu_gain_fcal = 0;
199 simu_gain_bcal = 0;
200 }
201 }
202
203 // runnumber = 30942;
204
205 if(use_rcdb == 1){
206 status = Read_RCDB(runnumber, print_messages);
207 if(print_messages) PrintTriggers();
208 }
209
210
211 if( (use_rcdb == 0) || (status > 0) || (triggers_enabled.size() == 0)){
212
213 // Simulate FCAL & BCAL main production trigger only
214
215 trigger_conf trig_tmp;
216 trig_tmp.bit = 0;
217 trig_tmp.gtp.fcal = FCAL_EN_SC;
218 trig_tmp.gtp.bcal = BCAL_EN_SC;
219 trig_tmp.gtp.en_thr = FCAL_BCAL_EN;
220 trig_tmp.gtp.fcal_min = 200;
221 triggers_enabled.push_back(trig_tmp);
222
223 if(print_messages)
224 cout << " Do not use RCDB for the trigger simulation. Default (spring 2017) trigger settings are used " << endl;
225 }
226
227
228 // extract the FCAL Geometry
229 vector<const DFCALGeometry*> fcalGeomVect;
230 eventLoop->Get( fcalGeomVect );
231 if (fcalGeomVect.size() < 1)
232 return OBJECT_NOT_AVAILABLE;
233 const DFCALGeometry& fcalGeom = *(fcalGeomVect[0]);
234
235 if(print_messages) jout << "In DL1MCTrigger_factory, loading constants..." << endl;
236
237 vector< double > fcal_gains_ch;
238 vector< double > fcal_pedestals_ch;
239
240 if (eventLoop->GetCalib("/FCAL/gains", fcal_gains_ch)){
241 jout << "DL1MCTrigger_factory: Error loading /FCAL/gains !" << endl;
242 // Load default values of gains if CCDB table is not found
243 for(int ii = 0; ii < DFCALGeometry::kBlocksTall; ii++){
244 for(int jj = 0; jj < DFCALGeometry::kBlocksWide; jj++){
245 fcal_gains[ii][jj] = 1.;
246 }
247 }
248 } else {
249 LoadFCALConst(fcal_gains, fcal_gains_ch, fcalGeom);
250
251 if(debug){
252 for(int ch = 0; ch < (int)fcal_gains_ch.size(); ch++){
253 int row = fcalGeom.row(ch);
254 int col = fcalGeom.column(ch);
255 // Sanity check for regular FCAL (row,col) ranges (anticipating
256 // future upgrade to FCAL to include insert)
257 if(fcalGeom.isBlockActive(row,col)&&row<DFCALGeometry::kBlocksTall
258 && col<DFCALGeometry::kBlocksWide){
259 hfcal_gains->Fill(fcal_gains[row][col]);
260 DVector2 aaa = fcalGeom.positionOnFace(row,col);
261 hfcal_gains2->Fill(float(aaa.X()), float(aaa.Y()), fcal_gains[row][col]);
262 // cout << aaa.X() << " " << aaa.Y() << endl;
263
264 }
265 }
266 }
267
268 }
269
270 if (eventLoop->GetCalib("/FCAL/pedestals", fcal_pedestals_ch)){
271 jout << "DL1MCTrigger_factory: Error loading /FCAL/pedestals !" << endl;
272 // Load default values of pedestals if CCDB table is not found
273 for(int ii = 0; ii < DFCALGeometry::kBlocksTall; ii++){
274 for(int jj = 0; jj < DFCALGeometry::kBlocksWide; jj++){
275 fcal_pedestals[ii][jj] = 100.;
276 }
277 }
278 } else {
279 LoadFCALConst(fcal_pedestals, fcal_pedestals_ch, fcalGeom);
280
281 if(debug){
282 for(int ch = 0; ch < (int)fcal_gains_ch.size(); ch++){
283 int row = fcalGeom.row(ch);
284 int col = fcalGeom.column(ch);
285 // Sanity check for regular FCAL (row,col) ranges (anticipating
286 // future upgrade to FCAL to include insert)
287 if(fcalGeom.isBlockActive(row,col)&&row<DFCALGeometry::kBlocksTall
288 && col<DFCALGeometry::kBlocksWide){
289 hfcal_ped->Fill(fcal_pedestals[row][col]);
290 }
291 }
292 }
293
294 }
295
296 if(!SIMU_BASELINE){
297 simu_baseline_fcal = 0;
298 simu_baseline_bcal = 0;
299 }
300
301 if(!SIMU_GAIN){
302 simu_gain_fcal = 0;
303 simu_gain_bcal = 0;
304 }
305
306 if(debug){
307 for(int ii = 0; ii < 100; ii++){
308 cout << " Channel = " << ii << " Value = " <<
309 fcal_gains_ch[ii] << endl;
310 }
311 }
312
313
314 return NOERROR;
315}
316
317//------------------
318// evnt
319//------------------
320jerror_t DL1MCTrigger_factory::evnt(JEventLoop *loop, uint64_t eventnumber){
321
322 if(BYPASS) {
323 DL1MCTrigger *trigger = new DL1MCTrigger;
324 trigger->trig_mask = 1;
325 _data.push_back(trigger);
326 return NOERROR;
327 }
328
329 int l1_found = 1;
330
331 int status = 0;
332
333 fcal_signal_hits.clear();
334 bcal_signal_hits.clear();
335
336 fcal_merged_hits.clear();
337 bcal_merged_hits.clear();
338
339 memset(fcal_ssp,0,sizeof(fcal_ssp));
340 memset(fcal_gtp,0,sizeof(fcal_gtp));
341
342 memset(bcal_ssp,0,sizeof(bcal_ssp));
343 memset(bcal_gtp,0,sizeof(bcal_gtp));
344
345
346 vector<const DFCALHit*> fcal_hits;
347 vector<const DBCALHit*> bcal_hits;
348
349 loop->Get(fcal_hits);
350 loop->Get(bcal_hits);
351
352 DRandom2 gDRandom(0); // declared extern in DRandom2.h
353
354 // This is temporary, to allow this simulation to be run on data
355 // to help out with trigger efficiency studies - sdobbs (Aug. 26, 2020)
356 if( loop->GetJEvent().GetStatusBit(kSTATUS_EVIO) ){
357 if(print_data_message) {
358 jout << "WARNING: Running L1 trigger simulation on EVIO data" << endl;
359 print_data_message = false;
360 }
361
362 // for data, don't add in baseline shifts, since they already exist
363 simu_baseline_fcal = 0;
364 simu_baseline_bcal = 0;
365 } else {
366 // Initialize random number generator
367 // Read seeds from hddm file
368 // Generate seeds according to the event number if they are not stored in hddm
369 // The proceedure is consistent with the mcsmear
370
371 UInt_t seed1 = 0;
372 UInt_t seed2 = 0;
373 UInt_t seed3 = 0;
374
375 GetSeeds(loop, eventnumber, seed1, seed2, seed3);
376
377 gDRandom.SetSeeds(seed1, seed2, seed3);
378 }
379
380 // cout << endl;
381 // cout << " Event = " << eventnumber << endl;
382 // cout << " Seed 1: " << seed1 << endl;
383 // cout << " Seed 2: " << seed2 << endl;
384 // cout << " Seed 3: " << seed3 << endl;
385 // cout << endl;
386
387
388 DL1MCTrigger *trigger = new DL1MCTrigger;
389
390 // FCAL energy sum
391 double fcal_hit_en = 0;
392
393 for (unsigned int ii = 0; ii < fcal_hits.size(); ii++){
394 int row = fcal_hits[ii]->row;
395 int col = fcal_hits[ii]->column;
396
397 // Shift time to simulate pile up hits
398 double time = fcal_hits[ii]->t + time_shift;
399 if((time < time_min) || (time > time_max)){
400 continue;
401 }
402
403 // Check channels masked for trigger
404 int ch_masked = 0;
405 for(unsigned int jj = 0; jj < fcal_trig_mask.size(); jj++){
406 if( (row == fcal_trig_mask[jj].row) && (col == fcal_trig_mask[jj].col)){
407 ch_masked = 1;
408 break;
409 }
410 }
411
412 if(ch_masked == 0){
413 fcal_hit_en += fcal_hits[ii]->E;
414
415 fcal_signal fcal_tmp;
416
417 fcal_tmp.row = row;
418 fcal_tmp.column = col;
419
420 fcal_tmp.energy = fcal_hits[ii]->E;
421 fcal_tmp.time = time;
422 memset(fcal_tmp.adc_amp,0,sizeof(fcal_tmp.adc_amp));
423 memset(fcal_tmp.adc_en, 0,sizeof(fcal_tmp.adc_en));
424
425 double fcal_adc_en = fcal_tmp.energy*FCAL_ADC_PER_MEV*1000;
426
427 // Account for gain fluctuations
428 if(simu_gain_fcal && row<DFCALGeometry::kBlocksTall
429 && col<DFCALGeometry::kBlocksWide){
430
431 double gain = fcal_gains[row][col];
432
433 fcal_adc_en *= gain;
434 }
435
436 status = SignalPulse(fcal_adc_en, fcal_tmp.time, fcal_tmp.adc_en, 1);
437 status = 0;
438
439 fcal_signal_hits.push_back(fcal_tmp);
440 }
441
442 }
443
444
445 // Merge FCAL hits
446 for(unsigned int ii = 0; ii < fcal_signal_hits.size(); ii++){
447
448 if(fcal_signal_hits[ii].merged == 1) continue;
449
450 fcal_signal fcal_tmp;
451 fcal_tmp.row = fcal_signal_hits[ii].row;
452 fcal_tmp.column = fcal_signal_hits[ii].column;
453 fcal_tmp.energy = 0.;
454 fcal_tmp.time = 0.;
455 for(int kk = 0; kk < sample; kk++)
456 fcal_tmp.adc_en[kk] = fcal_signal_hits[ii].adc_en[kk];
457
458 for(unsigned int jj = ii + 1; jj < fcal_signal_hits.size(); jj++){
459 if((fcal_signal_hits[ii].row == fcal_signal_hits[jj].row) &&
460 (fcal_signal_hits[ii].column == fcal_signal_hits[jj].column)){
461
462 fcal_signal_hits[jj].merged = 1;
463
464 for(int kk = 0; kk < sample; kk++)
465 fcal_tmp.adc_en[kk] += fcal_signal_hits[jj].adc_en[kk];
466 }
467 }
468
469 fcal_merged_hits.push_back(fcal_tmp);
470 }
471
472 // Add baseline fluctuations for channels with hits
473 if(simu_baseline_fcal){
474 for(unsigned int ii = 0; ii < fcal_merged_hits.size(); ii++){
475 int row = fcal_merged_hits[ii].row;
476 int column = fcal_merged_hits[ii].column;
477 double pedestal = 100.0;
478 if (row<DFCALGeometry::kBlocksTall
479 && column<DFCALGeometry::kBlocksWide){
480 pedestal=fcal_pedestals[row][column];
481 }
482 AddBaseline(fcal_merged_hits[ii].adc_en, pedestal, gDRandom);
483 }
484 }
485
486
487 // Digitize
488 for(unsigned int ii = 0; ii < fcal_merged_hits.size(); ii++){
489 Digitize(fcal_merged_hits[ii].adc_en,fcal_merged_hits[ii].adc_amp);
490 // cout << " Digitize " << fcal_merged_hits[ii].adc_en[sample - 3]
491 // << " " << fcal_merged_hits[ii].adc_amp[sample - 3] << endl;
492 }
493
494
495 int fcal_hit_adc_en = 0;
496
497 for(unsigned int ii = 0; ii < fcal_merged_hits.size(); ii++)
498 for(int jj = 0; jj < sample; jj++)
499 fcal_hit_adc_en += fcal_merged_hits[ii].adc_amp[jj];
500
501
502 status += FADC_SSP(fcal_merged_hits, 1);
503
504 status += GTP(1);
505
506
507
508 // BCAL energy sum
509 double bcal_hit_en = 0;
510
511 for (unsigned int ii = 0; ii < bcal_hits.size(); ii++){
512
513 // Shift time to simulate pile up hits
514 double time = bcal_hits[ii]->t + time_shift;
515 if((time < time_min) || (time > time_max)){
516 continue;
517 }
518
519 int module = bcal_hits[ii]->module;
520 int layer = bcal_hits[ii]->layer;
521 int sector = bcal_hits[ii]->sector;
522 int end = bcal_hits[ii]->end;
523
524 // Check trigger masks
525 int ch_masked = 0;
526
527 for(unsigned int jj = 0; jj < bcal_trig_mask.size(); jj++){
528 if( (module == bcal_trig_mask[jj].module) && (layer == bcal_trig_mask[jj].layer) &&
529 (sector == bcal_trig_mask[jj].sector) && (end == bcal_trig_mask[jj].end) ){
530 ch_masked = 1;
531 break;
532 }
533 }
534
535
536 if(ch_masked == 0){
537 bcal_hit_en += bcal_hits[ii]->E;
538
539 bcal_signal bcal_tmp;
540 bcal_tmp.module = module;
541 bcal_tmp.layer = layer;
542 bcal_tmp.sector = sector;
543 bcal_tmp.end = end;
544
545 bcal_tmp.energy = bcal_hits[ii]->E;
546 bcal_tmp.time = time;
547 memset(bcal_tmp.adc_amp,0,sizeof(bcal_tmp.adc_amp));
548 memset(bcal_tmp.adc_en, 0,sizeof(bcal_tmp.adc_en));
549
550 double bcal_adc_en = bcal_tmp.energy*BCAL_ADC_PER_MEV*1000;
551
552 status = SignalPulse(bcal_adc_en, bcal_tmp.time, bcal_tmp.adc_en, 2);
553 status = 0;
554
555 bcal_signal_hits.push_back(bcal_tmp);
556 }
557 }
558
559
560 // Merge BCAL hits
561 for(unsigned int ii = 0; ii < bcal_signal_hits.size(); ii++){
562
563 if(bcal_signal_hits[ii].merged == 1) continue;
564
565 bcal_signal bcal_tmp;
566 bcal_tmp.module = bcal_signal_hits[ii].module;
567 bcal_tmp.layer = bcal_signal_hits[ii].layer;
568 bcal_tmp.sector = bcal_signal_hits[ii].sector;
569 bcal_tmp.end = bcal_signal_hits[ii].end;
570
571 bcal_tmp.energy = 0.;
572 bcal_tmp.time = 0.;
573
574 for(int kk = 0; kk < sample; kk++)
575 bcal_tmp.adc_en[kk] = bcal_signal_hits[ii].adc_en[kk];
576
577 for(unsigned int jj = ii + 1; jj < bcal_signal_hits.size(); jj++){
578 if((bcal_signal_hits[ii].module == bcal_signal_hits[jj].module) &&
579 (bcal_signal_hits[ii].layer == bcal_signal_hits[jj].layer) &&
580 (bcal_signal_hits[ii].sector == bcal_signal_hits[jj].sector) &&
581 (bcal_signal_hits[ii].end == bcal_signal_hits[jj].end)){
582
583 bcal_signal_hits[jj].merged = 1;
584
585 for(int kk = 0; kk < sample; kk++)
586 bcal_tmp.adc_en[kk] += bcal_signal_hits[jj].adc_en[kk];
587 }
588 }
589
590 bcal_merged_hits.push_back(bcal_tmp);
591 }
592
593
594 // Add baseline fluctuations for channels with hits
595 if(simu_baseline_bcal){
596 for(unsigned int ii = 0; ii < bcal_merged_hits.size(); ii++){
597 // Assume that all BCAL pedestals are 100
598 double pedestal = TRIG_BASELINE;
599 AddBaseline(bcal_merged_hits[ii].adc_en, pedestal, gDRandom);
600 }
601 }
602
603
604 // Digitize
605 for(unsigned int ii = 0; ii < bcal_merged_hits.size(); ii++)
606 Digitize(bcal_merged_hits[ii].adc_en,bcal_merged_hits[ii].adc_amp);
607
608
609 int bcal_hit_adc_en = 0;
610 for(unsigned int ii = 0; ii < bcal_merged_hits.size(); ii++)
611 for(int jj = 0; jj < sample; jj++)
612 bcal_hit_adc_en += bcal_merged_hits[ii].adc_amp[jj];
613
614
615 status = FADC_SSP(bcal_merged_hits, 2);
616
617 status = GTP(2);
Value stored to 'status' is never read
618
619 // Search for triggers
620
621 l1_found = FindTriggers(trigger);
622
623 if(l1_found){
624
625 int fcal_gtp_max = 0;
626 int bcal_gtp_max = 0;
627
628 for(unsigned int ii = 0; ii < sample; ii++){
629 if(fcal_gtp[ii] > fcal_gtp_max) fcal_gtp_max = fcal_gtp[ii];
630 if(bcal_gtp[ii] > bcal_gtp_max) bcal_gtp_max = bcal_gtp[ii];
631 }
632
633 trigger->fcal_en = fcal_hit_en;
634 trigger->fcal_adc = fcal_hit_adc_en;
635 trigger->fcal_adc_en = fcal_hit_adc_en/FCAL_ADC_PER_MEV/1000.;
636 trigger->fcal_gtp = fcal_gtp_max;
637 trigger->fcal_gtp_en = fcal_gtp_max/FCAL_ADC_PER_MEV/1000.;
638
639 trigger->bcal_en = bcal_hit_en;
640 trigger->bcal_adc = bcal_hit_adc_en;
641 trigger->bcal_adc_en = bcal_hit_adc_en/BCAL_ADC_PER_MEV_CORRECT/2./1000.;
642 trigger->bcal_gtp = bcal_gtp_max;
643 trigger->bcal_gtp_en = bcal_gtp_max/BCAL_ADC_PER_MEV_CORRECT/2./1000.;
644
645 _data.push_back(trigger);
646
647 } else{
648 delete trigger;
649 }
650
651 return NOERROR;
652}
653
654//------------------
655// erun
656//------------------
657jerror_t DL1MCTrigger_factory::erun(void)
658{
659 return NOERROR;
660}
661
662//------------------
663// fini
664//------------------
665jerror_t DL1MCTrigger_factory::fini(void)
666{
667 return NOERROR;
668}
669
670//*********************
671// Read RCDB
672//*********************
673
674int DL1MCTrigger_factory::Read_RCDB(int32_t runnumber, bool print_messages)
675{
676
677#if HAVE_RCDB1
678
679 vector<const DTranslationTable*> ttab;
680 eventLoop->Get(ttab);
681
682 vector<string> SectionNames = {"TRIGGER", "GLOBAL", "FCAL", "BCAL", "TOF", "ST", "TAGH",
683 "TAGM", "PS", "PSC", "TPOL", "CDC", "FDC"};
684
685 string RCDB_CONNECTION;
686
687 if( getenv("RCDB_CONNECTION")!= NULL__null )
688 RCDB_CONNECTION = getenv("RCDB_CONNECTION");
689 else
690 RCDB_CONNECTION = "mysql://rcdb@hallddb.jlab.org/rcdb"; // default to outward-facing MySQL DB
691
692
693 rcdb::Connection connection(RCDB_CONNECTION);
694
695 auto rtvsCnd = connection.GetCondition(runnumber, "rtvs");
696
697 if( !rtvsCnd ) {
698 cout<<"'rtvs' condition is not set for run " << runnumber << endl;
699 return 2;
700 }
701
702
703 auto json = rtvsCnd->ToJsonDocument(); // The CODA rtvs is serialized as JSon dictionary.
704 string fileName(json["%(config)"].GetString()); // We need item with name '%(config)'
705
706 auto file = connection.GetFile(runnumber, fileName);
707
708 if(!file) { // If there is no such file, null is returned
709 cout<<"File with name: "<< fileName
710 <<" doesn't exist (not associated) with run: "<< runnumber << endl;
711 return 3;
712 }
713
714 string fileContent = file->GetContent(); // Get file content
715 auto result = rcdb::ConfigParser::Parse(fileContent, SectionNames); // Parse it!
716
717
718 // Pulse parameters for trigger
719
720 auto trig_thr = result.Sections["FCAL"].NameValues["FADC250_TRIG_THR"];
721 auto trig_nsb = result.Sections["FCAL"].NameValues["FADC250_TRIG_NSB"];
722 auto trig_nsa = result.Sections["FCAL"].NameValues["FADC250_TRIG_NSA"];
723
724 if(trig_thr.size() > 0){
725 FCAL_CELL_THR = stoi(trig_thr);
726 if(FCAL_CELL_THR < 0) FCAL_CELL_THR = 0;
727 }
728
729 if(trig_nsb.size() > 0)
730 FCAL_NSB = stoi(trig_nsb);
731
732 if(trig_thr.size() > 0)
733 FCAL_NSA = stoi(trig_nsa);
734
735 trig_thr = result.Sections["BCAL"].NameValues["FADC250_TRIG_THR"];
736 trig_nsb = result.Sections["BCAL"].NameValues["FADC250_TRIG_NSB"];
737 trig_nsa = result.Sections["BCAL"].NameValues["FADC250_TRIG_NSA"];
738
739 if(trig_thr.size() > 0){
740 BCAL_CELL_THR = stoi(trig_thr);
741 if(BCAL_CELL_THR < 0) BCAL_CELL_THR = 0;
742 }
743
744 if(trig_nsb.size() > 0)
745 BCAL_NSB = stoi(trig_nsb);
746
747 if(trig_thr.size() > 0)
748 BCAL_NSA = stoi(trig_nsa);
749
750
751 // List of enabled GTP equations and triggers
752 vector<vector<string>> triggerTypes;
753
754 for(auto row : result.Sections["TRIGGER"].Rows) {
755
756 if(row[0] == "TRIG_TYPE") {
757
758 if(row.size() >= 9){
759 triggerTypes.push_back(row); // The line in config file starts with TRIG_TYPE
760 } else {
761 cout << " Cannot parse TRIG_TYPE. Insufficient number of parameters " << row.size() << endl;
762 }
763 }
764
765
766 // Integration windows for BCAL and FCAL
767 if(row[0] == "TRIG_EQ") {
768 if(row.size() >= 5){
769 if(stoi(row[4]) == 1){ // Trigger enabled
770 if(row[1] == "FCAL")
771 FCAL_WINDOW = stoi(row[3]);
772 if(row[1] == "BCAL_E")
773 BCAL_WINDOW = stoi(row[3]);
774 }
775 }
776 }
777 }
778
779
780
781
782 for(int ii = 0; ii < 32; ii++){
783
784 int trig_found = 0;
785
786 trigger_conf trigger_tmp;
787
788 memset(&trigger_tmp,0,sizeof(trigger_tmp));
789
790 trigger_tmp.bit = ii;
791
792 for(unsigned int jj = 0; jj < triggerTypes.size(); jj++){
793
794 if(triggerTypes[jj].size() < 9) continue; // Minimum 9 parameters are required in the config file
795
796 int bit = stoi(triggerTypes[jj][8]); // Trigger bit
797
798 if( bit == ii){
799
800 // FCAL, BCAL triggers
801 if(triggerTypes[jj][1] == "BFCAL"){
802
803 int fcal = 0;
804 int bcal = 0;
805
806 if(triggerTypes[jj][4].size() > 0) fcal = stoi(triggerTypes[jj][4]);
807 if(triggerTypes[jj][5].size() > 0) bcal = stoi(triggerTypes[jj][5]);
808
809
810 if( (fcal > 0) && (bcal > 0)){ // FCAL & BCAL
811 trigger_tmp.type = 0x3;
812 trigger_tmp.gtp.fcal = fcal;
813 trigger_tmp.gtp.bcal = bcal;
814 if(triggerTypes[jj][6].size() > 0) trigger_tmp.gtp.en_thr = stoi(triggerTypes[jj][6]);
815
816 if(triggerTypes[jj].size() > 9){
817 if((triggerTypes[jj].size() >= 10) && (triggerTypes[jj][9].size() > 0)) trigger_tmp.gtp.fcal_min = stoi(triggerTypes[jj][9]);
818 if((triggerTypes[jj].size() >= 11) && (triggerTypes[jj][10].size() > 0)) trigger_tmp.gtp.fcal_max = stoi(triggerTypes[jj][10]);
819 if((triggerTypes[jj].size() >= 12) && (triggerTypes[jj][11].size() > 0)) trigger_tmp.gtp.bcal_min = stoi(triggerTypes[jj][11]);
820 if((triggerTypes[jj].size() >= 13) && (triggerTypes[jj][12].size() > 0)) trigger_tmp.gtp.bcal_max = stoi(triggerTypes[jj][12]);
821 }
822
823 trig_found++;
824 } else if ((fcal > 0) && (bcal == 0)){ // FCAL only
825 trigger_tmp.type = 0x1;
826 trigger_tmp.gtp.fcal = fcal;
827 if(triggerTypes[jj][6].size() > 0) trigger_tmp.gtp.en_thr = stoi(triggerTypes[jj][6]);
828 if(triggerTypes[jj].size() > 9){
829 if(triggerTypes[jj][9].size() > 0) trigger_tmp.gtp.fcal_min = stoi(triggerTypes[jj][9]);
830 if(triggerTypes[jj][10].size() > 0) trigger_tmp.gtp.fcal_max = stoi(triggerTypes[jj][10]);
831 }
832 trig_found++;
833 } else if ((bcal > 0) && (fcal == 0)){ // BCAL only
834 trigger_tmp.type = 0x2;
835 trigger_tmp.gtp.bcal = bcal;
836 if(triggerTypes[jj][6].size() > 0) trigger_tmp.gtp.en_thr = stoi(triggerTypes[jj][6]);
837 if((triggerTypes[jj].size() >= 12) && (triggerTypes[jj][11].size() > 0)) trigger_tmp.gtp.bcal_min = stoi(triggerTypes[jj][11]);
838 if((triggerTypes[jj].size() >= 13) && (triggerTypes[jj][12].size() > 0)) trigger_tmp.gtp.bcal_max = stoi(triggerTypes[jj][12]);
839
840
841 trig_found++;
842 } else {
843 cout << " Incorrect parameters for BFCAL trigger " << endl;
844 }
845
846 } else if(triggerTypes[jj][1] == "ST"){ // ST
847 trigger_tmp.type = 0x4;
848 if(triggerTypes[jj][7].size() > 0) trigger_tmp.gtp.st_nhit = stoi(triggerTypes[jj][7]);
849 if((triggerTypes[jj].size() >= 14) && (triggerTypes[jj][13].size() > 0)) trigger_tmp.gtp.st_pattern = stoi(triggerTypes[jj][13]);
850
851 trig_found++;
852 } else if(triggerTypes[jj][1] == "PS"){ // PS
853 trigger_tmp.type = 0x8;
854 if(triggerTypes[jj][7].size() > 0) trigger_tmp.gtp.ps_nhit = stoi(triggerTypes[jj][7]);
855 if((triggerTypes[jj].size() >= 14) && (triggerTypes[jj][13].size() > 0)) trigger_tmp.gtp.ps_pattern = stoi(triggerTypes[jj][13]);
856 trig_found++;
857 } else if(triggerTypes[jj][1] == "TAGH"){ // TAGH
858 trigger_tmp.type = 0x10;
859 if(triggerTypes[jj][7].size() > 0) trigger_tmp.gtp.tof_nhit = stoi(triggerTypes[jj][7]);
860 if((triggerTypes[jj].size() >= 14) && (triggerTypes[jj][13].size() > 0)){
861 trigger_tmp.gtp.tagh_pattern = stoi(triggerTypes[jj][13],nullptr,0);
862 }
863
864 trig_found++;
865 } else if(triggerTypes[jj][1] == "TOF"){ // TOF
866 trigger_tmp.type = 0x20;
867 if(triggerTypes[jj][7].size() > 0) trigger_tmp.gtp.ps_nhit = stoi(triggerTypes[jj][7]);
868 if((triggerTypes[jj].size() >= 14) && (triggerTypes[jj][13].size() > 0)) trigger_tmp.gtp.tof_pattern = stoi(triggerTypes[jj][13]);
869 trig_found++;
870 } else {
871 cout << " Incorrect Trigger type " << triggerTypes[jj][1] << endl;
872 }
873
874 if(trig_found > 0){
875 triggers_enabled.push_back(trigger_tmp);
876 }
877
878 } // Trigger lane is enabled
879 } // Trigger types
880
881 }
882
883
884
885 // Load FCAL Trigger Masks
886
887 string comDir = result.Sections["FCAL"].NameValues["FADC250_COM_DIR"];
888 string comVer = result.Sections["FCAL"].NameValues["FADC250_COM_VER"];
889 string userDir = result.Sections["FCAL"].NameValues["FADC250_USER_DIR"];
890 string userVer = result.Sections["FCAL"].NameValues["FADC250_USER_VER"];
891
892
893 for(int crate = 1; crate <= 12; crate++){
894
895 std::string s = std::to_string(crate);
896
897 string comFileName = comDir + "/rocfcal" +s + "_fadc250_" + comVer + ".cnf";
898 string userFileName = userDir + "/rocfcal" +s + "_" + userVer + ".cnf";
899
900 auto userFile = connection.GetFile(runnumber, userFileName);
901
902 if(!userFile) {
903 // cout<<" USER File with name: "<< userFileName
904 // <<" doesn't exist (not associated) with run: "<< runnumber <<endl;
905 continue;
906 }
907
908
909 auto userParseResult = rcdb::ConfigParser::ParseWithSlots(userFile->GetContent(), "FADC250_SLOTS");
910
911
912 for(unsigned int slot = 3; slot <= 21; slot++){
913
914 auto userValues = userParseResult.SectionsBySlotNumber[slot].NameVectors["FADC250_TRG_MASK"]; // Parse it and return
915
916 if(userValues.size() > 0){
917 for (unsigned int ch = 0; ch < userValues.size(); ++ch) {
918
919 if(userValues[ch].size() == 0) continue;
920
921 if(stoi(userValues[ch]) > 0){
922
923 uint32_t roc_id = 10 + crate;
924
925 DTranslationTable::csc_t daq_index = {roc_id, slot, ch };
926
927 DTranslationTable::DChannelInfo channel_info;
928
929 try {
930 channel_info = ttab[0]->GetDetectorIndex(daq_index);
931 }
932
933 catch(...){
934 if(print_messages) cout << "Exception: FCAL channel is not in the translation table " << " Crate = " << 10 + crate << " Slot = " << slot <<
935 " Channel = " << ch << endl;
936 continue;
937 }
938
939 fcal_mod tmp;
940
941 tmp.roc = crate;
942 tmp.slot = slot;
943 tmp.ch = ch;
944
945 tmp.row = channel_info.fcal.row;
946 tmp.col = channel_info.fcal.col;
947
948 fcal_trig_mask.push_back(tmp);
949 }
950
951 }
952
953 }
954
955 } // Loop over slots
956 } // Loop over crates
957
958
959
960
961 // Load BCAL Trigger Masks
962
963 comDir = result.Sections["BCAL"].NameValues["FADC250_COM_DIR"];
964 comVer = result.Sections["BCAL"].NameValues["FADC250_COM_VER"];
965 userDir = result.Sections["BCAL"].NameValues["FADC250_USER_DIR"];
966 userVer = result.Sections["BCAL"].NameValues["FADC250_USER_VER"];
967
968
969 for(int crate = 1; crate <= 12; crate++){
970
971 if( (crate == 3) || (crate == 6) || (crate == 9) || (crate == 12)) continue;
972
973 std::string s = std::to_string(crate);
974
975 string comFileName = comDir + "/rocbcal" +s + "_fadc250_" + comVer + ".cnf";
976 string userFileName = userDir + "/rocbcal" +s + "_" + userVer + ".cnf";
977
978 auto userFile = connection.GetFile(runnumber, userFileName);
979
980 if(!userFile){
981 // cout<<" USER File with name: "<< userFileName
982 // <<" doesn't exist (not associated) with run: "<< runnumber <<endl;
983 continue;
984 }
985
986
987 auto userParseResult = rcdb::ConfigParser::ParseWithSlots(userFile->GetContent(), "FADC250_SLOTS");
988
989 for(unsigned int slot = 3; slot <= 21; slot++){
990
991 auto userValues = userParseResult.SectionsBySlotNumber[slot].NameVectors["FADC250_TRG_MASK"]; // Parse it and return
992
993 if(userValues.size() > 0){
994
995 for (unsigned int ch = 0; ch < userValues.size(); ++ch) {
996
997 if(userValues[ch].size() == 0) continue;
998
999 if(stoi(userValues[ch]) > 0){
1000
1001 uint32_t roc_id = 30 + crate;
1002
1003 DTranslationTable::csc_t daq_index = {roc_id, slot, ch };
1004
1005 DTranslationTable::DChannelInfo channel_info;
1006
1007 try {
1008 channel_info = ttab[0]->GetDetectorIndex(daq_index);
1009 }
1010
1011 catch(...){
1012 cout << "Exception: BCAL channel is not in the translation table " << " Crate = " << 30 + crate << " Slot = " << slot <<
1013 " Channel = " << ch << endl;
1014 continue;
1015 }
1016
1017
1018 bcal_mod tmp;
1019
1020 tmp.roc = crate;
1021 tmp.slot = slot;
1022 tmp.ch = ch;
1023
1024 tmp.module = channel_info.bcal.module;
1025 tmp.layer = channel_info.bcal.layer;
1026 tmp.sector = channel_info.bcal.sector;
1027 tmp.end = channel_info.bcal.end;
1028
1029 bcal_trig_mask.push_back(tmp);
1030 }
1031
1032 }
1033
1034 }
1035
1036 } // Loop over slots
1037 } // Loop over crates
1038
1039 return 0;
1040
1041#else // HAVE_RCDB
1042
1043 return 10; // RCDB is not available
1044
1045#endif
1046
1047}
1048
1049
1050int DL1MCTrigger_factory::SignalPulse(double en, double time, double amp_array[sample], int type){
1051
1052
1053 // Parameterize and digitize pulse shapes. Sum up amplitudes
1054 // type = 1 - FCAL
1055 // type = 2 - BCAL
1056
1057 float exp_par = 0.358;
1058
1059 int pulse_length = 20;
1060
1061 if(type == 2) exp_par = 0.18;
1062
1063 int sample_first = (int)floor(time/time_stamp);
1064
1065 // digitization range
1066 int ind_min = sample_first + 1;
1067 int ind_max = ind_min + pulse_length + 1;
1068
1069 if( (ind_min > sample) || (ind_min < 0)){
1070 // cout << " SignalPulse() FATAL error: time out of range " << time << " " << ind_min << endl;
1071 return 1;
1072 }
1073
1074 if(ind_max > sample){
1075 // cout << " SignalPulse(): ind_max set to maximum" << time << " " << ind_max << endl;
1076 ind_max = sample - 1;
1077 }
1078
1079 for(int i = ind_min; i < ind_max; i++ ){
1080 double adc_t = time_stamp*i - time;
1081 double amp = exp_par*exp_par*exp(-adc_t*exp_par)*adc_t;
1082
1083 // amp_array[i] += (int)(amp*time_stamp*en + 0.5);
1084 // if(amp_array[i] > max_adc_bins){
1085 // amp_array[i] = max_adc_bins;
1086 // }
1087
1088 amp_array[i] += amp*time_stamp*en;
1089
1090 }
1091
1092 return 0;
1093}
1094
1095int DL1MCTrigger_factory::GTP(int detector){
1096
1097 // 1 - FCAL
1098 // 2 - BCAL
1099
1100 int INT_WINDOW = 20;
1101
1102 switch(detector){
1103 case 1:
1104 INT_WINDOW = FCAL_WINDOW;
1105 break;
1106 case 2:
1107 INT_WINDOW = BCAL_WINDOW;
1108 break;
1109 default:
1110 break;
1111 }
1112
1113 int index_min = 0;
1114 int index_max = 0;
1115
1116 for(unsigned int samp = 0; samp < sample; samp++){
1117
1118 index_max = samp;
1119 index_min = samp - INT_WINDOW;
1120
1121 if(index_min < 0) index_min = 0;
1122
1123 int energy_sum = 0;
1124
1125 for(int ii = index_min; ii <= index_max; ii++){
1126 if(detector == 1)
1127 energy_sum += fcal_ssp[ii];
1128 else
1129 energy_sum += bcal_ssp[ii];
1130 }
1131
1132 if(detector == 1)
1133 fcal_gtp[samp] = energy_sum;
1134 else
1135 bcal_gtp[samp] = energy_sum;
1136 }
1137
1138 return 0;
1139
1140}
1141
1142
1143template <typename T> int DL1MCTrigger_factory::FADC_SSP(vector<T> merged_hits, int detector){
1144
1145 // 1 - FCAL
1146 // 2 - BCAL
1147
1148 int EN_THR = 4096;
1149 int NSA = 10;
1150 int NSB = 3;;
1151
1152 switch(detector){
1153 case 1:
1154 EN_THR = FCAL_CELL_THR;
1155 NSA = FCAL_NSA;
1156 NSB = FCAL_NSB;
1157 break;
1158 case 2:
1159 EN_THR = BCAL_CELL_THR;
1160 NSA = BCAL_NSA;
1161 NSB = BCAL_NSB;
1162 break;
1163 default:
1164 break;
1165 }
1166
1167 for(unsigned int hit = 0; hit < merged_hits.size(); hit++){
1168
1169 int index_min = -10;
1170 int index_max = -10;
1171
1172 for(int ii = 0; ii < sample; ii++){
1173
1174 int pulse_found = 0;
1175
1176 if(merged_hits[hit].adc_amp[ii] >= EN_THR){
1177 pulse_found = 1;
1178 } else {
1179 continue;
1180 }
1181
1182 index_min = ii - NSB;
1183
1184 if(index_max > index_min) index_min = index_max + 1;
1185
1186 index_max = ii + NSA - 1;
1187
1188 if(index_min < 0) index_min = 0;
1189
1190 if(index_max >= sample){
1191 index_max = sample - 1;
1192 }
1193
1194 int extend_nsa = 1;
1195
1196 // Extend FADC range if needed
1197
1198 while(extend_nsa){
1199
1200 int index_tmp = index_max + 1;
1201
1202 if(index_tmp < sample){
1203 if(merged_hits[hit].adc_amp[index_tmp] >= EN_THR){
1204 index_max += NSA;
1205 } else extend_nsa = 0;
1206 } else extend_nsa = 0;
1207 }
1208
1209 if(index_max >= sample)
1210 index_max = sample - 1;
1211
1212 for(int kk = index_min; kk <= index_max; kk++){
1213 if(detector == 1){
1214 if((merged_hits[hit].adc_amp[kk] - 100) > 0)
1215 fcal_ssp[kk] += (merged_hits[hit].adc_amp[kk] - TRIG_BASELINE);
1216 }
1217 else if(detector == 2){
1218 if((merged_hits[hit].adc_amp[kk] - 100) > 0)
1219 bcal_ssp[kk] += (merged_hits[hit].adc_amp[kk] - TRIG_BASELINE);
1220 }
1221 }
1222
1223 if(pulse_found == 1){
1224 ii = index_max + 1;
1225 pulse_found = 0;
1226 }
1227
1228 }
1229
1230 }
1231
1232 return 0;
1233}
1234
1235void DL1MCTrigger_factory::PrintTriggers(){
1236
1237 string detector;
1238 int nhit = 0;
1239 unsigned int pattern = 0;
1240
1241 cout << endl << endl;
1242 cout << " ------------ Trigger Settings --------------- " << endl;
1243 cout << endl << endl;
1244
1245 cout << "----------- FCAL ----------- " << endl << endl;
1246
1247 cout << "FCAL_CELL_THR = " << setw(10) << FCAL_CELL_THR << endl;
1248 cout << "FCAL_NSA = " << setw(10) << FCAL_NSA << endl;
1249 cout << "FCAL_NSB = " << setw(10) << FCAL_NSB << endl;
1250 cout << "FCAL_WINDOW = " << setw(10) << FCAL_WINDOW << endl;
1251
1252 cout << endl;
1253
1254 cout << "----------- BCAL ----------- " << endl << endl;
1255
1256 cout << "BCAL_CELL_THR = " << setw(10) << BCAL_CELL_THR << endl;
1257 cout << "BCAL_NSA = " << setw(10) << BCAL_NSA << endl;
1258 cout << "BCAL_NSB = " << setw(10) << BCAL_NSB << endl;
1259 cout << "BCAL_WINDOW = " << setw(10) << BCAL_WINDOW << endl;
1260
1261 cout << endl << endl;
1262
1263
1264 if(triggers_enabled.size() > 0){
1265 cout << "TYPE " << "FCAL_E " << "BCAL_E " <<
1266 "EN_THR " << "NHIT " << "LANE " << "FCAL_EMIN " << "FCAL_EMAX " <<
1267 "BCAL_EMIN " << "BCAL_EMAX " << "PATTERN " << endl;
1268 }
1269
1270
1271 for(unsigned int ii = 0; ii < triggers_enabled.size(); ii++){
1272
1273 switch(triggers_enabled[ii].type){
1274 case 1: detector = "BFCAL ";
1275 break;
1276 case 2: detector = "BFCAL ";
1277 break;
1278 case 3: detector = "BFCAL ";
1279 break;
1280 case 4: detector = "ST ";
1281 nhit = triggers_enabled[ii].gtp.st_nhit;
1282 pattern = triggers_enabled[ii].gtp.st_pattern;
1283 break;
1284 case 8: detector = "PS ";
1285 nhit = triggers_enabled[ii].gtp.ps_nhit;
1286 pattern = triggers_enabled[ii].gtp.ps_pattern;
1287 break;
1288 case 16: detector = "TAGH ";
1289 nhit = triggers_enabled[ii].gtp.tagh_nhit;
1290 pattern = triggers_enabled[ii].gtp.tagh_pattern;
1291 break;
1292 case 32: detector = "TOF ";
1293 nhit = triggers_enabled[ii].gtp.tof_nhit;
1294 pattern = triggers_enabled[ii].gtp.tof_pattern;
1295 break;
1296 default: detector = "NONE ";
1297 cout << " Unknown detector ===== " << triggers_enabled[ii].type << endl;
1298 break;
1299 }
1300
1301 cout << detector << setw(6) <<
1302 triggers_enabled[ii].gtp.fcal << setw(9) <<
1303 triggers_enabled[ii].gtp.bcal << setw(11) <<
1304 triggers_enabled[ii].gtp.en_thr << setw(6) <<
1305 nhit << setw(9) <<
1306 triggers_enabled[ii].bit << setw(12) <<
1307 triggers_enabled[ii].gtp.fcal_min << setw(14) <<
1308 triggers_enabled[ii].gtp.fcal_max << setw(10) <<
1309 triggers_enabled[ii].gtp.bcal_min << setw(14) <<
1310 triggers_enabled[ii].gtp.bcal_max << setw(8) <<
1311 hex << uppercase << "0x" << pattern << nouppercase << dec << endl;
1312
1313 }
1314
1315 cout << endl << endl;
1316
1317}
1318
1319
1320
1321int DL1MCTrigger_factory::FindTriggers(DL1MCTrigger *trigger){
1322
1323 int trig_found = 0;
1324
1325 // Main production trigger
1326 for(unsigned int ii = 0; ii < triggers_enabled.size(); ii++){
1327
1328 if(triggers_enabled[ii].bit == 0){ // Main production trigger found
1329
1330 int gtp_energy = 0;
1331 int bcal_energy = 0;
1332
1333 for(unsigned int samp = 0; samp < sample; samp++){
1334
1335 int bcal_samp = samp - BCAL_OFFSET;
1336
1337 if(bcal_samp < 0){
1338 bcal_energy = 0;
1339 } else if(bcal_samp >= sample){
1340 bcal_energy = 0;
1341 } else{
1342 bcal_energy = bcal_gtp[bcal_samp];
1343 }
1344
1345
1346 gtp_energy = triggers_enabled[ii].gtp.fcal*fcal_gtp[samp] +
1347 triggers_enabled[ii].gtp.bcal*bcal_energy;
1348
1349 if(gtp_energy >= triggers_enabled[ii].gtp.en_thr){
1350
1351 if(triggers_enabled[ii].gtp.fcal_min > 0) { // FCAL > 0
1352 if(fcal_gtp[samp] > triggers_enabled[ii].gtp.fcal_min){
1353 trigger->trig_mask = (trigger->trig_mask | 0x1);
1354 trigger->trig_time[0] = samp - 25;
1355 trig_found++;
1356 }
1357 } else {
1358
1359 trigger->trig_mask = (trigger->trig_mask | 0x1);
1360 trigger->trig_time[0] = samp - 25;
1361 trig_found++;
1362 }
1363
1364 break;
1365
1366 } // Check energy threshold
1367 }
1368 } // Trigger Bit 0
1369
1370
1371 if(triggers_enabled[ii].bit == 2){ // BCAL trigger found
1372
1373 int gtp_energy = 0;
1374 int bcal_energy = 0;
1375
1376 for(unsigned int samp = 0; samp < sample; samp++){
1377
1378 int bcal_samp = samp - BCAL_OFFSET;
1379
1380 if(bcal_samp < 0){
1381 bcal_energy = 0;
1382 } else if(bcal_samp >= sample){
1383 bcal_energy = 0;
1384 } else{
1385 bcal_energy = bcal_gtp[bcal_samp];
1386 }
1387
1388 gtp_energy = triggers_enabled[ii].gtp.bcal*bcal_energy;
1389
1390 if(gtp_energy >= triggers_enabled[ii].gtp.en_thr){
1391
1392 trigger->trig_mask = (trigger->trig_mask | 0x4);
1393 trigger->trig_time[2] = samp - 25;
1394 trig_found++;
1395
1396 break;
1397
1398 } // Energy threshold
1399 }
1400 } // Trigger Bit 3
1401
1402 } // Loop over triggers found in the config file
1403
1404 return trig_found;
1405
1406}
1407
1408
1409// Fill fcal calibration tables similar to FCALHit factory
1410void DL1MCTrigger_factory::LoadFCALConst(fcal_constants_t &table, const vector<double> &fcal_const_ch,
1411 const DFCALGeometry &fcalGeom){
1412 for (int ch = 0; ch < static_cast<int>(fcal_const_ch.size()); ch++) {
1413 int row = fcalGeom.row(ch);
1414 int col = fcalGeom.column(ch);
1415 table[row][col] = fcal_const_ch[ch];
1416 }
1417
1418}
1419
1420void DL1MCTrigger_factory::Digitize(double adc_amp[sample], int adc_count[sample]){
1421
1422 for(int samp = 0; samp < sample; samp++ ){
1423
1424 adc_count[samp] += (int)(adc_amp[samp] + TRIG_BASELINE + 0.5);
1425
1426 if(adc_count[samp] > max_adc_bins)
1427 adc_count[samp] = max_adc_bins;
1428
1429 }
1430}
1431
1432
1433void DL1MCTrigger_factory::AddBaseline(double adc_amp[sample], double pedestal, DRandom2 &gDRandom){
1434
1435 double pedestal_correct = pedestal - TRIG_BASELINE;
1436
1437 for(int samp = 0; samp < sample; samp++ ){
1438 double tmp = gDRandom.Gaus(pedestal_correct, pedestal_sigma);
1439 adc_amp[samp] += tmp;
1440 // cout << " " << tmp;
1441 }
1442
1443 if(debug){
1444 cout << endl;
1445 cout << " Corrected pedestals = " << pedestal_correct << " " << adc_amp[sample - 2]
1446 << " " << pedestal_sigma << endl;
1447 }
1448
1449}
1450
1451
1452
1453void DL1MCTrigger_factory::GetSeeds(JEventLoop *loop, uint64_t eventnumber, UInt_t &seed1, UInt_t &seed2, UInt_t &seed3){
1454
1455 // Use seeds similar to mcsmear
1456
1457 JEvent& event = loop->GetJEvent();
1458
1459 JEventSource *source = event.GetJEventSource();
1460
1461 DEventSourceHDDM *hddm_source = dynamic_cast<DEventSourceHDDM*>(source);
1462
1463 if (!hddm_source) {
1464
1465 cerr << "DL1MCTrigger_factory: This program MUST be used with an HDDM file as input!" << endl;
1466 cerr << " Default seeds will be used for the random generator " << endl;
1467 seed1 = 259921049 + eventnumber;
1468 seed2 = 442249570 + eventnumber;
1469 seed3 = 709975946 + eventnumber;
1470 } else {
1471
1472 hddm_s::HDDM *record = (hddm_s::HDDM*)event.GetRef();
1473 if (!record){
1474 seed1 = 259921049 + eventnumber;
1475 seed2 = 442249570 + eventnumber;
1476 seed3 = 709975946 + eventnumber;
1477 } else {
1478
1479
1480 hddm_s::ReactionList::iterator reiter = record->getReactions().begin();
1481
1482 hddm_s::Random my_rand = reiter->getRandom();
1483
1484 // Copy seeds from event record to local variables
1485 seed1 = my_rand.getSeed1();
1486 seed2 = my_rand.getSeed2();
1487 seed3 = my_rand.getSeed3();
1488
1489 // If no seeds are stored in the hddm file, generate them in the same way
1490 // as in mcsmear
1491 if ((seed1 == 0) || (seed2 == 0) || (seed3 == 0)){
1492 uint64_t eventNo = record->getPhysicsEvent().getEventNo();
1493 seed1 = 259921049 + eventNo;
1494 seed2 = 442249570 + eventNo;
1495 seed3 = 709975946 + eventNo;
1496
1497 }
1498
1499 } // Record doesn't exist
1500 } // Not an HDDM file
1501}