CbmRoot
Loading...
Searching...
No Matches
CbmTrdModuleSimR.cxx
Go to the documentation of this file.
1/* Copyright (C) 2018-2020 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
2 SPDX-License-Identifier: GPL-3.0-only
3 Authors: Etienne Bechtel, Florian Uhlig [committer] */
4
5// Includes from TRD
6#include "CbmTrdModuleSimR.h"
7
8#include "CbmTrdAddress.h"
9#include "CbmTrdDigi.h"
10#include "CbmTrdDigitizer.h"
11#include "CbmTrdParModAsic.h"
12#include "CbmTrdParModDigi.h"
13#include "CbmTrdParSpadic.h"
14#include "CbmTrdPoint.h"
15#include "CbmTrdRadiator.h"
16
17// Includes from CbmRoot
18//#include "CbmDaqBuffer.h"
19#include "CbmMatch.h"
20
21// Includes from FairRoot
22#include <Logger.h>
23
24// Includes from Root
25#include <TClonesArray.h>
26#include <TGeoManager.h>
27#include <TMath.h>
28#include <TRandom3.h>
29#include <TStopwatch.h>
30#include <TVector3.h>
31
32// Includes from C++
33#include "CbmTrdRawToDigiR.h"
34
35#include <fstream>
36#include <iomanip>
37#include <iostream>
38
39#include <cmath>
40
41
42//_________________________________________________________________________________
43CbmTrdModuleSimR::CbmTrdModuleSimR(Int_t mod, Int_t ly, Int_t rot)
44 : CbmTrdModuleSim(mod, ly, rot)
45 , fSigma_noise_keV(0.2)
46 , fMinimumChargeTH(.5e-06)
47 , fCurrentTime(-1.)
48 , fAddress(-1.)
49 , fLastEventTime(-1)
50 , fEventTime(-1)
51 , fLastTime(-1)
52 , fCollectTime(2048)
53 , //in pulsed mode
54 fnClusterConst(0)
55 , fnScanRowConst(0)
56 ,
57
58 fPulseSwitch(true)
59 , fPrintPulse(false)
60 , //only for debugging
61 fTimeShift(true)
62 , //for pulsed mode
63 fAddCrosstalk(true)
64 , //for pulsed mode
65 fClipping(true)
66 , //for pulsed mode
67 fepoints(0)
68 , fAdcNoise(2)
69 , fDistributionMode(4)
70 , fCrosstalkLevel(0.01)
71 ,
72
73 nofElectrons(0)
74 , nofLatticeHits(0)
75 , nofPointsAboveThreshold(0)
76 ,
77
78 fAnalogBuffer()
79 , fPulseBuffer()
80 , fMultiBuffer()
81 , fTimeBuffer()
82 , fShiftQA()
83 , fLinkQA()
84 , fMCQA()
85 , fMCBuffer()
86
87{
88 FairRootManager* ioman = FairRootManager::Instance();
89 fTimeSlice = (CbmTimeSlice*) ioman->GetObject("TimeSlice.");
90
91 if (!fPulseSwitch && CbmTrdDigitizer::IsTimeBased()) fCollectTime = 200;
92 SetNameTitle(Form("TrdSimR%d", mod), "Simulator for rectangular read-out.");
93
95 TFile* oldFile = gFile;
96 TDirectory* oldDir = gDirectory;
97
98 TString dir = getenv("VMCWORKDIR");
99 TString filename = dir + "/parameters/trd/FeatureExtractionLookup.root";
100 TFile* f = new TFile(filename, "OPEN");
101 LOG_IF(fatal, !f->IsOpen()) << "parameter file " << filename << " does not exist!";
102 fDriftTime = f->Get<TH2D>("Drift");
103 LOG_IF(fatal, !fDriftTime) << "No histogram Drift founfd in file " << filename;
104 fDriftTime->SetDirectory(0);
105 f->Close();
106
108 gFile = oldFile;
109 gDirectory = oldDir;
110
111 oldFile = nullptr;
112 delete oldFile;
113
114 if (fPulseSwitch) {
115 SetSpadicResponse(fCalibration, fTau);
116 SetPulsePars(fRecoMode);
117 }
118}
119
120//_______________________________________________________________________________
121void CbmTrdModuleSimR::AddDigi(Int_t address, Double_t charge, Double_t /*chargeTR*/, Double_t time, Int_t trigger)
122{
123 Double_t weighting = charge;
124 if (CbmTrdDigitizer::UseWeightedDist()) {
125 TVector3 padPos, padPosErr;
126 fDigiPar->GetPadPosition(address, padPos, padPosErr);
127 Double_t distance = sqrt(pow(fXYZ[0] - padPos[0], 2) + pow(fXYZ[1] - padPos[1], 2));
128 weighting = 1. / distance;
129 }
130
131 // std::map<Int_t, std::pair<CbmTrdDigi*, CbmMatch*> >::iterator it = fDigiMap.find(address);
132 // std::map<Int_t, std::pair<CbmTrdDigi*, CbmMatch*> >::iterator it;
133
134 Int_t col = CbmTrdAddress::GetColumnId(address);
135 Int_t row = CbmTrdAddress::GetRowId(address);
136 Int_t sec = CbmTrdAddress::GetSectorId(address);
137 Int_t ncols = fDigiPar->GetNofColumns();
138 Int_t channel(0);
139 for (Int_t isec(0); isec < sec; isec++)
140 channel += ncols * fDigiPar->GetNofRowsInSector(isec);
141 channel += ncols * row + col;
142
143 std::map<Int_t, std::pair<CbmTrdDigi*, CbmMatch*>>::iterator it = fDigiMap.find(address);
144 if (it == fDigiMap.end()) { // Pixel not yet in map -> Add new pixel
145 CbmMatch* digiMatch = new CbmMatch();
146 digiMatch->AddLink(CbmLink(weighting, fPointId, fEventId, fInputId));
147 AddNoise(charge);
148
149 CbmTrdDigi::eTriggerType triggertype = CbmTrdDigi::eTriggerType::kNTrg;
150 if (trigger == 1) triggertype = CbmTrdDigi::eTriggerType::kSelf;
151 if (trigger == 2) triggertype = CbmTrdDigi::eTriggerType::kNeighbor;
152
153 CbmTrdDigi* digi = new CbmTrdDigi(channel, fModAddress, charge * 1e6, ULong64_t(time), triggertype, 0);
154
155 digi->SetFlag(0, true);
156 if (fDigiPar->GetPadSizeY(sec) == 1.5) digi->SetErrorClass(1);
157 if (fDigiPar->GetPadSizeY(sec) == 4.) digi->SetErrorClass(2);
158 if (fDigiPar->GetPadSizeY(sec) == 6.75) digi->SetErrorClass(3);
159 if (fDigiPar->GetPadSizeY(sec) == 12.) digi->SetErrorClass(4);
160
161 fDigiMap[address] = std::make_pair(digi, digiMatch);
162
163 it = fDigiMap.find(address);
164 // it->second.first->SetAddressModule(fModAddress);//module); <- now handled in the digi contructor
165 }
166 else {
167 it->second.first->AddCharge(charge * 1e6);
168 it->second.second->AddLink(CbmLink(weighting, fPointId, fEventId, fInputId));
169 }
170}
171
172//_______________________________________________________________________________
173void CbmTrdModuleSimR::ProcessBuffer(Int_t address)
174{
175
176 std::vector<std::pair<CbmTrdDigi*, CbmMatch*>> analog = fAnalogBuffer[address];
177 std::vector<std::pair<CbmTrdDigi*, CbmMatch*>>::iterator it;
178 CbmTrdDigi* digi = analog.begin()->first;
179 CbmMatch* digiMatch = new CbmMatch();
180 //printf("CbmTrdModuleSimR::ProcessBuffer(%10d)=%3d\n", address, digi->GetAddressChannel());
181
182 // Int_t trigger = 0;
183 Float_t digicharge = 0;
184 // Float_t digiTRcharge=0;
185 for (it = analog.begin(); it != analog.end(); it++) {
186 digicharge += it->first->GetCharge();
187 digiMatch->AddLink(it->second->GetLink(0));
188 //printf(" add charge[%f] trigger[%d]\n", it->first->GetCharge(), it->first->GetTriggerType());
189 }
190 digi->SetCharge(digicharge);
191 digi->SetTriggerType(fAnalogBuffer[address][0].first->GetTriggerType());
192
193 // std::cout<<digicharge<<std::endl;
194
195 // if(analog.size()>1) for (it=analog.begin()+1 ; it != analog.end(); it++) if(it->first) delete it->first;
196
197 fDigiMap[address] = std::make_pair(digi, digiMatch);
198
199 fAnalogBuffer.erase(address);
200}
201
202//_______________________________________________________________________________
203void CbmTrdModuleSimR::ProcessPulseBuffer(Int_t address, Bool_t FNcall, Bool_t MultiCall = false, Bool_t down = true,
204 Bool_t up = true)
205{
206
207 std::map<Int_t, std::pair<std::vector<Double_t>, CbmMatch*>>::iterator iBuff = fPulseBuffer.find(address);
208 std::map<Int_t, Double_t>::iterator tBuff = fTimeBuffer.find(address);
209
210 if (iBuff == fPulseBuffer.end() || tBuff == fTimeBuffer.end()) return;
211 Int_t trigger = CheckTrigger(fPulseBuffer[address].first);
212 if (fPulseBuffer[address].first.size() < 32) { return; }
213
214 if (trigger == 0 && !FNcall) { return; }
215 if (trigger == 1 && FNcall) { FNcall = false; }
216
217 Int_t col = CbmTrdAddress::GetColumnId(address);
218 Int_t row = CbmTrdAddress::GetRowId(address);
219 Int_t sec = CbmTrdAddress::GetSectorId(address);
220 Int_t ncols = fDigiPar->GetNofColumns();
221 Int_t channel(0);
222 for (Int_t isec(0); isec < sec; isec++)
223 channel += ncols * fDigiPar->GetNofRowsInSector(isec);
224 channel += ncols * row + col;
225
226
227 CbmMatch* digiMatch = new CbmMatch();
228 std::vector<Int_t> temp;
229 for (size_t i = 0; i < fPulseBuffer[address].first.size(); i++) {
230 Int_t noise = AddNoiseADC();
231 Int_t cross = AddCrosstalk(address, i, sec, row, col, ncols);
232
233 fPulseBuffer[address].first[i] = fPulseBuffer[address].first[i] + noise + cross;
234 if (fPulseBuffer[address].first[i] <= fClipLevel && fPulseBuffer[address].first[i] >= -12.)
235 temp.push_back((Int_t) fPulseBuffer[address].first[i]);
236 if (fPulseBuffer[address].first[i] > fClipLevel) temp.push_back(fClipLevel - 1);
237 if (fPulseBuffer[address].first[i] < -12.) temp.push_back(-12.);
238 if (fDebug) fQA->Fill("Pulse", i, fPulseBuffer[address].first[i]);
239 if (fDebug) fQA->CreateHist("Pulse self", 32, -0.5, 31.5, 512, -12., 500.);
240 if (fDebug) fQA->CreateHist("Pulse FN", 32, -0.5, 31.5, 512, -12., 500.);
241 if (fDebug && trigger == 1) fQA->Fill("Pulse self", i, fPulseBuffer[address].first[i]);
242 if (fDebug && trigger == 0) fQA->Fill("Pulse FN", i, fPulseBuffer[address].first[i]);
243 }
244
245 Float_t shift = fMessageConverter->GetTimeShift(temp);
246 Float_t corr = fMinDrift; //correction of average sampling to clock difference and systematic average drift time
247 if (!CbmTrdDigitizer::IsTimeBased()) corr = CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC);
248 //correction for EB case is done later, due to the event time 0 and the unsigned data member for the time in the digi
249
250 if (fTimeSlice) {
251 if (fTimeBuffer[address] + corr - shift < fTimeSlice->GetStartTime()) {
252 delete digiMatch;
253 delete fPulseBuffer[address].second;
254 fPulseBuffer.erase(address);
255 fTimeBuffer.erase(address);
256 return;
257 }
258 }
259
260 CbmTrdDigi* digi = nullptr;
261 if (fTimeBuffer[address] + corr - shift > 0.)
262 digi = fMessageConverter->MakeDigi(temp, channel, fModAddress, fTimeBuffer[address] + corr - shift, true);
263 else
264 digi = fMessageConverter->MakeDigi(temp, channel, fModAddress, fTimeBuffer[address] + corr, true);
265
266 if (fDigiPar->GetPadSizeY(sec) == 1.5) digi->SetErrorClass(1);
267 if (fDigiPar->GetPadSizeY(sec) == 4.) digi->SetErrorClass(2);
268 if (fDigiPar->GetPadSizeY(sec) == 6.75) digi->SetErrorClass(3);
269 if (fDigiPar->GetPadSizeY(sec) == 12.) digi->SetErrorClass(4);
270 if (!CbmTrdDigitizer::IsTimeBased()) digi->SetFlag(1, true);
271
272 Int_t shiftcut = fShiftQA[address] * 10;
273 Float_t timeshift = shiftcut / 10.0;
274 if (temp[fMinBin] == fClipLevel - 1 && temp[fMaxBin] == fClipLevel - 1) digi->SetCharge(35.);
275
276 std::vector<CbmLink> links = fPulseBuffer[address].second->GetLinks();
277 Int_t Links = 0.;
278 for (UInt_t iLinks = 0; iLinks < links.size(); iLinks++) {
279 digiMatch->AddLink(links[iLinks]);
280 Links++;
281 }
282
283 if (fDebug) {
284 fQA->CreateHist("MC", 200, 0., 50.);
285 fQA->Fill("MC", fMCQA[address]);
286 fQA->CreateHist("rec shift", 63, -0.5, 62.5);
287 fQA->Fill("rec shift", shift);
288 fQA->CreateHist("T res", 70, -35, 35);
289 fQA->Fill("T res", shift - timeshift);
290 fQA->CreateHist("Time", 100, -50, 50);
291 fQA->Fill("Time", digi->GetTime() - fLinkQA[address][0]["Time"]);
292
293
294 fQA->CreateProfile("MAX ADC", 63, -0.5, 62.5, 512, -12., 500.);
295 fQA->FillProfile("MAX ADC", timeshift, temp[fMaxBin], fMCQA[address]);
296 fQA->CreateProfile("ASYM MAP", 512, -12., 500., 512, -12., 500.);
297 fQA->FillProfile("ASYM MAP", temp[fMinBin], temp[fMaxBin], timeshift);
298
299
300 if (trigger == 1 && MultiCall && digi->GetCharge() > 0.) fQA->Fill("Multi Quote", 1);
301 else
302 fQA->Fill("Multi Quote", 0);
303
304 if (trigger == 1 && !MultiCall) {
305 fQA->CreateHist("E self", 200, 0., 50.0);
306 fQA->Fill("E self", digi->GetCharge());
307 fQA->CreateHist("E res", 400., -1., 1.);
308 fQA->Fill("E res", digi->GetCharge() - fMCQA[address]);
309 fQA->CreateHist("E res rel", 400., -1., 1.);
310 fQA->Fill("E res rel", (digi->GetCharge() - fMCQA[address]) / fMCQA[address]);
311 fQA->CreateHist("Charge Max", 200, 0., 50.0, 512, -12., 500.);
312 fQA->Fill("Charge Max", fMessageConverter->GetCharge(temp, timeshift), temp[fMaxBin]);
313 fQA->CreateHist("Links Res Self", 400., -1., 1., 5, -0.5, 4.5);
314 fQA->Fill("Links Res Self", (digi->GetCharge() - fMCQA[address]) / fMCQA[address], digiMatch->GetNofLinks());
315
316 if (Links == 1) {
317 fQA->CreateProfile("1 Link Prof", 32, 0., 32.);
318 for (size_t i = 0; i < fPulseBuffer[address].first.size(); i++) {
319 fQA->FillProfile("1 Link Prof", i, temp[i]);
320 }
321 }
322
323 if (fLinkQA[address].size() == 2 && fLinkQA[address][1]["Event"] != fLinkQA[address][0]["Event"]) {
324 fQA->CreateHist("Links QA time", 400., -1., 1., 500, 0., 2000.);
325 fQA->Fill("Links QA time", (digi->GetCharge() - fMCQA[address]) / fMCQA[address],
326 fLinkQA[address][1]["Time"] - fLinkQA[address][0]["Time"]);
327 fQA->CreateProfile("2 Link Prof", 32, 0., 32., 100, 0., 2000.);
328 for (size_t i = 0; i < fPulseBuffer[address].first.size(); i++) {
329 fQA->FillProfile("2 Link Prof", i, fLinkQA[address][1]["Time"] - fLinkQA[address][0]["Time"], temp[i]);
330 }
331 }
332
333 if (Links == 2 && links[0].GetEntry() == links[1].GetEntry()) {
334 fQA->CreateHist("In Event Res", 400., -1., 1.);
335 fQA->Fill("In Event Res", (digi->GetCharge() - fMCQA[address]) / fMCQA[address]);
336 }
337 if (Links == 2 && links[0].GetEntry() != links[1].GetEntry()) {
338 fQA->CreateHist("Inter Event Res", 400., -1., 1.);
339 fQA->Fill("Inter Event Res", (digi->GetCharge() - fMCQA[address]) / fMCQA[address]);
340 }
341 }
342 else {
343 fQA->CreateHist("E FN", 200, 0., 10.0);
344 fQA->Fill("E FN", digi->GetCharge());
345 fQA->CreateHist("E FN max", 512, -12., 500.);
346 fQA->Fill("E FN max", temp[fMaxBin]);
347 fQA->CreateHist("E FN MC max", 512, -12., 500., 200, 0., 50.);
348 fQA->Fill("E FN MC max", temp[fMaxBin], fMCQA[address]);
349 fQA->CreateHist("E FN res", 400., -1., 1.);
350 fQA->Fill("E FN res", digi->GetCharge() - fMCQA[address]);
351 fQA->CreateHist("E FN rel", 400., -1., 1.);
352 fQA->Fill("E FN rel", (digi->GetCharge() - fMCQA[address]) / fMCQA[address]);
353 fQA->CreateHist("Links Res FN", 400., -1., 1., 5, -0.5, 4.5);
354 fQA->Fill("Links Res FN", digi->GetCharge() - fMCQA[address], digiMatch->GetNofLinks());
355
356 fQA->CreateProfile("1 Link Prof FN", 32, 0., 32.);
357 for (size_t i = 0; i < fPulseBuffer[address].first.size(); i++) {
358 fQA->FillProfile("1 Link Prof FN", i, temp[i]);
359 }
360
361
362 if (fLinkQA[address].size() == 2 && fLinkQA[address][1]["Event"] != fLinkQA[address][0]["Event"]) {
363 fQA->CreateHist("Links QA time FN", 400., -1., 1., 500, 0., 2000.);
364 fQA->Fill("Links QA time FN", (digi->GetCharge() - fMCQA[address]) / fMCQA[address],
365 fLinkQA[address][1]["Time"] - fLinkQA[address][0]["Time"]);
366 fQA->CreateProfile("2 Link Prof FN", 32, 0., 32., 100, 0., 2000.);
367 for (size_t i = 0; i < fPulseBuffer[address].first.size(); i++) {
368 fQA->FillProfile("2 Link Prof FN", i, fLinkQA[address][1]["Time"] - fLinkQA[address][0]["Time"], temp[i]);
369 }
370 }
371 }
372 fShiftQA.erase(address);
373 fMCQA.erase(address);
374 fLinkQA.erase(address);
375 }
376
377
378 // digi->SetAddressModule(fModAddress); Not required anymore, now handled in the digi c'tor
379
380 if (trigger == 1) { digi->SetTriggerType(CbmTrdDigi::eTriggerType::kSelf); }
381 if (trigger == 0 && FNcall) { digi->SetTriggerType(CbmTrdDigi::eTriggerType::kNeighbor); }
382 if (trigger == 1 && MultiCall) { digi->SetTriggerType(CbmTrdDigi::eTriggerType::kMulti); }
383
384 //digi->SetMatch(digiMatch);
385 if (fDebug) {
386 fQA->CreateHist("Links", 10, -0.5, 9.5);
387 fQA->Fill("Links", digiMatch->GetNofLinks());
388 }
389 if (!FNcall && fPrintPulse)
390 std::cout << "main call charge: " << digi->GetCharge() << " col : " << col
391 << " lay : " << CbmTrdAddress::GetLayerId(address) << " trigger: " << trigger
392 << " time: " << digi->GetTime() << std::endl;
393 if (FNcall && fPrintPulse)
394 std::cout << "FN call charge: " << digi->GetCharge() << " col : " << col
395 << " lay : " << CbmTrdAddress::GetLayerId(address) << " trigger: " << trigger
396 << " time: " << digi->GetTime() << std::endl;
397
398 // if(!FNcall && MultiCall) std::cout<<"main call charge: "<<digi->GetCharge()<<" col : " << col<<" lay : " << CbmTrdAddress::GetLayerId(address)<<" trigger: " << trigger<<" time: " << digi->GetTime()<<std::endl;
399 // if(FNcall && MultiCall) std::cout<<"FN call charge: "<<digi->GetCharge()<<" col : " << col<<" lay : " << CbmTrdAddress::GetLayerId(address)<< " trigger: " << trigger<<" time: " << digi->GetTime()<<std::endl;
400
401 fDigiMap[address] = std::make_pair(digi, digiMatch);
402
403 fPulseBuffer.erase(address);
404
405 if (!FNcall && !MultiCall && trigger == 1) {
406 if (down) {
407 Int_t FNaddress = 0;
408 if (col >= 1)
409 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
410 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col - 1);
411 Double_t timediff = TMath::Abs(fTimeBuffer[address] - fTimeBuffer[FNaddress]);
412 if (FNaddress != 0 && timediff < CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC))
413 ProcessPulseBuffer(FNaddress, true, false, true, false);
414 }
415
416 if (up) {
417 Int_t FNaddress = 0;
418 if (col < (ncols - 1))
419 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(address),
420 sec, row, col + 1);
421 Double_t timediff = TMath::Abs(fTimeBuffer[address] - fTimeBuffer[FNaddress]);
422 if (FNaddress != 0 && timediff < CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC))
423 ProcessPulseBuffer(FNaddress, true, false, false, true);
424 }
425 }
426
427
428 if (!FNcall && MultiCall && trigger == 1) {
429 if (down) {
430 Int_t FNaddress = 0;
431 if (col >= 1)
432 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
433 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col - 1);
434 // Double_t timediff = TMath::Abs(fTimeBuffer[address]-fTimeBuffer[FNaddress]);
435 if (FNaddress != 0) ProcessPulseBuffer(FNaddress, true, true, true, false);
436 }
437
438 if (up) {
439 Int_t FNaddress = 0;
440 if (col < (ncols - 1))
441 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(address),
442 sec, row, col + 1);
443 // Double_t timediff = TMath::Abs(fTimeBuffer[address]-fTimeBuffer[FNaddress]);
444 if (FNaddress != 0) ProcessPulseBuffer(FNaddress, true, true, false, true);
445 }
446 }
447
448 fTimeBuffer.erase(address);
449}
450
451
452//_______________________________________________________________________________
453void CbmTrdModuleSimR::AddDigitoBuffer(Int_t address, Double_t charge, Double_t /*chargeTR*/, Double_t time,
454 Int_t trigger)
455{
456 Double_t weighting = charge;
457 if (CbmTrdDigitizer::UseWeightedDist()) {
458 TVector3 padPos, padPosErr;
459 fDigiPar->GetPadPosition(address, padPos, padPosErr);
460 Double_t distance = sqrt(pow(fXYZ[0] - padPos[0], 2) + pow(fXYZ[1] - padPos[1], 2));
461 weighting = 1. / distance;
462 }
463
464 //compare times of the buffer content with the actual time and process the buffer if collecttime is over
465 Bool_t eventtime = false;
466 if (time > 0.000) eventtime = true;
467 if (eventtime) CheckTime(address);
468
469 AddNoise(charge);
470
471 //fill digis in the buffer
472 CbmMatch* digiMatch = new CbmMatch();
473 digiMatch->AddLink(CbmLink(weighting, fPointId, fEventId, fInputId));
474
475 Int_t col = CbmTrdAddress::GetColumnId(address);
476 Int_t row = CbmTrdAddress::GetRowId(address);
477 Int_t sec = CbmTrdAddress::GetSectorId(address);
478 Int_t ncols = fDigiPar->GetNofColumns();
479 Int_t channel(0);
480 for (Int_t isec(0); isec < sec; isec++)
481 channel += ncols * fDigiPar->GetNofRowsInSector(isec);
482 channel += ncols * row + col;
483
484 auto triggertype = CbmTrdDigi::eTriggerType::kNTrg;
485 if (trigger == 1) triggertype = CbmTrdDigi::eTriggerType::kSelf;
486 if (trigger == 2) triggertype = CbmTrdDigi::eTriggerType::kNeighbor;
487 // std::cout<<charge*1e6<<" "<<fTimeBuffer[address]/CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)<<std::endl;
488 CbmTrdDigi* digi =
489 new CbmTrdDigi(channel, fModAddress, charge * 1e6,
490 ULong64_t(time / CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)), triggertype, 0);
491
492
493 //digi->SetMatch(digiMatch);
494 // printf("CbmTrdModuleSimR::AddDigitoBuffer(%10d)=%3d [%d] col[%3d] row[%d] sec[%d]\n", address, channel, fModAddress, col, row, sec);
495
496 fAnalogBuffer[address].push_back(std::make_pair(digi, digiMatch));
497 fTimeBuffer[address] = fCurrentTime;
498 if (!eventtime) ProcessBuffer(address);
499}
500
501//_______________________________________________________________________________
502void CbmTrdModuleSimR::AddDigitoPulseBuffer(Int_t address, Double_t reldrift, Double_t charge, Double_t /*chargeTR*/,
503 Double_t time, Int_t /*trigger*/, Int_t /*epoints*/, Int_t /*ipoint*/)
504{
505
506 Double_t weighting = charge * fepoints;
507 if (CbmTrdDigitizer::UseWeightedDist()) {
508 TVector3 padPos, padPosErr;
509 fDigiPar->GetPadPosition(address, padPos, padPosErr);
510 Double_t distance = sqrt(pow(fXYZ[0] - padPos[0], 2) + pow(fXYZ[1] - padPos[1], 2));
511 weighting = 1. / distance;
512 }
513
514 // if(!CbmTrdDigitizer::IsTimeBased() && fPointId-fLastPoint!=0) {CheckBuffer(true);CleanUp(true);}
515
516 if (CbmTrdDigitizer::IsTimeBased() && reldrift == 0.) {
517 Bool_t gotMulti = CheckMulti(address, fPulseBuffer[address].first);
518 fMultiBuffer.erase(address);
519 if (gotMulti) fMCBuffer[address].clear();
520 }
521
522 if (CbmTrdDigitizer::IsTimeBased() && reldrift == 0.) CheckTime(address);
523
524 fMultiBuffer[address] = std::make_pair(fMultiBuffer[address].first + charge, 0.);
525 if (fMCBuffer[address].size() > 0) {
526 fMCBuffer[address][0].push_back(fPointId);
527 fMCBuffer[address][1].push_back(fEventId);
528 fMCBuffer[address][2].push_back(fInputId);
529 }
530 else
531 for (auto ini = 0; ini < 3; ini++) {
532 std::vector<Int_t> v;
533 fMCBuffer[address].push_back(v);
534 }
535 std::vector<Double_t> pulse;
536 if (fTimeBuffer[address] > 0) {
537 pulse = fPulseBuffer[address].first;
538 if (pulse.size() < 32) return;
539 AddToPulse(address, charge, reldrift, pulse);
540 Bool_t found = false;
541 for (Int_t links = 0; links < fPulseBuffer[address].second->GetNofLinks(); links++)
542 if (fPulseBuffer[address].second->GetLink(links).GetIndex() == fPointId) found = true;
543 if (!found) {
544 fPulseBuffer[address].second->AddLink(CbmLink(weighting, fPointId, fEventId, fInputId));
545 if (fDebug) {
546 std::map<TString, Int_t> LinkQA;
547 LinkQA["Event"] = fEventId;
548 LinkQA["Point"] = fPointId;
549 LinkQA["Time"] = fEventTime;
550 LinkQA["Charge"] = charge * 1e6 * fepoints;
551 LinkQA["PosX"] = fQAPosition[0];
552 LinkQA["PosY"] = fQAPosition[1];
553 fLinkQA[address].push_back(LinkQA);
554 }
555 }
556 }
557 else {
558 fMCBuffer[address].clear();
559 pulse = MakePulse(charge, pulse, address);
560 fPulseBuffer[address].first = pulse;
561 CbmMatch* digiMatch = new CbmMatch();
562 digiMatch->AddLink(CbmLink(weighting, fPointId, fEventId, fInputId));
563
564 if (fDebug) {
565 std::vector<std::map<TString, Int_t>> vecLink;
566 std::map<TString, Int_t> LinkQA;
567 LinkQA["Event"] = fEventId;
568 LinkQA["Point"] = fPointId;
569 LinkQA["Time"] = fEventTime;
570 LinkQA["Charge"] = charge * 1e6 * fepoints;
571 LinkQA["PosX"] = fQAPosition[0];
572 LinkQA["PosY"] = fQAPosition[1];
573 vecLink.push_back(LinkQA);
574 fLinkQA[address] = vecLink;
575 }
576
577 fPulseBuffer[address].second = digiMatch;
578 fTimeBuffer[address] = time;
579 }
580
581 // if(!CbmTrdDigitizer::IsTimeBased() && finish) {CheckBuffer(true);CleanUp(true);}
582}
583
584
585std::vector<Double_t> CbmTrdModuleSimR::MakePulse(Double_t charge, std::vector<Double_t> pulse, Int_t address)
586{
587
588 Double_t sample[32];
589 for (Int_t i = 0; i < 32; i++)
590 sample[i] = 0;
591 // Double_t timeshift = gRandom->Uniform(0.,CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
592 Double_t timeshift =
593 ((Int_t)(fCurrentTime * 10) % (Int_t)(CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) * 10)) / 10.0;
594 if (fDebug) fQA->Fill("Shift", timeshift);
595 // Int_t shift=timeshift;
596 //fShiftQA[address]=shift;
597 fShiftQA[address] = timeshift;
598 for (Int_t i = 0; i < 32; i++) {
599 if (i < fPresamples) {
600 sample[i] = 0.;
601 continue;
602 }
603 if (fTimeShift)
604 sample[i] = fCalibration * charge * 1e6
605 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
606 if (!fTimeShift)
607 sample[i] = fCalibration * charge * 1e6
608 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
609 if (sample[i] > fClipLevel && fClipping) sample[i] = fClipLevel - 1; //clipping
610 }
611
612 for (Int_t i = 0; i < 32; i++) {
613 pulse.push_back(sample[i]);
614 }
615 AddCorrelatedNoise(pulse);
616
617 // if(fDebug && CheckTrigger(pulse) == 1) fMCQA[address]=charge*1e6;
618 if (fDebug) fMCQA[address] = charge * 1e6;
619 return pulse;
620}
621
622void CbmTrdModuleSimR::AddToPulse(Int_t address, Double_t charge, Double_t reldrift, std::vector<Double_t> pulse)
623{
624
625 Int_t comptrigger = CheckTrigger(pulse);
626 std::vector<Double_t> temppulse;
627 for (Int_t i = 0; i < 32; i++)
628 temppulse.push_back(pulse[i]);
629 Double_t dt = fCurrentTime - fTimeBuffer[address];
630 Int_t startbin = (dt + reldrift) / CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC);
631 Double_t timeshift =
632 ((Int_t)((dt + reldrift) * 10) % (Int_t)(CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) * 10)) / 10.0;
633 if (startbin > 31 || dt < 0.) return;
634 if (fDebug) fMCQA[address] += charge * 1e6;
635
636 for (Int_t i = 0; i < 32; i++) {
637 if (i < fPresamples + startbin) {
638 pulse[i] = temppulse[i];
639 continue;
640 }
641 Int_t addtime = i - startbin - fPresamples;
642 pulse[i] += fCalibration * charge * 1e6
643 * CalcResponse(addtime * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
644 if (pulse[i] > fClipLevel && fClipping) pulse[i] = fClipLevel - 1; //clipping
645 }
646
647 // std::vector<Int_t> newpulse;
648 // for(Int_t i=0;i<32;i++){
649 // if(i < fPresamples) continue;
650 // if(fTimeShift){
651 // Int_t sample = fCalibration*charge*1e6*CalcResponse((i-fPresamples)*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)+timeshift);
652 // if(sample > fClipLevel && fClipping) sample=fClipLevel-1; //clipping
653 // newpulse.push_back(sample);
654 // }
655 // if(!fTimeShift){
656 // Int_t sample = fCalibration*charge*1e6*CalcResponse((i-fPresamples)*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
657 // if(sample > fClipLevel && fClipping) sample=fClipLevel-1; //clipping
658 // newpulse.push_back(sample);
659 // }
660 // }
661
662 Int_t trigger = CheckTrigger(pulse);
663 if (comptrigger == 0 && trigger == 1) {
664 for (Int_t i = 0; i < 32; i++) {
665 if (i < fPresamples && startbin >= fPresamples) {
666 pulse[i] = temppulse[i + startbin - fPresamples];
667 continue;
668 }
669 if (i < fPresamples && startbin < fPresamples) {
670 pulse[i] = temppulse[i + startbin];
671 continue;
672 }
673 Int_t addtime = i - fPresamples;
674 Int_t shift = startbin + i;
675 if (fTimeShift) {
676 if (shift < 32)
677 pulse[i] = temppulse[shift]
678 + fCalibration * charge * 1e6
679 * CalcResponse(addtime * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
680 else
681 pulse[i] = fCalibration * charge * 1e6
682 * CalcResponse(addtime * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
683 if (pulse[i] > fClipLevel && fClipping) pulse[i] = fClipLevel - 1; //clipping
684 }
685 if (!fTimeShift) {
686 if (shift < 32)
687 pulse[i] = temppulse[shift]
688 + fCalibration * charge * 1e6
689 * CalcResponse(addtime * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
690 else
691 pulse[i] =
692 fCalibration * charge * 1e6 * CalcResponse(addtime * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
693 if (pulse[i] > fClipLevel && fClipping) pulse[i] = fClipLevel - 1; //clipping
694 }
695 }
696 fTimeBuffer[address] = fCurrentTime;
697 }
698
699 if (trigger == 2) { fMultiBuffer[address].second = fCurrentTime; }
700
701 fPulseBuffer[address].first = pulse;
702}
703
704Bool_t CbmTrdModuleSimR::CheckMulti(Int_t address, std::vector<Double_t> pulse)
705{
706
707 Bool_t processed = false;
708 Int_t trigger = CheckTrigger(pulse);
709 if (trigger == 2) {
710 processed = true;
711 Int_t maincol = CbmTrdAddress::GetColumnId(address);
712 Int_t row = CbmTrdAddress::GetRowId(address);
713 Int_t sec = CbmTrdAddress::GetSectorId(address);
714 Int_t shift = GetMultiBin(pulse);
715 Int_t ncols = fDigiPar->GetNofColumns();
716 // Double_t timeshift = gRandom->Uniform(0.,CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
717 Double_t timeshift =
718 ((Int_t)(fMultiBuffer[address].second * 10) % (Int_t)(CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) * 10))
719 / 10.0;
720
721 std::vector<Double_t> temppulse;
722 std::map<Int_t, std::vector<Double_t>> templow;
723 std::map<Int_t, std::vector<Double_t>> temphigh;
724
725 temppulse = pulse;
726 for (Int_t i = 0; i < 32; i++) {
727 if (fDebug) fQA->CreateHist("Multi self", 32, -0.5, 31.5, 512, -12., 500.);
728 if (fDebug) fQA->Fill("Multi self", i, pulse[i]);
729 if (i >= shift) pulse[i] = 0.;
730 }
731 fPulseBuffer[address].first = pulse;
732
733 // std::cout<<"multi time: " << fTimeBuffer[address]<<" col: "<<maincol<<std::endl;
734 // Double_t dt=TMath::Abs(fMultiBuffer[address].second - fTimeBuffer[address]);
735
736
737 Int_t col = maincol;
738 Int_t FNaddress = 0;
739 Double_t FNshift = 0;
740 std::vector<Double_t> FNpulse = fPulseBuffer[address].first;
741 if (col >= 1)
742 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(fModAddress),
743 sec, row, col - 1);
744 Int_t FNtrigger = 1;
745
746 while (FNtrigger == 1 && FNaddress != 0) {
747 if (col >= 1)
748 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
749 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col - 1);
750 else
751 break;
752 col--;
753 FNpulse = fPulseBuffer[FNaddress].first;
754 templow[FNaddress] = FNpulse;
755 FNshift =
756 (fTimeBuffer[address] - fTimeBuffer[FNaddress]) / CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + shift;
757 for (Int_t i = 0; i < 32; i++) {
758 if (i >= FNshift) FNpulse[i] = 0.;
759 }
760 FNtrigger = CheckTrigger(FNpulse);
761 fPulseBuffer[FNaddress].first = FNpulse;
762 if (col == 0) break;
763
764 // std::cout<<"col: "<<col<<" "<< fTimeBuffer[FNaddress]<<" FNaddress: " << FNaddress<<" FNtrigger: "<< FNtrigger<<" samples: " << fPulseBuffer[FNaddress].first.size()<<" time: " << fTimeBuffer[FNaddress]<<std::endl;
765 }
766
767 col = maincol;
768 FNaddress = 0;
769 if (col < (ncols - 1))
770 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(fModAddress),
771 sec, row, col + 1);
772 FNtrigger = 1;
773
774 while (FNtrigger == 1 && FNaddress != 0) {
775 if (col < (ncols - 1))
776 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
777 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col + 1);
778 else
779 break;
780 col++;
781 FNpulse = fPulseBuffer[FNaddress].first;
782 temphigh[FNaddress] = FNpulse;
783 FNshift =
784 (fTimeBuffer[address] - fTimeBuffer[FNaddress]) / CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + shift;
785 for (Int_t i = 0; i < 32; i++) {
786 if (i >= FNshift) FNpulse[i] = 0.;
787 }
788 FNtrigger = CheckTrigger(FNpulse);
789 fPulseBuffer[FNaddress].first = FNpulse;
790 if (col == ncols - 1) break;
791 }
792 ProcessPulseBuffer(address, false, true, true, true);
793
794 for (Int_t i = 0; i < 32; i++) {
795 if (i < fPresamples) {
796 pulse[i] = temppulse[shift + i - fPresamples];
797 continue;
798 }
799 if (shift + i - fPresamples < 32) {
800 if (fTimeShift)
801 pulse[i] =
802 temppulse[shift + i - fPresamples]
803 + fCalibration * fMultiBuffer[address].first * 1e6
804 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
805 if (!fTimeShift)
806 pulse[i] = temppulse[shift + i - fPresamples]
807 + fCalibration * fMultiBuffer[address].first * 1e6
808 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
809 }
810 else {
811 if (fTimeShift)
812 pulse[i] =
813 fCalibration * fMultiBuffer[address].first * 1e6
814 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
815 if (!fTimeShift)
816 pulse[i] = fCalibration * fMultiBuffer[address].first * 1e6
817 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
818 }
819
820 // if(shift+i < 32){
821 // if(fTimeShift) pulse[i]=temppulse[shift+i]+fCalibration*fMultiBuffer[address].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)+timeshift);
822 // if(!fTimeShift) pulse[i]=temppulse[shift+i]+fCalibration*fMultiBuffer[address].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
823 // }
824 // else{
825 // if(fTimeShift) pulse[i]=fCalibration*fMultiBuffer[address].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)+timeshift);
826 // if(!fTimeShift) pulse[i]=fCalibration*fMultiBuffer[address].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
827 // }
828 }
829 for (Int_t i = 0; i < 32; i++) {
830 if (fDebug) fQA->CreateHist("Multi self after", 32, -0.5, 31.5, 512, -12., 500.);
831 if (fDebug) fQA->Fill("Multi self after", i, pulse[i]);
832 }
833
834 fTimeBuffer[address] = fMultiBuffer[address].second;
835 fPulseBuffer[address].first = pulse;
836
837 if (col < ncols - 1)
838 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(fModAddress),
839 sec, row, col + 1);
840 FNtrigger = 1;
841
842 while (FNtrigger == 1 && FNaddress != 0) {
843 if (col < (ncols - 1))
844 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
845 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col + 1);
846 else
847 break;
848 col++;
849 FNpulse = fPulseBuffer[FNaddress].first;
850 for (Int_t i = 0; i < 32; i++) {
851 if (i < fPresamples && temphigh[FNaddress].size() > 0.) {
852 pulse[i] = temphigh[FNaddress][shift + i - fPresamples];
853 continue;
854 }
855 if (fTimeShift) {
856 if ((size_t) shift + i - fPresamples < temphigh[FNaddress].size())
857 FNpulse[i] =
858 temphigh[FNaddress][shift + i - fPresamples]
859 + fCalibration * fMultiBuffer[FNaddress].first * 1e6
860 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
861 else
862 FNpulse[i] =
863 fCalibration * fMultiBuffer[FNaddress].first * 1e6
864 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
865 if (FNpulse[i] > fClipLevel && fClipping) FNpulse[i] = fClipLevel - 1; //clipping
866 }
867 if (!fTimeShift) {
868 if ((size_t) shift + i - fPresamples < temphigh[FNaddress].size())
869 FNpulse[i] = temphigh[FNaddress][shift + i - fPresamples]
870 + fCalibration * fMultiBuffer[FNaddress].first * 1e6
871 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
872 else
873 FNpulse[i] = fCalibration * fMultiBuffer[FNaddress].first * 1e6
874 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
875 if (FNpulse[i] > fClipLevel && fClipping) FNpulse[i] = fClipLevel - 1; //clipping
876 }
877
878 // if(fTimeShift){
879 // if(shift+i<32 && temphigh[FNaddress].size()>0) FNpulse[i]=temphigh[FNaddress][shift+i]+fCalibration*fMultiBuffer[FNaddress].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)+timeshift);
880 // else FNpulse[i]=fCalibration*fMultiBuffer[FNaddress].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC)+timeshift);
881 // if(FNpulse[i] > fClipLevel && fClipping) FNpulse[i]=fClipLevel-1; //clipping
882 // }
883 // if(!fTimeShift){
884 // if(shift+i<32 && temphigh[FNaddress].size()>0) FNpulse[i]=temphigh[FNaddress][shift+i]+fCalibration*fMultiBuffer[FNaddress].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
885 // else FNpulse[i]=fCalibration*fMultiBuffer[FNaddress].first*1e6*CalcResponse(i*CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
886 // if(FNpulse[i] > fClipLevel && fClipping) FNpulse[i]=fClipLevel-1; //clipping
887 // }
888 }
889
890 fPulseBuffer[FNaddress].first = FNpulse;
891 CbmMatch* FNdigiMatch = new CbmMatch();
892 if (!fMCBuffer[FNaddress].empty())
893 for (UInt_t links = 0; links < fMCBuffer[FNaddress][0].size(); links++)
894 FNdigiMatch->AddLink(CbmLink(fMultiBuffer[FNaddress].first, fMCBuffer[FNaddress][0][links],
895 fMCBuffer[FNaddress][1][links], fMCBuffer[FNaddress][2][links]));
896 fPulseBuffer[FNaddress].second = FNdigiMatch;
897 if (fDebug) {
898 std::vector<std::map<TString, Int_t>> vecLink;
899 std::map<TString, Int_t> LinkQA;
900 LinkQA["Event"] = fEventId;
901 LinkQA["Point"] = fPointId;
902 LinkQA["Time"] = fEventTime;
903 LinkQA["Charge"] = fMultiBuffer[FNaddress].first * 1e6 * fepoints;
904 vecLink.push_back(LinkQA);
905 fLinkQA[FNaddress] = vecLink;
906 }
907
908 FNtrigger = CheckTrigger(FNpulse);
909
910 fTimeBuffer[FNaddress] = fMultiBuffer[address].second;
911 fMultiBuffer.erase(FNaddress);
912 if (col == ncols - 1) break;
913 }
914
915 if (col >= 1)
916 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(fModAddress),
917 sec, row, col - 1);
918 FNtrigger = 1;
919
920 while (FNtrigger == 1 && FNaddress != 0) {
921 if (col >= 1)
922 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address),
923 CbmTrdAddress::GetModuleId(fModAddress), sec, row, col - 1);
924 else
925 break;
926 col--;
927 FNpulse = fPulseBuffer[FNaddress].first;
928 for (Int_t i = 0; i < 32; i++) {
929 if (i < fPresamples && templow[FNaddress].size() > 0.) {
930 pulse[i] = templow[FNaddress][shift + i - fPresamples];
931 continue;
932 }
933 if (fTimeShift) {
934 if ((size_t) shift + i - fPresamples < templow[FNaddress].size())
935 FNpulse[i] =
936 templow[FNaddress][shift + i - fPresamples]
937 + fCalibration * fMultiBuffer[FNaddress].first * 1e6
938 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
939 else
940 FNpulse[i] =
941 fCalibration * fMultiBuffer[FNaddress].first * 1e6
942 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC) + timeshift);
943 if (FNpulse[i] > fClipLevel && fClipping) FNpulse[i] = fClipLevel - 1; //clipping
944 }
945 if (!fTimeShift) {
946 if ((size_t) shift + i - fPresamples < templow[FNaddress].size())
947 FNpulse[i] = templow[FNaddress][shift + i - fPresamples]
948 + fCalibration * fMultiBuffer[FNaddress].first * 1e6
949 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
950 else
951 FNpulse[i] = fCalibration * fMultiBuffer[FNaddress].first * 1e6
952 * CalcResponse((i - fPresamples) * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
953 if (FNpulse[i] > fClipLevel && fClipping) FNpulse[i] = fClipLevel - 1; //clipping
954 }
955 }
956 fPulseBuffer[FNaddress].first = FNpulse;
957 CbmMatch* FNdigiMatch = new CbmMatch();
958 if (!fMCBuffer[FNaddress].empty())
959 for (UInt_t links = 0; links < fMCBuffer[FNaddress][0].size(); links++)
960 FNdigiMatch->AddLink(CbmLink(fMultiBuffer[FNaddress].first, fMCBuffer[FNaddress][0][links],
961 fMCBuffer[FNaddress][1][links], fMCBuffer[FNaddress][2][links]));
962 fPulseBuffer[FNaddress].second = FNdigiMatch;
963 if (fDebug) {
964 std::vector<std::map<TString, Int_t>> vecLink;
965 std::map<TString, Int_t> LinkQA;
966 LinkQA["Event"] = fEventId;
967 LinkQA["Point"] = fPointId;
968 LinkQA["Time"] = fEventTime;
969 LinkQA["Charge"] = fMultiBuffer[FNaddress].first * 1e6 * fepoints;
970 vecLink.push_back(LinkQA);
971 fLinkQA[FNaddress] = vecLink;
972 }
973
974 FNtrigger = CheckTrigger(FNpulse);
975
976 fTimeBuffer[FNaddress] = fMultiBuffer[address].second;
977 fMultiBuffer.erase(FNaddress);
978 if (col == 0) break;
979 }
980
981 CbmMatch* digiMatch = new CbmMatch();
982 if (!fMCBuffer[address].empty())
983 for (UInt_t links = 0; links < fMCBuffer[address][0].size(); links++)
984 digiMatch->AddLink(CbmLink(fMultiBuffer[address].first, fMCBuffer[address][0][links],
985 fMCBuffer[address][1][links], fMCBuffer[address][2][links]));
986 fPulseBuffer[address].second = digiMatch;
987 if (fDebug) {
988 std::vector<std::map<TString, Int_t>> vecLink;
989 std::map<TString, Int_t> LinkQA;
990 LinkQA["Event"] = fEventId;
991 LinkQA["Point"] = fPointId;
992 LinkQA["Time"] = fEventTime;
993 LinkQA["Charge"] = fMultiBuffer[address].first * 1e6 * fepoints;
994 vecLink.push_back(LinkQA);
995 fLinkQA[address] = vecLink;
996 }
997 }
998
999 return processed;
1000}
1001
1002Int_t CbmTrdModuleSimR::CheckTrigger(std::vector<Double_t> pulse)
1003{
1004
1005 Int_t slope = 0;
1006 Bool_t trigger = false;
1007 Bool_t falling = false;
1008 Bool_t multihit = false;
1009 for (size_t i = 0; i < pulse.size(); i++) {
1010 if (i < pulse.size() - 1) slope = pulse[i + 1] - pulse[i];
1011 if (slope >= fTriggerSlope && !trigger) trigger = true;
1012 if (slope < 0 && trigger) falling = true;
1013 if (slope >= fTriggerSlope && trigger && falling && (Int_t) i > fMaxBin) multihit = true;
1014 }
1015
1016 if (!trigger && !multihit) return 0;
1017 if (trigger && !multihit) return 1;
1018 if (trigger && multihit) return 2;
1019
1020 return 0;
1021}
1022
1023Int_t CbmTrdModuleSimR::GetMultiBin(std::vector<Double_t> pulse)
1024{
1025
1026 Int_t slope = 0;
1027 Int_t startbin = 0;
1028 Bool_t trigger = false;
1029 Bool_t falling = false;
1030 for (size_t i = 0; i < pulse.size(); i++) {
1031 if (i < 31) slope = pulse[i + 1] - pulse[i];
1032 if (slope >= fTriggerSlope && !trigger) trigger = true;
1033 if (slope < 0 && trigger) falling = true;
1034 if (slope >= fTriggerSlope && trigger && falling && (Int_t) i > fMaxBin) startbin = i;
1035 }
1036
1037 return startbin;
1038}
1039
1040
1041//_______________________________________________________________________________
1042Double_t CbmTrdModuleSimR::CalcPRF(Double_t x, Double_t W, Double_t h)
1043{
1044 Float_t K3 = 0.525;
1045 Double_t SqrtK3 = sqrt(K3);
1046
1047 return std::fabs(-1. / (2. * atan(SqrtK3))
1048 * (atan(SqrtK3 * tanh(TMath::Pi() * (-2. + SqrtK3) * (W + 2. * x) / (8. * h)))
1049 + atan(SqrtK3 * tanh(TMath::Pi() * (-2. + SqrtK3) * (W - 2. * x) / (8. * h)))));
1050}
1051
1052//_______________________________________________________________________________
1053Double_t CbmTrdModuleSimR::CalcResponse(Double_t t)
1054{
1055
1056 if (fShapingOrder == 1) return (t / fTau) * TMath::Exp(-(t / fTau));
1057 else
1058 return (t / fTau) * (t / fTau) * TMath::Exp(-(t / fTau));
1059}
1060
1061//_______________________________________________________________________________
1062Bool_t CbmTrdModuleSimR::DistributeCharge(Double_t pointin[3], Double_t pointout[3], Double_t delta[3], Double_t pos[3],
1063 Int_t ipoints)
1064{
1065 if (fDistributionMode == 1) {
1066 for (Int_t i = 0; i < 3; i++) {
1067 pos[i] = pointin[i] + (0.01) * delta[i] + 0.95 * delta[i] / fepoints * ipoints;
1068 }
1069 }
1070
1071 if (fDistributionMode == 5) {
1072 for (Int_t i = 0; i < 3; i++) {
1073 pos[i] = pointin[i] + 0.5 * delta[i];
1074 }
1075 }
1076
1077
1078 //in development
1079 if (fDistributionMode == 2) {
1080 Double_t lastpos[3] = {pointin[0], pointin[1], pointin[2]};
1081 Double_t dist_gas = TMath::Sqrt(delta[0] * delta[0] + delta[1] * delta[1] + delta[2] * delta[2]);
1082 if (ipoints > 0)
1083 for (Int_t i = 0; i < 3; i++)
1084 lastpos[i] = pos[i];
1085 Double_t roll = gRandom->Integer(100);
1086 Double_t s = (GetStep(dist_gas, roll) / dist_gas) / 2;
1087
1088 if (TMath::Abs(lastpos[0] + s * delta[0]) > fDigiPar->GetSizeX()
1089 || TMath::Abs(lastpos[1] + s * delta[1]) > fDigiPar->GetSizeY() || (lastpos[2] + s * delta[2]) >= pointout[2]) {
1090 for (Int_t i = 0; i < 3; i++) {
1091 pos[i] = lastpos[i];
1092 }
1093 return true;
1094 }
1095 for (Int_t i = 0; i < 3; i++)
1096 pos[i] = lastpos[i] + s * delta[i];
1097
1098 return true;
1099 }
1100
1101 if (fDistributionMode == 3) {
1102 for (Int_t i = 0; i < 3; i++) {
1103 pos[i] = pointin[i] + (0.5 + ipoints) * delta[i];
1104 }
1105 }
1106
1107
1108 return true;
1109}
1110
1111//_______________________________________________________________________________
1112Bool_t CbmTrdModuleSimR::MakeDigi(CbmTrdPoint* point, Double_t time, Bool_t TR)
1113{
1114 // if(!CbmTrdDigitizer::IsTimeBased()) fPulseSwitch = false;
1115 if (fPulseSwitch && !fMessageConverter->GetSetter()) {
1116 if (fDebug) fQA->CreateHist("Pulse", 32, -0.5, 31.5, 512, -12., 500.);
1117 if (fDebug) fQA->CreateHist("Shift", 63, -0.5, 62.5);
1118 if (fDebug) fQA->CreateHist("Multi Quote", 4, -0.5, 3.5);
1119 // if(fDebug) fMessageConverter->SetDebug(true);
1120 if (fDebug) fMessageConverter->SetQA(fQA);
1121 std::vector<Int_t> recomask;
1122 for (Int_t i = frecostart; i <= frecostop; i++)
1123 recomask.push_back(i);
1124 fMessageConverter->SetRecoMask(recomask);
1125 fMessageConverter->SetCalibration(fCalibration);
1126 fMessageConverter->SetShapingOrder(fShapingOrder);
1127 fMessageConverter->SetTau(fTau);
1128 fMessageConverter->SetPresamples(fPresamples);
1129 fMessageConverter->SetMaxBin(fMaxBin);
1130 fMessageConverter->SetMinBin(fMinBin);
1131 //fMessageConverter->SetLookup(1);
1132 TString dir = getenv("VMCWORKDIR");
1133 TString filename = dir + "/parameters/trd/FeatureExtractionLookup.root";
1134 // if(fDistributionMode == 1 && fepoints == 3) filename= dir + "/parameters/trd/Feature_mode1_3points.root";
1135 // if(fDistributionMode == 1 && fepoints == 5) filename= dir + "/parameters/trd/Feature_mode1_5points.root";
1136 // if(fDistributionMode == 4) filename= dir + "/parameters/trd/Feature_mode4.root";
1137 fMessageConverter->SetReadFile(filename.Data());
1138 fMessageConverter->Init();
1139 fMessageConverter->SetSetter(true);
1140 }
1141
1142
1143 // calculate current physical time
1144 fCurrentTime = time + point->GetTime(); //+ AddDrifttime(gRandom->Integer(240))*1000; //convert to ns;
1145 fEventTime = time;
1146
1147 const Double_t nClusterPerCm = 1.0;
1148 Double_t point_in[3] = {point->GetXIn(), point->GetYIn(), point->GetZIn()};
1149 Double_t point_out[3] = {point->GetXOut(), point->GetYOut(), point->GetZOut()};
1150 Double_t local_point_out[3]; // exit point coordinates in local module cs
1151 Double_t local_point_in[3]; // entrace point coordinates in local module cs
1152 gGeoManager->cd(GetPath());
1153 gGeoManager->MasterToLocal(point_in, local_point_in);
1154 gGeoManager->MasterToLocal(point_out, local_point_out);
1155 SetPositionMC(local_point_out);
1156
1157 for (Int_t i = 0; i < 3; i++)
1158 fQAPosition[i] = point_in[i];
1159 for (Int_t i = 0; i < 3; i++)
1160 fQAPos_out[i] = point_out[i];
1161
1162 //fCurrentTime -= 1e9 * (point_in[2] / 100)
1163 // / 2.99e8; // subtract time of flight to the detector layer
1164
1165 // General processing on the MC point
1166 Double_t ELoss(0.), ELossTR(0.), ELossdEdX(point->GetEnergyLoss());
1167 if (fRadiator && TR) {
1168 nofElectrons++;
1169 if (
1170 fRadiator->LatticeHit(
1171 point)) { // electron has passed lattice grid (or frame material) befor reaching the gas volume -> TR-photons have been absorbed by the lattice grid
1172 nofLatticeHits++;
1173 }
1174 else if (point_out[2] >= point_in[2]) { //electron has passed the radiator
1175 TVector3 mom;
1176 point->Momentum(mom);
1177 ELossTR = fRadiator->GetTR(mom);
1178 }
1179 }
1180 ELoss = ELossTR + ELossdEdX;
1181
1182 if (fDebug) fQA->CreateHist("E MC", 200, 0., 50.0);
1183 if (fDebug) fQA->Fill("E MC", ELoss * 1e6);
1184
1185 if (ELoss > fMinimumChargeTH) nofPointsAboveThreshold++;
1186
1187 Double_t cluster_pos[3]; // cluster position in local module coordinate system
1188 Double_t cluster_delta
1189 [3]; // vector pointing in MC-track direction with length of one path slice within chamber volume to creat n cluster
1190
1191 Double_t trackLength = 0;
1192
1193 for (Int_t i = 0; i < 3; i++) {
1194 cluster_delta[i] = (local_point_out[i] - local_point_in[i]);
1195 trackLength += cluster_delta[i] * cluster_delta[i];
1196 }
1197 trackLength = TMath::Sqrt(trackLength);
1198 Int_t nCluster =
1199 trackLength / nClusterPerCm + 0.9; // Track length threshold of minimum 0.1cm track length in gas volume
1200
1201 if (fnClusterConst > 0) {
1202 nCluster = fnClusterConst; // Set number of cluster to constant value
1203 }
1204
1205 if (nCluster < 1) { return kFALSE; }
1206 nCluster = 1;
1207
1208 for (Int_t i = 0; i < 3; i++) {
1209 cluster_delta[i] /= Double_t(nCluster);
1210 }
1211
1212 Double_t clusterELoss = ELoss / Double_t(nCluster);
1213 Double_t clusterELossTR = ELossTR / Double_t(nCluster);
1214
1215
1216 //to change the number of ionization points in the gas
1217 Int_t epoints = fepoints;
1218
1219 if (fDistributionMode == 3) { epoints = nCluster; }
1220 if (fDistributionMode == 5) { epoints = 1; }
1221
1222 //in development
1223 std::vector<Double_t> vec;
1224 std::pair<Int_t, std::vector<Double_t>> steps = std::make_pair(0, vec);
1225 if (fDistributionMode == 4) {
1226 Double_t dist_gas = TMath::Sqrt(cluster_delta[0] * cluster_delta[0] + cluster_delta[1] * cluster_delta[1]
1227 + cluster_delta[2] * cluster_delta[2]);
1228 steps = (GetTotalSteps(local_point_in, local_point_out, dist_gas));
1229 epoints = steps.first;
1230 fepoints = steps.first;
1231 if (fDebug) fQA->CreateHist("Ionization", 63, -0.5, 62.5);
1232 if (fDebug) fQA->Fill("Ionization", steps.first);
1233 if (fDebug) fQA->CreateHist("Dist Ionization", 63, -0.5, 62.5, 20, 0., 2.);
1234 if (fDebug) fQA->Fill("Dist Ionization", steps.first, dist_gas);
1235 }
1236
1237 if (epoints != 1) {
1238 clusterELoss = ELoss / epoints;
1239 clusterELossTR = ELossTR / epoints;
1240 }
1241
1242 if (!fPulseSwitch) {
1243 for (Int_t ipoints = 0; ipoints < epoints; ipoints++) {
1244 Bool_t dist = DistributeCharge(local_point_in, local_point_out, cluster_delta, cluster_pos, ipoints);
1245 if (!dist) break;
1246 if (fDistributionMode == 4)
1247 for (Int_t i = 0; i < 3; i++)
1248 cluster_pos[i] = steps.second[i + ipoints * 3];
1249
1250 if (fDigiPar->GetSizeX() < std::fabs(cluster_pos[0]) || fDigiPar->GetSizeY() < std::fabs(cluster_pos[1])) {
1251 printf("-> nC %i/%i x: %7.3f y: %7.3f \n", ipoints, nCluster - 1, cluster_pos[0], cluster_pos[1]);
1252 for (Int_t i = 0; i < 3; i++)
1253 printf(" (%i) | in: %7.3f + delta: %7.3f * cluster: %i/%i = "
1254 "cluster_pos: %7.3f out: %7.3f g_in:%f g_out:%f\n",
1255 i, local_point_in[i], cluster_delta[i], ipoints, (Int_t) nCluster, cluster_pos[i], local_point_out[i],
1256 point_in[i], point_out[i]);
1257 }
1258
1259 if (CbmTrdDigitizer::AddNoise()) {
1260 Int_t noiserate = gRandom->Uniform(0, 3); //still in development
1261 Double_t simtime = fCurrentTime;
1262 for (Int_t ndigi = 0; ndigi < noiserate; ndigi++) {
1263 NoiseTime(time);
1264 // ScanPadPlane(cluster_pos, gRandom->Gaus(0, fSigma_noise_keV * 1.E-6), 0,epoints,ipoints);
1265 }
1266 fCurrentTime = simtime;
1267 }
1268 ScanPadPlane(cluster_pos, 0., clusterELoss, clusterELossTR, epoints, ipoints);
1269 }
1270 }
1271
1272 Double_t driftcomp = 10000;
1273 Int_t start = -1;
1274 Double_t Ionizations[epoints][3];
1275 if (fPulseSwitch) {
1276 for (Int_t ipoints = 0; ipoints < epoints; ipoints++) {
1277 Bool_t dist = DistributeCharge(local_point_in, local_point_out, cluster_delta, cluster_pos, ipoints);
1278 if (!dist) break;
1279 if (fDistributionMode == 4)
1280 for (Int_t i = 0; i < 3; i++)
1281 cluster_pos[i] = steps.second[i + ipoints * 3];
1282 for (Int_t i = 0; i < 3; i++)
1283 Ionizations[ipoints][i] = cluster_pos[i];
1284
1285 if (fDigiPar->GetSizeX() < std::fabs(cluster_pos[0]) || fDigiPar->GetSizeY() < std::fabs(cluster_pos[1])) {
1286 printf("-> nC %i/%i x: %7.3f y: %7.3f \n", ipoints, nCluster - 1, cluster_pos[0], cluster_pos[1]);
1287 for (Int_t i = 0; i < 3; i++)
1288 printf(" (%i) | in: %7.3f + delta: %7.3f * cluster: %i/%i = "
1289 "cluster_pos: %7.3f out: %7.3f g_in:%f g_out:%f\n",
1290 i, local_point_in[i], cluster_delta[i], ipoints, (Int_t) nCluster, cluster_pos[i], local_point_out[i],
1291 point_in[i], point_out[i]);
1292 }
1293
1294 fDigiPar->ProjectPositionToNextAnodeWire(cluster_pos);
1295 Int_t relz = 239 - (cluster_pos[2] - local_point_in[2]) / 0.005;
1296 if (relz > 239 || relz < 0) relz = 239;
1297
1298 Double_t Drift_x =
1299 TMath::Abs(double(int(cluster_pos[0] * 1000000) % int(fDigiPar->GetAnodeWireSpacing() * 100)) / 100);
1300 // std::cout<<" pos: " << Drift_x<< " " << int(cluster_pos[0]*100000)<< " " << int(fDigiPar->GetAnodeWireSpacing()*100)<<" " << (int(cluster_pos[0]*100000) % int(fDigiPar->GetAnodeWireSpacing()*100))<<std::endl;
1301 if (TMath::Abs(AddDrifttime(Drift_x, cluster_pos[2])) < driftcomp
1302 && TMath::Abs(AddDrifttime(Drift_x, cluster_pos[2])) > 0.) {
1303 driftcomp = TMath::Abs(AddDrifttime(Drift_x, cluster_pos[2]));
1304 start = ipoints;
1305 }
1306 }
1307
1308 if (start == -1) return false;
1309 for (Int_t i = 0; i < 3; i++)
1310 cluster_pos[i] = Ionizations[start][i];
1311
1312 Int_t relz = 239 - (cluster_pos[2] - local_point_in[2]) / 0.005;
1313 if (relz > 239 || relz < 0) relz = 239;
1314 Double_t Drift_x =
1315 TMath::Abs(double(int(cluster_pos[0] * 1000000) % int(fDigiPar->GetAnodeWireSpacing() * 100)) / 100);
1316 Double_t reldrift = TMath::Abs(AddDrifttime(Drift_x, cluster_pos[2]) - driftcomp) * 1000;
1317 fDriftStart = driftcomp * 1000;
1318 fCurrentTime += fDriftStart;
1319 if (reldrift < 250.) ScanPadPlane(cluster_pos, 0., clusterELoss, clusterELossTR, epoints, start);
1320 if (fPrintPulse) std::cout << std::endl;
1321
1322 for (Int_t ipoints = 0; ipoints < epoints; ipoints++) {
1323 if (ipoints == start) continue;
1324 for (Int_t i = 0; i < 3; i++)
1325 cluster_pos[i] = Ionizations[ipoints][i];
1326
1327 relz = 239 - (cluster_pos[2] - local_point_in[2]) / 0.005;
1328 if (relz > 239 || relz < 0) relz = 239;
1329 Drift_x = TMath::Abs(double(int(cluster_pos[0] * 1000000) % int(fDigiPar->GetAnodeWireSpacing() * 100)) / 100);
1330 reldrift = TMath::Abs(AddDrifttime(Drift_x, cluster_pos[2]) - driftcomp) * 1000;
1331 if (reldrift < 250.) ScanPadPlane(cluster_pos, reldrift, clusterELoss, clusterELossTR, epoints, ipoints);
1332 if (fPrintPulse) std::cout << std::endl;
1333 }
1334 }
1335
1336 fLastEventTime = time;
1337 fLastPoint = fPointId;
1338 fLastEvent = fEventId;
1339 fLastTime = fCurrentTime;
1340 return true;
1341}
1342
1343
1344//_______________________________________________________________________________
1345void CbmTrdModuleSimR::ScanPadPlane(const Double_t* local_point, Double_t reldrift, Double_t clusterELoss,
1346 Double_t clusterELossTR, Int_t epoints, Int_t ipoint)
1347{
1348 Int_t sectorId(-1), columnId(-1), rowId(-1);
1349 fDigiPar->GetPadInfo(local_point, sectorId, columnId, rowId);
1350 if (sectorId < 0 && columnId < 0 && rowId < 0) { return; }
1351 else {
1352 for (Int_t i = 0; i < sectorId; i++) {
1353 rowId += fDigiPar->GetNofRowsInSector(i); // local -> global row
1354 }
1355
1356 // for (Int_t i = 0; i < 3; i++) fQAPosition[i]=local_point[i];
1357
1358 Double_t displacement_x(0), displacement_y(0); //mm
1359 Double_t h = fDigiPar->GetAnodeWireToPadPlaneDistance();
1360 Double_t W(fDigiPar->GetPadSizeX(sectorId)), H(fDigiPar->GetPadSizeY(sectorId));
1361 fDigiPar->TransformToLocalPad(local_point, displacement_x, displacement_y);
1362
1363 Int_t maxRow(6);
1364 if (fnScanRowConst > 0) maxRow = fnScanRowConst;
1365
1366 Int_t startRow(rowId - maxRow / 2);
1367 Int_t secRow(-1), targCol(-1), targRow(-1), targSec(-1), address(-1), fnRow(fDigiPar->GetNofRows()),
1368 fnCol(fDigiPar->GetNofColumns());
1369
1370 for (Int_t iRow = startRow; iRow <= rowId + maxRow / 2; iRow++) {
1371 Int_t iCol = columnId;
1372 if (((iCol >= 0) && (iCol <= fnCol - 1)) && ((iRow >= 0) && (iRow <= fnRow - 1))) { // real adress
1373 targSec = fDigiPar->GetSector(iRow, secRow);
1374 address = CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), targSec, secRow, iCol);
1375 }
1376 else {
1377 targRow = iRow;
1378 targCol = iCol;
1379 if (iCol < 0) { targCol = 0; }
1380 else if (iCol > fnCol - 1) {
1381 targCol = fnCol - 1;
1382 }
1383 if (iRow < 0) { targRow = 0; }
1384 else if (iRow > fnRow - 1) {
1385 targRow = fnRow - 1;
1386 }
1387
1388 targSec = fDigiPar->GetSector(targRow, secRow);
1389 address =
1390 CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), targSec, secRow, targCol);
1391 }
1392
1393 Bool_t print = false;
1394
1395 //distribute the mc charge fraction over the channels wit the PRF
1396 Float_t chargeFraction = 0;
1397 Float_t ch = 0;
1398 Float_t tr = 0;
1399
1400 // std::cout<<" prf half: " << CalcPRF(0 * W , W, h)<<std::endl;
1401 // std::cout<<" prf half -1 : " << CalcPRF(-1 * W , W, h)<<std::endl;
1402 // std::cout<<" prf half +1: " << CalcPRF(1 * W , W, h)<<std::endl;
1403 chargeFraction =
1404 CalcPRF((iCol - columnId) * W - displacement_x, W, h) * CalcPRF((iRow - rowId) * H - displacement_y, H, h);
1405
1406 ch = chargeFraction * clusterELoss;
1407 tr = chargeFraction * clusterELossTR;
1408
1409 Bool_t lowerend = false;
1410 Bool_t upperend = false;
1411 Int_t collow = 1;
1412 Int_t colhigh = 1;
1413
1414 if (fDebug) fQA->CreateHist("E self MC", 200, 0., 50.0);
1415 if (fDebug) fQA->CreateHist("E FN MC", 200, 0., 10.0);
1416
1417 if (ch >= (fMinimumChargeTH / epoints)) {
1418 if (!CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && !print)
1419 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1420 if (!CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && print) {
1421 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1422 std::cout << " time: " << fCurrentTime << " col: " << iCol << " row: " << iRow - rowId
1423 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1424 }
1425 if (CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch) AddDigitoBuffer(address, ch, tr, fCurrentTime, Int_t(1));
1426 if (fPulseSwitch) {
1427 if (fDebug) fQA->Fill("E self MC", ch * epoints * 1e6);
1428 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1429 }
1430 if (fPulseSwitch && print) {
1431 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1432 std::cout << " time: " << fCurrentTime << " col: " << iCol << " row: " << iRow - rowId
1433 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1434 }
1435
1436
1437 while (!lowerend) {
1438 if ((((iCol - collow) >= 0) && ((iCol - collow) <= fnCol - 1))
1439 && ((iRow >= 0) && (iRow <= fnRow - 1))) { // real adress
1440 targSec = fDigiPar->GetSector(iRow, secRow);
1441 address = CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), targSec, secRow,
1442 iCol - collow);
1443 }
1444 else {
1445 break;
1446 }
1447
1448 chargeFraction = CalcPRF(((iCol - collow) - columnId) * W - displacement_x, W, h)
1449 * CalcPRF((iRow - rowId) * H - displacement_y, H, h);
1450 ch = chargeFraction * clusterELoss;
1451 tr = chargeFraction * clusterELossTR;
1452
1453
1454 if (ch >= (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && !print) {
1455 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1456 collow++;
1457 }
1458 if (ch < (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && !print) {
1459 AddDigi(address, ch, tr, fCurrentTime, Int_t(2));
1460 lowerend = true;
1461 }
1462 if (ch >= (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && print) {
1463 std::cout << " time: " << fCurrentTime << " col: " << iCol - collow << " row: " << iRow - rowId
1464 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1465 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1466 collow++;
1467 }
1468 if (ch < (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && print) {
1469 std::cout << " time: " << fCurrentTime << " col: " << iCol - collow << " row: " << iRow - rowId
1470 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 0 " << std::endl;
1471 AddDigi(address, ch, tr, fCurrentTime, Int_t(2));
1472 lowerend = true;
1473 }
1474 if (ch >= (fMinimumChargeTH / epoints) && CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch) {
1475 AddDigitoBuffer(address, ch, tr, fCurrentTime, Int_t(1));
1476 collow++;
1477 }
1478 if (ch < (fMinimumChargeTH / epoints) && CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch) {
1479 AddDigitoBuffer(address, ch, tr, fCurrentTime, Int_t(2));
1480 lowerend = true;
1481 }
1482 if (ch >= (fMinimumChargeTH / epoints) && fPulseSwitch && !print) {
1483 if (fDebug) fQA->Fill("E self MC", ch * epoints * 1e6);
1484 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1485 collow++;
1486 }
1487 if (ch < (fMinimumChargeTH / epoints) && fPulseSwitch && !print) {
1488 if (fDebug) fQA->Fill("E FN MC", ch * epoints * 1e6);
1489 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1490 lowerend = true;
1491 }
1492
1493 if (ch >= (fMinimumChargeTH / epoints) && fPulseSwitch && print) {
1494 std::cout << " time: " << fCurrentTime << " col: " << iCol - collow << " row: " << iRow - rowId
1495 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1496 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1497 collow++;
1498 }
1499 if (ch < (fMinimumChargeTH / epoints) && fPulseSwitch && print) {
1500 std::cout << " time: " << fCurrentTime << " col: " << iCol - collow << " row: " << iRow - rowId
1501 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 0 " << std::endl;
1502 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1503 lowerend = true;
1504 }
1505 }
1506
1507 while (!upperend) {
1508
1509 if ((((iCol + colhigh) >= 0) && ((iCol + colhigh) <= fnCol - 1))
1510 && ((iRow >= 0) && (iRow <= fnRow - 1))) { // real adress
1511 targSec = fDigiPar->GetSector(iRow, secRow);
1512 address = CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), targSec, secRow,
1513 iCol + colhigh);
1514 }
1515 else {
1516 break;
1517 }
1518
1519
1520 chargeFraction = CalcPRF(((iCol + colhigh) - columnId) * W - displacement_x, W, h)
1521 * CalcPRF((iRow - rowId) * H - displacement_y, H, h);
1522 ch = chargeFraction * clusterELoss;
1523 tr = chargeFraction * clusterELossTR;
1524
1525 if (ch >= (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && !print) {
1526 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1527 colhigh++;
1528 }
1529 if (ch < (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && !print) {
1530 AddDigi(address, ch, tr, fCurrentTime, Int_t(2));
1531 upperend = true;
1532 }
1533 if (ch >= (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && print) {
1534 std::cout << " time: " << fCurrentTime << " col: " << iCol + colhigh << " row: " << iRow - rowId
1535 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1536 AddDigi(address, ch, tr, fCurrentTime, Int_t(1));
1537 colhigh++;
1538 }
1539 if (ch < (fMinimumChargeTH / epoints) && !CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch && print) {
1540 std::cout << " time: " << fCurrentTime << " col: " << iCol + colhigh << " row: " << iRow - rowId
1541 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 0 " << std::endl;
1542 AddDigi(address, ch, tr, fCurrentTime, Int_t(2));
1543 upperend = true;
1544 }
1545 if (ch >= (fMinimumChargeTH / epoints) && CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch) {
1546 AddDigitoBuffer(address, ch, tr, fCurrentTime, Int_t(1));
1547 colhigh++;
1548 }
1549 if (ch < (fMinimumChargeTH / epoints) && CbmTrdDigitizer::IsTimeBased() && !fPulseSwitch) {
1550 AddDigitoBuffer(address, ch, tr, fCurrentTime, Int_t(2));
1551 upperend = true;
1552 }
1553 if (ch >= (fMinimumChargeTH / epoints) && fPulseSwitch && !print) {
1554 if (fDebug) fQA->Fill("E self MC", ch * epoints * 1e6);
1555 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1556 colhigh++;
1557 }
1558 if (ch < (fMinimumChargeTH / epoints) && fPulseSwitch && !print) {
1559 if (fDebug) fQA->Fill("E FN MC", ch * epoints * 1e6);
1560 if (ipoint == epoints - 1 && epoints > 1)
1561 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1562 if (ipoint != epoints - 1 && epoints > 1)
1563 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1564 if (epoints == 1) AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1565 upperend = true;
1566 }
1567
1568 if (ch >= (fMinimumChargeTH / epoints) && fPulseSwitch && print) {
1569 std::cout << " time: " << fCurrentTime << " col: " << iCol + colhigh << " row: " << iRow - rowId
1570 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 1 " << std::endl;
1571 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(1), epoints, ipoint);
1572 colhigh++;
1573 }
1574 if (ch < (fMinimumChargeTH / epoints) && fPulseSwitch && print) {
1575 std::cout << " time: " << fCurrentTime << " col: " << iCol + colhigh << " row: " << iRow - rowId
1576 << " secrow: " << secRow << " charge: " << ch * 1e6 << " 0 " << std::endl;
1577 if (ipoint == epoints - 1 && epoints > 1)
1578 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1579 if (ipoint != epoints - 1 && epoints > 1)
1580 AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1581 if (epoints == 1) AddDigitoPulseBuffer(address, reldrift, ch, tr, fCurrentTime, Int_t(2), epoints, ipoint);
1582 upperend = true;
1583 }
1584 }
1585
1586 if (print) std::cout << std::endl;
1587
1588 } //if charge > trigger
1589 } //for rows
1590 }
1591}
1592
1593
1594//_______________________________________________________________________________
1595void CbmTrdModuleSimR::SetAsicPar(CbmTrdParModAsic* /*p*/)
1596{
1600 if (!fDigiPar) {
1601 LOG(warn) << GetName() << "::SetAsicPar : No Digi params for module " << fModAddress
1602 << ". Try calling first CbmTrdModSim::SetDigiPar.";
1603 return;
1604 }
1605
1606 if (fAsicPar) {
1607 LOG(warn) << GetName() << "::SetAsicPar : The list for module " << fModAddress << " already initialized.";
1608 return;
1609 }
1610 fAsicPar = new CbmTrdParModAsic();
1611 CbmTrdParSpadic* asic(NULL);
1612
1613 Int_t iFebGroup = 0; // 1; 2; // normal, super, ultimate
1614 Int_t gRow[3] = {2, 2, 2}; // re-ordering on the feb -> same mapping for normal and super
1615 Int_t gCol[3] = {16, 8, 4}; // re-ordering on the feb -> same mapping for normal and super
1616 Double_t xAsic = 0; // x position of Asic
1617 Double_t yAsic = 0; // y position of Asic
1618
1619 Int_t rowId(0), isecId(0), irowId(0), iAsic(0);
1620 for (Int_t s = 0; s < fDigiPar->GetNofSectors(); s++) {
1621 for (Int_t r = 0; r < fDigiPar->GetNofRowsInSector(s); r++) {
1622 for (Int_t c = 0; c < fDigiPar->GetNofColumnsInSector(s); c++) {
1623 // ultimate density 6 rows, 5 pads
1624 // super density 4 rows, 8 pads
1625 // normal density 2 rows, 16 pads
1626 if ((rowId % gRow[iFebGroup]) == 0) {
1627 if ((c % gCol[iFebGroup]) == 0) {
1628 xAsic = c + gCol[iFebGroup] / 2.;
1629 yAsic = r + gRow[iFebGroup] / 2.;
1630
1631 Double_t local_point[3];
1632 Double_t padsizex = fDigiPar->GetPadSizeX(s);
1633 Double_t padsizey = fDigiPar->GetPadSizeY(s);
1634
1635 // calculate position in sector coordinate system
1636 // with the origin in the lower left corner (looking upstream)
1637 local_point[0] = ((Int_t)(xAsic + 0.5) * padsizex);
1638 local_point[1] = ((Int_t)(yAsic + 0.5) * padsizey);
1639
1640 // calculate position in module coordinate system
1641 // with the origin in the lower left corner (looking upstream)
1642 local_point[0] += fDigiPar->GetSectorBeginX(s);
1643 local_point[1] += fDigiPar->GetSectorBeginY(s);
1644
1645 // local_point[i] must be >= 0 at this point Double_t local_point[3];
1646 Double_t fDx(GetDx()), fDy(GetDy());
1647 asic = new CbmTrdParSpadic(iAsic, iFebGroup, local_point[0] - fDx, local_point[1] - fDy);
1648 fAsicPar->SetAsicPar(asic);
1649 if (local_point[0] > 2 * fDx) {
1650 LOG(error) << "CbmTrdModuleSimR::SetAsicPar: asic position x=" << local_point[0]
1651 << " is out of bounds [0," << 2 * fDx << "]!";
1652 fDigiPar->Print("all");
1653 }
1654 if (local_point[1] > 2 * fDy) {
1655 LOG(error) << "CbmTrdModuleSimR::SetAsicPar: asic position y=" << local_point[1]
1656 << " is out of bounds [0," << 2 * fDy << "]!";
1657 fDigiPar->Print("all");
1658 }
1659 for (Int_t ir = rowId; ir < rowId + gRow[iFebGroup]; ir++) {
1660 for (Int_t ic = c; ic < c + gCol[iFebGroup]; ic++) {
1661 if (ir >= fDigiPar->GetNofRows())
1662 LOG(error) << "CbmTrdModuleSimR::SetAsicPar: ir " << ir << " is out of bounds!";
1663 if (ic >= fDigiPar->GetNofColumns())
1664 LOG(error) << "CbmTrdModuleSimR::SetAsicPar: ic " << ic << " is out of bounds!";
1665 isecId = fDigiPar->GetSector((Int_t) ir, irowId);
1666 asic->SetChannelAddress(
1667 CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), isecId, irowId, ic));
1668 //s, ir, ic));//new
1669 if (false)
1670 printf(" M:%10i(%4i) s: %i irowId: %4i ic: "
1671 "%4i r: %4i c: %4i address:%10i\n",
1672 fModAddress, CbmTrdAddress::GetModuleId(fModAddress), isecId, irowId, ic, r, c,
1673 CbmTrdAddress::GetAddress(fLayerId, fModAddress, isecId, irowId, ic));
1674 }
1675 }
1676 iAsic++; // next Asic
1677 }
1678 }
1679 }
1680 rowId++;
1681 }
1682 }
1683
1684 // Self Test
1685 for (Int_t s = 0; s < fDigiPar->GetNofSectors(); s++) {
1686 const Int_t nRow = fDigiPar->GetNofRowsInSector(s);
1687 const Int_t nCol = fDigiPar->GetNofColumnsInSector(s);
1688 for (Int_t r = 0; r < nRow; r++) {
1689 for (Int_t c = 0; c < nCol; c++) {
1690 Int_t channelAddress = CbmTrdAddress::GetAddress(fLayerId, CbmTrdAddress::GetModuleId(fModAddress), s, r, c);
1691 if (fAsicPar->GetAsicAddress(channelAddress) == -1)
1692 LOG(error) << "CbmTrdModuleSimR::SetAsicPar: Channel address:" << channelAddress
1693 << " is not or multible initialized in module " << fModAddress
1694 << "(ID:" << CbmTrdAddress::GetModuleId(fModAddress) << ")"
1695 << "(s:" << s << ", r:" << r << ", c:" << c << ")";
1696 }
1697 }
1698 }
1699}
1700
1701
1702//_______________________________________________________________________________
1703void CbmTrdModuleSimR::SetNoiseLevel(Double_t sigma_keV) { fSigma_noise_keV = sigma_keV; }
1704
1705//_______________________________________________________________________________
1706void CbmTrdModuleSimR::SetDistributionPoints(Int_t points) { fepoints = points; }
1707
1708//_______________________________________________________________________________
1709void CbmTrdModuleSimR::SetSpadicResponse(Double_t calibration, Double_t tau)
1710{
1711
1712 fCalibration = calibration;
1713 fTau = tau;
1714 Double_t sum = 0;
1715 for (auto i = frecostart; i <= frecostop; i++)
1716 sum += fCalibration * CalcResponse(i * CbmTrdDigi::Clk(CbmTrdDigi::eCbmTrdAsicType::kSPADIC));
1717 fEReco = sum;
1718}
1719
1720//_______________________________________________________________________________
1721void CbmTrdModuleSimR::SetPulsePars(Int_t mode)
1722{
1723
1724 if (mode == 1) {
1725 frecostart = 2;
1726 frecostop = 8;
1727 }
1728 if (mode == 2) {
1729 frecostart = 2;
1730 frecostop = 6;
1731 }
1732}
1733
1734
1735//_______________________________________________________________________________
1736void CbmTrdModuleSimR::SetPulseMode(Bool_t pulsed = true) { fPulseSwitch = pulsed; }
1737
1738
1739//_______________________________________________________________________________
1740void CbmTrdModuleSimR::SetPadPlaneScanArea(Int_t row)
1741{
1742 if (row % 2 == 0) row += 1;
1743 fnScanRowConst = row;
1744}
1745
1746
1747//_______________________________________________________________________________
1748Double_t CbmTrdModuleSimR::AddNoise(Double_t charge)
1749{
1750
1751 if (fSigma_noise_keV > 0.0 && CbmTrdDigitizer::AddNoise() && !fPulseSwitch) {
1752 Double_t noise = gRandom->Gaus(0,
1753 fSigma_noise_keV * 1e-6); // keV->GeV // add only once per digi and event noise !!!
1754 charge += noise; // resulting charge can be < 0 -> possible problems with position reconstruction
1755 return charge;
1756 }
1757 else
1758 return 0.;
1759}
1760
1761//_______________________________________________________________________________
1762Int_t CbmTrdModuleSimR::AddNoiseADC()
1763{
1764
1765
1766 if (CbmTrdDigitizer::AddNoise() && fPulseSwitch) {
1767 if (fAdcNoise == 0) return 0.;
1768 Int_t noise = gRandom->Gaus(0, fAdcNoise);
1769 return noise;
1770 // return 0;
1771 }
1772 else
1773 return 0.;
1774}
1775
1776std::vector<Double_t> CbmTrdModuleSimR::AddCorrelatedNoise(std::vector<Double_t> pulse)
1777{
1778
1779 //dummy for now
1780 return pulse;
1781
1782 // for(size_t i=0;i<pulse.size();i++){
1783 // pulse[i] += TMath::Sin(i)*5;
1784 // }
1785
1786 // return pulse;
1787}
1788
1789
1790//_______________________________________________________________________________
1791Int_t CbmTrdModuleSimR::AddCrosstalk(Double_t address, Int_t i, Int_t sec, Int_t row, Int_t col, Int_t ncols)
1792{
1793
1794 Double_t cross = 0.;
1795 if (fAddCrosstalk && fPulseSwitch) {
1796 Int_t FNaddress = 0;
1797 if (col >= 1)
1798 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(fModAddress),
1799 sec, row, col - 1);
1800 if (fPulseBuffer[FNaddress].first.size() >= (size_t) i + 1)
1801 cross += fPulseBuffer[FNaddress].first[i] * fCrosstalkLevel;
1802
1803 FNaddress = 0;
1804 if (col < (ncols - 1))
1805 FNaddress = CbmTrdAddress::GetAddress(CbmTrdAddress::GetLayerId(address), CbmTrdAddress::GetModuleId(address),
1806 sec, row, col + 1);
1807 if (fPulseBuffer[FNaddress].first.size() >= (size_t) i + 1)
1808 cross += fPulseBuffer[FNaddress].first[i] * fCrosstalkLevel;
1809 }
1810
1811 return cross;
1812}
1813
1814//_______________________________________________________________________________
1815void CbmTrdModuleSimR::CheckBuffer(Bool_t EB = false)
1816{
1817
1818
1819 std::map<Int_t, Double_t>::iterator timeit;
1820 std::vector<Int_t> toBeErased;
1821
1822 Bool_t done = false;
1823
1824 while (!done) {
1825 done = true;
1826 for (timeit = fTimeBuffer.begin(); timeit != fTimeBuffer.end(); timeit++) {
1827 Int_t add = timeit->first;
1828 if (fCurrentTime < fTimeBuffer[add]) continue;
1829 Double_t dt = fCurrentTime - fTimeBuffer[add];
1830 if ((dt < fCollectTime || dt == fCurrentTime) && !EB) continue;
1831 // if(!fPulseSwitch) {ProcessBuffer(add);fTimeBuffer.erase(add);}
1832 if (!fPulseSwitch) {
1833 ProcessBuffer(add);
1834 toBeErased.push_back(add);
1835 }
1836 if (fPulseSwitch) {
1837 std::vector<Double_t> pulse;
1838 pulse = fPulseBuffer[add].first;
1839
1840 if (CheckTrigger(pulse) == 1 && EB) {
1841 ProcessPulseBuffer(add, false, false, true, true);
1842 break;
1843 }
1844
1845
1846 if (CheckTrigger(pulse) == 1 && !EB) {
1847 ProcessPulseBuffer(add, false, false, true, true);
1848 done = false;
1849 break;
1850 }
1851
1852 if (fPrintPulse) std::cout << std::endl;
1853 }
1854 }
1855 }
1856
1857 for (auto& address : toBeErased) {
1858 fTimeBuffer.erase(address);
1859 }
1860}
1861
1862//_______________________________________________________________________________
1863Int_t CbmTrdModuleSimR::FlushBuffer(ULong64_t time)
1864{
1865 Bool_t closeTS(kFALSE);
1866 if (fTimeSlice && (time - fTimeSlice->GetEndTime()) > -1000.) closeTS = true;
1867
1868 // if(!CbmTrdDigitizer::IsTimeBased()) return 0;
1869 //process channels before timeslice ending and release memory
1870 // if(closeTS && CbmTrdDigitizer::IsTimeBased()){
1871 if (closeTS || time == 0) {
1872 std::map<Int_t, Double_t>::iterator timeit;
1873 Bool_t done = false;
1874
1875 while (!done) {
1876 done = true;
1877 for (timeit = fTimeBuffer.begin(); timeit != fTimeBuffer.end(); timeit++) {
1878 Int_t add = timeit->first;
1879 if (!fPulseSwitch) { ProcessBuffer(add); }
1880 if (fPulseSwitch) {
1881 std::vector<Double_t> pulse;
1882 pulse = fPulseBuffer[add].first;
1883
1884 if (CheckTrigger(pulse) == 1) {
1885 ProcessPulseBuffer(add, false, false, true, true);
1886 done = false;
1887 break;
1888 }
1889
1890 if (fPrintPulse) std::cout << std::endl;
1891 }
1892 }
1893 }
1894 std::map<Int_t, std::pair<std::vector<Double_t>, CbmMatch*>>::iterator itBuffer;
1895 for (itBuffer = fPulseBuffer.begin(); itBuffer != fPulseBuffer.end(); itBuffer++) {
1896 if (fPulseBuffer[itBuffer->first].second) delete fPulseBuffer[itBuffer->first].second;
1897 }
1898 fPulseBuffer.clear();
1899 fTimeBuffer.clear();
1900 fMultiBuffer.clear();
1901 fMCBuffer.clear();
1902 }
1903 return 0;
1904}
1905
1906
1907//_______________________________________________________________________________
1908void CbmTrdModuleSimR::CleanUp(Bool_t EB = false)
1909{
1910
1911 std::map<Int_t, Double_t>::iterator timeit;
1912 // clean up
1913 std::vector<Int_t> erase_list;
1914
1915 if (fPulseSwitch) {
1916 for (timeit = fTimeBuffer.begin(); timeit != fTimeBuffer.end(); timeit++) {
1917 Int_t add = timeit->first;
1918 Double_t dt = fCurrentTime - fTimeBuffer[add];
1919 if (fTimeSlice) {
1920 if (fTimeBuffer[add] < fTimeSlice->GetStartTime()) {
1921 erase_list.push_back(add);
1922 continue;
1923 }
1924 }
1925 if ((dt < fCollectTime || dt == fCurrentTime) && !EB) continue;
1926 erase_list.push_back(add);
1927 }
1928 }
1929
1930 // //release memory for all non-triggered channels after signal collection time
1931 for (UInt_t i = 0; i < erase_list.size(); i++) {
1932 if (fPulseBuffer[erase_list[i]].second) delete fPulseBuffer[erase_list[i]].second;
1933 fPulseBuffer.erase(erase_list[i]);
1934 fTimeBuffer.erase(erase_list[i]);
1935 fMultiBuffer.erase(erase_list[i]);
1936 fMCBuffer.erase(erase_list[i]);
1937 }
1938}
1939
1940//_______________________________________________________________________________
1941void CbmTrdModuleSimR::CheckTime(Int_t address)
1942{
1943
1944 //compare last entry in the actual channel with the current time
1945 std::map<Int_t, Double_t>::iterator timeit;
1946 Double_t dt = fCurrentTime - fTimeBuffer[address];
1947 // std::cout<<" dt: " << dt<<std::endl;
1948 Bool_t go = false;
1949 if (fCurrentTime > fTimeBuffer[address] && dt > 0.0000000) {
1950 if (dt > fCollectTime && dt != fCurrentTime && !fPulseSwitch) {
1951 ProcessBuffer(address);
1952 fTimeBuffer.erase(address);
1953 go = true;
1954 }
1955 // if(dt>fCollectTime && dt!=fCurrentTime && fPulseSwitch) {ProcessPulseBuffer(address,false,false);std::cout<<" ------ " <<std::endl;go=true;}
1956 if (dt > fCollectTime && dt != fCurrentTime && fPulseSwitch) {
1957 //ProcessPulseBuffer(address,false,false,true,true);
1958 go = true;
1959 if (fPrintPulse) std::cout << std::endl;
1960 }
1961 }
1962
1963 if (go && fPulseSwitch) {
1964 CheckBuffer(false);
1965 CleanUp(false);
1966 }
1967 if (go && !fPulseSwitch) CheckBuffer();
1968}
1969
1970//_______________________________________________________________________________
1971void CbmTrdModuleSimR::NoiseTime(ULong64_t eventTime) { fCurrentTime = gRandom->Uniform(fLastEventTime, eventTime); }
1972
1973//_______________________________________________________________________________
1974Double_t CbmTrdModuleSimR::AddDrifttime(Double_t x, Double_t z)
1975{
1976
1977 if (fDebug) fQA->CreateHist("All drift", 300, 0., 300.);
1978 if (fDebug) fQA->Fill("All drift", fDriftTime->GetBinContent(fDriftTime->FindBin(x, z)) * 1000);
1979
1980 return fDriftTime->GetBinContent(fDriftTime->FindBin(x, z));
1981}
1982
1983//_______________________________________________________________________________
1984Double_t CbmTrdModuleSimR::AddDrifttime(Int_t x)
1985{
1986 Double_t drifttime[241] = {
1987 0.11829, 0.11689, 0.11549, 0.11409, 0.11268, 0.11128, 0.10988, 0.10847, 0.10707, 0.10567, 0.10427,
1988 0.10287, 0.10146, 0.10006, 0.09866, 0.09726, 0.095859, 0.094459, 0.09306, 0.091661, 0.090262, 0.088865,
1989 0.087467, 0.086072, 0.084677, 0.083283, 0.08189, 0.080499, 0.07911, 0.077722, 0.076337, 0.074954, 0.073574,
1990 0.072197, 0.070824, 0.069455, 0.06809, 0.066731, 0.065379, 0.064035, 0.0627, 0.061376, 0.060063, 0.058764,
1991 0.05748, 0.056214, 0.054967, 0.053743, 0.052544, 0.051374, 0.05024, 0.049149, 0.048106, 0.047119, 0.046195,
1992 0.045345, 0.044583, 0.043925, 0.043403, 0.043043, 0.042872, 0.042932, 0.043291, 0.044029, 0.045101, 0.04658,
1993 0.048452, 0.050507, 0.052293, 0.053458, 0.054021, 0.053378, 0.052139, 0.53458, 0.050477, 0.048788, 0.047383,
1994 0.046341, 0.045631, 0.045178, 0.045022, 0.045112, 0.045395, 0.045833, 0.046402, 0.047084, 0.047865, 0.048726,
1995 0.049651, 0.050629, 0.051654, 0.052718, 0.053816, 0.054944, 0.056098, 0.057274, 0.058469, 0.059682, 0.060909,
1996 0.062149, 0.0634, 0.064661, 0.06593, 0.067207, 0.06849, 0.069778, 0.07107, 0.072367, 0.073666, 0.074968,
1997 0.076272, 0.077577, 0.078883, 0.080189, 0.081495, 0.082801, 0.084104, 0.085407, 0.086707, 0.088004, 0.089297,
1998 0.090585, 0.091867, 0.093142, 0.094408, 0.095664, 0.096907, 0.098134, 0.099336, 0.10051, 0.10164, 0.10273,
1999 0.10375, 0.10468, 0.10548, 0.10611, 0.10649, 0.10655, 0.10608, 0.10566, 0.1072, 0.10799, 0.10875,
2000 0.11103, 0.11491, 0.11819, 0.12051, 0.12211, 0.12339, 0.12449, 0.12556, 0.12663, 0.12771, 0.12881,
2001 0.12995, 0.13111, 0.13229, 0.13348, 0.13468, 0.13589, 0.13711, 0.13834, 0.13957, 0.1408, 0.14204,
2002 0.14328, 0.14452, 0.14576, 0.14701, 0.14825, 0.1495, 0.15075, 0.152, 0.15325, 0.1545, 0.15576,
2003 0.15701, 0.15826, 0.15952, 0.16077, 0.16203, 0.16328, 0.16454, 0.16579, 0.16705, 0.1683, 0.16956,
2004 0.17082, 0.17207, 0.17333, 0.17458, 0.17584, 0.1771, 0.17835, 0.17961, 0.18087, 0.18212, 0.18338,
2005 0.18464, 0.18589, 0.18715, 0.18841, 0.18966, 0.19092, 0.19218, 0.19343, 0.19469, 0.19595, 0.19721,
2006 0.19846, 0.19972, 0.20098, 0.20223, 0.20349, 0.20475, 0.20601, 0.20726, 0.20852, 0.20978, 0.21103,
2007 0.21229, 0.21355, 0.2148, 0.21606, 0.21732, 0.21857, 0.21983, 0.22109, 0.22234, 0.2236, 0.22486,
2008 0.22612, 0.22737, 0.22863, 0.22989, 0.23114, 0.2324, 0.23366, 0.23491, 0.23617, 0.2363};
2009
2010
2011 // Int_t xindex=0;
2012
2013 return drifttime[Int_t(x)];
2014}
2015
2016
2017//_______________________________________________________________________________
2018Double_t CbmTrdModuleSimR::GetStep(Double_t dist, Int_t roll)
2019{
2020 Double_t prob = 0.;
2021 Int_t steps = 1000;
2022 Double_t CalcGamma = 0.;
2023
2024 std::pair<Double_t, Double_t> bethe[12] = {
2025 std::make_pair(1.5, 1.5), std::make_pair(2, 1.1), std::make_pair(3, 1.025), std::make_pair(4, 1),
2026 std::make_pair(10, 1.1), std::make_pair(20, 1.2), std::make_pair(100, 1.5), std::make_pair(200, 1.6),
2027 std::make_pair(300, 1.65), std::make_pair(400, 1.675), std::make_pair(500, 1.7), std::make_pair(1000, 1.725)};
2028
2029 for (Int_t n = 0; n < 12; n++) {
2030 if (fGamma < bethe[0].first) {
2031 CalcGamma = bethe[0].second;
2032 break;
2033 }
2034 if (n == 11) {
2035 CalcGamma = bethe[11].second;
2036 break;
2037 }
2038
2039 if (fGamma >= bethe[n].first && fGamma <= bethe[n + 1].first) {
2040 Double_t dx = bethe[n + 1].first - bethe[n].first;
2041 Double_t dy = bethe[n + 1].second - bethe[n].second;
2042 Double_t slope = dy / dx;
2043 CalcGamma = (fGamma - bethe[n].first) * slope + bethe[n].second;
2044 break;
2045 }
2046 }
2047
2048 Double_t s = 0;
2049 Double_t D = 1 / (20.5 * CalcGamma);
2050 for (Int_t i = 1; i < steps; i++) {
2051 s = (dist / steps) * i;
2052 prob = (1 - TMath::Exp(-s / D)) * 100;
2053 if (prob >= roll) return s;
2054 }
2055
2056 return s;
2057}
2058
2059std::pair<Int_t, std::vector<Double_t>> CbmTrdModuleSimR::GetTotalSteps(Double_t In[3], Double_t Out[3], Double_t dist)
2060{
2061 Double_t prob = 0.;
2062 Int_t steps = 1000;
2063 Double_t CalcGamma = 0.;
2064 Double_t roll = gRandom->Integer(100);
2065
2066 std::pair<Double_t, Double_t> bethe[12] = {
2067 std::make_pair(1.5, 1.5), std::make_pair(2, 1.1), std::make_pair(3, 1.025), std::make_pair(4, 1),
2068 std::make_pair(10, 1.1), std::make_pair(20, 1.2), std::make_pair(100, 1.5), std::make_pair(200, 1.6),
2069 std::make_pair(300, 1.65), std::make_pair(400, 1.675), std::make_pair(500, 1.7), std::make_pair(1000, 1.725)};
2070
2071 for (Int_t n = 0; n < 12; n++) {
2072 if (fGamma < bethe[0].first) {
2073 CalcGamma = bethe[0].second;
2074 break;
2075 }
2076 if (n == 11) {
2077 CalcGamma = bethe[11].second;
2078 break;
2079 }
2080
2081 if (fGamma >= bethe[n].first && fGamma <= bethe[n + 1].first) {
2082 Double_t dx = bethe[n + 1].first - bethe[n].first;
2083 Double_t dy = bethe[n + 1].second - bethe[n].second;
2084 Double_t slope = dy / dx;
2085 CalcGamma = (fGamma - bethe[n].first) * slope + bethe[n].second;
2086 break;
2087 }
2088 }
2089
2090 Double_t pos[3] = {In[0], In[1], In[2]};
2091 Double_t lastpos[3] = {In[0], In[1], In[2]};
2092 Int_t pointcount = 0;
2093 std::vector<Double_t> posvec;
2094 Double_t D = 1 / (20.5 * CalcGamma);
2095 Double_t delta[3];
2096 Double_t s = 0;
2097 for (Int_t i = 0; i < 3; i++)
2098 delta[i] = (Out[i] - In[i]);
2099
2100 roll = gRandom->Integer(100);
2101 for (Int_t i = 1; i < steps; i++) {
2102 s = (dist / steps) * i;
2103 prob = (1 - TMath::Exp(-s / D)) * 100;
2104 if (prob >= roll) {
2105 Double_t move = 2 * (s / dist);
2106 for (Int_t n = 0; n < 3; n++)
2107 pos[n] = lastpos[n] + move * delta[n];
2108 for (Int_t n = 0; n < 3; n++)
2109 lastpos[n] = pos[n];
2110 // Double_t r = TMath::Sqrt((In[0]-pos[0])*(In[0]-pos[0])+(In[1]-pos[1])*(In[1]-pos[1]));
2111 if (TMath::Abs(pos[0]) < fDigiPar->GetSizeX() && TMath::Abs(pos[1]) < fDigiPar->GetSizeY() && (pos[2]) < Out[2]) {
2112 posvec.push_back(pos[0]);
2113 posvec.push_back(pos[1]);
2114 posvec.push_back(pos[2]);
2115 pointcount++;
2116 i = 1;
2117 roll = gRandom->Integer(100);
2118 continue;
2119 }
2120 }
2121 if (TMath::Abs(pos[0]) < fDigiPar->GetSizeX() && TMath::Abs(pos[1]) < fDigiPar->GetSizeY() && (pos[2]) < Out[2]) {
2122 continue;
2123 }
2124 break;
2125 }
2126
2127 return std::make_pair(pointcount, posvec);
2128}
2129
2130
2131ClassImp(CbmTrdModuleSimR)
points
TClonesArray * points
Definition Analyze_matching.h:22
ClassImp
ClassImp(CbmConverterManager)
CbmTrdAddress.h
Helper class to convert unique channel ID back and forth.
CbmTrdDigitizer.h
TRD digitizer. Updated 24/04/2013 by Andrey Lebedev andrey.lebedev@gsi.de Updated 4/06/2018 by Alex B...
sqrt
friend fvec sqrt(const fvec &a)
Definition KfSimdPseudo.h:118
size
static constexpr size_t size()
Definition KfSimdPseudo.h:2
v
fscal v[fmask::Size]
Definition KfSimdPseudo.h:4
first
bool first
Definition LKFMinuit.cxx:149
h
Generates beam ions for transport simulation.
CbmMatch::GetNofLinks
int32_t GetNofLinks() const
Definition CbmMatch.h:42
CbmMatch::AddLink
void AddLink(const CbmLink &newLink)
Definition CbmMatch.cxx:47
CbmTimeSlice
Bookkeeping of time-slice content.
Definition CbmTimeSlice.h:33
CbmTimeSlice::GetEndTime
double GetEndTime() const
Definition CbmTimeSlice.cxx:101
CbmTimeSlice::GetStartTime
double GetStartTime() const
Definition CbmTimeSlice.h:114
CbmTrdAddress::GetModuleId
static uint32_t GetModuleId(uint32_t address)
Return module ID from address.
Definition CbmTrdAddress.h:79
CbmTrdAddress::GetSectorId
static uint32_t GetSectorId(uint32_t address)
Return sector ID from address.
Definition CbmTrdAddress.h:89
CbmTrdAddress::GetColumnId
static uint32_t GetColumnId(uint32_t address)
Return column ID from address.
Definition CbmTrdAddress.h:106
CbmTrdAddress::GetLayerId
static uint32_t GetLayerId(uint32_t address)
Return layer ID from address.
Definition CbmTrdAddress.h:69
CbmTrdAddress::GetAddress
static uint32_t GetAddress(int32_t layerId, int32_t moduleId, int32_t sectorId, int32_t rowId, int32_t columnId)
Return address from system ID, layer, module, sector, column and row IDs.
Definition CbmTrdAddress.h:44
CbmTrdAddress::GetRowId
static uint32_t GetRowId(uint32_t address)
Return row ID from address.
Definition CbmTrdAddress.h:99
CbmTrdCheckUtil::CreateHist
void CreateHist(std::string name, Int_t xbins, Double_t xlow, Double_t xhigh, Int_t ybins=0, Double_t ylow=1., Double_t yhigh=1.)
Definition CbmTrdCheckUtil.cxx:30
CbmTrdCheckUtil::Fill
void Fill(std::string name, Double_t x, Double_t y=9999.)
Definition CbmTrdCheckUtil.cxx:65
CbmTrdCheckUtil::FillProfile
void FillProfile(std::string name, Double_t x, Double_t y, Double_t z=9999.)
Definition CbmTrdCheckUtil.cxx:88
CbmTrdCheckUtil::CreateProfile
void CreateProfile(std::string name, Int_t xbins, Double_t xlow, Double_t xhigh, Int_t ybins=0, Double_t ylow=1., Double_t yhigh=1.)
Definition CbmTrdCheckUtil.cxx:47
CbmTrdDigi::SetTriggerType
void SetTriggerType(const eTriggerType triggerType)
Set digi trigger type.
Definition CbmTrdDigi.cxx:280
CbmTrdDigi::SetFlag
void SetFlag(const int32_t iflag, bool set=true)
Generic flag status setter.
Definition CbmTrdDigi.cxx:246
CbmTrdDigi::Clk
static float Clk(eCbmTrdAsicType ty)
DAQ clock accessor for each ASIC.
Definition CbmTrdDigi.h:109
CbmTrdDigi::SetErrorClass
void SetErrorClass(const int32_t n)
Set digi error class (SPADIC only)
Definition CbmTrdDigi.h:231
CbmTrdDigi::GetTime
double GetTime() const
Getter for physical time [ns]. Accounts for clock representation of each ASIC. In SPADIC case physica...
Definition CbmTrdDigi.h:153
CbmTrdDigi::GetCharge
double GetCharge() const
Common purpose charge getter.
Definition CbmTrdDigi.cxx:150
CbmTrdDigi::SetCharge
void SetCharge(float c)
Charge setter for SPADIC ASIC.
Definition CbmTrdDigi.cxx:239
CbmTrdDigitizer::AddNoise
static Bool_t AddNoise()
Definition CbmTrdDigitizer.h:69
CbmTrdDigitizer::UseWeightedDist
static Bool_t UseWeightedDist()
Definition CbmTrdDigitizer.h:71
CbmTrdDigitizer::IsTimeBased
static Bool_t IsTimeBased()
Definition CbmTrdDigitizer.h:70
CbmTrdModuleAbstract::fLayerId
Char_t fLayerId
layer identifier
Definition CbmTrdModuleAbstract.h:89
CbmTrdModuleAbstract::fDigiPar
const CbmTrdParModDigi * fDigiPar
read-out description of module
Definition CbmTrdModuleAbstract.h:93
CbmTrdModuleAbstract::GetDx
virtual Double_t GetDx() const
Shortcut getter size x/2 [cm].
Definition CbmTrdModuleAbstract.h:44
CbmTrdModuleAbstract::GetPath
virtual const Char_t * GetPath() const
Definition CbmTrdModuleAbstract.h:72
CbmTrdModuleAbstract::fModAddress
UShort_t fModAddress
unique identifier for current module
Definition CbmTrdModuleAbstract.h:88
CbmTrdModuleAbstract::fAsicPar
CbmTrdParModAsic * fAsicPar
the set of ASIC operating on the module (owned)
Definition CbmTrdModuleAbstract.h:95
CbmTrdModuleAbstract::GetDy
virtual Double_t GetDy() const
Shortcut getter size y/2 [cm].
Definition CbmTrdModuleAbstract.h:46
CbmTrdModuleSimR
Simulation module implementation for rectangular pad geometry.
Definition CbmTrdModuleSimR.h:23
CbmTrdModuleSimR::SetPulsePars
void SetPulsePars(Int_t mode)
Definition CbmTrdModuleSimR.cxx:1721
CbmTrdModuleSimR::GetTotalSteps
std::pair< Int_t, std::vector< Double_t > > GetTotalSteps(Double_t In[3], Double_t Out[3], Double_t dist)
Definition CbmTrdModuleSimR.cxx:2059
CbmTrdModuleSimR::FlushBuffer
Int_t FlushBuffer(ULong64_t time=0)
Flush local digi buffer.
Definition CbmTrdModuleSimR.cxx:1863
CbmTrdModuleSimR::SetSpadicResponse
void SetSpadicResponse(Double_t calibration, Double_t tau)
Definition CbmTrdModuleSimR.cxx:1709
CbmTrdModuleSimR::SetNoiseLevel
void SetNoiseLevel(Double_t sigma_keV)
Definition CbmTrdModuleSimR.cxx:1703
CbmTrdModuleSimR::MakeDigi
Bool_t MakeDigi(CbmTrdPoint *p, Double_t time, Bool_t TR)
Steering routine for converting MC point to digits.
Definition CbmTrdModuleSimR.cxx:1112
CbmTrdModuleSimR::AddDigi
void AddDigi(Int_t address, Double_t charge, Double_t chargeTR, Double_t time, Int_t trigger)
Definition CbmTrdModuleSimR.cxx:121
CbmTrdModuleSimR::AddNoise
Double_t AddNoise(Double_t charge)
Definition CbmTrdModuleSimR.cxx:1748
CbmTrdModuleSimR::CleanUp
void CleanUp(Bool_t EB)
Definition CbmTrdModuleSimR.cxx:1908
CbmTrdModuleSimR::CalcPRF
Double_t CalcPRF(Double_t x, Double_t W, Double_t h)
Definition CbmTrdModuleSimR.cxx:1042
CbmTrdModuleSimR::CalcResponse
Double_t CalcResponse(Double_t t)
Definition CbmTrdModuleSimR.cxx:1053
CbmTrdModuleSimR::AddDigitoBuffer
void AddDigitoBuffer(Int_t address, Double_t charge, Double_t chargeTR, Double_t time, Int_t trigger)
Definition CbmTrdModuleSimR.cxx:453
CbmTrdModuleSimR::ProcessBuffer
void ProcessBuffer(Int_t address)
Definition CbmTrdModuleSimR.cxx:173
CbmTrdModuleSimR::fMCBuffer
std::map< Int_t, std::vector< std::vector< Int_t > > > fMCBuffer
Definition CbmTrdModuleSimR.h:163
CbmTrdModuleSimR::NoiseTime
void NoiseTime(ULong64_t eventTime)
Definition CbmTrdModuleSimR.cxx:1971
CbmTrdModuleSimR::fQA
CbmTrdCheckUtil * fQA
Definition CbmTrdModuleSimR.h:169
CbmTrdModuleSimR::ScanPadPlane
void ScanPadPlane(const Double_t *local_point, Double_t reldrift, Double_t clusterELoss, Double_t clusterELossTR, Int_t epoints, Int_t ipoint)
Definition CbmTrdModuleSimR.cxx:1345
CbmTrdModuleSimR::AddCorrelatedNoise
std::vector< Double_t > AddCorrelatedNoise(std::vector< Double_t > pulse)
Definition CbmTrdModuleSimR.cxx:1776
CbmTrdModuleSimR::fMessageConverter
CbmTrdRawToDigiR * fMessageConverter
Definition CbmTrdModuleSimR.h:167
CbmTrdModuleSimR::CbmTrdModuleSimR
CbmTrdModuleSimR(Int_t mod, Int_t ly, Int_t rot)
Definition CbmTrdModuleSimR.cxx:43
CbmTrdModuleSimR::fAnalogBuffer
std::map< Int_t, std::vector< std::pair< CbmTrdDigi *, CbmMatch * > > > fAnalogBuffer
Definition CbmTrdModuleSimR.h:156
CbmTrdModuleSimR::fTimeSlice
CbmTimeSlice * fTimeSlice
link to CBM time slice
Definition CbmTrdModuleSimR.h:146
CbmTrdModuleSimR::CheckTrigger
Int_t CheckTrigger(std::vector< Double_t > pulse)
Definition CbmTrdModuleSimR.cxx:1002
CbmTrdModuleSimR::AddDigitoPulseBuffer
void AddDigitoPulseBuffer(Int_t address, Double_t reldrift, Double_t charge, Double_t chargeTR, Double_t time, Int_t trigger, Int_t epoints, Int_t ipoint)
Definition CbmTrdModuleSimR.cxx:502
CbmTrdModuleSimR::fShiftQA
std::map< Int_t, Double_t > fShiftQA
Definition CbmTrdModuleSimR.h:160
CbmTrdModuleSimR::fMCQA
std::map< Int_t, Double_t > fMCQA
Definition CbmTrdModuleSimR.h:162
CbmTrdModuleSimR::AddCrosstalk
Int_t AddCrosstalk(Double_t address, Int_t i, Int_t sec, Int_t row, Int_t col, Int_t ncols)
Definition CbmTrdModuleSimR.cxx:1791
CbmTrdModuleSimR::SetAsicPar
void SetAsicPar(CbmTrdParModAsic *p=NULL)
Definition CbmTrdModuleSimR.cxx:1595
CbmTrdModuleSimR::MakePulse
std::vector< Double_t > MakePulse(Double_t charge, std::vector< Double_t > pulse, Int_t address)
Definition CbmTrdModuleSimR.cxx:585
CbmTrdModuleSimR::AddDrifttime
Double_t AddDrifttime(Double_t x, Double_t z)
Definition CbmTrdModuleSimR.cxx:1974
CbmTrdModuleSimR::fLinkQA
std::map< Int_t, std::vector< std::map< TString, Int_t > > > fLinkQA
Definition CbmTrdModuleSimR.h:161
CbmTrdModuleSimR::CheckTime
void CheckTime(Int_t address)
Definition CbmTrdModuleSimR.cxx:1941
CbmTrdModuleSimR::GetMultiBin
Int_t GetMultiBin(std::vector< Double_t > pulse)
Definition CbmTrdModuleSimR.cxx:1023
CbmTrdModuleSimR::CheckMulti
Bool_t CheckMulti(Int_t address, std::vector< Double_t > pulse)
Definition CbmTrdModuleSimR.cxx:704
CbmTrdModuleSimR::fMultiBuffer
std::map< Int_t, std::pair< Double_t, Int_t > > fMultiBuffer
Definition CbmTrdModuleSimR.h:158
CbmTrdModuleSimR::fPulseBuffer
std::map< Int_t, std::pair< std::vector< Double_t >, CbmMatch * > > fPulseBuffer
Definition CbmTrdModuleSimR.h:157
CbmTrdModuleSimR::fTimeBuffer
std::map< Int_t, Double_t > fTimeBuffer
Definition CbmTrdModuleSimR.h:159
CbmTrdModuleSimR::ProcessPulseBuffer
void ProcessPulseBuffer(Int_t address, Bool_t FNcall, Bool_t MultiCall, Bool_t down, Bool_t up)
Definition CbmTrdModuleSimR.cxx:203
CbmTrdModuleSimR::SetPadPlaneScanArea
void SetPadPlaneScanArea(Int_t row)
Definition CbmTrdModuleSimR.cxx:1740
CbmTrdModuleSimR::SetDistributionPoints
void SetDistributionPoints(Int_t points)
Definition CbmTrdModuleSimR.cxx:1706
CbmTrdModuleSimR::SetPulseMode
void SetPulseMode(Bool_t pulsed)
Definition CbmTrdModuleSimR.cxx:1736
CbmTrdModuleSimR::GetStep
Double_t GetStep(Double_t dist, Int_t roll)
Definition CbmTrdModuleSimR.cxx:2018
CbmTrdModuleSimR::AddToPulse
void AddToPulse(Int_t address, Double_t charge, Double_t reldrift, std::vector< Double_t > pulse)
Definition CbmTrdModuleSimR.cxx:622
CbmTrdModuleSimR::CheckBuffer
void CheckBuffer(Bool_t EB)
Definition CbmTrdModuleSimR.cxx:1815
CbmTrdModuleSimR::DistributeCharge
Bool_t DistributeCharge(Double_t pointin[3], Double_t pointout[3], Double_t delta[3], Double_t pos[3], Int_t ipoints)
Definition CbmTrdModuleSimR.cxx:1062
CbmTrdModuleSim
Abstract class for module wise digitization and raw format producing.
Definition CbmTrdModuleSim.h:24
CbmTrdModuleSim::fRadiator
std::shared_ptr< CbmTrdRadiator > fRadiator
Pointer to digitizer.
Definition CbmTrdModuleSim.h:102
CbmTrdModuleSim::fInputId
Int_t fInputId
MC input file number.
Definition CbmTrdModuleSim.h:95
CbmTrdModuleSim::fPointId
Int_t fPointId
MC point id being processed.
Definition CbmTrdModuleSim.h:93
CbmTrdModuleSim::fDigiMap
std::map< Int_t, std::pair< CbmTrdDigi *, CbmMatch * > > fDigiMap
Temporary storage for complete digis for each CBM address.
Definition CbmTrdModuleSim.h:105
CbmTrdModuleSim::SetPositionMC
virtual void SetPositionMC(Double_t pos[3])
Definition CbmTrdModuleSim.h:79
CbmTrdModuleSim::fEventId
Int_t fEventId
MC event id being processed.
Definition CbmTrdModuleSim.h:94
CbmTrdModuleSim::fXYZ
Double_t fXYZ[3]
MC position of the point in module coordinates.
Definition CbmTrdModuleSim.h:97
CbmTrdParAsic::SetChannelAddress
virtual void SetChannelAddress(Int_t address)
Definition CbmTrdParAsic.cxx:41
CbmTrdParModAsic
Describe TRD module ASIC settings (electronic gain, delays, etc)
Definition CbmTrdParModAsic.h:30
CbmTrdParModAsic::GetAsicAddress
virtual Int_t GetAsicAddress(Int_t chAddress) const
Look for the ASIC which operates on a specific channel It applies to the list of ASICs.
Definition CbmTrdParModAsic.cxx:68
CbmTrdParModAsic::SetAsicPar
virtual void SetAsicPar(CbmTrdParAsic *p)
Initialize the ASIC parameters for DAQ id It applies to the list of ASICs.
Definition CbmTrdParModAsic.cxx:169
CbmTrdParModDigi::GetSectorBeginY
Double_t GetSectorBeginY(Int_t i) const
Definition CbmTrdParModDigi.h:35
CbmTrdParModDigi::Print
void Print(Option_t *opt="") const
Definition CbmTrdParModDigi.cxx:151
CbmTrdParModDigi::GetPadInfo
Bool_t GetPadInfo(const Double_t *local_point, Int_t &sectorId, Int_t &columnId, Int_t &rowId) const
Definition CbmTrdParModDigi.cxx:412
CbmTrdParModDigi::GetNofSectors
Int_t GetNofSectors() const
Definition CbmTrdParModDigi.h:43
CbmTrdParModDigi::GetSizeY
Double_t GetSizeY() const
Definition CbmTrdParModDigi.h:84
CbmTrdParModDigi::GetNofColumnsInSector
Int_t GetNofColumnsInSector(Int_t i) const
Definition CbmTrdParModDigi.cxx:336
CbmTrdParModDigi::GetNofRows
Int_t GetNofRows() const
Definition CbmTrdParModDigi.cxx:317
CbmTrdParModDigi::GetPadPosition
void GetPadPosition(const Int_t sector, const Int_t col, const Int_t row, TVector3 &padPos, TVector3 &padPosErr) const
Definition CbmTrdParModDigi.cxx:672
CbmTrdParModDigi::GetAnodeWireSpacing
Double_t GetAnodeWireSpacing() const
Definition CbmTrdParModDigi.h:41
CbmTrdParModDigi::GetPadSizeY
Double_t GetPadSizeY(Int_t i) const
Definition CbmTrdParModDigi.h:33
CbmTrdParModDigi::GetNofRowsInSector
Int_t GetNofRowsInSector(Int_t i) const
Definition CbmTrdParModDigi.cxx:343
CbmTrdParModDigi::GetSectorBeginX
Double_t GetSectorBeginX(Int_t i) const
Definition CbmTrdParModDigi.h:34
CbmTrdParModDigi::GetSector
Int_t GetSector(Int_t npady, Int_t &rowId) const
Definition CbmTrdParModDigi.cxx:257
CbmTrdParModDigi::ProjectPositionToNextAnodeWire
void ProjectPositionToNextAnodeWire(Double_t *local_point) const
Definition CbmTrdParModDigi.cxx:176
CbmTrdParModDigi::GetPadSizeX
Double_t GetPadSizeX(Int_t i) const
Definition CbmTrdParModDigi.h:32
CbmTrdParModDigi::GetAnodeWireToPadPlaneDistance
Double_t GetAnodeWireToPadPlaneDistance() const
Definition CbmTrdParModDigi.h:39
CbmTrdParModDigi::GetNofColumns
Int_t GetNofColumns() const
Definition CbmTrdParModDigi.cxx:298
CbmTrdParModDigi::GetSizeX
Double_t GetSizeX() const
Definition CbmTrdParModDigi.h:83
CbmTrdParModDigi::TransformToLocalPad
void TransformToLocalPad(const Double_t *local_point, Double_t &posX, Double_t &posY) const
Definition CbmTrdParModDigi.cxx:590
CbmTrdParSpadic
Definition of SPADIC parameters.
Definition CbmTrdParSpadic.h:22
CbmTrdPoint::GetYOut
double GetYOut() const
Definition CbmTrdPoint.h:69
CbmTrdPoint::GetXIn
double GetXIn() const
Definition CbmTrdPoint.h:65
CbmTrdPoint::GetZIn
double GetZIn() const
Definition CbmTrdPoint.h:67
CbmTrdPoint::GetXOut
double GetXOut() const
Definition CbmTrdPoint.h:68
CbmTrdPoint::GetYIn
double GetYIn() const
Definition CbmTrdPoint.h:66
CbmTrdPoint::GetZOut
double GetZOut() const
Definition CbmTrdPoint.h:70
CbmTrdRawToDigiR::MakeDigi
CbmTrdDigi * MakeDigi(std::vector< Int_t > samples, Int_t channel, Int_t uniqueModuleId, ULong64_t time, Bool_t FN=false)
Definition CbmTrdRawToDigiR.cxx:412
CbmTrdRawToDigiR::SetPresamples
void SetPresamples(Int_t pre)
Definition CbmTrdRawToDigiR.h:60
CbmTrdRawToDigiR::SetMinBin
void SetMinBin(Int_t bin)
Definition CbmTrdRawToDigiR.h:59
CbmTrdRawToDigiR::SetQA
void SetQA(CbmTrdCheckUtil *qa)
Definition CbmTrdRawToDigiR.h:68
CbmTrdRawToDigiR::SetReadFile
void SetReadFile(std::string file)
Definition CbmTrdRawToDigiR.h:56
CbmTrdRawToDigiR::SetRecoMask
void SetRecoMask(std::vector< Int_t > mask)
Definition CbmTrdRawToDigiR.h:54
CbmTrdRawToDigiR::SetCalibration
void SetCalibration(Double_t cal)
Definition CbmTrdRawToDigiR.h:50
CbmTrdRawToDigiR::SetSetter
void SetSetter(Bool_t set)
Definition CbmTrdRawToDigiR.h:55
CbmTrdRawToDigiR::SetMaxBin
void SetMaxBin(Int_t bin)
Definition CbmTrdRawToDigiR.h:58
CbmTrdRawToDigiR::GetCharge
Double_t GetCharge(std::vector< Int_t > samples, Int_t shift=-1)
Definition CbmTrdRawToDigiR.cxx:462
CbmTrdRawToDigiR::SetTau
void SetTau(Double_t tau)
Definition CbmTrdRawToDigiR.h:51
CbmTrdRawToDigiR::GetTimeShift
Float_t GetTimeShift(std::vector< Int_t > samples)
Definition CbmTrdRawToDigiR.cxx:444
CbmTrdRawToDigiR::SetShapingOrder
void SetShapingOrder(Int_t order)
Definition CbmTrdRawToDigiR.h:53