146{
147
148
149 constexpr bool use3D = true;
150
151 cellIntersections.reserve(extension.caloLayerIntersections().size() * 1.3);
152
156
160
162 for (const auto* cell : cells) {
163
165
168 continue;
169
170
171
172
173
175 if (sample == CaloSampling::PreSamplerB) {
179 }
180 }
183
184 double drFix =
cell->caloDDE()->dr();
185 double dzFix =
cell->caloDDE()->dz();
186
187
188 int isample =
cell->caloDDE()->getSampling();
192 if (sample == CaloSampling::PreSamplerB || sample == CaloSampling::EMB1 ||
193 sample == CaloSampling::EMB2 || sample == CaloSampling::EMB3)
195
196 double drTG = fabs((
pos->second.first - pos2->second.second).perp());
197 double dzTG = fabs((
pos->second.first - pos2->second.second).z());
198
201 <<
cell->caloDDE()->getSampling() <<
" dr "
202 <<
cell->caloDDE()->dr() <<
" drTG " << drTG);
205 <<
cell->caloDDE()->getSampling() <<
" dz "
206 <<
cell->caloDDE()->dz() <<
" dzTG " << dzTG);
207
208 if (drFix == 0.) {
209
210
211
212
213 if (
cell->caloDDE()->deta() > 0) {
215 double dtheta =
220 drFix = fabs(dr);
221 double detaCheck =
224 <<
cell->caloDDE()->getSampling() <<
" deta "
225 <<
cell->caloDDE()->deta() <<
" detaCheck "
226 << detaCheck << " drFix " << drFix);
227 } else {
229 <<
theta <<
" dtheta " << dtheta <<
" sum/pi "
230 << (
theta + dtheta) * M_1_PI <<
" deta "
231 <<
cell->caloDDE()->deta());
232 }
233
234
235
236 } else {
237 double drMin = 100000.;
238 int dscut = 1;
240 dscut = 0;
241 const CaloCell* cellFound = nullptr;
242 for (const auto* celln : cells) {
243 if (cell == celln)
244 continue;
245 if (
cell->caloDDE()->getSubCalo() == celln->caloDDE()->getSubCalo()) {
246 int dsample = isample - celln->caloDDE()->getSampling();
247 if (abs(dsample) == dscut) {
248 double drNew = fabs(
cell->caloDDE()->r() - celln->caloDDE()->r());
249 if (drNew < 1)
250 continue;
251 if (drNew < drMin) {
252 drMin = drNew;
253 cellFound = celln;
254 }
255 }
256 }
257 }
258 drFix = drMin;
260 <<
cell->caloDDE()->getSampling() <<
" x "
261 <<
cell->caloDDE()->x() <<
" y " <<
cell->caloDDE()->y()
262 <<
" z " <<
cell->caloDDE()->z() <<
" dr "
263 <<
cell->caloDDE()->dr() <<
" drFix " << drFix
264 << " drTG " << drTG);
265 if (cellFound)
268 << cellFound->
caloDDE()->
x() <<
" y "
269 << cellFound->
caloDDE()->
y() <<
" z "
270 << cellFound->
caloDDE()->
z() <<
" dr "
271 << cellFound->
caloDDE()->
dr() <<
" dscut " << dscut
272 << " drFix " << drFix);
273 }
274 }
275
276 if (dzFix == 0.) {
277
278
279
280 if (
cell->caloDDE()->deta() > 0) {
282 double dtheta =
285 double dz = fabs(
cell->caloDDE()->r() /
tan(
theta + dtheta) -
287 dzFix = dz;
288 } else {
290 <<
theta <<
" dtheta " << dtheta <<
" sum/pi "
291 << (
theta + dtheta) * M_1_PI <<
" deta "
292 <<
cell->caloDDE()->deta());
293 }
294 double detaCheck =
297 <<
cell->caloDDE()->getSampling() <<
" deta "
298 <<
cell->caloDDE()->deta() <<
" detaCheck "
299 << detaCheck << " dtheta " << dtheta << " dzFix "
300 << dzFix);
301 } else {
302 double dzMin = 100000.;
303 int dscut = 1;
305 dscut = 0;
306 const CaloCell* cellFound = nullptr;
307 for (const auto* celln : cells) {
308 if (cell == celln)
309 continue;
310 if (
cell->caloDDE()->getSubCalo() == celln->caloDDE()->getSubCalo()) {
311 int isample2 = celln->caloDDE()->getSampling();
312 if (abs(isample - isample2) == dscut) {
313 double dzNew = fabs(
cell->caloDDE()->z() - celln->caloDDE()->z());
314 if (dzNew < 1)
315 continue;
316 if (dzNew < dzMin) {
317 dzMin = dzNew;
318 cellFound = celln;
319 }
320 }
321 }
322 }
323 dzFix = dzMin;
325 <<
cell->caloDDE()->getSampling() <<
" x "
326 <<
cell->caloDDE()->x() <<
" y " <<
cell->caloDDE()->y()
327 <<
" z " <<
cell->caloDDE()->z() <<
" dz "
328 <<
cell->caloDDE()->dz() <<
" dzFix " << dzFix
329 << " dzTG " << dzTG);
330 if (cellFound)
333 << cellFound->
caloDDE()->
x() <<
" y "
334 << cellFound->
caloDDE()->
y() <<
" z "
335 << cellFound->
caloDDE()->
z() <<
" dz "
336 << cellFound->
caloDDE()->
dz() <<
" dscut " << dscut
337 << " dzFix " << dzFix);
338 }
339 }
340
341
342
343
344
345 if (
cell->energy() > 50.)
347 <<
cell->caloDDE()->getSampling() <<
" cell energy "
348 <<
cell->energy() <<
" dzFix " << dzFix <<
" dzTG " << dzTG
349 << " drFix " << drFix << " drTG " << drTG << " barrel "
351
353 dzFix = dzTG;
355 drFix = drTG;
356
357 if (use3D) {
358
360 *cell,
pos->second.first, pos2->second.second, drFix, dzFix);
361 double totpath = (
pos->second.first - pos2->second.second).
mag();
362 path = totpath != 0 ? pathInMM / totpath : 0.;
363 if (path > 0 ||
cell->energy() > 50.) {
365 <<
cell->caloDDE()->getSampling() <<
" cell energy "
366 <<
cell->energy() <<
" drFix " << drFix <<
" dzFix "
367 << dzFix << " path " << path << " length TG " << totpath);
369 <<
cell->caloDDE()->dz() <<
" deta "
370 <<
cell->caloDDE()->deta());
371 }
372 }
373
375 double path2 = 0.;
376
377 if (!use3D)
379
380 if (path2 <= 0. && path <= 0.)
381 continue;
382
383
385 double eLoss = 0.;
386
387
388
390 eLoss = eLossPair->second;
391 if (sample == CaloSampling::PreSamplerB) {
394 eLoss = 0.5 * (eLossPair->second) + 0.5 * (eLossPair2->second);
395 }
396 } else if (sample == CaloSampling::EMB1) {
397 auto eLossPair2 =
eLossLayerMap.find(CaloSampling::PreSamplerB);
399 eLoss = 0.5 * (eLossPair->second) + 0.5 * (eLossPair2->second);
400 }
401 }
402 }
403
405 << path << " PATH2D = " << path2 << " eLoss " << eLoss
406 <<
" cell energy " << (cell)->
energy() <<
" radius "
407 <<
cell->caloDDE()->r() <<
" phi " <<
cell->caloDDE()->phi()
408 <<
" dr " <<
cell->caloDDE()->dr() <<
" dphi "
409 <<
cell->caloDDE()->dphi() <<
" x " <<
cell->caloDDE()->x()
410 <<
" y " <<
cell->caloDDE()->y() <<
" z "
411 <<
cell->caloDDE()->z() <<
" dx " <<
cell->caloDDE()->dx()
412 <<
" dy " <<
cell->caloDDE()->dy() <<
" dz "
413 <<
cell->caloDDE()->dz() <<
" volume "
414 <<
cell->caloDDE()->volume());
415
416 cellIntersections.emplace_back(
417 cell, new ParticleCellIntersection(*cell, eLoss, use3D ? path : path2));
418 }
419 ATH_MSG_DEBUG(
" added cell intersections " << cellIntersections.size());
420}
Scalar theta() const
theta method
Scalar mag() const
mag method
#define ATH_MSG_WARNING(x,...)
#define ATH_MSG_VERBOSE(x,...)
double pathInsideCell(const CaloCell &cell, const Amg::Vector3D &entry, const Amg::Vector3D &exit)
Return the % of the path crossed inside the cell, given the parameters for the extrapolation at entra...
const CaloDetDescrElement * caloDDE() const
get pointer to CaloDetDescrElement (data member)
CaloCell_ID::CaloSample getSampling() const
cell sampling
static double get3DPathLength(const CaloCell &cell, const Amg::Vector3D &entry, const Amg::Vector3D &exit, double drFix, double dzFix)
void entryExitLayerMap(const Trk::CaloExtension &extension, EntryExitLayerMap &result, const LayersToSelect *selection=nullptr)
std::map< CaloSampling::CaloSample, double > ScalarLayerMap
void pathLenLayerMap(const Trk::CaloExtension &extension, ScalarLayerMap &result)
void eLossLayerMap(const Trk::CaloExtension &extension, ScalarLayerMap &result)
std::map< CaloSampling::CaloSample, std::pair< Amg::Vector3D, Amg::Vector3D > > EntryExitLayerMap
path
python interpreter configuration --------------------------------------—