49 {
50
51
52 SG::ReadHandle<LArDigitContainer> inputContainer(
m_digitKey,ctx);
53
54
55 auto outputContainer = std::make_unique<LArRawChannelContainer>();
56
57
59 const ILArPedestal* peds=*pedHdl;
60
61 SG::ReadCondHandle<LArADC2MeV> adc2mevHdl(
m_adc2MeVKey,ctx);
62 const LArADC2MeV* adc2MeVs=*adc2mevHdl;
63
64 SG::ReadCondHandle<ILArOFC> ofcHdl(
m_ofcKey,ctx);
65 const ILArOFC* ofcs=*ofcHdl;
66
67 SG::ReadCondHandle<ILArShape> shapeHdl(
m_shapeKey,ctx);
68 const ILArShape* shapes=*shapeHdl;
69
71
72 std::unique_ptr<LArDSPThresholdsFlat> run2DSPThresh;
73 const LArDSPThresholdsComplete* run1DSPThresh = nullptr;
77 run2DSPThresh = std::make_unique<LArDSPThresholdsFlat>(*dspThrshAttr);
80 << ". Aborting." );
81 return StatusCode::FAILURE;
82 }
83 }
86 run1DSPThresh = dspThresh.cptr();
87 }
88 else {
90 return StatusCode::FAILURE;
91 }
92 }
93
94
95 for (const LArDigit* digit : *inputContainer) {
96
98
99 const HWIdentifier
id=
digit->hardwareID();
100 const bool connected=(*cabling)->isOnlineConnected(id);
101
102
104
105 const std::vector<short>& samples=
digit->samples();
108
109
110
111
112 const auto& ofca=ofcs->
OFC_a(
id,gain);
113 const auto& adc2mev=adc2MeVs->
ADC2MEV(
id,gain);
114
115 const size_t nOFC=ofca.size();
116
117
118
119 const size_t nSamples=samples.size()-firstSample;
120 if (nSamples<nOFC) {
121 ATH_MSG_ERROR(
"effective sample size: "<< nSamples <<
", must be >= OFC_a size: " << ofca.size());
122 return StatusCode::FAILURE;
123 }
124
125
127 if (!connected) continue;
129 << " gain " << gain);
130 return StatusCode::FAILURE;
131 }
132
134 if (!connected) continue;
136 << " gain " << gain);
137 return StatusCode::FAILURE;
138 }
139
140
141
142
143
146
147
148 std::vector<float> samp_no_ped(nOFC,0.0);
149 for (
size_t i=0;
i<nOFC;++
i) {
150 if (samples[i+firstSample]==4096 || samples[i+firstSample]==0)
saturated=
true;
151 samp_no_ped[
i]=samples[
i+firstSample]-
p;
152 }
153 for (
size_t i=0;
i<nOFC;++
i) {
154 A+=
static_cast<double>(samp_no_ped[
i])*ofca[i];
155 }
156
157
158 const float E=adc2mev[0]+
A*adc2mev[1];
159
161 float tau=0;
162
163
166
168 float ecut(0.);
170 if (run2DSPThresh) {
171 ecut = run2DSPThresh->tQThr(id);
172 }
173 else if (run1DSPThresh) {
174 ecut = run1DSPThresh->
tQThr(
id);
175 }
176 else {
178 return StatusCode::FAILURE;
179 }
180 }
181 else {
183 }
184 if (E1 > ecut) {
185 ATH_MSG_VERBOSE(
"Channel " <<
m_onlineId->channel_name(
id) <<
" gain " << gain <<
" above threshold for tQ computation");
187
188
189 const auto& ofcb=ofcs->
OFC_b(
id,gain);
190 double At=0;
191 for (
size_t i=0;
i<nOFC;++
i) {
192 At+=
static_cast<double>(samp_no_ped[
i])*ofcb[i];
193 }
194
195 tau=(std::fabs(A)>0.1) ? At/A : 0.0;
196 const auto& fullShape=shapes->
Shape(
id,gain);
197
198
199
200 const size_t nSamples=samples.size();
201 if (fullShape.size()>nSamples && nSamples==4 &&
m_firstSample==0) {
203 firstSample=1;
204 }
205 }
206
208 if (!connected) continue;
210 << " gain " << gain);
211 ATH_MSG_ERROR(
"Got size " << fullShape.size() <<
", expected at least " << nSamples+firstSample);
212 return StatusCode::FAILURE;
213 }
214
215 const float* shape=&*fullShape.begin()+firstSample;
216
219 const auto& fullshapeDer=shapes->
ShapeDer(
id,gain);
220 if (
ATH_UNLIKELY(fullshapeDer.size()<nOFC+firstSample)) {
222 << " gain " << gain);
223 ATH_MSG_ERROR(
"Got size " << fullshapeDer.size() <<
", expected at least " << nOFC+firstSample);
224 return StatusCode::FAILURE;
225 }
226
227 const float* shapeDer=&*fullshapeDer.begin()+firstSample;
228 for (
size_t i=0;
i<nOFC;++
i) {
229 q += std::pow((A*(shape[i]-tau*shapeDer[i])-(samp_no_ped[i])),2);
230 }
231 }
232 else {
233
234 for (
size_t i=0;
i<nOFC;++
i) {
235 q += std::pow((A*shape[i]-(samp_no_ped[i])),2);
236 }
237 }
238
239 int iqua =
static_cast<int>(
q);
240 if (iqua > 0xFFFF) iqua=0xFFFF;
241 iquaShort =
static_cast<uint16_t>(iqua & 0xFFFF);
242
244 tau*=(Gaudi::Units::nanosecond/Gaudi::Units::picosecond);
245 }
246
247 outputContainer->emplace_back(id,static_cast<int>(std::floor(E+0.5)),
248 static_cast<int>(std::floor(tau+0.5)),
250 }
251
252
253 SG::WriteHandle<LArRawChannelContainer>outputHandle(
m_rawChannelKey,ctx);
254 ATH_CHECK(outputHandle.record(std::move(outputContainer) ) );
255
256 return StatusCode::SUCCESS;
257}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_VERBOSE(x)
virtual OFCRef_t OFC_b(const HWIdentifier &id, int gain, int tbin=0) const =0
virtual OFCRef_t OFC_a(const HWIdentifier &id, int gain, int tbin=0) const =0
access to OFCs by online ID, gain, and tbin (!=0 for testbeam)
virtual float timeOffset(const HWIdentifier &CellID, int gain) const =0
virtual float pedestal(const HWIdentifier &id, int gain) const =0
virtual ShapeRef_t Shape(const HWIdentifier &id, int gain, int tbin=0, int mode=0) const =0
virtual ShapeRef_t ShapeDer(const HWIdentifier &id, int gain, int tbin=0, int mode=0) const =0
const LArVectorProxy ADC2MEV(const HWIdentifier &id, int gain) const
float tQThr(const HWIdentifier chid) const
Gaudi::Property< int > m_firstSample
SG::ReadCondHandleKey< ILArPedestal > m_pedestalKey
SG::ReadHandleKey< LArDigitContainer > m_digitKey
Gaudi::Property< bool > m_useShapeDer
SG::ReadCondHandleKey< LArOnOffIdMapping > m_cablingKey
SG::ReadCondHandleKey< AthenaAttributeList > m_run2DSPThresholdsKey
SG::ReadCondHandleKey< LArDSPThresholdsComplete > m_run1DSPThresholdsKey
const LArOnlineID * m_onlineId
SG::WriteHandleKey< LArRawChannelContainer > m_rawChannelKey
Gaudi::Property< bool > m_useDBFortQ
Gaudi::Property< bool > m_absECutFortQ
SG::ReadCondHandleKey< LArADC2MeV > m_adc2MeVKey
SG::ReadCondHandleKey< ILArShape > m_shapeKey
SG::ReadCondHandleKey< ILArOFC > m_ofcKey
Gaudi::Property< float > m_eCutFortQ
setScaleOne setStatusOne saturated