242 {
243
244 std::unique_ptr<CaloEnergy> caloParamMean{
m_caloParamTool->meanParametrizedEnergy(trackMomentum,
eta,
phi)};
245
246 std::unique_ptr<CaloEnergy> caloParamMop{
m_caloParamTool->mopParametrizedEnergy(trackMomentum,
eta,
phi)};
247
248 std::unique_ptr<CaloEnergy> caloParamMip{
m_caloParamTool->meanParametrizedEnergy(10. * Units::GeV,
eta,
phi)};
249
251
252
253
254
255
256
257 double MopLoss = std::abs(caloParamMop->deltaE());
258
259 double MopError = mopPeak->sigmaDeltaE();
260
262
263
265
267
269
270 double MaterialCorrection =
InertMaterial * MopLossCorrected;
271
273
274 const double TileMeasurementMaterial = caloMeas.Tile_SamplingFraction();
275
276 const double LArHECMeasurementMaterial = caloMeas.LArHEC_SamplingFraction();
277
278 const double LArEmMeasurementMaterial = caloMeas.LArEM_SamplingFraction();
279
280 const double TileEnergy = caloMeas.Tile_EnergyMeasured();
281
282 const double EmEnergy = caloMeas.LArEM_EnergyMeasured();
283
284
285 double ForwardHECCorrection = 0.;
286 if (std::abs(
eta) > 2. && caloMeas.LArHEC_EnergyMeasured() > 100.)
287 ForwardHECCorrection = (1. - LArHECMeasurementMaterial) * HECMaterial * MopLossCorrected;
288 const double LArHECEnergy = caloMeas.LArHEC_EnergyMeasured() + ForwardHECCorrection;
289
290 double TotalMeasuredEnergy = TileEnergy + EmEnergy + LArHECEnergy;
291
292 ATH_MSG_VERBOSE(
"Energy Deposition:Tile= " << TileEnergy / Units::GeV <<
",LArHEC= " << LArHECEnergy / Units::GeV
293 << ",EM= " << EmEnergy / Units::GeV << " ForwardHECCorrection "
294 << ForwardHECCorrection / Units::GeV << " HECMaterial " << HECMaterial
295 << " MopLossCorrected " << MopLossCorrected / Units::GeV);
296
297 bool bHEC = false;
298 bool bEM = false;
299
300
301 if (std::abs(
eta) < 1.4) {
302 if (LArHECEnergy + TileEnergy > 0.1 * MopLoss * HECMaterial) bHEC = true;
303 }
else if (std::abs(
eta) > 1.8) {
304 if (LArHECEnergy + TileEnergy > 0.2 * MopLoss * HECMaterial) bHEC = true;
305 } else {
306 if (LArHECEnergy + TileEnergy > 0.25 * MopLoss * HECMaterial) bHEC = true;
307 }
308 if (EmEnergy > 0.5 * MopLoss * EmMaterial) bEM = true;
309
310 if (!bHEC) {
311
312
313 }
314 if (!bEM) {
315
316
317 }
318
319 double MeasCorrected = TotalMeasuredEnergy + MaterialCorrection;
320
321
322
323
324 const double IonizationLoss = (1. / 1.4) * std::abs(caloParamMip->deltaE());
325
326 double eOverMipCorrectionEm = 0.;
327 double eOverMipCorrectionHEC = 0.;
328
329
330
331 if (bEM) {
332 const double emipEM = 0.78;
333 eOverMipCorrectionEm = -(1. / emipEM - 1.) * IonizationLoss * EmMaterial * LArEmMeasurementMaterial;
334 if (EmEnergy + eOverMipCorrectionEm < 0.) eOverMipCorrectionEm = 0.;
335 }
336 if (bHEC) {
337 const double emipTile = 0.86;
338 const double emipLAr = 0.94;
339 const double HECEnergy = TileEnergy + LArHECEnergy;
340 const double eOverMipCorrectionTile =
341 -(1. / emipTile - 1.) * TileEnergy / HECEnergy * IonizationLoss * HECMaterial * TileMeasurementMaterial;
342 const double eOverMipCorrectionLAr =
343 -(1. / emipLAr - 1.) * LArHECEnergy / HECEnergy * IonizationLoss * HECMaterial * LArHECMeasurementMaterial;
344 eOverMipCorrectionHEC = eOverMipCorrectionTile + eOverMipCorrectionLAr;
345 if (LArHECEnergy + TileEnergy + eOverMipCorrectionHEC < 0.) eOverMipCorrectionHEC = 0.;
346 }
347 const double eOverMipCorrection = eOverMipCorrectionEm + eOverMipCorrectionHEC;
348
349
350
351 constexpr double fix1FromPeter[26] = {0.424104, 0.479637, 0.483419, 0.490242, 0.52806, 0.573582, 0.822098,
352 0.767301, 0.809919, 0.658745, 0.157187, 0.413214, 0.771074, 0.61815,
353 0.350113, 0.322785, 0.479294, 0.806183, 0.822161, 0.757731, -0.0857186,
354 -0.0992693, -0.0492252, 0.0650174, 0.261538, 0.360413};
355
356 constexpr double fix2FromPeter[26] = {-0.647703, -0.303498, -0.268645, -0.261292, -0.260152, -0.269253, -0.266212,
357 -0.240837, -0.130172, -0.111638, -0.329423, -0.321011, -0.346050, -0.305592,
358 -0.313293, -0.317111, -0.428393, -0.524839, -0.599547, -0.464013, -0.159663,
359 -0.140879, -0.0975618, 0.0225352, 0.0701925, -0.24778};
360 int ieta =
static_cast<int>(std::abs(
eta) / 0.10);
361 if (ieta > 25) ieta = 25;
362 double FinalMeasuredEnergy = MeasCorrected + eOverMipCorrection + (fix1FromPeter[ieta] + fix2FromPeter[ieta]) * Units::GeV;
363
364 ATH_MSG_VERBOSE(
"Sum of cells " << (TileEnergy + EmEnergy + LArHECEnergy) / Units::GeV <<
" Total energy deposition "
365 << TotalMeasuredEnergy / Units::GeV << " corrected energy deposition " << MeasCorrected / Units::GeV
366 << " e/mip corre " << FinalMeasuredEnergy / Units::GeV << std::endl
367 << "\nFinal Energy Measurement = "
368 << FinalMeasuredEnergy / Units::GeV
369
370
371 << "\nMop Energy Deposition = " << MopLoss / Units::GeV << " +- "
372 << MopError / Units::GeV
373
374
375
376
377 << std::endl
378 << "Final Meas = " << FinalMeasuredEnergy / Units::GeV << " Mop Dep = " << MopLoss / Units::GeV
379 << " +- " << MopError / Units::GeV);
380
381 const double HECIso = caloMeas.Tile_Isolation() + caloMeas.LArHEC_Isolation();
382 const double EmIso = caloMeas.LArEM_Isolation();
383 const double sinTheta = 1. / std::cosh(
eta);
385
386 constexpr double OneOver85 = 1. / 85;
387 constexpr double FifteenGeV = 15. * Units::GeV;
391 bool PassCut = true;
393 PassCut = false;
394
395 int nTracks = 0;
396
399 nTracks = inner.first;
400
401 if (pT < 100. * Units::GeV && nTracks >
m_maxNTracksIso) PassCut =
false;
402 }
403
404 ATH_MSG_VERBOSE(
"pT= " << pT / Units::GeV <<
",HECIso= " << HECIso / Units::GeV <<
",EmIso= " << EmIso / Units::GeV
405 << ", nTracks= " << nTracks << ",PassCut= " << PassCut);
406
408 std::unique_ptr<CaloEnergy> caloEnergy;
409
410
411
412 if (FinalMeasuredEnergy < MopLoss + 2. * MopError && FinalMeasuredEnergy >
m_minFinalEnergy) {
414 caloEnergy.swap(mopPeak);
415 } else if (PassCut) {
416
417
418 double F1 = 0.;
419 if (caloMeas.LArEM_EnergyMeasured() >
m_emEtCut) F1 = caloMeas.LArEM_FirstCompartmentEnergy() / caloMeas.LArEM_EnergyMeasured();
420 ATH_MSG_VERBOSE(
" start Tail and FSR treatment: Et in e.m. " << EmEnergy * sinTheta / Units::GeV <<
" F1 ratio " << F1);
423 double FinalEnergyErrorPlus = 0.50 * std::sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
424 double FinalEnergyErrorMinus = FinalEnergyErrorPlus;
425
426
427 if (LArHECEnergy > 1. * Units::GeV) {
428 FinalEnergyErrorMinus = FinalEnergyErrorPlus = 0.50 * std::sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
429 }
430 double FinalEnergyError = 0.5 * (FinalEnergyErrorMinus + FinalEnergyErrorPlus);
432 caloEnergy = std::make_unique<CaloEnergy>(FinalMeasuredEnergy, FinalEnergyError, FinalEnergyErrorMinus,
433 FinalEnergyErrorPlus, lossType);
435 } else {
436
438 double FinalMeasuredEnergyNoEm = FinalMeasuredEnergy - EmEnergy + MopLossEm;
439
440
441 ATH_MSG_VERBOSE(
" CaloEnergy FSR : EmEnergy " << EmEnergy <<
" FinalMeasuredEnergy " << FinalMeasuredEnergy
442 << " MopLossEm " << MopLossEm << " MopErrorEm " << MopErrorEm
443 << " FinalMeasuredEnergyNoEm " << FinalMeasuredEnergyNoEm);
444 if (FinalMeasuredEnergyNoEm < MopLoss + 2. * MopError) {
445
449 caloEnergy =
450 std::make_unique<CaloEnergy>(caloParamMop->deltaE(), caloParamMop->sigmaDeltaE(),
451 caloParamMop->sigmaMinusDeltaE(), caloParamMop->sigmaPlusDeltaE(), lossType);
452 } else {
454 caloEnergy =
455 std::make_unique<CaloEnergy>(caloParamMean->deltaE(), caloParamMean->sigmaDeltaE(),
456 caloParamMean->sigmaMinusDeltaE(), caloParamMean->sigmaPlusDeltaE(), lossType);
457 }
458 ATH_MSG_VERBOSE(
" CaloEnergy FSR : Small deposit, FinalMeasuredEnergyNoEm " << FinalMeasuredEnergyNoEm
459 << " using Eloss " << caloEnergy->deltaE());
460 } else {
461
464 double FinalEnergyErrorNoEm = 0.50 * std::sqrt(FinalMeasuredEnergyNoEm / Units::GeV) * Units::GeV;
465 FinalEnergyErrorNoEm = std::hypot(FinalEnergyErrorNoEm, MopErrorEm);
466 caloEnergy = std::make_unique<CaloEnergy>(FinalMeasuredEnergyNoEm, FinalEnergyErrorNoEm, FinalEnergyErrorNoEm,
467 FinalEnergyErrorNoEm, lossType);
468 ATH_MSG_VERBOSE(
" CaloEnergy FSR : Large deposit, FinalMeasuredEnergyNoEm " << FinalMeasuredEnergyNoEm);
469 }
470 }
471 } else {
472
475 caloEnergy = std::make_unique<CaloEnergy>(caloParamMop->deltaE(), caloParamMop->sigmaDeltaE(),
476 caloParamMop->sigmaMinusDeltaE(), caloParamMop->sigmaPlusDeltaE(), lossType);
477 } else {
479 caloEnergy = std::make_unique<CaloEnergy>(caloParamMean->deltaE(), caloParamMean->sigmaDeltaE(),
480 caloParamMean->sigmaMinusDeltaE(), caloParamMean->sigmaPlusDeltaE(), lossType);
481 }
482 }
483
484 caloEnergy->set_fsrCandidateEnergy(MopLossEm);
485 return caloEnergy;
486 }
EnergyLossType
Calo Energy Loss Type Parametrized : reconstruction configured to use the parametrization w/o looking...
const std::string trackMomentum