1
/* This file is part of MAUS: http://micewww.pp.rl.ac.uk:8080/projects/maus
3
* MAUS is free software: you can redistribute it and/or modify
4
* it under the terms of the GNU General Public License as published by
5
* the Free Software Foundation, either version 3 of the License, or
6
* (at your option) any later version.
8
* MAUS is distributed in the hope that it will be useful,
9
* but WITHOUT ANY WARRANTY; without even the implied warranty of
10
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11
* GNU General Public License for more details.
13
* You should have received a copy of the GNU General Public License
14
* along with MAUS. If not, see <http://www.gnu.org/licenses/>.
18
#include "Utils/CppErrorHandler.hh"
19
#include "Utils/JsonWrapper.hh"
20
#include "Interface/dataCards.hh"
21
#include "DataStructure/MCEvent.hh"
22
#include "DataStructure/Spill.hh"
23
#include "API/PyWrapMapBase.hh"
24
#include "Converter/DataConverters/CppJsonSpillConverter.hh"
25
#include "Converter/DataConverters/JsonCppSpillConverter.hh"
27
#include "src/map/MapCppEMRMCDigitization/MapCppEMRMCDigitization.hh"
31
PyMODINIT_FUNC init_MapCppEMRMCDigitization(void) {
32
PyWrapMapBase<MAUS::MapCppEMRMCDigitization>::PyWrapMapBaseModInit
33
("MapCppEMRMCDigitization", "", "", "", "");
36
MapCppEMRMCDigitization::MapCppEMRMCDigitization()
37
: MapBase<MAUS::Data>("MapCppEMRMCDigitization") {
40
void MapCppEMRMCDigitization::_birth(const std::string& argJsonConfigDocument) {
41
_classname = "MapCppEMRMCDigitization";
42
char* pMAUS_ROOT_DIR = getenv("MAUS_ROOT_DIR");
44
if (!pMAUS_ROOT_DIR) {
45
throw MAUS::Exception(Exception::recoverable,
46
"Could not resolve ${MAUS_ROOT_DIR} environment variable",
47
"MapCppEMRMCDigitization::birth");
51
Json::Value configJSON;
52
configJSON = JsonWrapper::StringToJson(argJsonConfigDocument);
55
_number_of_planes = configJSON["EMRnumberOfPlanes"].asInt();
56
_number_of_bars = configJSON["EMRnumberOfBars"].asInt();
58
_tot_noise_up = configJSON["EMRtotNoiseUp"].asInt();
59
_tot_noise_low = configJSON["EMRtotNoiseLow"].asInt();
61
_do_sampling = configJSON["EMRdoSampling"].asInt();
62
_nph_per_MeV = configJSON["EMRnphPerMeV"].asInt();
63
_seed = configJSON["EMRseed"].asInt();
64
_trap_eff = configJSON["EMRtrapEff"].asDouble();
65
_QE_MAPMT = configJSON["EMRqeMAPMT"].asDouble();
66
_QE_SAPMT = configJSON["EMRqeSAPMT"].asDouble();
67
_nADC_per_pe_MAPMT = configJSON["EMRnadcPerPeMAPMT"].asDouble();
68
_nADC_per_pe_SAPMT = configJSON["EMRnadcPerPeSAPMT"].asDouble();
69
_electronics_response_spread_MAPMT = configJSON["EMRelectronicsResponseSpreadMAPMT"].asDouble();
70
_electronics_response_spread_SAPMT = configJSON["EMRelectronicsResponseSpreadSAPMT"].asDouble();
71
_dbb_count = configJSON["EMRdbbCount"].asDouble();
72
_fadc_count = configJSON["EMRfadcCount"].asDouble();
73
_time_response_spread = configJSON["EMRtimeResponseSpread"].asInt();
74
_tot_func_p1 = configJSON["EMRtotFuncP1"].asDouble();
75
_tot_func_p2 = configJSON["EMRtotFuncP2"].asDouble();
76
_tot_func_p3 = configJSON["EMRtotFuncP3"].asDouble();
77
_tot_func_p4 = configJSON["EMRtotFuncP4"].asDouble();
78
_acquisition_window = configJSON["EMRacquisitionWindow"].asInt();
79
_signal_integration_window = configJSON["EMRsignalIntegrationWindow"].asInt();
80
_arrival_time_shift = configJSON["EMRarrivalTimeShift"].asInt();
81
_arrival_time_gaus_width = configJSON["EMRarrivalTimeGausWidth"].asDouble();
82
_arrival_time_uniform_width = configJSON["EMRarrivalTimeUniformWidth"].asDouble();
83
_pulse_shape_landau_width = configJSON["EMRpulseShapeLandauWidth"].asDouble();
84
_fom = configJSON["EMRfom"].asString();
85
_birks_constant = configJSON["EMRbirksConstant"].asDouble();
86
_average_path_length = configJSON["EMRaveragePathLength"].asDouble();
87
_signal_energy_threshold = configJSON["EMRsignalEnergyThreshold"].asDouble();
88
_baseline_spread = configJSON["EMRbaselineSpread"].asInt();
89
_maximum_noise_level = configJSON["EMRmaximumNoiseLevel"].asInt();
90
_deltat_shift = configJSON["EMRdeltatShift"].asDouble();
91
_baseline_position = configJSON["EMRbaselinePosition"].asDouble();
92
_arrival_time_spread = configJSON["EMRarrivalTimeSpread"].asDouble();
94
// Generate random noise in each SAPMT
95
_rand = new TRandom3(_seed);
97
for (int iPlane = 0; iPlane < _number_of_planes; iPlane++) {
98
_baseline[iPlane] = _baseline_position
99
+ static_cast<int>(_rand->Uniform(-_baseline_spread, _baseline_spread));
100
_noise_level[iPlane] = static_cast<int>(_rand->Uniform(0, _maximum_noise_level));
101
_noise_position[iPlane] = static_cast<int>(_rand->Uniform(0, 1.9999999));
104
// Load the EMR calibration map
105
bool loaded = _calibMap.InitializeFromCards(configJSON);
107
throw(Exception(Exception::recoverable,
108
"Could not find EMR calibration map",
109
"MapCppEMRMCDigitizer::birth"));
111
// Load the EMR attenuation map
112
loaded = _attenMap.InitializeFromCards(configJSON);
114
throw(Exception(Exception::recoverable,
115
"Could not find EMR attenuation map",
116
"MapCppEMRMCDigitizer::birth"));
119
void MapCppEMRMCDigitization::_death() {
122
void MapCppEMRMCDigitization::_process(Data *data) const {
124
// Get spill, break if there's no DAQ data
125
Spill *spill = data->GetSpill();
127
if (spill->GetMCEvents() == NULL)
130
int nPartEvents = spill->GetMCEventSize();
132
// Check the Recon event array is initialised, and if not make it so
133
if (!spill->GetReconEvents()) {
134
ReconEventPArray* revts = new ReconEventPArray();
135
spill->SetReconEvents(revts);
138
// Set primary hits/noise time windows
139
std::vector<double> delta_t_array;
140
delta_t_array.resize(0);
142
for (int iPe = 0; iPe < nPartEvents; iPe++) {
143
EMRHitArray *hits = spill->GetAnMCEvent(iPe).GetEMRHits();
146
int nHits = hits->size();
148
Primary *primary = spill->GetAnMCEvent(iPe).GetPrimary();
149
double pTime = primary->GetTime(); // ns
151
for (int iHit = 0; iHit < nHits; iHit++) {
152
EMRHit hit = hits->at(iHit);
154
int time = static_cast<int>(hit.GetTime()); // ns
155
double delta_t = time - pTime;
157
delta_t_array.push_back(delta_t);
162
if (!delta_t_array.size()) return;
164
int lt = static_cast<int>(*std::min_element(delta_t_array.begin(), delta_t_array.end())
166
int deltat_limits[4] = {lt, lt+20, lt+20, lt+40};
168
// Create the DBB and fADC arrays with n events (1 per trigger)
169
EMRDBBEventVector emr_dbb_events_tmp = get_dbb_data_tmp(nPartEvents);
170
EMRfADCEventVector emr_fadc_events_tmp = get_fadc_data_tmp(nPartEvents);
172
// Fill the Recon event array with G4 Spill information
173
for (int iPe = 0; iPe < nPartEvents; iPe++) {
174
EMRHitArray *EMRhits = spill->GetAnMCEvent(iPe).GetEMRHits();
177
Primary *primary = spill->GetAnMCEvent(iPe).GetPrimary();
178
double pTime = primary->GetTime(); // ns
180
processDBB(EMRhits, iPe, pTime, emr_dbb_events_tmp, emr_fadc_events_tmp);
181
processFADC(EMRhits, iPe, emr_fadc_events_tmp);
185
fill(spill, nPartEvents, emr_dbb_events_tmp, emr_fadc_events_tmp);
187
// Reset/Resize DBB and fADC arrays with n+2 events (1 per trigger + noise + decays)
188
emr_dbb_events_tmp = get_dbb_data_tmp(nPartEvents+2);
189
emr_fadc_events_tmp = get_fadc_data_tmp(nPartEvents+2);
191
// Digitize the Recon event array
192
for (int iPe = 0; iPe < nPartEvents; iPe++) {
193
EMREvent *EMRevt = spill->GetReconEvents()->at(iPe)->GetEMREvent();
196
digitize(EMRevt, iPe, nPartEvents, deltat_limits, emr_dbb_events_tmp, emr_fadc_events_tmp);
200
fill(spill, nPartEvents+2, emr_dbb_events_tmp, emr_fadc_events_tmp);
203
void MapCppEMRMCDigitization::processDBB(MAUS::EMRHitArray *EMRhits,
206
EMRDBBEventVector& emr_dbb_events_tmp,
207
EMRfADCEventVector& emr_fadc_events_tmp) const {
209
int nHits = EMRhits->size();
211
for (int iHit = 0; iHit < nHits; iHit++) {
213
EMRHit hit = EMRhits->at(iHit);
216
EMRChannelId *id = hit.GetChannelId();
218
int g4barid = id->GetBar();
219
int planeid = g4barid / 59;
220
int barid = g4barid % 59 + 1;
223
double edep = hit.GetEnergyDeposited(); /*MeV*/
226
ThreeVector pos = hit.GetPosition();
227
double pozx = pos.x(); /*mm*/
228
double pozy = pos.y(); /*mm*/
229
double pozz = pos.z(); /*mm*/
232
double time = hit.GetTime(); /*ns*/
234
int lt = static_cast<int>(time);/*ns*/
235
int tot = static_cast<int>(edep*1000000);/*eV*/
237
int xPlane = planeid;
239
int delta_t = lt-pTime;
248
bHit.SetDeltaT(delta_t);
249
emr_dbb_events_tmp[xPe][xPlane][xBar].push_back(bHit);
253
void MapCppEMRMCDigitization::processFADC(MAUS::EMRHitArray *EMRhits,
255
EMRfADCEventVector& emr_fadc_events_tmp) const {
257
int nHits = EMRhits->size();
259
for (int iHit = 0; iHit < nHits; iHit++) {
261
EMRHit hit = EMRhits->at(iHit);
264
EMRChannelId *id = hit.GetChannelId();
266
int g4barid = id->GetBar();
267
int planeid = g4barid / 59;
270
double edep = hit.GetEnergyDeposited(); /*MeV*/
272
double xCharge = edep*1000000;/*eV*/
273
int xPlane = planeid;
275
xCharge += emr_fadc_events_tmp[xPe][xPlane]._charge;
277
data._orientation = planeid % 2;
278
data._charge = xCharge;
280
emr_fadc_events_tmp[xPe][xPlane] = data;
284
void MapCppEMRMCDigitization::digitize(MAUS::EMREvent *EMRevt,
288
EMRDBBEventVector& emr_dbb_events_tmp,
289
EMRfADCEventVector& emr_fadc_events_tmp) const {
291
int nPlHits = EMRevt->GetEMRPlaneHitArray().size();
293
for (int iPlane = 0; iPlane < nPlHits; iPlane++) {
294
EMRPlaneHit *plHit = EMRevt->GetEMRPlaneHitArray().at(iPlane);
295
int xPlane = plHit->GetPlane();
297
// Find PMT gain correction coefficient
298
EMRChannelKey xAverageKey(xPlane, xPlane%2, -1, "emr");
299
double epsilon_MA = _calibMap.Eps(xAverageKey, "MA");
300
double epsilon_SA = _calibMap.Eps(xAverageKey, "SA");
302
int nBars = plHit->GetEMRBarArray().size();
305
for (int iBar = 0; iBar < nBars; iBar++) {
306
EMRBar *bar = plHit->GetEMRBarArray().at(iBar);
307
int xBar = bar->GetBar();
308
int nBarHits = bar->GetEMRBarHitArray().size();
309
bool matched = false;
311
for (int iBarHit = 0; iBarHit < nBarHits; iBarHit++) {
312
EMRBarHit barHit = bar->GetEMRBarHitArray().at(iBarHit);
313
EMRChannelKey xKey(xPlane, xPlane%2, xBar, "emr");
314
double x = barHit.GetX();/*mm*/
315
double y = barHit.GetY();/*mm*/
316
int xTot = barHit.GetTot();/*eV*/
317
int xDeltaT = barHit.GetDeltaT();/*ns*/
318
int xHitTime = barHit.GetHitTime();/*ns*/
320
/*--------- MAPMT signal simulation -------*/
322
// simulate electronics noise and PMT dark current
325
// set g4 information
326
double energy = xTot/1000000.0;/*MeV*/
327
int time = xDeltaT;/*ns*/
329
if (energy < _signal_energy_threshold) continue;
331
// convert energy to the number of scintillation photons (nsph) according to Birks's law
332
// !!! TO DO !!! calculate path length
333
int nsph = static_cast<int>(static_cast<double>(_nph_per_MeV*energy) // MeV
334
/(1.0+_birks_constant // mm/MeV
335
/_average_path_length*energy)); // MeV
337
// sample nsph with Poisson
338
if (_do_sampling) nsph = _rand->Poisson(nsph);
340
// covert nsph to the number of trapped photons (ntph)
341
int ntph = static_cast<int>(static_cast<double>(nsph*_trap_eff));
343
// sample ntph with Poisson
344
if (_do_sampling) ntph = _rand->Poisson(ntph);
346
// reduce ntph according to fibre length (naph)
348
int naph_MAPMT = static_cast<int>(static_cast<double>(ntph)
349
/2.0*_attenMap.fibreAtten(xKey, x, y, "MA"));
352
int naph_SAPMT_hit = static_cast<int>(static_cast<double>(ntph)
353
/2.0*_attenMap.fibreAtten(xKey, x, y, "SA"));
355
// apply channel attenuation map
357
naph_MAPMT = static_cast<int>(static_cast<double>(naph_MAPMT)
358
*_attenMap.connectorAtten(xKey, "MA"));
361
naph_SAPMT_hit = static_cast<int>(static_cast<double>(naph_SAPMT_hit)
362
*_attenMap.connectorAtten(xKey, "SA"));
364
if (time > deltat_limits[0] && time < deltat_limits[1]) naph_SAPMT += naph_SAPMT_hit;
366
// sample naph with Poisson
367
if (_do_sampling) naph_MAPMT = _rand->Poisson(naph_MAPMT);
369
// simulate cross-talk and misalignment
372
// convert naph to the number of photoelectrons (npe)
373
int npe = static_cast<int>(static_cast<double>(naph_MAPMT)*_QE_MAPMT);
375
// sample npe with Poisson
376
if (_do_sampling) npe = _rand->Poisson(npe);
378
// correct npe for the photocathode non-uniformity
381
// correct npe for pmt gain difference
382
npe = static_cast<int>(static_cast<double>(npe)*epsilon_MA);
384
if (npe == 0) continue;
386
// convert npe to the number of ADC counts
387
int nADC = static_cast<int>(static_cast<double>(npe)*_nADC_per_pe_MAPMT);
389
if (nADC == 0) continue;
391
// simulate electronics response
392
if (_do_sampling) nADC = static_cast<int>(_rand->Gaus(nADC, _electronics_response_spread_MAPMT));
394
// convert nADC to Tot
395
int xTotDigi = static_cast<int>(_tot_func_p1+_tot_func_p2
396
*log(static_cast<double>(nADC)
397
/_tot_func_p4+_tot_func_p3));
399
if (xTotDigi <= 0) continue;
401
// convert g4 time to deltaT time in ADC counts
402
int xDeltaTDigi = static_cast<int>(static_cast<double>(time)/_dbb_count);
405
if (_do_sampling) xDeltaTDigi = static_cast<int>(_rand->Gaus(xDeltaTDigi,
406
_time_response_spread));
407
if (xDeltaTDigi < 0) xDeltaTDigi = 0;
409
// correct deltaT with fibre length
410
xDeltaTDigi += static_cast<int>(static_cast<double>(_attenMap.fibreDelay(xKey, x, y, "MA"))
414
int xHitTimeDigi = static_cast<int>(static_cast<double>(xHitTime)/_dbb_count);
418
bHit.SetTot(xTotDigi);
419
bHit.SetDeltaT(xDeltaTDigi);
420
bHit.SetHitTime(xHitTimeDigi);
422
// discriminate noise and decays from signal events
423
if (xDeltaTDigi >= deltat_limits[0] && xDeltaTDigi < deltat_limits[1]) {
424
emr_dbb_events_tmp[xPe][xPlane][xBar].push_back(bHit);
426
} else if (xDeltaTDigi >= deltat_limits[2] && xDeltaTDigi < deltat_limits[3] &&
427
xTotDigi > _tot_noise_low && xTotDigi < _tot_noise_up ) {
428
emr_dbb_events_tmp[nPartEvents][xPlane][xBar].push_back(bHit);
430
} else if (iBarHit == nBarHits-1 && !matched) {
432
emr_dbb_events_tmp[nPartEvents+1][xPlane][xBar].push_back(bHit);
437
int xOri = plHit->GetOrientation();
439
/*------------ SAPMT signal simulation -------*/
441
// simulate PMT dark current
444
// sample naph with Poisson
445
if (_do_sampling) naph_SAPMT = _rand->Poisson(naph_SAPMT);
447
// convert naph to the number of photoelectrons (npe)
448
int npe = static_cast<int>(static_cast<double>(naph_SAPMT)*_QE_SAPMT);
450
// sample npe with Poisson
451
if (_do_sampling) npe = _rand->Poisson(npe);
453
// correct npe for pmt gain difference
454
npe = static_cast<int>(static_cast<double>(npe)*epsilon_SA);
456
// convert npe to the number of ADC counts
457
int nADC = static_cast<int>(static_cast<double>(npe)*_nADC_per_pe_SAPMT);
459
// simulate electronics response
460
if (_do_sampling) nADC = static_cast<int>(_rand->Gaus(nADC,
461
_electronics_response_spread_SAPMT));
463
// correct non physical values
464
if (nADC < 0) nADC = 0;
466
// set signal baseline (8bit ADC)
467
int baseline = _baseline[xPlane];
469
// set electronics noise level - number of fluctuations within acquisition window
470
int noise_level = _noise_level[xPlane];
472
// set electronics noise position (upwards/downwards fluctuations)
473
int noise_position = _noise_position[xPlane];
475
// simulate negative voltage pulse with random electronics noise
476
TH1F *pulse_shape = new TH1F("pulse_shape", "pulse_shape",
477
_acquisition_window, 0, _acquisition_window);
479
double arrival_time = _signal_integration_window*2 + _arrival_time_shift;
481
if (static_cast<int>(_rand->Uniform(0, 1.9999999))) arrival_time += _arrival_time_spread;
483
arrival_time += _rand->Gaus(0, _arrival_time_gaus_width)
484
+_rand->Uniform(-_arrival_time_uniform_width, _arrival_time_uniform_width);
486
for (int i = 0; i < noise_level; i++) {
487
double noise = _rand->Uniform(0, _acquisition_window);
488
int bin = static_cast<Int_t>(noise+1);
489
if (pulse_shape->GetBinContent(bin) != 1) pulse_shape->Fill(noise);
492
if (noise_position) for (Int_t i = 1; i <= _acquisition_window; i++) {
493
pulse_shape->SetBinContent(i, -pulse_shape->GetBinContent(i));}
495
for (int i = 0; i < nADC; i++) {
496
pulse_shape->Fill(_rand->Landau(arrival_time, _pulse_shape_landau_width));}
498
for (int i = 1; i <= _acquisition_window; i++) {
499
pulse_shape->SetBinContent(i, baseline-pulse_shape->GetBinContent(i));}
501
double pedestal_baseline = pulse_shape->Integral(1, _signal_integration_window)
502
/static_cast<double>(_signal_integration_window);
503
double pedestal_area = pedestal_baseline*static_cast<double>(_signal_integration_window)
504
- pulse_shape->Integral(_signal_integration_window+1,
505
_signal_integration_window*2);
506
double signal_area = pedestal_baseline*static_cast<double>(_signal_integration_window)
507
- pulse_shape->Integral(arrival_time-10+1,
508
arrival_time+_signal_integration_window-10);
511
double xAreaDigi = signal_area;
513
// correct non physical values
514
if (nADC == 0) xAreaDigi = 0;
515
if (xAreaDigi < 0) xAreaDigi = 0;
518
double xPedestalAreaDigi = pedestal_area;
520
// set hit arrival time
521
int xArrivalTimeDigi = arrival_time;
523
// simulate time delay in cables
527
std::vector<int> xSamplesDigi;
528
for (Int_t i = 1; i <= _acquisition_window; i++) {
529
xSamplesDigi.push_back(pulse_shape->GetBinContent(i));}
532
data._orientation = xOri;
533
data._charge = xAreaDigi;/*ADC*/
534
data._pedestal_area = xPedestalAreaDigi;/*ADC*/
535
data._time = xArrivalTimeDigi;/*ADC*/
536
data._samples = xSamplesDigi;/*#*/
537
emr_fadc_events_tmp[xPe][xPlane] = data;
543
void MapCppEMRMCDigitization::fill(MAUS::Spill *spill,
545
EMRDBBEventVector emr_dbb_events_tmp,
546
EMRfADCEventVector emr_fadc_events_tmp) const {
548
int recPartEvents = spill->GetReconEventSize();
549
ReconEventPArray *recEvts = spill->GetReconEvents();
550
int xSpill = spill->GetSpillNumber();
552
// Resize the recon event to harbour all the EMR noise+decays
553
if (recPartEvents == 0) { // No recEvts yet
554
for (int iPe = 0; iPe < nPartEvents; iPe++) {
555
recEvts->push_back(new ReconEvent);
559
if (nPartEvents == recPartEvents+2) { // Regular sized recEvts already created
560
for (int iPe = 0; iPe < 2; iPe++) {
561
recEvts->push_back(new ReconEvent);
565
// Fill it with DBB and fADC arrays
566
for (int iPe = 0; iPe < nPartEvents; iPe++) {
567
EMREvent *evt = new EMREvent;
568
EMRPlaneHitArray plArray;
570
for (int iPlane = 0; iPlane < _number_of_planes; iPlane++) {
571
int xOri = emr_fadc_events_tmp[iPe][iPlane]._orientation;
572
double xCharge = emr_fadc_events_tmp[iPe][iPlane]._charge;
573
int xPos = emr_fadc_events_tmp[iPe][iPlane]._time;
574
double xPedestalArea = emr_fadc_events_tmp[iPe][iPlane]._pedestal_area;
575
std::vector<int> xSamples = emr_fadc_events_tmp[iPe][iPlane]._samples;
577
EMRPlaneHit *plHit = new EMRPlaneHit;
578
plHit->SetPlane(iPlane);
579
plHit->SetTrigger(iPe);
581
plHit->SetOrientation(xOri);
582
plHit->SetCharge(xCharge);
583
plHit->SetDeltaT(xPos);
584
plHit->SetSpill(xSpill);
585
plHit->SetPedestalArea(xPedestalArea);
586
plHit->SetSamples(xSamples);
588
EMRBarArray barArray;
589
for (int iBar = 0; iBar < _number_of_bars; iBar++) {
590
int nHits = emr_dbb_events_tmp[iPe][iPlane][iBar].size();
592
EMRBar *bar = new EMRBar;
594
bar->SetEMRBarHitArray(emr_dbb_events_tmp[iPe][iPlane][iBar]);
595
barArray.push_back(bar);
599
plHit->SetEMRBarArray(barArray);
600
if (barArray.size() || xCharge) {
601
plArray.push_back(plHit);
605
if (iPe < 1) {evt->SetInitialTrigger(true);
606
} else {evt->SetInitialTrigger(false);}
608
evt->SetEMRPlaneHitArray(plArray);
609
recEvts->at(iPe)->SetEMREvent(evt);
610
recEvts->at(iPe)->SetPartEventNumber(iPe);
613
spill->SetReconEvents(recEvts);
616
EMRDBBEventVector MapCppEMRMCDigitization::get_dbb_data_tmp(int nPartEvts) const {
617
EMRDBBEventVector emr_dbb_events_tmp;
618
emr_dbb_events_tmp.resize(nPartEvts);
619
for (int iPe = 0; iPe < nPartEvts; iPe++) {
620
emr_dbb_events_tmp[iPe].resize(_number_of_planes); // number of planes
621
for (int iPlane = 0; iPlane < _number_of_planes; iPlane++) {
622
emr_dbb_events_tmp[iPe][iPlane].resize(_number_of_bars); // number of bars in a plane
625
return emr_dbb_events_tmp;
628
EMRfADCEventVector MapCppEMRMCDigitization::get_fadc_data_tmp(int nPartEvts) const {
629
EMRfADCEventVector emr_fadc_events_tmp;
630
emr_fadc_events_tmp.resize(nPartEvts);
631
for (int iPe = 0; iPe < nPartEvts ;iPe++) {
632
emr_fadc_events_tmp[iPe].resize(_number_of_planes);
633
for (int iPlane = 0; iPlane < _number_of_planes; iPlane++) {
635
data._orientation = iPlane%2;
637
data._pedestal_area = 0.0;
639
std::vector<int> xSamples;
640
data._samples = xSamples;
641
emr_fadc_events_tmp[iPe][iPlane] = data;
644
return emr_fadc_events_tmp;