LDMX Software
HgcrocEmulator.cxx
1
2#include "Tools/HgcrocEmulator.h"
3
4namespace ldmx {
5
7 // settings of readout chip that are the same for all chips
8 // used in actual digitization
9 noise_ = ps.getParameter<bool>("noise");
10 timingJitter_ = ps.getParameter<double>("timingJitter");
11 rateUpSlope_ = ps.getParameter<double>("rateUpSlope");
12 timeUpSlope_ = ps.getParameter<double>("timeUpSlope");
13 rateDnSlope_ = ps.getParameter<double>("rateDnSlope");
14 timeDnSlope_ = ps.getParameter<double>("timeDnSlope");
15 timePeak_ = ps.getParameter<double>("timePeak");
16 clockCycle_ = ps.getParameter<double>("clockCycle");
17 nADCs_ = ps.getParameter<int>("nADCs");
18 iSOI_ = ps.getParameter<int>("iSOI");
19
20 // Time -> clock counts conversion
21 // time [ns] * ( 2^10 / max time in ns ) = clock counts
22 ns_ = 1024. / clockCycle_;
23
24 hit_merge_ns_ = 0.05; // combine at 50 ps level
25
26 // Configure the pulse shape function
28 TF1("pulseFunc",
29 "[0]*((1.0+exp([1]*(-[2]+[3])))*(1.0+exp([5]*(-[6]+[3]))))/"
30 "((1.0+exp([1]*(x-[2]+[3]-[4])))*(1.0+exp([5]*(x-[6]+[3]-[4]))))",
31 0.0, (double)nADCs_ * clockCycle_);
32 pulseFunc_.FixParameter(0, 1.0); // amplitude is set externally
33 pulseFunc_.FixParameter(1, rateUpSlope_);
34 pulseFunc_.FixParameter(2, timeUpSlope_);
35 pulseFunc_.FixParameter(3, timePeak_);
36 pulseFunc_.FixParameter(4, 0); // not using time offset in this way
37 pulseFunc_.FixParameter(5, rateDnSlope_);
38 pulseFunc_.FixParameter(6, timeDnSlope_);
39}
40
41void HgcrocEmulator::seedGenerator(uint64_t seed) {
42 noiseInjector_ = std::make_unique<TRandom3>(seed);
43}
44
46 const int &channelID,
47 std::vector<std::pair<double, double>> &arriving_pulses,
48 std::vector<ldmx::HgcrocDigiCollection::Sample> &digiToAdd) const {
49 // step 0: prepare ourselves for emulation
50 digiToAdd.clear(); // make sure it is clean
51
52 // Configure chip settings based off of table (that may have been passed)
53 double totMax = getCondition(channelID, "TOT_MAX");
54 double padCapacitance = getCondition(channelID, "PAD_CAPACITANCE");
55 double gain = this->gain(channelID);
56 double pedestal = this->pedestal(channelID);
57 double toaThreshold = getCondition(channelID, "TOA_THRESHOLD");
58 double totThreshold = getCondition(channelID, "TOT_THRESHOLD");
59 // measTime defines the point in the BX where an in-time
60 // (time=0 in times vector) hit would arrive.
61 // Used to determine BX boundaries and TOA behavior.
62 double measTime = getCondition(channelID, "MEAS_TIME");
63 double drainRate = getCondition(channelID, "DRAIN_RATE");
64 double readoutThresholdFloat = this->readoutThreshold(channelID);
65 int readoutThreshold = int(readoutThresholdFloat);
66
67 // sort by amplitude
68 // ==> makes sure that puleses are merged towards higher ones
69 std::sort(
70 arriving_pulses.begin(), arriving_pulses.end(),
71 [](const std::pair<double, double> &a,
72 const std::pair<double, double> &b) { return a.first > b.first; });
73
74 // step 1: gather voltages into groups separated by (programmable) ns, single
75 // pass
77
78 for (auto hit : arriving_pulses) pulse.addOrMerge(hit, hit_merge_ns_);
79
80 // TODO step 2: add timing jitter
81 // if (noise_) pulse.jitter();
82
84
85 // step 3: go through each BX sample one by one
86 bool wasTOA = false;
87 for (int iADC = 0; iADC < nADCs_; iADC++) {
88 double startBX = (iADC - iSOI_) * clockCycle_ - measTime;
89 ldmx_log(trace) << " iADC = " << iADC << " at startBX = " << startBX;
90
91 // step 3b: check each merged hit to see if it peaks in this BX. If so,
92 // check its peak time to see if it's over TOT or TOA.
93 bool startTOT = false;
94 bool overTOA = false;
95 double toverTOA = -1;
96 double toverTOT = -1;
97 for (auto hit : pulse.hits()) {
98 int hitBX = int((hit.second + measTime) / clockCycle_ + iSOI_);
99 // if this hit wasn't in the current BX, continue...
100 if (hitBX != iADC) {
101 continue;
102 }
103
104 double vpeak = pulse(hit.second);
105
106 if (vpeak > totThreshold) {
107 startTOT = true;
108 // use the latest time in the window
109 if (toverTOT < hit.second) {
110 toverTOT = hit.second;
111 }
112 }
113
114 if (vpeak > toaThreshold) {
115 if (!overTOA || hit.second < toverTOA) toverTOA = hit.second;
116 overTOA = true;
117 }
118
119 } // loop over sim hits
120
121 // check for the case of a TOA even though the peak is in the next BX
122 if (!overTOA && pulse(startBX + clockCycle_) > toaThreshold) {
123 if (pulse(startBX) < toaThreshold) {
124 // pulse crossed TOA threshold somewhere between the start of this
125 // basket and the end
126 overTOA = true;
127 toverTOA = startBX + clockCycle_;
128 }
129 }
130
131 if (startTOT) {
132 // above TOT threshold -> do TOT readout mode
133
134 // @TODO NO NOISE
135 // CompositePulse includes pedestal, we need to remove it
136 // when calculating the charge deposited.
137 double charge_deposited =
138 (pulse(toverTOT) - gain * pedestal) * padCapacitance;
139
140 // Measure Time Over Threshold (TOT) by using the drain rate.
141 // 1. Use drain rate to see how long it takes for the charge to drain off
142 // 2. Translate this into DIGI samples
143
144 // Assume linear drain with slope drain rate:
145 // y-intercept = pulse amplitude
146 // slope = drain rate
147 // ==> x-intercept = amplitude / rate
148 // actual time over threshold using the real signal voltage amplitude
149 double tot = charge_deposited / drainRate;
150 ldmx_log(trace) << " we are in TOT read-out mode, TOT = " << tot;
151
152 // calculate the TDC counts for this tot measurement
153 // internally, the chip uses 12 bits (2^12 = 4096)
154 // to measure a maximum of tot Max [ns]
155 int tdc_counts = int(tot * 4096 / totMax) + pedestal;
156
157 // were we already over TOA? TOT is reported in BX where TOA went over
158 // threshold...
159 int toa{0};
160 if (wasTOA) {
161 // TOA was in the past
162 toa = digiToAdd.back().toa();
163 } else {
164 // TOA is here and we need to find it
165 double timecross = pulse.findCrossing(startBX, toverTOT, toaThreshold);
166 toa = int((timecross - startBX) * ns_);
167 // keep inside valid limits
168 if (toa == 0) toa = 1;
169 if (toa > 1023) toa = 1023;
170 }
171 ldmx_log(trace) << " Adding TOT hit with toa = " << toa
172 << ", tdc_counts = " << tdc_counts
173 << " adc_t at prev iADC = "
174 << digiToAdd.at(iADC - 1).adc_t();
175 // ADC at t-1
176 auto adc_at_tminus1 =
177 (iADC > 0) ? digiToAdd.at(iADC - 1).adc_t() : pedestal;
178 auto i_tot_sample = digiToAdd.size();
179 // mark as a TOT measurement with 2nd boolean as true
180 digiToAdd.emplace_back(false, true, adc_at_tminus1, tdc_counts, toa);
181
182 // TODO: properly handle saturation and recovery, eventually.
183 // Now just kill everything...
184 ldmx_log(trace) << " Adding further hits with ADC [t-1] = 0x3FF, toa = "
185 "0x3FF, until digiToAdd.size() = "
186 << digiToAdd.size() << " < nADCs_(" << nADCs_ << ")";
187 while (digiToAdd.size() < nADCs_) {
188 // flags to mark type of sample
189 digiToAdd.emplace_back(true, false, 0x3FF, 0x3FF, 0);
190 }
191 // Read out if the toa is within one Bx after nominal
192 return (i_tot_sample <= iSOI_ + 1);
193 } else {
194 // determine the voltage at the sampling time
195 double bxvolts = pulse((iADC - iSOI_) * clockCycle_);
196 // add noise if requested
197 if (noise_) bxvolts += noise(channelID);
198 // convert to integer and keep in range (handle low and high saturation)
199 int adc = bxvolts / gain;
200 ldmx_log(trace) << " we are in ADC read-out mode, adc = " << adc;
201 if (adc < 0) adc = 0;
202 if (adc > 1023) adc = 1023;
203
204 // check for TOA
205 int toa(0);
206 if (pulse(startBX) < toaThreshold && overTOA) {
207 double timecross = pulse.findCrossing(startBX, toverTOA, toaThreshold);
208 toa = int((timecross - startBX) * ns_);
209 // keep inside valid limits
210 if (toa == 0) toa = 1;
211 if (toa > 1023) toa = 1023;
212 wasTOA = true;
213 } else {
214 wasTOA = false;
215 }
216 // ADC at t-1
217 auto adc_t_minus1 =
218 (iADC > 0) ? digiToAdd.at(iADC - 1).adc_t() : pedestal;
219
220 digiToAdd.emplace_back(false, false, adc_t_minus1, adc, toa);
221 } // TOT or ADC Mode
222 } // sampling baskets
223
224 // we only get here if we never went into TOT mode
225 // check the SOI to see if we should read out
226 ldmx_log(trace) << " we are adding the hit IFF iSOI= " << iSOI_
227 << "'s adc_t = " << digiToAdd.at(iSOI_).adc_t()
228 << " >= thresh (" << readoutThreshold << ")";
229 return digiToAdd.at(iSOI_).adc_t() >= readoutThreshold;
230} // HgcrocEmulator::digitize
231
232std::vector<ldmx::HgcrocDigiCollection::Sample> HgcrocEmulator::noiseDigi(
233 const int &channel, const double &soi_amplitude) const {
234 // get chip conditions from emulator
235 double pedestal{this->pedestal(channel)};
236 double gain{this->gain(channel)};
237 // fill a digi with noise samples
238 std::vector<ldmx::HgcrocDigiCollection::Sample> noise_digi;
239 for (int iADC{0}; iADC < nADCs_; iADC++) {
240 // gen noise for ADC samples
241 // ADC at t-1
242 int adc_tm1{static_cast<int>(pedestal)};
243 if (iADC > 0) {
244 adc_tm1 = noise_digi.at(iADC - 1).adc_t();
245 } else {
246 adc_tm1 += noise(channel) / gain;
247 }
248 // ADC at t
249 int adc_t{static_cast<int>(pedestal + noise(channel) / gain)};
250
251 if (iADC == iSOI_) adc_t += soi_amplitude / gain;
252
253 // set toa to 0 (not determined)
254 // put new sample into noise digi
255 noise_digi.emplace_back(false, false, adc_tm1, adc_t, 0);
256 } // samples in noise digi
257 return noise_digi;
258}
259
260} // namespace ldmx
Class encapsulating parameters for configuring a processor.
Definition Parameters.h:29
const std::vector< std::pair< double, double > > & hits() const
Get list of individual pulses that are entering the chip.
double findCrossing(double low, double high, double level, double prec=0.01)
Find the time at which we cross the input level.
void addOrMerge(const std::pair< double, double > &hit, double hit_merge_ns)
Put another hit into this composite pulse.
double readoutThreshold(const int &id) const
Readout Threshold (ADC Counts)
double rateUpSlope_
Rate of Up Slope in Pulse Shape [1/ns].
void seedGenerator(uint64_t seed)
Seed the emulator for random number generation.
double noise(const int &channelID) const
Get random noise amplitdue for input channel [mV].
double pedestal(const int &id) const
Pedestal [ADC Counts] for input channel.
double getCondition(int id, const std::string &name) const
Get condition for input chip ID, condition name, and default value.
double hit_merge_ns_
Hit merging time [ns].
double timeUpSlope_
Time of Up Slope relative to Pulse Shape Fit [ns].
double timePeak_
Time of Peak relative to pulse shape fit [ns].
int iSOI_
Index for the Sample Of Interest in the list of digi samples.
bool digitize(const int &channelID, std::vector< std::pair< double, double > > &arriving_pulses, std::vector< ldmx::HgcrocDigiCollection::Sample > &digiToAdd) const
Digitize the signals from the simulated hits.
std::unique_ptr< TRandom3 > noiseInjector_
Generates Gaussian noise on top of real hits.
double gain(const int &channelID) const
Gain for input channel.
double timingJitter_
Jitter of timing mechanism in the chip [ns].
int nADCs_
Depth of ADC buffer.
bool noise_
Put noise in channels, only configure to false if testing.
double timeDnSlope_
Time of Down Slope relative to Pulse Shape Fit [ns].
TF1 pulseFunc_
Functional shape of signal pulse in time.
double clockCycle_
Time interval for chip clock [ns].
std::vector< ldmx::HgcrocDigiCollection::Sample > noiseDigi(const int &channel, const double &soi_amplitude=0) const
Generate a digi of pure noise.
double rateDnSlope_
Rate of Down Slope in Pulse Shape [1/ns].
HgcrocEmulator(const framework::config::Parameters &ps)
Constructor.
double ns_
Conversion from time [ns] to counts.