diff --git a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyMeas.h b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyMeas.h
index 70420ae31177469eb50c0e494db0a3ef80d1808c..e4806513800fd34531e11ed7ae99c1cda4efec81 100755
--- a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyMeas.h
+++ b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyMeas.h
@@ -1,5 +1,5 @@
/*
- Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+ Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
*/
///////////////////////////////////////////////////////////////////
@@ -7,7 +7,7 @@
AlgTool to access the energy deposited by a muon in the
calorimeters. The energy deposit is estimated according to
the active measurements in the traversed calorimeter cells.
-
+
@author Konstantinos.Nikolopoulos@cern.ch, Alan.Poppleton@cern.ch
(c) ATLAS Combined Muon software
*/
@@ -18,14 +18,14 @@
//<<<<<< INCLUDES >>>>>>
+#include <atomic>
+
#include "AthenaBaseComps/AthAlgTool.h"
+#include "CaloConditions/CaloNoise.h"
+#include "CaloEvent/CaloCellContainer.h"
#include "GaudiKernel/ToolHandle.h"
#include "MuidInterfaces/IMuidCaloEnergyMeas.h"
#include "StoreGate/ReadHandleKey.h"
-#include "CaloEvent/CaloCellContainer.h"
-#include "CaloConditions/CaloNoise.h"
-
-#include <atomic>
//<<<<<< CLASS DECLARATIONS >>>>>>
@@ -35,91 +35,76 @@ class LArEM_ID;
class LArHEC_ID;
class TileID;
-namespace Rec
-{
+namespace Rec {
- class caloMeas;
- class IMuidCaloEnergyParam;
+class caloMeas;
+class IMuidCaloEnergyParam;
-class MuidCaloEnergyMeas: public AthAlgTool,
- virtual public IMuidCaloEnergyMeas
-{
-public:
- MuidCaloEnergyMeas(const std::string& type,
- const std::string& name,
- const IInterface* parent);
- ~MuidCaloEnergyMeas(void); // destructor
+class MuidCaloEnergyMeas : public AthAlgTool, virtual public IMuidCaloEnergyMeas {
+ public:
+ MuidCaloEnergyMeas(const std::string& type, const std::string& name, const IInterface* parent);
+ ~MuidCaloEnergyMeas(void); // destructor
+
+ StatusCode initialize();
+ StatusCode finalize();
- StatusCode initialize();
- StatusCode finalize();
-
/**IMuidCaloEnergyMeas interface:
get the muon energy loss measurement from the calorimeter,
knowing the track intersection at the em and had cals*/
- CaloMeas* energyMeasurement(double etaEM,
- double phiEM,
- double etaHad,
- double phiHad) const;
-private:
+ CaloMeas* energyMeasurement(double etaEM, double phiEM, double etaHad, double phiHad) const;
+
+ private:
// private methods
- void energyInCalo(CaloMeas& caloMeas,
- const CaloCellContainer* cellContainer,
- double eta,
- double phi,
- int iSubCalo) const;
- void isolationEnergy(CaloMeas& caloMeas,
- const CaloCellContainer* cellContainer,
- double eta,
- double phi,
- int iSubCalo) const;
- double energyInTile(const CaloCellContainer* cellContainer,
- double eta,
- double phi,
- int,
- int) const;
- double energyInLArHEC(const CaloCellContainer* cellContainer,
- double eta,
- double phi,
- int,
- int) const;
- double energyInLArEM(const CaloCellContainer* cellContainer,
- double eta,
- double phi,
- int,
- int) const;
- int cellCounting(const CaloCellContainer* cellContainer,
- double eta,
- double phi) const;
-
+ void energyInCalo(CaloMeas& caloMeas, const CaloCellContainer* cellContainer, double eta, double phi,
+ int iSubCalo) const;
+ void isolationEnergy(CaloMeas& caloMeas, const CaloCellContainer* cellContainer, double eta, double phi,
+ int iSubCalo) const;
+ double energyInTile(const CaloCellContainer* cellContainer, double eta, double phi, int, int) const;
+ double energyInLArHEC(const CaloCellContainer* cellContainer, double eta, double phi, int, int) const;
+ double energyInLArEM(const CaloCellContainer* cellContainer, double eta, double phi, int, int) const;
+ int cellCounting(const CaloCellContainer* cellContainer, double eta, double phi) const;
+
// helpers, managers, tools
- SG::ReadCondHandleKey<CaloNoise> m_noiseCDOKey{this,"CaloNoiseKey","totalNoise","SG Key of CaloNoise data object"};
- ToolHandle<IMuidCaloEnergyParam> m_caloParamTool;
-
- const TileID* m_tileID;
- const LArEM_ID* m_emID;
- const LArHEC_ID* m_hecID;
- SG::ReadHandleKey<CaloCellContainer> m_cellContainerLocation{this,"CellContainerLocation","AllCalo","calo cell container location"};
-
- double m_measurementConeTile;
- double m_measurementConeLArHEC;
- double m_measurementConeLArEM;
- double m_isolationConeTile;
- double m_isolationConeLArHEC;
- double m_isolationConeLArEM;
-
- double m_sigmasAboveNoise; // The minimum sigmas above the noise tool rms
- double m_sigmasAboveNoiseCore; // The minimum sigmas above the noise tool rms
-
- mutable std::atomic_int m_totalCoreCellsEM;
- mutable std::atomic_int m_totalCoreCellsHEC;
- mutable std::atomic_int m_totalCoreCellsTile;
- mutable std::atomic_int m_totalSelectedEM;
- mutable std::atomic_int m_totalSelectedHEC;
- mutable std::atomic_int m_totalSelectedTile;
+ SG::ReadCondHandleKey<CaloNoise> m_noiseCDOKey{
+ this,
+ "CaloNoiseKey",
+ "totalNoise",
+ "SG Key of CaloNoise data object",
+ };
+ ToolHandle<IMuidCaloEnergyParam> m_caloParamTool{
+ this,
+ "CaloParamTool",
+ "",
+ };
+
+ const TileID* m_tileID;
+ const LArEM_ID* m_emID;
+ const LArHEC_ID* m_hecID;
+ SG::ReadHandleKey<CaloCellContainer> m_cellContainerLocation{
+ this,
+ "CellContainerLocation",
+ "AllCalo",
+ "calo cell container location",
+ };
+
+ double m_measurementConeTile;
+ double m_measurementConeLArHEC;
+ double m_measurementConeLArEM;
+ double m_isolationConeTile;
+ double m_isolationConeLArHEC;
+ double m_isolationConeLArEM;
+
+ double m_sigmasAboveNoise; // The minimum sigmas above the noise tool rms
+ double m_sigmasAboveNoiseCore; // The minimum sigmas above the noise tool rms
+
+ mutable std::atomic_int m_totalCoreCellsEM;
+ mutable std::atomic_int m_totalCoreCellsHEC;
+ mutable std::atomic_int m_totalCoreCellsTile;
+ mutable std::atomic_int m_totalSelectedEM;
+ mutable std::atomic_int m_totalSelectedHEC;
+ mutable std::atomic_int m_totalSelectedTile;
};
-} // end of namespace
-
-#endif // MUIDCALOENERGYTOOLS_MUIDCALOENERGYMEAS_H
-
+} // namespace Rec
+#endif // MUIDCALOENERGYTOOLS_MUIDCALOENERGYMEAS_H
diff --git a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyTool.h b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyTool.h
index 16b83ab5884ef9ff28ee4a9eeef229f8f3f91e07..a4ea560f9deeb0bc7978960e25da8db926f354dd 100755
--- a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyTool.h
+++ b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/MuidCaloEnergyTools/MuidCaloEnergyTool.h
@@ -1,5 +1,5 @@
/*
- Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+ Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
*/
//////////////////////////////////////////////////////////////////////////////
@@ -13,7 +13,7 @@
K. Nikolopoulos, D. Fassouliotis, C. Kourkoumelis, and A. Poppleton,
"Event-by-Event Estimate of Muon Energy Loss in ATLAS",
IEEE Trans. Nucl. Sci., vol. 54, no. 5, pt. 2, pp. 1792-1796, Oct. 2007.
-
+
@author Konstantinos.Nikolopoulos@cern.ch, Alan.Poppleton@cern.ch
(c) ATLAS Combined Muon software
*/
@@ -24,95 +24,89 @@
//<<<<<< INCLUDES >>>>>>
+#include <atomic>
+
#include "AthenaBaseComps/AthAlgTool.h"
#include "GaudiKernel/ToolHandle.h"
#include "MuidInterfaces/IMuidCaloEnergy.h"
-
-#include <atomic>
+#include "MuidInterfaces/IMuidCaloEnergyParam.h"
+#include "MuidInterfaces/IMuidTrackIsolation.h"
//<<<<<< CLASS DECLARATIONS >>>>>>
-namespace Rec
-{
+namespace Rec {
class CaloMeas;
class IMuidCaloEnergyMeas;
class IMuidCaloEnergyParam;
class IMuidTrackIsolation;
-class MuidCaloEnergyTool: public AthAlgTool,
- virtual public IMuidCaloEnergy
-{
+class MuidCaloEnergyTool : public AthAlgTool, virtual public IMuidCaloEnergy {
+
+ public:
+ MuidCaloEnergyTool(const std::string& type, const std::string& name, const IInterface* parent);
+ ~MuidCaloEnergyTool(void); // destructor
-public:
- MuidCaloEnergyTool(const std::string& type,
- const std::string& name,
- const IInterface* parent);
- ~MuidCaloEnergyTool(void); // destructor
-
StatusCode initialize();
StatusCode finalize();
/**IMuidCaloEnergy interface:
to get the total energyLoss in the calorimeters */
- CaloEnergy* energyLoss(double trackMomentum,
- double eta,
- double phi) const;
-
+ CaloEnergy* energyLoss(double trackMomentum, double eta, double phi) const;
+
/**IMuidCaloEnergy interface:
TrackStateOnSurface for parameters and energyLoss at the calorimeter mid-surface */
- const Trk::TrackStateOnSurface* trackStateOnSurface(
- const Trk::TrackParameters& middleParameters,
- const Trk::TrackParameters* innerParameters,
- const Trk::TrackParameters* outerParameters) const;
-
-private:
+ const Trk::TrackStateOnSurface* trackStateOnSurface(const Trk::TrackParameters& middleParameters,
+ const Trk::TrackParameters* innerParameters,
+ const Trk::TrackParameters* outerParameters) const;
+
+ private:
// private methods
- CaloEnergy* measurement(double trackMomentum,
- double eta,
- double phi,
- CaloMeas* caloMeas) const;
- double muSpecResolParam(double trackMomentum,
- double eta) const;
- double paramCorrection(double trackMomentum,
- double eta,
- double MopLoss,
- double MopSigma) const;
- double landau(double x,
- double mpv,
- double sigma,
- bool norm) const;
+ CaloEnergy* measurement(double trackMomentum, double eta, double phi, CaloMeas* caloMeas) const;
+ double muSpecResolParam(double trackMomentum, double eta) const;
+ double paramCorrection(double trackMomentum, double eta, double MopLoss, double MopSigma) const;
+ double landau(double x, double mpv, double sigma, bool norm) const;
// helpers, managers, tools
- ToolHandle<IMuidCaloEnergyMeas> m_caloMeasTool;
- ToolHandle<IMuidCaloEnergyParam> m_caloParamTool;
- ToolHandle<IMuidTrackIsolation> m_trackIsolationTool;
-
+ ToolHandle<IMuidCaloEnergyMeas> m_caloMeasTool{
+ this,
+ "CaloMeasTool",
+ "Rec::MuidCaloEnergyMeas/MuidCaloEnergyMeas",
+ };
+ ToolHandle<IMuidCaloEnergyParam> m_caloParamTool{
+ this,
+ "CaloParamTool",
+ "Rec::MuidCaloEnergyParam/MuidCaloEnergyParam",
+ };
+ ToolHandle<IMuidTrackIsolation> m_trackIsolationTool{
+ this,
+ "TrackIsolationTool",
+ "Rec::MuidTrackIsolation/MuidTrackIsolation",
+ };
+
// configurable options
- bool m_cosmics;
- bool m_energyLossMeasurement;
- bool m_forceIsolationFailure;
- bool m_FSRtreatment;
- bool m_MOPparametrization;
- bool m_trackIsolation;
+ bool m_cosmics;
+ bool m_energyLossMeasurement;
+ bool m_forceIsolationFailure;
+ bool m_FSRtreatment;
+ bool m_MOPparametrization;
+ bool m_trackIsolation;
// thresholds for use of energy measurement
- double m_emEtCut; // minimum Et in em for FSR treatment
- double m_emF1Cut; // minimum F1 in em for FSR treatment
- double m_emMinEnergy; // mininum energy in the EM
- double m_hecMinEnergy; // minimum energy in the HEC
- int m_maxNTracksIso; // max #tracks in the isolation cone
- double m_minFinalEnergy; // minimum measured final energy
- double m_minMuonPt; // minimum pt of the muon
+ double m_emEtCut; // minimum Et in em for FSR treatment
+ double m_emF1Cut; // minimum F1 in em for FSR treatment
+ double m_emMinEnergy; // mininum energy in the EM
+ double m_hecMinEnergy; // minimum energy in the HEC
+ int m_maxNTracksIso; // max #tracks in the isolation cone
+ double m_minFinalEnergy; // minimum measured final energy
+ double m_minMuonPt; // minimum pt of the muon
// counters (for finalize)
- mutable std::atomic_int m_countMean; //number of tracks using mean
- mutable std::atomic_int m_countMeasurement; //number of tracks using measurement
- mutable std::atomic_int m_countMop; //number of tracks using mop
+ mutable std::atomic_int m_countMean; // number of tracks using mean
+ mutable std::atomic_int m_countMeasurement; // number of tracks using measurement
+ mutable std::atomic_int m_countMop; // number of tracks using mop
};
-} // end of namespace
-
-#endif // MUIDCALOENERGYTOOLS_MUIDCALOENERGYTOOL_H
-
+} // namespace Rec
+#endif // MUIDCALOENERGYTOOLS_MUIDCALOENERGYTOOL_H
diff --git a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyMeas.cxx b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyMeas.cxx
index ccc58ff0315acacddd10ff9daec40f58e1275d31..e1798e456c84eb5f9c4cf597e0e27fe9c33ae111 100755
--- a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyMeas.cxx
+++ b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyMeas.cxx
@@ -1,5 +1,5 @@
/*
- Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+ Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
*/
///////////////////////////////////////////////////////////////////
@@ -13,58 +13,53 @@
//<<<<<< INCLUDES >>>>>>
+#include "MuidCaloEnergyTools/MuidCaloEnergyMeas.h"
+
+#include "AthenaKernel/Units.h"
#include "CLHEP/Geometry/Vector3D.h"
#include "CaloEvent/CaloCell.h"
#include "CaloIdentifier/CaloCell_ID.h"
-#include "CaloIdentifier/TileID.h"
+#include "CaloIdentifier/CaloID.h"
#include "CaloIdentifier/LArEM_ID.h"
#include "CaloIdentifier/LArHEC_ID.h"
-#include "CaloIdentifier/CaloID.h"
+#include "CaloIdentifier/TileID.h"
#include "CaloUtils/CaloCellList.h"
-#include "MuidCaloEnergyTools/MuidCaloEnergyMeas.h"
#include "MuidEvent/CaloMeas.h"
#include "MuidInterfaces/IMuidCaloEnergyParam.h"
-#include "AthenaKernel/Units.h"
//<<<<<< CLASS STRUCTURE INITIALIZATION >>>>>>
namespace Units = Athena::Units;
-namespace Rec
-{
-
-MuidCaloEnergyMeas::MuidCaloEnergyMeas (const std::string& type,
- const std::string& name,
- const IInterface* parent)
- : AthAlgTool (type, name, parent),
- m_caloParamTool ("", this),
- m_tileID (0),
- m_emID (0),
- m_hecID (0),
- m_sigmasAboveNoise (4.),
- m_sigmasAboveNoiseCore (1.5),
- m_totalCoreCellsEM (0),
- m_totalCoreCellsHEC (0),
- m_totalCoreCellsTile (0),
- m_totalSelectedEM (0),
- m_totalSelectedHEC (0),
- m_totalSelectedTile (0)
+namespace Rec {
+
+MuidCaloEnergyMeas::MuidCaloEnergyMeas(const std::string& type, const std::string& name, const IInterface* parent)
+ : AthAlgTool(type, name, parent),
+ m_tileID(0),
+ m_emID(0),
+ m_hecID(0),
+ m_sigmasAboveNoise(4.),
+ m_sigmasAboveNoiseCore(1.5),
+ m_totalCoreCellsEM(0),
+ m_totalCoreCellsHEC(0),
+ m_totalCoreCellsTile(0),
+ m_totalSelectedEM(0),
+ m_totalSelectedHEC(0),
+ m_totalSelectedTile(0)
{
declareInterface<IMuidCaloEnergyMeas>(this);
- declareProperty ("CaloParamTool", m_caloParamTool);
- declareProperty ("NoiseThresInSigmas", m_sigmasAboveNoise);
- declareProperty ("NoiseThresInSigmasCore", m_sigmasAboveNoiseCore);
-
- m_measurementConeTile = 0.15;
- m_measurementConeLArHEC = 0.15;
- m_measurementConeLArEM = 0.075;
- m_isolationConeTile = 0.3;
- m_isolationConeLArHEC = 0.3;
- m_isolationConeLArEM = 0.15;
+ declareProperty("NoiseThresInSigmas", m_sigmasAboveNoise);
+ declareProperty("NoiseThresInSigmasCore", m_sigmasAboveNoiseCore);
+
+ m_measurementConeTile = 0.15;
+ m_measurementConeLArHEC = 0.15;
+ m_measurementConeLArEM = 0.075;
+ m_isolationConeTile = 0.3;
+ m_isolationConeLArHEC = 0.3;
+ m_isolationConeLArEM = 0.15;
}
-MuidCaloEnergyMeas::~MuidCaloEnergyMeas (void)
-{}
+MuidCaloEnergyMeas::~MuidCaloEnergyMeas(void) {}
//<<<<<< PUBLIC MEMBER FUNCTION DEFINITIONS >>>>>>
@@ -72,34 +67,31 @@ MuidCaloEnergyMeas::~MuidCaloEnergyMeas (void)
StatusCode
MuidCaloEnergyMeas::initialize()
{
- ATH_MSG_INFO( "Initializing MuidCaloEnergyMeas AlgTool"
- << " - package version " << PACKAGE_VERSION );
+ ATH_MSG_INFO("Initializing MuidCaloEnergyMeas AlgTool"
+ << " - package version " << PACKAGE_VERSION);
// retrieve TileID helper and TileIfno from det store
- if (StatusCode::SUCCESS != detStore()->retrieve(m_tileID))
- {
- ATH_MSG_FATAL( "Unable to retrieve TileID helper from DetectorStore" );
- return StatusCode::FAILURE;
+ if (StatusCode::SUCCESS != detStore()->retrieve(m_tileID)) {
+ ATH_MSG_FATAL("Unable to retrieve TileID helper from DetectorStore");
+ return StatusCode::FAILURE;
}
- ATH_MSG_VERBOSE( "Accessed TileID helper" );
-
- // retrieve LArEM helper
- if (StatusCode::SUCCESS != detStore()->retrieve(m_emID,"LArEM_ID"))
- {
- ATH_MSG_FATAL( "Unable to retrieve LArEM_ID helper from DetectorStore" );
- return StatusCode::FAILURE;
+ ATH_MSG_VERBOSE("Accessed TileID helper");
+
+ // retrieve LArEM helper
+ if (StatusCode::SUCCESS != detStore()->retrieve(m_emID, "LArEM_ID")) {
+ ATH_MSG_FATAL("Unable to retrieve LArEM_ID helper from DetectorStore");
+ return StatusCode::FAILURE;
}
- ATH_MSG_VERBOSE( "Accessed LArEM helper" );
+ ATH_MSG_VERBOSE("Accessed LArEM helper");
// retrieve LArHEC helper
- if (StatusCode::SUCCESS != detStore()->retrieve(m_hecID,"LArHEC_ID"))
- {
- ATH_MSG_FATAL( "Unable to retrieve LArHEC_ID helper from DetectorStore" );
- return StatusCode::FAILURE;
+ if (StatusCode::SUCCESS != detStore()->retrieve(m_hecID, "LArHEC_ID")) {
+ ATH_MSG_FATAL("Unable to retrieve LArHEC_ID helper from DetectorStore");
+ return StatusCode::FAILURE;
}
- ATH_MSG_VERBOSE( "Accessed LArHEC helper" );
-
+ ATH_MSG_VERBOSE("Accessed LArHEC helper");
+ ATH_CHECK(m_caloParamTool.retrieve());
ATH_CHECK(m_noiseCDOKey.initialize());
ATH_CHECK(m_cellContainerLocation.initialize());
@@ -109,122 +101,99 @@ MuidCaloEnergyMeas::initialize()
StatusCode
MuidCaloEnergyMeas::finalize()
{
- ATH_MSG_INFO( "Finalizing MuidCaloEnergyMeas Tool" );
- ATH_MSG_INFO( " EM: selected " << m_totalSelectedEM << " from "
- << m_totalCoreCellsEM << " cells in core" );
-
- ATH_MSG_INFO( " Tile: selected " << m_totalSelectedTile << " from "
- << m_totalCoreCellsTile << " cells in core" );
-
- ATH_MSG_INFO( " HEC: selected " << m_totalSelectedHEC << " from "
- << m_totalCoreCellsHEC << " cells in core" );
+ ATH_MSG_INFO("Finalizing MuidCaloEnergyMeas Tool");
+ ATH_MSG_INFO(" EM: selected " << m_totalSelectedEM << " from " << m_totalCoreCellsEM << " cells in core");
+
+ ATH_MSG_INFO(" Tile: selected " << m_totalSelectedTile << " from " << m_totalCoreCellsTile << " cells in core");
+
+ ATH_MSG_INFO(" HEC: selected " << m_totalSelectedHEC << " from " << m_totalCoreCellsHEC << " cells in core");
return StatusCode::SUCCESS;
}
CaloMeas*
-MuidCaloEnergyMeas::energyMeasurement (double etaEM,
- double phiEM,
- double etaHad,
- double phiHad) const
+MuidCaloEnergyMeas::energyMeasurement(double etaEM, double phiEM, double etaHad, double phiHad) const
{
- SG::ReadHandle<CaloCellContainer> cellContainer(m_cellContainerLocation);
- if(!cellContainer.isPresent()){
- ATH_MSG_DEBUG("No calo cell container "<<m_cellContainerLocation.key()<<", energy measurement is 0");
- return 0;
- }
- if(!cellContainer.isValid()){
- ATH_MSG_WARNING("Calo cell container "<<m_cellContainerLocation.key()<<" not valid!");
- return 0;
- }
+ SG::ReadHandle<CaloCellContainer> cellContainer(m_cellContainerLocation);
+ if (!cellContainer.isPresent()) {
+ ATH_MSG_DEBUG("No calo cell container " << m_cellContainerLocation.key() << ", energy measurement is 0");
+ return 0;
+ }
+ if (!cellContainer.isValid()) {
+ ATH_MSG_WARNING("Calo cell container " << m_cellContainerLocation.key() << " not valid!");
+ return 0;
+ }
// set measured tile energy, measured sampling fraction and isolation energy into CaloMeas
- CaloMeas* caloMeas = new CaloMeas();
- energyInCalo(*caloMeas,cellContainer.cptr(),etaHad,phiHad,0);
- isolationEnergy(*caloMeas,cellContainer.cptr(),etaHad,phiHad,0);
+ CaloMeas* caloMeas = new CaloMeas();
+ energyInCalo(*caloMeas, cellContainer.cptr(), etaHad, phiHad, 0);
+ isolationEnergy(*caloMeas, cellContainer.cptr(), etaHad, phiHad, 0);
// similar for LArHEC
- energyInCalo(*caloMeas,cellContainer.cptr(),etaHad,phiHad,1);
- isolationEnergy(*caloMeas,cellContainer.cptr(),etaHad,phiHad,1);
+ energyInCalo(*caloMeas, cellContainer.cptr(), etaHad, phiHad, 1);
+ isolationEnergy(*caloMeas, cellContainer.cptr(), etaHad, phiHad, 1);
// and for the em calo
- energyInCalo(*caloMeas,cellContainer.cptr(),etaEM,phiEM,2);
- isolationEnergy(*caloMeas,cellContainer.cptr(),etaEM,phiEM,2);
-
- if (msgLvl(MSG::DEBUG))
- {
- // int cellsOverThreshold = 0;
- // cellsOverThreshold = cellCounting(cellContainer,etaEM,phiEM);
- ATH_MSG_DEBUG( std::setiosflags(std::ios::fixed)
- << " Tile energy (GeV) :" << std::setw(8) << std::setprecision(3)
- << caloMeas->Tile_EnergyMeasured()/Units::GeV
- << " Material percent:" << std::setw(4) << std::setprecision(0)
- << 100.*caloMeas->Tile_SamplingFraction()
- << " ISO :" << std::setw(8) << std::setprecision(3)
- << caloMeas->Tile_Isolation()/Units::GeV
- << endmsg
- << " LArHEC energy (GeV):" << std::setw(8) << std::setprecision(3)
- << caloMeas->LArHEC_EnergyMeasured()/Units::GeV
- << " Material percent:" << std::setw(4) << std::setprecision(0)
- << 100.*caloMeas->LArHEC_SamplingFraction()
- << " ISO :" << std::setw(8) << std::setprecision(3)
- << caloMeas->LArHEC_Isolation()/Units::GeV
- << endmsg
- << " EM energy :" << std::setw(8) << std::setprecision(3)
- << caloMeas->LArEM_EnergyMeasured()/Units::GeV
- << " first compartment:" << std::setw(8) << std::setprecision(3)
- << caloMeas->LArEM_FirstCompartmentEnergy()/Units::GeV
- << " Material percent:" << std::setw(4) << std::setprecision(0)
- << 100.*caloMeas->LArEM_SamplingFraction()
- << " ISO :" << std::setw(8) << std::setprecision(3)
- << caloMeas->LArEM_Isolation()/Units::GeV);
+ energyInCalo(*caloMeas, cellContainer.cptr(), etaEM, phiEM, 2);
+ isolationEnergy(*caloMeas, cellContainer.cptr(), etaEM, phiEM, 2);
+
+ if (msgLvl(MSG::DEBUG)) {
+ // int cellsOverThreshold = 0;
+ // cellsOverThreshold = cellCounting(cellContainer,etaEM,phiEM);
+ ATH_MSG_DEBUG(std::setiosflags(std::ios::fixed)
+ << " Tile energy (GeV) :" << std::setw(8) << std::setprecision(3)
+ << caloMeas->Tile_EnergyMeasured() / Units::GeV << " Material percent:" << std::setw(4)
+ << std::setprecision(0) << 100. * caloMeas->Tile_SamplingFraction() << " ISO :" << std::setw(8)
+ << std::setprecision(3) << caloMeas->Tile_Isolation() / Units::GeV << endmsg
+ << " LArHEC energy (GeV):" << std::setw(8) << std::setprecision(3)
+ << caloMeas->LArHEC_EnergyMeasured() / Units::GeV << " Material percent:" << std::setw(4)
+ << std::setprecision(0) << 100. * caloMeas->LArHEC_SamplingFraction() << " ISO :" << std::setw(8)
+ << std::setprecision(3) << caloMeas->LArHEC_Isolation() / Units::GeV << endmsg
+ << " EM energy :" << std::setw(8) << std::setprecision(3)
+ << caloMeas->LArEM_EnergyMeasured() / Units::GeV << " first compartment:" << std::setw(8)
+ << std::setprecision(3) << caloMeas->LArEM_FirstCompartmentEnergy() / Units::GeV
+ << " Material percent:" << std::setw(4) << std::setprecision(0)
+ << 100. * caloMeas->LArEM_SamplingFraction() << " ISO :" << std::setw(8) << std::setprecision(3)
+ << caloMeas->LArEM_Isolation() / Units::GeV);
}
-
+
return caloMeas;
}
//<<<<<< PRIVATE MEMBER FUNCTION DEFINITIONS >>>>>>
int
-MuidCaloEnergyMeas::cellCounting(const CaloCellContainer* cellContainer,
- double mu_eta,
- double mu_phi) const
+MuidCaloEnergyMeas::cellCounting(const CaloCellContainer* cellContainer, double mu_eta, double mu_phi) const
{
- //int isubcalo = 2;
+ // int isubcalo = 2;
constexpr double lowest_threshold = 4 * 50.;
-
+
SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
+ const CaloNoise* noiseCDO = *noiseHdl;
CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
- CaloCellList myList(cellContainer,iCalo); // Construct the list
- myList.select(mu_eta,mu_phi,0.2,0.2);
+ CaloCellList myList(cellContainer, iCalo); // Construct the list
+ myList.select(mu_eta, mu_phi, 0.2, 0.2);
- int count = 0;
+ int count = 0;
CaloCellList::list_iterator ilistfirst = myList.begin();
CaloCellList::list_iterator ilistlast = myList.end();
- for(;ilistfirst!=ilistlast;ilistfirst++)
- {
- const double cellEnergy= (*ilistfirst)->energy();
- const double noise_rms =noiseCDO->getNoise((*ilistfirst)->ID(),(*ilistfirst)->gain());
-
- if( cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise )
- {
- count+=1;
- }
+ for (; ilistfirst != ilistlast; ilistfirst++) {
+ const double cellEnergy = (*ilistfirst)->energy();
+ const double noise_rms = noiseCDO->getNoise((*ilistfirst)->ID(), (*ilistfirst)->gain());
+
+ if (cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise) {
+ count += 1;
+ }
}
- ATH_MSG_DEBUG( " counted "<< count
- << " cells over threshold out of a total of " << myList.ncells() << " cells" );
+ ATH_MSG_DEBUG(" counted " << count << " cells over threshold out of a total of " << myList.ncells() << " cells");
return count;
}
void
-MuidCaloEnergyMeas::energyInCalo (CaloMeas& caloMeas,
- const CaloCellContainer* cellContainer,
- double muEta,
- double muPhi,
- int isubcalo) const
+MuidCaloEnergyMeas::energyInCalo(CaloMeas& caloMeas, const CaloCellContainer* cellContainer, double muEta, double muPhi,
+ int isubcalo) const
{
/* -------------------------------------------
Tile Cal
@@ -241,409 +210,288 @@ MuidCaloEnergyMeas::energyInCalo (CaloMeas& caloMeas,
2 --> Sampling 3
3 --> Sampling 4
===========================================
- LarEM calorimeter
- sample_number
+ LarEM calorimeter
+ sample_number
0 --> Presampler
1 --> Sampling 1
2 --> Sampling 2
3 --> Sampling 3
leadingEnergy is contribution from presampler and first compartment
- -------------------------------------------*/
+ -------------------------------------------*/
SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
-
- double totalEnergy = 0.;
- double leadingEnergy = 0.;
- CaloCellList* myList = 0;
-
- if (isubcalo == 0)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta,muPhi,m_measurementConeTile,m_measurementConeTile);
- }
- else if (isubcalo == 1)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LARHEC;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta,muPhi,m_measurementConeLArHEC,m_measurementConeLArHEC);
- }
- else if (isubcalo == 2)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta,muPhi,m_measurementConeLArEM,m_measurementConeLArEM);
+ const CaloNoise* noiseCDO = *noiseHdl;
+
+ double totalEnergy = 0.;
+ double leadingEnergy = 0.;
+ CaloCellList* myList = 0;
+
+ if (isubcalo == 0) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_measurementConeTile, m_measurementConeTile);
+ } else if (isubcalo == 1) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LARHEC;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_measurementConeLArHEC, m_measurementConeLArHEC);
+ } else if (isubcalo == 2) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_measurementConeLArEM, m_measurementConeLArEM);
}
- int countCoreCells = 0;
- int countSelected = 0;
- bool compartment[4] = {false,false,false,false};
- double measuredSamplings = 0.;
-
- if (myList)
- {
- // get last cell (as CaloCellList->back() method doesn't work)
- const CaloCell* lastCell= 0;
- CaloCellList::list_iterator cell = myList->begin();
- for ( ; cell != myList->end(); ++cell) lastCell = *cell;
-
- // flag core in each sampling
- for (int coreSampling = 0; coreSampling != 4; ++coreSampling)
- {
- const CaloCell* coreCell = 0;
- double coreRadius = 0.;
- for (cell = myList->begin(); cell != myList->end(); ++cell)
- {
- int sampling = -1;
- if (isubcalo == 0)
- {
- sampling = m_tileID->sample((**cell).ID());
- }
- else if (isubcalo == 1)
- {
- sampling = m_hecID->sampling((**cell).ID());
- }
- else if (isubcalo == 2)
- {
- sampling = m_emID->sampling((**cell).ID());
- }
- if (*cell != lastCell && sampling != coreSampling) continue;
-
- double deltaEta = (**cell).eta() - muEta;
- double deltaPhi = (**cell).phi() - muPhi;
- double radius = deltaEta*deltaEta + deltaPhi*deltaPhi;
-
- if (sampling == coreSampling)
- {
- if (! coreCell || radius < coreRadius)
- {
- coreCell = *cell;
- coreRadius = radius;
- }
- }
-
- if (*cell != lastCell || ! coreCell) continue;
-
- for (CaloCellList::list_iterator cell2 = myList->begin();
- cell2 != myList->end();
- ++cell2)
- {
- if (isubcalo == 0)
- {
- sampling = m_tileID->sample((**cell2).ID());
- }
- else if (isubcalo == 1)
- {
- sampling = m_hecID->sampling((**cell2).ID());
- }
- else if (isubcalo == 2)
- {
- sampling = m_emID->sampling((**cell2).ID());
- }
- if (sampling != coreSampling) continue;
-
- double cellEnergy = (**cell2).energy();
- double noiseRms = noiseCDO->getNoise((**cell2).ID(),(**cell2).gain());
-
- // looser selection for core cell where at least mip is expected
- bool cellSelected = cellEnergy > m_sigmasAboveNoise * noiseRms;
- if (*cell2 == coreCell
- && cellEnergy > m_sigmasAboveNoiseCore * noiseRms) cellSelected = true;
-
- if (cellSelected)
- {
- ++countSelected;
- totalEnergy += cellEnergy;
- compartment[coreSampling] = true;
- if (coreSampling < 2)
- leadingEnergy += cellEnergy;
- }
- if (*cell2 == coreCell)
- {
- ++countCoreCells;
- if (isubcalo == 0)
- {
- ++m_totalCoreCellsTile;
- if (cellSelected) ++m_totalSelectedTile;
- }
- else if (isubcalo == 1)
- {
- ++m_totalCoreCellsHEC;
- if (cellSelected) ++m_totalSelectedHEC;
- }
- else if (isubcalo == 2)
- {
- ++m_totalCoreCellsEM;
- if (cellSelected) ++m_totalSelectedEM;
- }
- }
- if (msgLvl(MSG::DEBUG))
- {
- std::string info = "";
- std::string type = " Tile ";
- if (isubcalo == 1)
- {
- type = " LArHEC";
- }
- else if (isubcalo == 2)
- {
- type = " EM ";
- }
- if (cellSelected && *cell2 == coreCell)
- {
- info = " selected core cell# ";
- }
- else if (cellSelected)
- {
- info = " selected cell# ";
- }
- else if (*cell2 == coreCell)
- {
- info = " cell in core NOT selected";
- }
-
- if (info == "")
- {
- ATH_MSG_VERBOSE( std::setiosflags(std::ios::fixed) << type
- << " Sampling: "<< std::setw(1)
- << coreSampling
- << " Radius :"
- << std::setw(6) << std::setprecision(0)
- << (**cell2).caloDDE()->r()
- << " Eta : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).eta()
- << " Phi : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).phi()
- << " Noise level : "
- << std::setw(6) << std::setprecision(0)
- << noiseRms
- << " Energy : "
- << std::setw(7) << std::setprecision(0)
- << (**cell2).energy() << info );
- }
- else if (cellSelected)
- {
- ATH_MSG_DEBUG( std::setiosflags(std::ios::fixed) << type
- << " Sampling: "<< std::setw(1)
- << coreSampling
- << " Radius :"
- << std::setw(6) << std::setprecision(0)
- << (**cell2).caloDDE()->r()
- << " Eta : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).eta()
- << " Phi : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).phi()
- << " Noise level : "
- << std::setw(6) << std::setprecision(0)
- << noiseRms
- << " Energy : "
- << std::setw(7) << std::setprecision(0)
- << (**cell2).energy() << info << countSelected);
- }
- else
- {
- ATH_MSG_DEBUG( std::setiosflags(std::ios::fixed) << type
- << " Sampling: "<< std::setw(1)
- << coreSampling
- << " Radius :"
- << std::setw(6) << std::setprecision(0)
- << (**cell2).caloDDE()->r()
- << " Eta : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).eta()
- << " Phi : "
- << std::setw(6) << std::setprecision(2)
- << (**cell2).phi()
- << " Noise level : "
- << std::setw(6) << std::setprecision(0)
- << noiseRms
- << " Energy : "
- << std::setw(7) << std::setprecision(0)
- << (**cell2).energy() << info );
- }
- }
- }
- }
- }
-
- if (msgLvl(MSG::DEBUG))
- {
- std::string info = "";
- std::string type = " Tile ";
- if (isubcalo == 1)
- {
- type = " LArHEC";
- }
- else if (isubcalo == 2)
- {
- type = " EM ";
- }
-
- ATH_MSG_DEBUG ( type << " at eta = " << muEta << " : selected "
- << countSelected << " from measured cone with "
- << countCoreCells << " cells forming core" );
- }
-
- delete myList;
-
- for (int i = 0; i < 4 ; ++i)
- if (compartment[i]) measuredSamplings += m_caloParamTool->caloCompartmentDepth(isubcalo,i);
-
- // store result in caloMeas
- if (isubcalo == 0)
- {
- caloMeas.Tile_EnergyMeasured(totalEnergy);
- caloMeas.Tile_SamplingFraction(measuredSamplings);
- }
- else if (isubcalo == 1)
- {
- caloMeas.LArHEC_EnergyMeasured(totalEnergy);
- caloMeas.LArHEC_SamplingFraction(measuredSamplings);
- }
- else if (isubcalo == 2)
- {
- caloMeas.LArEM_EnergyMeasured(totalEnergy);
- caloMeas.LArEM_FirstCompartmentEnergy(leadingEnergy);
- caloMeas.LArEM_SamplingFraction(measuredSamplings);
- }
+ int countCoreCells = 0;
+ int countSelected = 0;
+ bool compartment[4] = {false, false, false, false};
+ double measuredSamplings = 0.;
+
+ if (myList) {
+ // get last cell (as CaloCellList->back() method doesn't work)
+ const CaloCell* lastCell = 0;
+ CaloCellList::list_iterator cell = myList->begin();
+ for (; cell != myList->end(); ++cell) lastCell = *cell;
+
+ // flag core in each sampling
+ for (int coreSampling = 0; coreSampling != 4; ++coreSampling) {
+ const CaloCell* coreCell = 0;
+ double coreRadius = 0.;
+ for (cell = myList->begin(); cell != myList->end(); ++cell) {
+ int sampling = -1;
+ if (isubcalo == 0) {
+ sampling = m_tileID->sample((**cell).ID());
+ } else if (isubcalo == 1) {
+ sampling = m_hecID->sampling((**cell).ID());
+ } else if (isubcalo == 2) {
+ sampling = m_emID->sampling((**cell).ID());
+ }
+ if (*cell != lastCell && sampling != coreSampling) continue;
+
+ double deltaEta = (**cell).eta() - muEta;
+ double deltaPhi = (**cell).phi() - muPhi;
+ double radius = deltaEta * deltaEta + deltaPhi * deltaPhi;
+
+ if (sampling == coreSampling) {
+ if (!coreCell || radius < coreRadius) {
+ coreCell = *cell;
+ coreRadius = radius;
+ }
+ }
+
+ if (*cell != lastCell || !coreCell) continue;
+
+ for (CaloCellList::list_iterator cell2 = myList->begin(); cell2 != myList->end(); ++cell2) {
+ if (isubcalo == 0) {
+ sampling = m_tileID->sample((**cell2).ID());
+ } else if (isubcalo == 1) {
+ sampling = m_hecID->sampling((**cell2).ID());
+ } else if (isubcalo == 2) {
+ sampling = m_emID->sampling((**cell2).ID());
+ }
+ if (sampling != coreSampling) continue;
+
+ double cellEnergy = (**cell2).energy();
+ double noiseRms = noiseCDO->getNoise((**cell2).ID(), (**cell2).gain());
+
+ // looser selection for core cell where at least mip is expected
+ bool cellSelected = cellEnergy > m_sigmasAboveNoise * noiseRms;
+ if (*cell2 == coreCell && cellEnergy > m_sigmasAboveNoiseCore * noiseRms) cellSelected = true;
+
+ if (cellSelected) {
+ ++countSelected;
+ totalEnergy += cellEnergy;
+ compartment[coreSampling] = true;
+ if (coreSampling < 2) leadingEnergy += cellEnergy;
+ }
+ if (*cell2 == coreCell) {
+ ++countCoreCells;
+ if (isubcalo == 0) {
+ ++m_totalCoreCellsTile;
+ if (cellSelected) ++m_totalSelectedTile;
+ } else if (isubcalo == 1) {
+ ++m_totalCoreCellsHEC;
+ if (cellSelected) ++m_totalSelectedHEC;
+ } else if (isubcalo == 2) {
+ ++m_totalCoreCellsEM;
+ if (cellSelected) ++m_totalSelectedEM;
+ }
+ }
+ if (msgLvl(MSG::DEBUG)) {
+ std::string info = "";
+ std::string type = " Tile ";
+ if (isubcalo == 1) {
+ type = " LArHEC";
+ } else if (isubcalo == 2) {
+ type = " EM ";
+ }
+ if (cellSelected && *cell2 == coreCell) {
+ info = " selected core cell# ";
+ } else if (cellSelected) {
+ info = " selected cell# ";
+ } else if (*cell2 == coreCell) {
+ info = " cell in core NOT selected";
+ }
+
+ if (info == "") {
+ ATH_MSG_VERBOSE(std::setiosflags(std::ios::fixed)
+ << type << " Sampling: " << std::setw(1) << coreSampling << " Radius :"
+ << std::setw(6) << std::setprecision(0) << (**cell2).caloDDE()->r()
+ << " Eta : " << std::setw(6) << std::setprecision(2) << (**cell2).eta()
+ << " Phi : " << std::setw(6) << std::setprecision(2) << (**cell2).phi()
+ << " Noise level : " << std::setw(6) << std::setprecision(0) << noiseRms
+ << " Energy : " << std::setw(7) << std::setprecision(0)
+ << (**cell2).energy() << info);
+ } else if (cellSelected) {
+ ATH_MSG_DEBUG(std::setiosflags(std::ios::fixed)
+ << type << " Sampling: " << std::setw(1) << coreSampling << " Radius :"
+ << std::setw(6) << std::setprecision(0) << (**cell2).caloDDE()->r()
+ << " Eta : " << std::setw(6) << std::setprecision(2) << (**cell2).eta()
+ << " Phi : " << std::setw(6) << std::setprecision(2) << (**cell2).phi()
+ << " Noise level : " << std::setw(6) << std::setprecision(0) << noiseRms
+ << " Energy : " << std::setw(7) << std::setprecision(0) << (**cell2).energy()
+ << info << countSelected);
+ } else {
+ ATH_MSG_DEBUG(std::setiosflags(std::ios::fixed)
+ << type << " Sampling: " << std::setw(1) << coreSampling << " Radius :"
+ << std::setw(6) << std::setprecision(0) << (**cell2).caloDDE()->r()
+ << " Eta : " << std::setw(6) << std::setprecision(2) << (**cell2).eta()
+ << " Phi : " << std::setw(6) << std::setprecision(2) << (**cell2).phi()
+ << " Noise level : " << std::setw(6) << std::setprecision(0) << noiseRms
+ << " Energy : " << std::setw(7) << std::setprecision(0) << (**cell2).energy()
+ << info);
+ }
+ }
+ }
+ }
+ }
+
+ if (msgLvl(MSG::DEBUG)) {
+ std::string info = "";
+ std::string type = " Tile ";
+ if (isubcalo == 1) {
+ type = " LArHEC";
+ } else if (isubcalo == 2) {
+ type = " EM ";
+ }
+
+ ATH_MSG_DEBUG(type << " at eta = " << muEta << " : selected " << countSelected
+ << " from measured cone with " << countCoreCells << " cells forming core");
+ }
+
+ delete myList;
+
+ for (int i = 0; i < 4; ++i)
+ if (compartment[i]) measuredSamplings += m_caloParamTool->caloCompartmentDepth(isubcalo, i);
+
+ // store result in caloMeas
+ if (isubcalo == 0) {
+ caloMeas.Tile_EnergyMeasured(totalEnergy);
+ caloMeas.Tile_SamplingFraction(measuredSamplings);
+ } else if (isubcalo == 1) {
+ caloMeas.LArHEC_EnergyMeasured(totalEnergy);
+ caloMeas.LArHEC_SamplingFraction(measuredSamplings);
+ } else if (isubcalo == 2) {
+ caloMeas.LArEM_EnergyMeasured(totalEnergy);
+ caloMeas.LArEM_FirstCompartmentEnergy(leadingEnergy);
+ caloMeas.LArEM_SamplingFraction(measuredSamplings);
+ }
}
}
void
-MuidCaloEnergyMeas::isolationEnergy (CaloMeas& caloMeas,
- const CaloCellContainer* cellContainer,
- double muEta,
- double muPhi,
- int isubcalo) const
+MuidCaloEnergyMeas::isolationEnergy(CaloMeas& caloMeas, const CaloCellContainer* cellContainer, double muEta,
+ double muPhi, int isubcalo) const
{
- double totalEnergy = 0.;
- CaloCellList* myList = 0;
+ double totalEnergy = 0.;
+ CaloCellList* myList = 0;
- SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
-
- if (isubcalo == 0)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta, muPhi, m_isolationConeTile, m_isolationConeTile);
- }
- else if (isubcalo == 1)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LARHEC;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta, muPhi, m_isolationConeLArHEC, m_isolationConeLArHEC);
- }
- else if (isubcalo == 2)
- {
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
- myList = new CaloCellList(cellContainer,iCalo);
- myList->select(muEta, muPhi, m_isolationConeLArEM, m_isolationConeLArEM);
+ SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
+ const CaloNoise* noiseCDO = *noiseHdl;
+
+ if (isubcalo == 0) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_isolationConeTile, m_isolationConeTile);
+ } else if (isubcalo == 1) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LARHEC;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_isolationConeLArHEC, m_isolationConeLArHEC);
+ } else if (isubcalo == 2) {
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
+ myList = new CaloCellList(cellContainer, iCalo);
+ myList->select(muEta, muPhi, m_isolationConeLArEM, m_isolationConeLArEM);
}
- if (myList)
- {
- // get last cell (as CaloCellList->back() method doesn't work)
- const CaloCell* lastCell= 0;
- CaloCellList::list_iterator cell = myList->begin();
- for ( ; cell != myList->end(); ++cell) lastCell = *cell;
-
- // flag core in each sampling
- for (int coreSampling = 0; coreSampling != 4; ++coreSampling)
- {
- const CaloCell* coreCell = 0;
- double coreRadius = 0.;
- for (cell = myList->begin(); cell != myList->end(); ++cell)
- {
- int sampling = -1;
- if (isubcalo == 0)
- {
- sampling = m_tileID->sample((**cell).ID());
- }
- else if (isubcalo == 1)
- {
- sampling = m_hecID->sampling((**cell).ID());
- }
- else if (isubcalo == 2)
- {
- sampling = m_emID->sampling((**cell).ID());
- }
- if (*cell != lastCell && sampling != coreSampling) continue;
-
- double deltaEta = (**cell).eta() - muEta;
- double deltaPhi = (**cell).phi() - muPhi;
- double radius = deltaEta*deltaEta + deltaPhi*deltaPhi;
-
- if (sampling == coreSampling)
- {
- if (! coreCell || radius < coreRadius)
- {
- coreCell = *cell;
- coreRadius = radius;
- }
- }
-
- if (*cell != lastCell || ! coreCell) continue;
-
- for (CaloCellList::list_iterator cell2 = myList->begin();
- cell2 != myList->end();
- ++cell2)
- {
- if (isubcalo == 0)
- {
- sampling = m_tileID->sample((**cell2).ID());
- }
- else if (isubcalo == 1)
- {
- sampling = m_hecID->sampling((**cell2).ID());
- }
- else if (isubcalo == 2)
- {
- sampling = m_emID->sampling((**cell2).ID());
- }
- if (sampling != coreSampling) continue;
-
- double cellEnergy = (**cell2).energy();
- double noiseRms = noiseCDO->getNoise((**cell2).ID(),(**cell2).gain());
- // looser selection for core cell where at least mip is expected
- bool cellSelected = cellEnergy > m_sigmasAboveNoise * noiseRms;
- if (*cell2 == coreCell
- && cellEnergy > m_sigmasAboveNoiseCore * noiseRms) cellSelected = true;
-
- if (cellSelected) totalEnergy += cellEnergy;
- }
- }
- }
-
- // store result in caloMeas
- if (isubcalo == 0)
- {
- caloMeas.Tile_Isolation(totalEnergy);
- }
- else if (isubcalo == 1)
- {
- caloMeas.LArHEC_Isolation(totalEnergy);
- }
- else if (isubcalo == 2)
- {
- caloMeas.LArEM_Isolation(totalEnergy);
- }
+ if (myList) {
+ // get last cell (as CaloCellList->back() method doesn't work)
+ const CaloCell* lastCell = 0;
+ CaloCellList::list_iterator cell = myList->begin();
+ for (; cell != myList->end(); ++cell) lastCell = *cell;
+
+ // flag core in each sampling
+ for (int coreSampling = 0; coreSampling != 4; ++coreSampling) {
+ const CaloCell* coreCell = 0;
+ double coreRadius = 0.;
+ for (cell = myList->begin(); cell != myList->end(); ++cell) {
+ int sampling = -1;
+ if (isubcalo == 0) {
+ sampling = m_tileID->sample((**cell).ID());
+ } else if (isubcalo == 1) {
+ sampling = m_hecID->sampling((**cell).ID());
+ } else if (isubcalo == 2) {
+ sampling = m_emID->sampling((**cell).ID());
+ }
+ if (*cell != lastCell && sampling != coreSampling) continue;
+
+ double deltaEta = (**cell).eta() - muEta;
+ double deltaPhi = (**cell).phi() - muPhi;
+ double radius = deltaEta * deltaEta + deltaPhi * deltaPhi;
+
+ if (sampling == coreSampling) {
+ if (!coreCell || radius < coreRadius) {
+ coreCell = *cell;
+ coreRadius = radius;
+ }
+ }
+
+ if (*cell != lastCell || !coreCell) continue;
+
+ for (CaloCellList::list_iterator cell2 = myList->begin(); cell2 != myList->end(); ++cell2) {
+ if (isubcalo == 0) {
+ sampling = m_tileID->sample((**cell2).ID());
+ } else if (isubcalo == 1) {
+ sampling = m_hecID->sampling((**cell2).ID());
+ } else if (isubcalo == 2) {
+ sampling = m_emID->sampling((**cell2).ID());
+ }
+ if (sampling != coreSampling) continue;
+
+ double cellEnergy = (**cell2).energy();
+ double noiseRms = noiseCDO->getNoise((**cell2).ID(), (**cell2).gain());
+ // looser selection for core cell where at least mip is expected
+ bool cellSelected = cellEnergy > m_sigmasAboveNoise * noiseRms;
+ if (*cell2 == coreCell && cellEnergy > m_sigmasAboveNoiseCore * noiseRms) cellSelected = true;
+
+ if (cellSelected) totalEnergy += cellEnergy;
+ }
+ }
+ }
+
+ // store result in caloMeas
+ if (isubcalo == 0) {
+ caloMeas.Tile_Isolation(totalEnergy);
+ } else if (isubcalo == 1) {
+ caloMeas.LArHEC_Isolation(totalEnergy);
+ } else if (isubcalo == 2) {
+ caloMeas.LArEM_Isolation(totalEnergy);
+ }
}
delete myList;
}
double
-MuidCaloEnergyMeas::energyInTile(const CaloCellContainer* cellContainer,
- double mu_eta,
- double mu_phi,
- int sample,
- int cone) const
+MuidCaloEnergyMeas::energyInTile(const CaloCellContainer* cellContainer, double mu_eta, double mu_phi, int sample,
+ int cone) const
{
// Tile Cal
// sample_number
@@ -651,72 +499,59 @@ MuidCaloEnergyMeas::energyInTile(const CaloCellContainer* cellContainer,
// 1 --> Sampling 2
// 2 --> Sampling 3
// 3 --> ITC
-
- //int i,j,k;
- SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
+ // int i,j,k;
- double tileTotalEnergy=0.;
- double tileTestEnergy=0.;
-
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
- CaloCellList myList(cellContainer,iCalo); // Construct the list
- if(cone == 1)
- {
- myList.select(mu_eta,mu_phi,0.15,0.15);
- }
- else if(cone == 2)
- {
- myList.select(mu_eta,mu_phi,0.3,0.3);
- }
- else
- {
- myList.select(mu_eta,mu_phi,0.,0.);
+ SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
+ const CaloNoise* noiseCDO = *noiseHdl;
+
+ double tileTotalEnergy = 0.;
+ double tileTestEnergy = 0.;
+
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::TILE;
+ CaloCellList myList(cellContainer, iCalo); // Construct the list
+ if (cone == 1) {
+ myList.select(mu_eta, mu_phi, 0.15, 0.15);
+ } else if (cone == 2) {
+ myList.select(mu_eta, mu_phi, 0.3, 0.3);
+ } else {
+ myList.select(mu_eta, mu_phi, 0., 0.);
}
-
+
CaloCellList::list_iterator ilistfirst = myList.begin();
CaloCellList::list_iterator ilistlast = myList.end();
-
- int count=0;
- double lowest_threshold=0.;
- for(;ilistfirst!=ilistlast;ilistfirst++)
- {
- double cellEnergy= (*ilistfirst)->energy();
- double noise_rms=noiseCDO->getNoise((*ilistfirst)->ID(),(*ilistfirst)->gain());
-
- if( cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise
- && m_tileID->sample((*ilistfirst)->ID()) == sample)
- {
- count+=1;
-
- ATH_MSG_DEBUG( "Energy : " << (*ilistfirst)->energy()
- << " Sampling: " << m_tileID->sample((*ilistfirst)->ID())
- << " Radius :" << (*ilistfirst)->caloDDE()->r()
- << " Eta : " << (*ilistfirst)->eta()
- << " Phi : " << (*ilistfirst)->phi()
- << " Noise Level : " << noise_rms );
-
- tileTotalEnergy += cellEnergy;
- }
- else
- {
- tileTestEnergy += cellEnergy;
- }
+
+ int count = 0;
+ double lowest_threshold = 0.;
+ for (; ilistfirst != ilistlast; ilistfirst++) {
+ double cellEnergy = (*ilistfirst)->energy();
+ double noise_rms = noiseCDO->getNoise((*ilistfirst)->ID(), (*ilistfirst)->gain());
+
+ if (cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise
+ && m_tileID->sample((*ilistfirst)->ID()) == sample)
+ {
+ count += 1;
+
+ ATH_MSG_DEBUG("Energy : " << (*ilistfirst)->energy()
+ << " Sampling: " << m_tileID->sample((*ilistfirst)->ID()) << " Radius :"
+ << (*ilistfirst)->caloDDE()->r() << " Eta : " << (*ilistfirst)->eta()
+ << " Phi : " << (*ilistfirst)->phi() << " Noise Level : " << noise_rms);
+
+ tileTotalEnergy += cellEnergy;
+ } else {
+ tileTestEnergy += cellEnergy;
+ }
}
-
- ATH_MSG_DEBUG( " counted " << count << " test energy " << tileTestEnergy );
+
+ ATH_MSG_DEBUG(" counted " << count << " test energy " << tileTestEnergy);
return tileTotalEnergy;
}
double
-MuidCaloEnergyMeas::energyInLArHEC(const CaloCellContainer* cellContainer,
- double mu_eta,
- double mu_phi,
- int sample,
- int cone) const
-{
+MuidCaloEnergyMeas::energyInLArHEC(const CaloCellContainer* cellContainer, double mu_eta, double mu_phi, int sample,
+ int cone) const
+{
// Look in the LarHEC calorimeter
// i.e. loop over its cells
/* sample_number
@@ -726,134 +561,113 @@ MuidCaloEnergyMeas::energyInLArHEC(const CaloCellContainer* cellContainer,
3 --> Sampling 4
*/
- //int i,j,k;
+ // int i,j,k;
SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
- double larhecTotal=0.;
-
+ const CaloNoise* noiseCDO = *noiseHdl;
+ double larhecTotal = 0.;
+
CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LARHEC;
- CaloCellList myList(cellContainer,iCalo); // Construct the list
- if(cone == 1)
- {
- myList.select(mu_eta,mu_phi,0.15,0.15);
- }
- else if(cone == 2)
- {
- myList.select(mu_eta,mu_phi,0.3,0.3);
+ CaloCellList myList(cellContainer, iCalo); // Construct the list
+ if (cone == 1) {
+ myList.select(mu_eta, mu_phi, 0.15, 0.15);
+ } else if (cone == 2) {
+ myList.select(mu_eta, mu_phi, 0.3, 0.3);
+ } else {
+ myList.select(mu_eta, mu_phi, 0., 0.);
}
- else
- {
- myList.select(mu_eta,mu_phi,0.,0.);
- }
-
+
CaloCellList::list_iterator ilistfirst = myList.begin();
CaloCellList::list_iterator ilistlast = myList.end();
- //std::vector<const CaloCell*> new_cell_list;
- int count=0;
- constexpr double lowest_threshold = 4.*150.;
-
- for(;ilistfirst!=ilistlast;ilistfirst++)
- {
- double cellEnergy= (*ilistfirst)->energy();
- const double noise_rms =noiseCDO->getNoise((*ilistfirst)->ID(),(*ilistfirst)->gain());
- if (cellEnergy > lowest_threshold
- && cellEnergy > noise_rms * m_sigmasAboveNoise
- && m_hecID->sampling((*ilistfirst)->ID()) == sample)
- {
- count+=1;
- //new_cell_list.push_back(*ilistfirst);
-
- ATH_MSG_DEBUG( "Energy : " << (*ilistfirst)->energy()
- << " Sampling: " << m_hecID->sampling((*ilistfirst)->ID())
- << " z :" << (*ilistfirst)->caloDDE()->z()
- << " Eta : " << (*ilistfirst)->eta()
- << " Phi : " << (*ilistfirst)->phi()
- << " Noise Level : " << noise_rms);
-
- larhecTotal += (*ilistfirst)->energy();
- }
+ // std::vector<const CaloCell*> new_cell_list;
+ int count = 0;
+ constexpr double lowest_threshold = 4. * 150.;
+
+ for (; ilistfirst != ilistlast; ilistfirst++) {
+ double cellEnergy = (*ilistfirst)->energy();
+ const double noise_rms = noiseCDO->getNoise((*ilistfirst)->ID(), (*ilistfirst)->gain());
+ if (cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise
+ && m_hecID->sampling((*ilistfirst)->ID()) == sample)
+ {
+ count += 1;
+ // new_cell_list.push_back(*ilistfirst);
+
+ ATH_MSG_DEBUG("Energy : " << (*ilistfirst)->energy()
+ << " Sampling: " << m_hecID->sampling((*ilistfirst)->ID())
+ << " z :" << (*ilistfirst)->caloDDE()->z() << " Eta : " << (*ilistfirst)->eta()
+ << " Phi : " << (*ilistfirst)->phi() << " Noise Level : " << noise_rms);
+
+ larhecTotal += (*ilistfirst)->energy();
+ }
}
-
- ATH_MSG_DEBUG( "larhec counted " << count );
-
+
+ ATH_MSG_DEBUG("larhec counted " << count);
+
return larhecTotal;
}
double
-MuidCaloEnergyMeas::energyInLArEM(const CaloCellContainer* cellContainer,
- double mu_eta,
- double mu_phi,
- int sample,
- int cone) const
-{
- // Look in the LarEM calorimeter
- // i.e. loop over its cells
- /*
- sample_number
+MuidCaloEnergyMeas::energyInLArEM(const CaloCellContainer* cellContainer, double mu_eta, double mu_phi, int sample,
+ int cone) const
+{
+ // Look in the LarEM calorimeter
+ // i.e. loop over its cells
+ /*
+ sample_number
0 --> Presampler
1 --> Sampling 1
2 --> Sampling 2
- 3 --> Sampling 3
+ 3 --> Sampling 3
*/
- //int i,j,k;
+ // int i,j,k;
SG::ReadCondHandle<CaloNoise> noiseHdl{m_noiseCDOKey};
- const CaloNoise* noiseCDO=*noiseHdl;
+ const CaloNoise* noiseCDO = *noiseHdl;
- double emTotalEnergy=0.;
+ double emTotalEnergy = 0.;
- CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
- CaloCellList myList(cellContainer,iCalo); // Construct the list
- if(cone == 1)
- {
- myList.select(mu_eta,mu_phi,0.075,0.075); // 0.1 0.1
- }
- else if(cone == 2)
- {
- myList.select(mu_eta,mu_phi,0.15,0.15);
- }
- else
- {
- myList.select(mu_eta,mu_phi,0.,0.);
+ CaloCell_ID::SUBCALO iCalo = CaloCell_ID::LAREM;
+ CaloCellList myList(cellContainer, iCalo); // Construct the list
+ if (cone == 1) {
+ myList.select(mu_eta, mu_phi, 0.075, 0.075); // 0.1 0.1
+ } else if (cone == 2) {
+ myList.select(mu_eta, mu_phi, 0.15, 0.15);
+ } else {
+ myList.select(mu_eta, mu_phi, 0., 0.);
}
CaloCellList::list_iterator ilistfirst = myList.begin();
CaloCellList::list_iterator ilistlast = myList.end();
-
- //std::vector<const CaloCell*> new_cell_list;
- int count=0;
-
- double lowest_threshold = 4.* 50.;
-
- for(;ilistfirst!=ilistlast;ilistfirst++)
- {
- double cellEnergy= (*ilistfirst)->energy();
- const double noise_rms =noiseCDO->getNoise((*ilistfirst)->ID(),(*ilistfirst)->gain());
-
- if( cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise
- && m_emID->sampling((*ilistfirst)->ID()) == sample)
- {
- count+=1;
- //new_cell_list.push_back(*ilistfirst);
-
- ATH_MSG_DEBUG( "Energy : " << (*ilistfirst)->energy()
- << " Sampling: " << m_emID->sampling((*ilistfirst)->ID())
- << " Radius :" << (*ilistfirst)->caloDDE()->r()
- << " z :" << (*ilistfirst)->caloDDE()->z()
- << " Eta : " << (*ilistfirst)->eta()
- << " Phi : " << (*ilistfirst)->phi()
- << " Noise Level : " << noise_rms);
-
- emTotalEnergy+=(*ilistfirst)->energy();
- }
+
+ // std::vector<const CaloCell*> new_cell_list;
+ int count = 0;
+
+ double lowest_threshold = 4. * 50.;
+
+ for (; ilistfirst != ilistlast; ilistfirst++) {
+ double cellEnergy = (*ilistfirst)->energy();
+ const double noise_rms = noiseCDO->getNoise((*ilistfirst)->ID(), (*ilistfirst)->gain());
+
+ if (cellEnergy > lowest_threshold && cellEnergy > noise_rms * m_sigmasAboveNoise
+ && m_emID->sampling((*ilistfirst)->ID()) == sample)
+ {
+ count += 1;
+ // new_cell_list.push_back(*ilistfirst);
+
+ ATH_MSG_DEBUG("Energy : " << (*ilistfirst)->energy()
+ << " Sampling: " << m_emID->sampling((*ilistfirst)->ID())
+ << " Radius :" << (*ilistfirst)->caloDDE()->r()
+ << " z :" << (*ilistfirst)->caloDDE()->z() << " Eta : " << (*ilistfirst)->eta()
+ << " Phi : " << (*ilistfirst)->phi() << " Noise Level : " << noise_rms);
+
+ emTotalEnergy += (*ilistfirst)->energy();
+ }
}
- ATH_MSG_DEBUG( "larem counted " << count );
-
+ ATH_MSG_DEBUG("larem counted " << count);
+
return emTotalEnergy;
}
-} // end of namespace
-
+} // namespace Rec
diff --git a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyTool.cxx b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyTool.cxx
index 98bb45031a0d4a6752aa53768e85bbb8b83d3bb3..c146b92a59f7464df60be2ac97c5868f345d95dd 100755
--- a/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyTool.cxx
+++ b/Reconstruction/MuonIdentification/MuidCaloEnergyTools/src/MuidCaloEnergyTool.cxx
@@ -1,5 +1,5 @@
/*
- Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+ Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
*/
//////////////////////////////////////////////////////////////////////////////
@@ -19,8 +19,11 @@
//<<<<<< INCLUDES >>>>>>
-#include <cmath>
#include "MuidCaloEnergyTools/MuidCaloEnergyTool.h"
+
+#include <cmath>
+
+#include "AthenaKernel/Units.h"
#include "MuidEvent/CaloMeas.h"
#include "MuidInterfaces/IMuidCaloEnergyMeas.h"
#include "MuidInterfaces/IMuidCaloEnergyParam.h"
@@ -29,122 +32,81 @@
#include "TrkMaterialOnTrack/MaterialEffectsOnTrack.h"
#include "TrkTrack/TrackStateOnSurface.h"
#include "muonEvent/CaloEnergy.h"
-#include "AthenaKernel/Units.h"
namespace Units = Athena::Units;
//<<<<<< CLASS STRUCTURE INITIALIZATION >>>>>>
-namespace Rec
-{
-
-MuidCaloEnergyTool::MuidCaloEnergyTool (const std::string&type,
- const std::string&name,
- const IInterface*parent)
- : AthAlgTool (type, name, parent),
- m_caloMeasTool ("Rec::MuidCaloEnergyMeas/MuidCaloEnergyMeas", this),
- m_caloParamTool ("Rec::MuidCaloEnergyParam/MuidCaloEnergyParam", this),
- m_trackIsolationTool ("Rec::MuidTrackIsolation/MuidTrackIsolation", this),
- m_cosmics (false),
- m_energyLossMeasurement (true),
- m_forceIsolationFailure (false),
- m_FSRtreatment (true),
- m_MOPparametrization (true),
- m_trackIsolation (false),
- m_emEtCut (2.5*Units::GeV),
- m_emF1Cut (0.15),
- m_emMinEnergy ( 2.*Units::GeV),
- m_hecMinEnergy (10.*Units::GeV),
- m_maxNTracksIso (2),
- m_minFinalEnergy ( 0.*Units::GeV),
- m_minMuonPt (15.*Units::GeV),
- m_countMean (0),
- m_countMeasurement (0),
- m_countMop (0)
+namespace Rec {
+
+MuidCaloEnergyTool::MuidCaloEnergyTool(const std::string& type, const std::string& name, const IInterface* parent)
+ : AthAlgTool(type, name, parent),
+ m_cosmics(false),
+ m_energyLossMeasurement(true),
+ m_forceIsolationFailure(false),
+ m_FSRtreatment(true),
+ m_MOPparametrization(true),
+ m_trackIsolation(false),
+ m_emEtCut(2.5 * Units::GeV),
+ m_emF1Cut(0.15),
+ m_emMinEnergy(2. * Units::GeV),
+ m_hecMinEnergy(10. * Units::GeV),
+ m_maxNTracksIso(2),
+ m_minFinalEnergy(0. * Units::GeV),
+ m_minMuonPt(15. * Units::GeV),
+ m_countMean(0),
+ m_countMeasurement(0),
+ m_countMop(0)
{
declareInterface<IMuidCaloEnergy>(this);
- declareProperty ("CaloMeasTool", m_caloMeasTool);
- declareProperty ("CaloParamTool", m_caloParamTool);
- declareProperty ("TrackIsolationTool", m_trackIsolationTool);
- declareProperty ("Cosmics", m_cosmics);
- declareProperty ("EnergyLossMeasurement", m_energyLossMeasurement);
- declareProperty ("ForceIsolationFailure", m_forceIsolationFailure);
- declareProperty ("FSRtreatment", m_FSRtreatment);
- declareProperty ("MopParametrization", m_MOPparametrization);
- declareProperty ("TrackIsolation", m_trackIsolation);
- declareProperty ("EmEtCut", m_emEtCut);
- declareProperty ("EmF1Cut", m_emF1Cut);
- declareProperty ("EmMinEnergy", m_emMinEnergy);
- declareProperty ("HecMinEnergy", m_hecMinEnergy);
- declareProperty ("MaxNTracksIso", m_maxNTracksIso);
- declareProperty ("MinFinalEnergy", m_minFinalEnergy);
- declareProperty ("MinMuonPt", m_minMuonPt);
+ declareProperty("Cosmics", m_cosmics);
+ declareProperty("EnergyLossMeasurement", m_energyLossMeasurement);
+ declareProperty("ForceIsolationFailure", m_forceIsolationFailure);
+ declareProperty("FSRtreatment", m_FSRtreatment);
+ declareProperty("MopParametrization", m_MOPparametrization);
+ declareProperty("TrackIsolation", m_trackIsolation);
+ declareProperty("EmEtCut", m_emEtCut);
+ declareProperty("EmF1Cut", m_emF1Cut);
+ declareProperty("EmMinEnergy", m_emMinEnergy);
+ declareProperty("HecMinEnergy", m_hecMinEnergy);
+ declareProperty("MaxNTracksIso", m_maxNTracksIso);
+ declareProperty("MinFinalEnergy", m_minFinalEnergy);
+ declareProperty("MinMuonPt", m_minMuonPt);
}
-MuidCaloEnergyTool::~MuidCaloEnergyTool (void)
-{}
+MuidCaloEnergyTool::~MuidCaloEnergyTool(void) {}
//<<<<<< PUBLIC MEMBER FUNCTION DEFINITIONS >>>>>>
-StatusCode MuidCaloEnergyTool::initialize()
+StatusCode
+MuidCaloEnergyTool::initialize()
{
- if (m_energyLossMeasurement)
- {
- ATH_MSG_INFO( "Initializing MuidCaloEnergyTool"
- << " - package version " << PACKAGE_VERSION
- << " - measured calorimeter energy deposition " );
- }
- else
- {
- ATH_MSG_INFO( "Initializing MuidCaloEnergyTool"
- << " - package version " << PACKAGE_VERSION
- << " - parametrized calorimeter energy deposition " );
+ if (m_energyLossMeasurement) {
+ ATH_MSG_INFO("Initializing MuidCaloEnergyTool"
+ << " - package version " << PACKAGE_VERSION << " - measured calorimeter energy deposition ");
+ } else {
+ ATH_MSG_INFO("Initializing MuidCaloEnergyTool"
+ << " - package version " << PACKAGE_VERSION << " - parametrized calorimeter energy deposition ");
}
// get the Tools
- if (m_caloParamTool.retrieve().isFailure())
- {
- ATH_MSG_FATAL( "Failed to retrieve tool " << m_caloParamTool );
- return StatusCode::FAILURE;
- }
- else
- {
- ATH_MSG_INFO( "Retrieved tool " << m_caloParamTool );
- }
-
- if (m_energyLossMeasurement)
- {
- if (m_caloMeasTool.retrieve().isFailure())
- {
- ATH_MSG_FATAL( "Failed to retrieve tool " << m_caloMeasTool );
- return StatusCode::FAILURE;
- }
- else
- {
- ATH_MSG_INFO( "Retrieved tool " << m_caloMeasTool );
- }
-
- if (m_trackIsolation)
- {
- if (m_trackIsolationTool.retrieve().isFailure())
- {
- ATH_MSG_FATAL( "Failed to retrieve tool " << m_trackIsolationTool );
- return StatusCode::FAILURE;
- }
- else
- {
- ATH_MSG_INFO( "Retrieved tool " << m_trackIsolationTool );
- }
- }
- else
- {
- ATH_MSG_WARNING( " Using energy measurement without trackIsolation " );
- m_trackIsolationTool.disable();
- }
- }
- else{
- m_trackIsolationTool.disable();
- m_caloMeasTool.disable();
+ ATH_CHECK(m_caloParamTool.retrieve());
+ ATH_MSG_INFO("Retrieved tool " << m_caloParamTool);
+
+ if (m_energyLossMeasurement) {
+ ATH_CHECK(m_caloMeasTool.retrieve());
+ ATH_MSG_INFO("Retrieved tool " << m_caloMeasTool);
+
+ if (m_trackIsolation) {
+ ATH_CHECK(m_trackIsolationTool.retrieve());
+ ATH_MSG_INFO("Retrieved tool " << m_trackIsolationTool);
+ } else {
+ ATH_MSG_WARNING(" Using energy measurement without trackIsolation ");
+ m_trackIsolationTool.disable();
+ }
+ } else {
+ m_trackIsolationTool.disable();
+ m_caloMeasTool.disable();
}
return StatusCode::SUCCESS;
@@ -153,89 +115,73 @@ StatusCode MuidCaloEnergyTool::initialize()
StatusCode
MuidCaloEnergyTool::finalize()
{
- ATH_MSG_INFO( "Finalizing MuidCaloEnergyTool."
- << " Tracks used mean: " << m_countMean
- << ", tracks used mop: " << m_countMop
- << ", tracks used measurement: " << m_countMeasurement );
-
+ ATH_MSG_INFO("Finalizing MuidCaloEnergyTool."
+ << " Tracks used mean: " << m_countMean << ", tracks used mop: " << m_countMop
+ << ", tracks used measurement: " << m_countMeasurement);
+
return StatusCode::SUCCESS;
}
CaloEnergy*
-MuidCaloEnergyTool::energyLoss(double trackMomentum,
- double eta,
- double phi) const
+MuidCaloEnergyTool::energyLoss(double trackMomentum, double eta, double phi) const
{
// debug
- ATH_MSG_VERBOSE( "Muon with : p = " << trackMomentum/Units::GeV
- << " Phi = " << phi
- << " Eta = " << eta );
+ ATH_MSG_VERBOSE("Muon with : p = " << trackMomentum / Units::GeV << " Phi = " << phi << " Eta = " << eta);
+
-
CaloEnergy* caloEnergy = 0;
- if (m_energyLossMeasurement)
- {
- // Energy Dep/Iso from calorimeters (projective assumption)
- CaloMeas* caloMeas = m_caloMeasTool->energyMeasurement(eta,phi,eta,phi);
- if (caloMeas)
- {
- caloEnergy = measurement(trackMomentum,eta,phi,caloMeas);
- delete caloMeas;
- }
+ if (m_energyLossMeasurement) {
+ // Energy Dep/Iso from calorimeters (projective assumption)
+ CaloMeas* caloMeas = m_caloMeasTool->energyMeasurement(eta, phi, eta, phi);
+ if (caloMeas) {
+ caloEnergy = measurement(trackMomentum, eta, phi, caloMeas);
+ delete caloMeas;
+ }
}
- if (! caloEnergy)
- {
- if (m_MOPparametrization)
- {
- ++m_countMop;
- caloEnergy = m_caloParamTool->mopParametrizedEnergy(trackMomentum,eta,phi);
- ATH_MSG_VERBOSE( "Selected the Mop Parametrization value! ==> " );
- }
- else
- {
- ++m_countMean;
- caloEnergy = m_caloParamTool->meanParametrizedEnergy(trackMomentum,eta,phi);
- ATH_MSG_VERBOSE( "Selected the Mean Parametrization value! ==> " );
- }
+ if (!caloEnergy) {
+ if (m_MOPparametrization) {
+ ++m_countMop;
+ caloEnergy = m_caloParamTool->mopParametrizedEnergy(trackMomentum, eta, phi);
+ ATH_MSG_VERBOSE("Selected the Mop Parametrization value! ==> ");
+ } else {
+ ++m_countMean;
+ caloEnergy = m_caloParamTool->meanParametrizedEnergy(trackMomentum, eta, phi);
+ ATH_MSG_VERBOSE("Selected the Mean Parametrization value! ==> ");
+ }
}
-
- if (msgLvl(MSG::DEBUG))
- {
- std::string eLossType = " no Calo !!";
- switch (caloEnergy->energyLossType())
- {
- case CaloEnergy::Parametrized:
- eLossType = "Parametrized";
- break;
- case CaloEnergy::NotIsolated:
- eLossType = "NotIsolated ";
- break;
- case CaloEnergy::MOP:
- eLossType = "MOP ";
- break;
- case CaloEnergy::Tail:
- eLossType = "Tail ";
- break;
- case CaloEnergy::FSRcandidate:
- eLossType = "FSRcandidate";
- break;
- default:
- break;
- };
- ATH_MSG_DEBUG( std::setiosflags(std::ios::fixed) << " energyLoss with"
- << " momentum =" << std::setw(6) << std::setprecision(1) << trackMomentum/Units::GeV
- << " phi =" << std::setw(6) << std::setprecision(3) << phi
- << " eta =" << std::setw(6) << std::setprecision(3) << eta
- << ". CaloEnergy: deltaE = " << std::setw(8) << std::setprecision(3)
- << caloEnergy->deltaE()/Units::GeV
- << " +" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaPlusDeltaE()/Units::GeV
- << " -" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaMinusDeltaE()/Units::GeV
- << " (" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaDeltaE()/Units::GeV
- << ") GeV, CaloEnergy::Type " << eLossType );
+
+ if (msgLvl(MSG::DEBUG)) {
+ std::string eLossType = " no Calo !!";
+ switch (caloEnergy->energyLossType()) {
+ case CaloEnergy::Parametrized:
+ eLossType = "Parametrized";
+ break;
+ case CaloEnergy::NotIsolated:
+ eLossType = "NotIsolated ";
+ break;
+ case CaloEnergy::MOP:
+ eLossType = "MOP ";
+ break;
+ case CaloEnergy::Tail:
+ eLossType = "Tail ";
+ break;
+ case CaloEnergy::FSRcandidate:
+ eLossType = "FSRcandidate";
+ break;
+ default:
+ break;
+ };
+ ATH_MSG_DEBUG(std::setiosflags(std::ios::fixed)
+ << " energyLoss with"
+ << " momentum =" << std::setw(6) << std::setprecision(1) << trackMomentum / Units::GeV
+ << " phi =" << std::setw(6) << std::setprecision(3) << phi << " eta =" << std::setw(6)
+ << std::setprecision(3) << eta << ". CaloEnergy: deltaE = " << std::setw(8)
+ << std::setprecision(3) << caloEnergy->deltaE() / Units::GeV << " +" << std::setw(5)
+ << std::setprecision(3) << caloEnergy->sigmaPlusDeltaE() / Units::GeV << " -" << std::setw(5)
+ << std::setprecision(3) << caloEnergy->sigmaMinusDeltaE() / Units::GeV << " (" << std::setw(5)
+ << std::setprecision(3) << caloEnergy->sigmaDeltaE() / Units::GeV << ") GeV, CaloEnergy::Type "
+ << eLossType);
}
return caloEnergy;
@@ -243,216 +189,182 @@ MuidCaloEnergyTool::energyLoss(double trackMomentum,
const Trk::TrackStateOnSurface*
MuidCaloEnergyTool::trackStateOnSurface(const Trk::TrackParameters& middleParameters,
- const Trk::TrackParameters* innerParameters,
- const Trk::TrackParameters* outerParameters) const
+ const Trk::TrackParameters* innerParameters,
+ const Trk::TrackParameters* outerParameters) const
{
- ATH_MSG_VERBOSE( "Muon with : p = " << middleParameters.momentum().mag()/Units::GeV
- << " Phi = " << middleParameters.position().phi()
- << " Eta = " << middleParameters.position().eta() );
+ ATH_MSG_VERBOSE("Muon with : p = " << middleParameters.momentum().mag() / Units::GeV
+ << " Phi = " << middleParameters.position().phi()
+ << " Eta = " << middleParameters.position().eta());
CaloEnergy* caloEnergy = 0;
- if (m_energyLossMeasurement)
- {
- // energy deposition according to the calo measurement
- double eta = middleParameters.position().eta();
- double phi = middleParameters.position().phi();
- double etaEM = eta;
- double phiEM = phi;
- double etaHad = eta;
- double phiHad = phi;
- if (innerParameters)
- {
- etaEM = innerParameters->position().eta();
- phiEM = innerParameters->position().phi();
- }
- if (outerParameters)
- {
- etaHad = outerParameters->position().eta();
- phiHad = outerParameters->position().phi();
- }
- CaloMeas* caloMeas = m_caloMeasTool->energyMeasurement(etaEM,phiEM,etaHad,phiHad);
- if (caloMeas)
- {
- caloEnergy = measurement(middleParameters.momentum().mag(),
- eta,
- phi,
- caloMeas);
- delete caloMeas;
- }
+ if (m_energyLossMeasurement) {
+ // energy deposition according to the calo measurement
+ double eta = middleParameters.position().eta();
+ double phi = middleParameters.position().phi();
+ double etaEM = eta;
+ double phiEM = phi;
+ double etaHad = eta;
+ double phiHad = phi;
+ if (innerParameters) {
+ etaEM = innerParameters->position().eta();
+ phiEM = innerParameters->position().phi();
+ }
+ if (outerParameters) {
+ etaHad = outerParameters->position().eta();
+ phiHad = outerParameters->position().phi();
+ }
+ CaloMeas* caloMeas = m_caloMeasTool->energyMeasurement(etaEM, phiEM, etaHad, phiHad);
+ if (caloMeas) {
+ caloEnergy = measurement(middleParameters.momentum().mag(), eta, phi, caloMeas);
+ delete caloMeas;
+ }
}
-
- if (! caloEnergy)
- {
- // parametrized energy deposition
- caloEnergy = energyLoss(middleParameters.momentum().mag(),
- middleParameters.position().eta(),
- middleParameters.position().phi());
-
- // WARN in case of anomalously high loss
- if (caloEnergy->deltaE() > 8.*Units::GeV && middleParameters.momentum().mag() < 500.*Units::GeV)
- ATH_MSG_WARNING(" high parametrized energy loss: "
- << caloEnergy->deltaE()/Units::GeV << " GeV"
- << " for p " << middleParameters.momentum().mag()/Units::GeV << " GeV"
- << " eta " << middleParameters.position().eta() );
+
+ if (!caloEnergy) {
+ // parametrized energy deposition
+ caloEnergy = energyLoss(middleParameters.momentum().mag(), middleParameters.position().eta(),
+ middleParameters.position().phi());
+
+ // WARN in case of anomalously high loss
+ if (caloEnergy->deltaE() > 8. * Units::GeV && middleParameters.momentum().mag() < 500. * Units::GeV)
+ ATH_MSG_WARNING(" high parametrized energy loss: "
+ << caloEnergy->deltaE() / Units::GeV << " GeV"
+ << " for p " << middleParameters.momentum().mag() / Units::GeV << " GeV"
+ << " eta " << middleParameters.position().eta());
}
// create MEOT owning CaloEnergy
std::bitset<Trk::MaterialEffectsBase::NumberOfMaterialEffectsTypes> typePattern(0);
typePattern.set(Trk::MaterialEffectsBase::EnergyLossEffects);
- const Trk::MaterialEffectsOnTrack* materialEffects =
- new const Trk::MaterialEffectsOnTrack(0.,
- caloEnergy,
- middleParameters.associatedSurface(),
- typePattern);
-
+ const Trk::MaterialEffectsOnTrack* materialEffects =
+ new const Trk::MaterialEffectsOnTrack(0., caloEnergy, middleParameters.associatedSurface(), typePattern);
+
// create TSOS
- const Trk::FitQualityOnSurface* fitQoS = 0;
- const Trk::MeasurementBase* measurementBase = 0;
+ const Trk::FitQualityOnSurface* fitQoS = 0;
+ const Trk::MeasurementBase* measurementBase = 0;
std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> pattern(0);
pattern.set(Trk::TrackStateOnSurface::CaloDeposit);
// debugging
- if (msgLvl(MSG::DEBUG))
- {
- std::string eLossType = " no Calo !!";
- switch (caloEnergy->energyLossType())
- {
- case CaloEnergy::Parametrized:
- eLossType = "Parametrized";
- break;
- case CaloEnergy::NotIsolated:
- eLossType = "NotIsolated ";
- break;
- case CaloEnergy::MOP:
- eLossType = "MOP ";
- break;
- case CaloEnergy::Tail:
- eLossType = "Tail ";
- break;
- case CaloEnergy::FSRcandidate:
- eLossType = "FSRcandidate";
- break;
- default:
- break;
- };
-
- ATH_MSG_DEBUG( std::setiosflags(std::ios::fixed) << " trackStateOnSurface with"
- << " momentum =" << std::setw(6) << std::setprecision(1)
- << middleParameters.momentum().mag()/Units::GeV
- << " phi =" << std::setw(6) << std::setprecision(3)
- << middleParameters.position().phi()
- << " eta =" << std::setw(6) << std::setprecision(3)
- << middleParameters.position().eta()
- << ". CaloEnergy: deltaE = " << std::setw(8) << std::setprecision(3)
- << caloEnergy->deltaE()/Units::GeV
- << " +" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaPlusDeltaE()/Units::GeV
- << " -" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaMinusDeltaE()/Units::GeV
- << " (" << std::setw(5) << std::setprecision(3)
- << caloEnergy->sigmaDeltaE()/Units::GeV
- << ") GeV, CaloEnergy::Type " << eLossType );
- }
-
- return new const Trk::TrackStateOnSurface(measurementBase,
- &middleParameters,
- fitQoS,
- materialEffects,
- pattern);
+ if (msgLvl(MSG::DEBUG)) {
+ std::string eLossType = " no Calo !!";
+ switch (caloEnergy->energyLossType()) {
+ case CaloEnergy::Parametrized:
+ eLossType = "Parametrized";
+ break;
+ case CaloEnergy::NotIsolated:
+ eLossType = "NotIsolated ";
+ break;
+ case CaloEnergy::MOP:
+ eLossType = "MOP ";
+ break;
+ case CaloEnergy::Tail:
+ eLossType = "Tail ";
+ break;
+ case CaloEnergy::FSRcandidate:
+ eLossType = "FSRcandidate";
+ break;
+ default:
+ break;
+ };
+
+ ATH_MSG_DEBUG(std::setiosflags(std::ios::fixed)
+ << " trackStateOnSurface with"
+ << " momentum =" << std::setw(6) << std::setprecision(1)
+ << middleParameters.momentum().mag() / Units::GeV << " phi =" << std::setw(6)
+ << std::setprecision(3) << middleParameters.position().phi() << " eta =" << std::setw(6)
+ << std::setprecision(3) << middleParameters.position().eta()
+ << ". CaloEnergy: deltaE = " << std::setw(8) << std::setprecision(3)
+ << caloEnergy->deltaE() / Units::GeV << " +" << std::setw(5) << std::setprecision(3)
+ << caloEnergy->sigmaPlusDeltaE() / Units::GeV << " -" << std::setw(5) << std::setprecision(3)
+ << caloEnergy->sigmaMinusDeltaE() / Units::GeV << " (" << std::setw(5) << std::setprecision(3)
+ << caloEnergy->sigmaDeltaE() / Units::GeV << ") GeV, CaloEnergy::Type " << eLossType);
+ }
+
+ return new const Trk::TrackStateOnSurface(measurementBase, &middleParameters, fitQoS, materialEffects, pattern);
}
//<<<<<< PRIVATE MEMBER FUNCTION DEFINITIONS >>>>>>
CaloEnergy*
-MuidCaloEnergyTool::measurement(double trackMomentum,
- double eta,
- double phi,
- CaloMeas* caloMeas) const
+MuidCaloEnergyTool::measurement(double trackMomentum, double eta, double phi, CaloMeas* caloMeas) const
{
// Mean Energy Loss parametrization
- CaloEnergy* caloParamMean = m_caloParamTool->meanParametrizedEnergy(trackMomentum,eta,phi);
+ CaloEnergy* caloParamMean = m_caloParamTool->meanParametrizedEnergy(trackMomentum, eta, phi);
// Mop Energy Loss parametrization
- CaloEnergy* caloParamMop = m_caloParamTool->mopParametrizedEnergy(trackMomentum,eta,phi);
+ CaloEnergy* caloParamMop = m_caloParamTool->mopParametrizedEnergy(trackMomentum, eta, phi);
// Mip Energy Loss
- CaloEnergy* caloParamMip = m_caloParamTool->meanParametrizedEnergy(10.*Units::GeV,eta,phi);
+ CaloEnergy* caloParamMip = m_caloParamTool->meanParametrizedEnergy(10. * Units::GeV, eta, phi);
// Mop Energy Loss peak
- CaloEnergy* mopPeak = m_caloParamTool->mopPeakEnergy(trackMomentum,eta,phi);
-
+ CaloEnergy* mopPeak = m_caloParamTool->mopPeakEnergy(trackMomentum, eta, phi);
+
// // flag upper hemisphere cosmic
// bool cosmic = false;
// if (caloParamMop->deltaE() < 0.) cosmic = true;
// mop energy deposition
- double MopLoss = fabs(caloParamMop->deltaE());
+ double MopLoss = fabs(caloParamMop->deltaE());
// mop energy deposition uncertainty
- double MopError = mopPeak->sigmaDeltaE();
+ double MopError = mopPeak->sigmaDeltaE();
// mop energy deposition corrected
- double MopLossCorrected = paramCorrection(trackMomentum,eta,MopLoss,MopError);
-
+ double MopLossCorrected = paramCorrection(trackMomentum, eta, MopLoss, MopError);
+
// percentage of inert material
- const double InertMaterial = m_caloParamTool->x0mapInertMaterial(eta);
+ const double InertMaterial = m_caloParamTool->x0mapInertMaterial(eta);
// percentage of em calorimeter material
- const double EmMaterial = m_caloParamTool->x0mapEmMaterial(eta);
+ const double EmMaterial = m_caloParamTool->x0mapEmMaterial(eta);
// percentage of hec calorimeter material
- const double HECMaterial = m_caloParamTool->x0mapHecMaterial(eta);
+ const double HECMaterial = m_caloParamTool->x0mapHecMaterial(eta);
// correction for the inert material
- double MaterialCorrection = InertMaterial * MopLossCorrected;
+ double MaterialCorrection = InertMaterial * MopLossCorrected;
// scale to get mop loss in em calo
- double MopLossEm = MopLoss*m_caloParamTool->emMopFraction(eta);
+ double MopLossEm = MopLoss * m_caloParamTool->emMopFraction(eta);
// fraction of Tile used for the measurement
- const double TileMeasurementMaterial = caloMeas->Tile_SamplingFraction();
+ const double TileMeasurementMaterial = caloMeas->Tile_SamplingFraction();
// fraction of LArHEC used for the measurement
- const double LArHECMeasurementMaterial = caloMeas->LArHEC_SamplingFraction();
+ const double LArHECMeasurementMaterial = caloMeas->LArHEC_SamplingFraction();
// fraction of LArEM used for the measurement
- const double LArEmMeasurementMaterial = caloMeas->LArEM_SamplingFraction();
+ const double LArEmMeasurementMaterial = caloMeas->LArEM_SamplingFraction();
// Measured energy deposition in Tile
- const double TileEnergy = caloMeas->Tile_EnergyMeasured();
+ const double TileEnergy = caloMeas->Tile_EnergyMeasured();
// Measured energy deposition in E/M
- const double EmEnergy = caloMeas->LArEM_EnergyMeasured();
+ const double EmEnergy = caloMeas->LArEM_EnergyMeasured();
// Correction for forward calorimetry
- double ForwardHECCorrection = 0.;
- if (fabs(eta)>2. && caloMeas->LArHEC_EnergyMeasured()>100.)
- ForwardHECCorrection = (1. - LArHECMeasurementMaterial) * HECMaterial * MopLossCorrected;
- const double LArHECEnergy = caloMeas->LArHEC_EnergyMeasured() + ForwardHECCorrection; // Measured energy deposition in LArHEC
-
- double TotalMeasuredEnergy = TileEnergy + EmEnergy + LArHECEnergy;
-
-
- ATH_MSG_VERBOSE( "Energy Deposition:Tile= " << TileEnergy/Units::GeV
- << ",LArHEC= " << LArHECEnergy/Units::GeV
- << ",EM= " << EmEnergy/Units::GeV
- << " ForwardHECCorrection " << ForwardHECCorrection/Units::GeV
- << " HECMaterial " << HECMaterial
- << " MopLossCorrected " << MopLossCorrected/Units::GeV );
-
- bool bHEC = false; // performed HEC measurement?
- bool bEM = false; // performed Em measurement?
-
+ double ForwardHECCorrection = 0.;
+ if (fabs(eta) > 2. && caloMeas->LArHEC_EnergyMeasured() > 100.)
+ ForwardHECCorrection = (1. - LArHECMeasurementMaterial) * HECMaterial * MopLossCorrected;
+ const double LArHECEnergy =
+ caloMeas->LArHEC_EnergyMeasured() + ForwardHECCorrection; // Measured energy deposition in LArHEC
+
+ double TotalMeasuredEnergy = TileEnergy + EmEnergy + LArHECEnergy;
+
+
+ ATH_MSG_VERBOSE("Energy Deposition:Tile= " << TileEnergy / Units::GeV << ",LArHEC= " << LArHECEnergy / Units::GeV
+ << ",EM= " << EmEnergy / Units::GeV << " ForwardHECCorrection "
+ << ForwardHECCorrection / Units::GeV << " HECMaterial " << HECMaterial
+ << " MopLossCorrected " << MopLossCorrected / Units::GeV);
+
+ bool bHEC = false; // performed HEC measurement?
+ bool bEM = false; // performed Em measurement?
+
// If muon isolated, and no significant measurement is made then use the mop parameterization, else the mean
- if (fabs(eta)<1.4)
- {
- if (LArHECEnergy + TileEnergy > 0.1 * MopLoss * HECMaterial) bHEC= true;
- }
- else if (fabs(eta)>1.8)
- {
- if (LArHECEnergy + TileEnergy > 0.2 * MopLoss * HECMaterial) bHEC= true;
- }
- else
- {
- if (LArHECEnergy + TileEnergy > 0.25 * MopLoss * HECMaterial) bHEC= true;
+ if (fabs(eta) < 1.4) {
+ if (LArHECEnergy + TileEnergy > 0.1 * MopLoss * HECMaterial) bHEC = true;
+ } else if (fabs(eta) > 1.8) {
+ if (LArHECEnergy + TileEnergy > 0.2 * MopLoss * HECMaterial) bHEC = true;
+ } else {
+ if (LArHECEnergy + TileEnergy > 0.25 * MopLoss * HECMaterial) bHEC = true;
}
- if (EmEnergy > 0.5 * MopLoss * EmMaterial) bEM = true;
+ if (EmEnergy > 0.5 * MopLoss * EmMaterial) bEM = true;
- if (!bHEC)
- {
-// TotalMeasuredEnergy -= TileEnergy + LArHECEnergy;
-// MaterialCorrection += HECMaterial * MopLossCorrected;
+ if (!bHEC) {
+ // TotalMeasuredEnergy -= TileEnergy + LArHECEnergy;
+ // MaterialCorrection += HECMaterial * MopLossCorrected;
}
- if (!bEM)
- {
-// TotalMeasuredEnergy -= EmEnergy;
-// MaterialCorrection += EmMaterial * MopLossCorrected;
+ if (!bEM) {
+ // TotalMeasuredEnergy -= EmEnergy;
+ // MaterialCorrection += EmMaterial * MopLossCorrected;
}
double MeasCorrected = TotalMeasuredEnergy + MaterialCorrection;
@@ -460,221 +372,173 @@ MuidCaloEnergyTool::measurement(double trackMomentum,
// Muons of 10 GeV are already in the relativistic rise region
// in order to obtain the mip deposition from the mean energy deposition of 10 GeV muons
// should divide by approximately 1.4 (Review of Particle Physics Figure 27.3 p.243)
- const double IonizationLoss = (1./1.4) * fabs(caloParamMip->deltaE());
-
- double eOverMipCorrectionEm = 0.;
- double eOverMipCorrectionHEC = 0.;
+ const double IonizationLoss = (1. / 1.4) * fabs(caloParamMip->deltaE());
+
+ double eOverMipCorrectionEm = 0.;
+ double eOverMipCorrectionHEC = 0.;
// Etrue = emip * Emeas
// -DE = Emeas - Etrue = Etrue ( 1./emip -1.)
- if (bEM)
- {
- const double emipEM = 0.78;
- eOverMipCorrectionEm = - (1./emipEM-1.) * IonizationLoss * EmMaterial * LArEmMeasurementMaterial;
- if (EmEnergy + eOverMipCorrectionEm<0.)eOverMipCorrectionEm=0.;
+ if (bEM) {
+ const double emipEM = 0.78;
+ eOverMipCorrectionEm = -(1. / emipEM - 1.) * IonizationLoss * EmMaterial * LArEmMeasurementMaterial;
+ if (EmEnergy + eOverMipCorrectionEm < 0.) eOverMipCorrectionEm = 0.;
}
- if (bHEC)
- {
- const double emipTile = 0.86;
- const double emipLAr = 0.94;
- const double HECEnergy = TileEnergy + LArHECEnergy;
- const double eOverMipCorrectionTile = - (1./emipTile-1.) * TileEnergy / HECEnergy * IonizationLoss * HECMaterial * TileMeasurementMaterial;
- const double eOverMipCorrectionLAr = - (1./emipLAr-1.) * LArHECEnergy / HECEnergy * IonizationLoss * HECMaterial * LArHECMeasurementMaterial;
- eOverMipCorrectionHEC = eOverMipCorrectionTile + eOverMipCorrectionLAr;
- if (LArHECEnergy + TileEnergy + eOverMipCorrectionHEC<0.)eOverMipCorrectionHEC=0.;
+ if (bHEC) {
+ const double emipTile = 0.86;
+ const double emipLAr = 0.94;
+ const double HECEnergy = TileEnergy + LArHECEnergy;
+ const double eOverMipCorrectionTile =
+ -(1. / emipTile - 1.) * TileEnergy / HECEnergy * IonizationLoss * HECMaterial * TileMeasurementMaterial;
+ const double eOverMipCorrectionLAr =
+ -(1. / emipLAr - 1.) * LArHECEnergy / HECEnergy * IonizationLoss * HECMaterial * LArHECMeasurementMaterial;
+ eOverMipCorrectionHEC = eOverMipCorrectionTile + eOverMipCorrectionLAr;
+ if (LArHECEnergy + TileEnergy + eOverMipCorrectionHEC < 0.) eOverMipCorrectionHEC = 0.;
}
- const double eOverMipCorrection = eOverMipCorrectionEm + eOverMipCorrectionHEC;
+ const double eOverMipCorrection = eOverMipCorrectionEm + eOverMipCorrectionHEC;
+
-
// additional offset from high-statistics Z->mumu MC (measured by Peter K 30/11/2011)
- double fix1FromPeter[26] = { 0.424104 , 0.479637 , 0.483419 , 0.490242 , 0.52806 ,
- 0.573582 , 0.822098 , 0.767301 , 0.809919 , 0.658745 ,
- 0.157187 , 0.413214 , 0.771074 , 0.61815 , 0.350113 ,
- 0.322785 , 0.479294 , 0.806183 , 0.822161 , 0.757731 ,
- -0.0857186, -0.0992693, -0.0492252, 0.0650174, 0.261538 ,
- 0.360413 };
+ double fix1FromPeter[26] = {0.424104, 0.479637, 0.483419, 0.490242, 0.52806, 0.573582, 0.822098,
+ 0.767301, 0.809919, 0.658745, 0.157187, 0.413214, 0.771074, 0.61815,
+ 0.350113, 0.322785, 0.479294, 0.806183, 0.822161, 0.757731, -0.0857186,
+ -0.0992693, -0.0492252, 0.0650174, 0.261538, 0.360413};
// (update from Peter K 09/12/2011)
- double fix2FromPeter[26] = { -0.647703 , -0.303498 , -0.268645 , -0.261292 , -0.260152 ,
- -0.269253 , -0.266212 , -0.240837 , -0.130172 , -0.111638 ,
- -0.329423 , -0.321011 , -0.346050 , -0.305592 , -0.313293 ,
- -0.317111 , -0.428393 , -0.524839 , -0.599547 , -0.464013 ,
- -0.159663 , -0.140879 , -0.0975618, 0.0225352, 0.0701925,
- -0.24778 };
- int ieta = static_cast<int> (fabs(eta)/0.10);
- if (ieta > 25) ieta = 25;
- double FinalMeasuredEnergy = MeasCorrected + eOverMipCorrection +
- (fix1FromPeter[ieta] + fix2FromPeter[ieta])*Units::GeV;
-
- ATH_MSG_VERBOSE( "Sum of cells " << (TileEnergy + EmEnergy + LArHECEnergy)/Units::GeV
- << " Total energy deposition " << TotalMeasuredEnergy/Units::GeV
- << " corrected energy deposition " << MeasCorrected/Units::GeV
- << " e/mip corre " << FinalMeasuredEnergy/Units::GeV << std::endl
- << "\nFinal Energy Measurement = " << FinalMeasuredEnergy /Units::GeV
- //<< "\nMean Energy Deposition = " << MeanLoss/Units::GeV
- //<< " - " << MeanErrorLeft/Units::GeV << " + "<< MeanErrorRight/Units::GeV
- << "\nMop Energy Deposition = " << MopLoss/Units::GeV << " +- " << MopError/Units::GeV
- //<< "\nOld parametrization energy= " << m_caloParamOld->deltaE()/Units::GeV
- //<< " +- " << m_caloParamOld->sigmaDeltaE()/Units::GeV
- //<< "\nTrack Momentum = " << trackMomentum/Units::GeV
- //<< " Eta= " << eta << " Phi= " << phi
- << std::endl
- << "Final Meas = " << FinalMeasuredEnergy / Units::GeV
- << " Mop Dep = " << MopLoss/Units::GeV << " +- " << MopError/Units::GeV );
-
- const double HECIso = caloMeas->Tile_Isolation() + caloMeas->LArHEC_Isolation();
- const double EmIso = caloMeas->LArEM_Isolation();
-
- const double Theta = 2.*atan(exp(-eta));
- const double pT = trackMomentum*sin(Theta)*Units::MeV;
- const double EmCut = m_emMinEnergy + (3.-2.)/(100.-15.)*(pT/Units::GeV-15.)*Units::GeV;
- const double HECCut = m_hecMinEnergy;
- const double pTCut = m_minMuonPt;
- bool PassCut = true;
- if (m_forceIsolationFailure
- || EmIso > EmCut
- || HECIso > HECCut
- || pT < pTCut
- || FinalMeasuredEnergy < m_minFinalEnergy ) PassCut = false;
-
- int nTracks = 0;
+ double fix2FromPeter[26] = {-0.647703, -0.303498, -0.268645, -0.261292, -0.260152, -0.269253, -0.266212,
+ -0.240837, -0.130172, -0.111638, -0.329423, -0.321011, -0.346050, -0.305592,
+ -0.313293, -0.317111, -0.428393, -0.524839, -0.599547, -0.464013, -0.159663,
+ -0.140879, -0.0975618, 0.0225352, 0.0701925, -0.24778};
+ int ieta = static_cast<int>(fabs(eta) / 0.10);
+ if (ieta > 25) ieta = 25;
+ double FinalMeasuredEnergy =
+ MeasCorrected + eOverMipCorrection + (fix1FromPeter[ieta] + fix2FromPeter[ieta]) * Units::GeV;
+
+ ATH_MSG_VERBOSE("Sum of cells " << (TileEnergy + EmEnergy + LArHECEnergy) / Units::GeV
+ << " Total energy deposition " << TotalMeasuredEnergy / Units::GeV
+ << " corrected energy deposition " << MeasCorrected / Units::GeV << " e/mip corre "
+ << FinalMeasuredEnergy / Units::GeV << std::endl
+ << "\nFinal Energy Measurement = "
+ << FinalMeasuredEnergy / Units::GeV
+ //<< "\nMean Energy Deposition = " << MeanLoss/Units::GeV
+ //<< " - " << MeanErrorLeft/Units::GeV << " + "<< MeanErrorRight/Units::GeV
+ << "\nMop Energy Deposition = " << MopLoss / Units::GeV << " +- "
+ << MopError / Units::GeV
+ //<< "\nOld parametrization energy= " << m_caloParamOld->deltaE()/Units::GeV
+ //<< " +- " << m_caloParamOld->sigmaDeltaE()/Units::GeV
+ //<< "\nTrack Momentum = " << trackMomentum/Units::GeV
+ //<< " Eta= " << eta << " Phi= " << phi
+ << std::endl
+ << "Final Meas = " << FinalMeasuredEnergy / Units::GeV
+ << " Mop Dep = " << MopLoss / Units::GeV << " +- " << MopError / Units::GeV);
+
+ const double HECIso = caloMeas->Tile_Isolation() + caloMeas->LArHEC_Isolation();
+ const double EmIso = caloMeas->LArEM_Isolation();
+
+ const double Theta = 2. * atan(exp(-eta));
+ const double pT = trackMomentum * sin(Theta) * Units::MeV;
+ const double EmCut = m_emMinEnergy + (3. - 2.) / (100. - 15.) * (pT / Units::GeV - 15.) * Units::GeV;
+ const double HECCut = m_hecMinEnergy;
+ const double pTCut = m_minMuonPt;
+ bool PassCut = true;
+ if (m_forceIsolationFailure || EmIso > EmCut || HECIso > HECCut || pT < pTCut
+ || FinalMeasuredEnergy < m_minFinalEnergy)
+ PassCut = false;
+
+ int nTracks = 0;
// double tracksEnergy = 0.;
- if (m_trackIsolation)
- {
- std::pair<int,double> inner = m_trackIsolationTool->trackIsolation(eta,phi);
- // double maxP = m_trackIsolationTool->maxP();
- nTracks = inner.first;
- // tracksEnergy = inner.second - maxP;
- if (pT < 100.*Units::GeV && nTracks > m_maxNTracksIso) PassCut = false;
+ if (m_trackIsolation) {
+ std::pair<int, double> inner = m_trackIsolationTool->trackIsolation(eta, phi);
+ // double maxP = m_trackIsolationTool->maxP();
+ nTracks = inner.first;
+ // tracksEnergy = inner.second - maxP;
+ if (pT < 100. * Units::GeV && nTracks > m_maxNTracksIso) PassCut = false;
}
- ATH_MSG_VERBOSE( "pT= " << pT/Units::GeV << ",HECIso= " << HECIso/Units::GeV
- << ",EmIso= " <<EmIso/Units::GeV << ", nTracks= "<< nTracks
- << ",PassCut= " << PassCut );
+ ATH_MSG_VERBOSE("pT= " << pT / Units::GeV << ",HECIso= " << HECIso / Units::GeV << ",EmIso= " << EmIso / Units::GeV
+ << ", nTracks= " << nTracks << ",PassCut= " << PassCut);
+
+ CaloEnergy::EnergyLossType lossType = CaloEnergy::NotIsolated;
+ CaloEnergy* caloEnergy = 0;
- CaloEnergy::EnergyLossType lossType = CaloEnergy::NotIsolated;
- CaloEnergy* caloEnergy = 0;
-
// choose between lossTypes MOP, Tail, FSR and NotIsolated according
// to measured energy, isolation cut and Et in em
- if (FinalMeasuredEnergy < MopLoss + 2. * MopError
- && FinalMeasuredEnergy > m_minFinalEnergy)
- {
- ++m_countMop;
- caloEnergy = mopPeak;
- mopPeak = 0;
- }
- else if (PassCut)
- {
- // tail offset from high-statistics Z->mumu MC (measured by Peter K 09/12/2011),
- // but next we try to separate any FSR contribution from the Landau tail
- double F1 = 0.;
- if (caloMeas->LArEM_EnergyMeasured() > m_emEtCut)
- F1 = caloMeas->LArEM_FirstCompartmentEnergy()/caloMeas->LArEM_EnergyMeasured();
- ATH_MSG_VERBOSE( " start Tail and FSR treatment: Et in e.m. " << EmEnergy*sin(Theta)/Units::GeV
- << " F1 ratio " << F1);
- if (! m_FSRtreatment
- || EmEnergy*sin(Theta) < m_emEtCut
- || F1 < m_emF1Cut)
- {
- ++m_countMeasurement;
- double FinalEnergyErrorMinus= 0.50 * sqrt(FinalMeasuredEnergy/Units::GeV) * Units::GeV;
- double FinalEnergyErrorPlus = 0.50 * sqrt(FinalMeasuredEnergy/Units::GeV) * Units::GeV;
-
- // overall also have 50% resolution in EC rather than the 70% naively expected from LArHEC
- if (LArHECEnergy > 1.*Units::GeV)
- {
- FinalEnergyErrorMinus = 0.50 * sqrt(FinalMeasuredEnergy/Units::GeV) * Units::GeV;
- FinalEnergyErrorPlus = 0.50 * sqrt(FinalMeasuredEnergy/Units::GeV) * Units::GeV;
- }
- double FinalEnergyError = 0.5*(FinalEnergyErrorMinus + FinalEnergyErrorPlus);
- lossType = CaloEnergy::Tail;
- caloEnergy = new CaloEnergy(FinalMeasuredEnergy,
- FinalEnergyError,
- FinalEnergyErrorMinus,
- FinalEnergyErrorPlus,
- lossType);
- ATH_MSG_VERBOSE( " CaloEnergy Tail : " << FinalMeasuredEnergy );
+ if (FinalMeasuredEnergy < MopLoss + 2. * MopError && FinalMeasuredEnergy > m_minFinalEnergy) {
+ ++m_countMop;
+ caloEnergy = mopPeak;
+ mopPeak = 0;
+ } else if (PassCut) {
+ // tail offset from high-statistics Z->mumu MC (measured by Peter K 09/12/2011),
+ // but next we try to separate any FSR contribution from the Landau tail
+ double F1 = 0.;
+ if (caloMeas->LArEM_EnergyMeasured() > m_emEtCut)
+ F1 = caloMeas->LArEM_FirstCompartmentEnergy() / caloMeas->LArEM_EnergyMeasured();
+ ATH_MSG_VERBOSE(" start Tail and FSR treatment: Et in e.m. " << EmEnergy * sin(Theta) / Units::GeV
+ << " F1 ratio " << F1);
+ if (!m_FSRtreatment || EmEnergy * sin(Theta) < m_emEtCut || F1 < m_emF1Cut) {
+ ++m_countMeasurement;
+ double FinalEnergyErrorMinus = 0.50 * sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
+ double FinalEnergyErrorPlus = 0.50 * sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
+
+ // overall also have 50% resolution in EC rather than the 70% naively expected from LArHEC
+ if (LArHECEnergy > 1. * Units::GeV) {
+ FinalEnergyErrorMinus = 0.50 * sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
+ FinalEnergyErrorPlus = 0.50 * sqrt(FinalMeasuredEnergy / Units::GeV) * Units::GeV;
+ }
+ double FinalEnergyError = 0.5 * (FinalEnergyErrorMinus + FinalEnergyErrorPlus);
+ lossType = CaloEnergy::Tail;
+ caloEnergy = new CaloEnergy(FinalMeasuredEnergy, FinalEnergyError, FinalEnergyErrorMinus,
+ FinalEnergyErrorPlus, lossType);
+ ATH_MSG_VERBOSE(" CaloEnergy Tail : " << FinalMeasuredEnergy);
+ } else {
+ // significant e.m. deposit and high F1
+ double MopErrorEm = MopError * m_caloParamTool->emMopFraction(eta);
+ double FinalMeasuredEnergyNoEm = FinalMeasuredEnergy - EmEnergy + MopLossEm;
+ // lossType = CaloEnergy::FSRcandidate;
+
+ ATH_MSG_VERBOSE(" CaloEnergy FSR : EmEnergy "
+ << EmEnergy << " FinalMeasuredEnergy " << FinalMeasuredEnergy << " MopLossEm "
+ << MopLossEm << " MopErrorEm " << MopErrorEm << " FinalMeasuredEnergyNoEm "
+ << FinalMeasuredEnergyNoEm);
+ if (FinalMeasuredEnergyNoEm < MopLoss + 2. * MopError) {
+ // small hadronic energy deposit: like NotIsolated (MOP or Mean according to configuration)
+ lossType = CaloEnergy::NotIsolated;
+ if (m_MOPparametrization) {
+ ++m_countMop;
+ caloEnergy =
+ new CaloEnergy(caloParamMop->deltaE(), caloParamMop->sigmaDeltaE(),
+ caloParamMop->sigmaMinusDeltaE(), caloParamMop->sigmaPlusDeltaE(), lossType);
+ } else {
+ ++m_countMean;
+ caloEnergy =
+ new CaloEnergy(caloParamMean->deltaE(), caloParamMean->sigmaDeltaE(),
+ caloParamMean->sigmaMinusDeltaE(), caloParamMean->sigmaPlusDeltaE(), lossType);
+ }
+ ATH_MSG_VERBOSE(" CaloEnergy FSR : Small deposit, FinalMeasuredEnergyNoEm "
+ << FinalMeasuredEnergyNoEm << " using Eloss " << caloEnergy->deltaE());
+ } else {
+ // significant hadronic energy deposit
+ ++m_countMeasurement;
+ lossType = CaloEnergy::FSRcandidate;
+ double FinalEnergyErrorNoEm = 0.50 * sqrt(FinalMeasuredEnergyNoEm / Units::GeV) * Units::GeV;
+ FinalEnergyErrorNoEm = sqrt(FinalEnergyErrorNoEm * FinalEnergyErrorNoEm + MopErrorEm * MopErrorEm);
+ caloEnergy = new CaloEnergy(FinalMeasuredEnergyNoEm, FinalEnergyErrorNoEm, FinalEnergyErrorNoEm,
+ FinalEnergyErrorNoEm, lossType);
+ ATH_MSG_VERBOSE(" CaloEnergy FSR : Large deposit, FinalMeasuredEnergyNoEm " << FinalMeasuredEnergyNoEm);
+ }
+ }
+ } else {
+ // lossType NotIsolated: MOP or Mean according to configuration
+ if (m_MOPparametrization) {
+ ++m_countMop;
+ caloEnergy = new CaloEnergy(caloParamMop->deltaE(), caloParamMop->sigmaDeltaE(),
+ caloParamMop->sigmaMinusDeltaE(), caloParamMop->sigmaPlusDeltaE(), lossType);
+ } else {
+ ++m_countMean;
+ caloEnergy = new CaloEnergy(caloParamMean->deltaE(), caloParamMean->sigmaDeltaE(),
+ caloParamMean->sigmaMinusDeltaE(), caloParamMean->sigmaPlusDeltaE(), lossType);
}
- else
- {
- // significant e.m. deposit and high F1
- double MopErrorEm = MopError*m_caloParamTool->emMopFraction(eta);
- double FinalMeasuredEnergyNoEm = FinalMeasuredEnergy - EmEnergy + MopLossEm;
- //lossType = CaloEnergy::FSRcandidate;
-
- ATH_MSG_VERBOSE( " CaloEnergy FSR : EmEnergy " << EmEnergy
- << " FinalMeasuredEnergy " << FinalMeasuredEnergy
- << " MopLossEm " << MopLossEm
- << " MopErrorEm " << MopErrorEm
- << " FinalMeasuredEnergyNoEm " << FinalMeasuredEnergyNoEm );
- if (FinalMeasuredEnergyNoEm < MopLoss + 2.*MopError)
- {
- // small hadronic energy deposit: like NotIsolated (MOP or Mean according to configuration)
- lossType = CaloEnergy::NotIsolated;
- if (m_MOPparametrization)
- {
- ++m_countMop;
- caloEnergy = new CaloEnergy(caloParamMop->deltaE(),
- caloParamMop->sigmaDeltaE(),
- caloParamMop->sigmaMinusDeltaE(),
- caloParamMop->sigmaPlusDeltaE(),
- lossType);
- }
- else
- {
- ++m_countMean;
- caloEnergy = new CaloEnergy(caloParamMean->deltaE(),
- caloParamMean->sigmaDeltaE(),
- caloParamMean->sigmaMinusDeltaE(),
- caloParamMean->sigmaPlusDeltaE(),
- lossType);
- }
- ATH_MSG_VERBOSE( " CaloEnergy FSR : Small deposit, FinalMeasuredEnergyNoEm "
- << FinalMeasuredEnergyNoEm
- << " using Eloss " << caloEnergy->deltaE());
- }
- else
- {
- // significant hadronic energy deposit
- ++m_countMeasurement;
- lossType = CaloEnergy::FSRcandidate;
- double FinalEnergyErrorNoEm = 0.50 * sqrt(FinalMeasuredEnergyNoEm/Units::GeV) * Units::GeV;
- FinalEnergyErrorNoEm = sqrt(FinalEnergyErrorNoEm*FinalEnergyErrorNoEm +
- MopErrorEm*MopErrorEm);
- caloEnergy = new CaloEnergy(FinalMeasuredEnergyNoEm,
- FinalEnergyErrorNoEm,
- FinalEnergyErrorNoEm,
- FinalEnergyErrorNoEm,
- lossType);
- ATH_MSG_VERBOSE( " CaloEnergy FSR : Large deposit, FinalMeasuredEnergyNoEm "
- << FinalMeasuredEnergyNoEm );
-
- }
- }
- }
- else
- {
- // lossType NotIsolated: MOP or Mean according to configuration
- if (m_MOPparametrization)
- {
- ++m_countMop;
- caloEnergy = new CaloEnergy(caloParamMop->deltaE(),
- caloParamMop->sigmaDeltaE(),
- caloParamMop->sigmaMinusDeltaE(),
- caloParamMop->sigmaPlusDeltaE(),
- lossType);
- }
- else
- {
- ++m_countMean;
- caloEnergy = new CaloEnergy(caloParamMean->deltaE(),
- caloParamMean->sigmaDeltaE(),
- caloParamMean->sigmaMinusDeltaE(),
- caloParamMean->sigmaPlusDeltaE(),
- lossType);
- }
}
-
+
// delete the various parametrized CaloEnergy's before return
delete caloParamMean;
delete caloParamMop;
@@ -684,57 +548,47 @@ MuidCaloEnergyTool::measurement(double trackMomentum,
return caloEnergy;
}
-double
-MuidCaloEnergyTool::muSpecResolParam(double trackMomentum,
- double eta) const
+double
+MuidCaloEnergyTool::muSpecResolParam(double trackMomentum, double eta) const
{
- const double Theta = 2.*atan(exp(-eta));
- const double pT = trackMomentum*sin(Theta)/Units::GeV; // pt in GeV
- double a = 0.;
- double b = 0.;
- if (fabs(eta)<1.)
- {
- a = 0.02255;
- b = 7.708e-5;
- }
- else if (fabs(eta)>1. && fabs(eta)<2.)
- {
- a = 0.04198;
- b = 8.912e-5;
- }
- else
- {
- a = 0.02181;
- b = 7.282e-5;
- }
- return sqrt(a*a+ (b*pT)*(b*pT));
+ const double Theta = 2. * atan(exp(-eta));
+ const double pT = trackMomentum * sin(Theta) / Units::GeV; // pt in GeV
+ double a = 0.;
+ double b = 0.;
+ if (fabs(eta) < 1.) {
+ a = 0.02255;
+ b = 7.708e-5;
+ } else if (fabs(eta) > 1. && fabs(eta) < 2.) {
+ a = 0.04198;
+ b = 8.912e-5;
+ } else {
+ a = 0.02181;
+ b = 7.282e-5;
+ }
+ return sqrt(a * a + (b * pT) * (b * pT));
}
-double
-MuidCaloEnergyTool::paramCorrection(double trackMomentum,
- double eta,
- double MopLoss,
- double MopLossSigma) const
+double
+MuidCaloEnergyTool::paramCorrection(double trackMomentum, double eta, double MopLoss, double MopLossSigma) const
{
- const double Nsigma = 2.;
- double MSres = muSpecResolParam(trackMomentum,eta);
- double MSsigma = trackMomentum*MSres;
- double sigma = sqrt(MSsigma*MSsigma + MopLossSigma*MopLossSigma);
- double sum = 0.;
- double weight = 0.;
- double xlow = MopLoss - Nsigma * sigma;
- if (xlow<0.) xlow = 0.;
- double xup = MopLoss + Nsigma * sigma;
- int Na = 50;
- const double inv_Na = 1. / static_cast<double> (Na);
- for (int j = 0; j < Na; ++j)
- {
- double x = xlow + j*(xup-xlow)*inv_Na;
- double w = landau(x,MopLoss,MopLossSigma,true);
- sum += x*w;
- weight += w;
- }
- double MopStat = sum/weight;
+ const double Nsigma = 2.;
+ double MSres = muSpecResolParam(trackMomentum, eta);
+ double MSsigma = trackMomentum * MSres;
+ double sigma = sqrt(MSsigma * MSsigma + MopLossSigma * MopLossSigma);
+ double sum = 0.;
+ double weight = 0.;
+ double xlow = MopLoss - Nsigma * sigma;
+ if (xlow < 0.) xlow = 0.;
+ double xup = MopLoss + Nsigma * sigma;
+ int Na = 50;
+ const double inv_Na = 1. / static_cast<double>(Na);
+ for (int j = 0; j < Na; ++j) {
+ double x = xlow + j * (xup - xlow) * inv_Na;
+ double w = landau(x, MopLoss, MopLossSigma, true);
+ sum += x * w;
+ weight += w;
+ }
+ double MopStat = sum / weight;
return MopStat;
}
@@ -745,81 +599,65 @@ MuidCaloEnergyTool::landau(double x, double mpv, double sigma, bool norm) const
// This function has been adapted from the CERNLIB routine G110 denlan.
// If norm=kTRUE (default is kFALSE) the result is divided by sigma
- double p1[5] = {0.4259894875,-0.1249762550, 0.03984243700, -0.006298287635, 0.001511162253};
- double q1[5] = {1.0 ,-0.3388260629, 0.09594393323, -0.01608042283, 0.003778942063};
+ double p1[5] = {0.4259894875, -0.1249762550, 0.03984243700, -0.006298287635, 0.001511162253};
+ double q1[5] = {1.0, -0.3388260629, 0.09594393323, -0.01608042283, 0.003778942063};
- double p2[5] = {0.1788541609, 0.1173957403, 0.01488850518, -0.001394989411, 0.0001283617211};
- double q2[5] = {1.0 , 0.7428795082, 0.3153932961, 0.06694219548, 0.008790609714};
+ double p2[5] = {0.1788541609, 0.1173957403, 0.01488850518, -0.001394989411, 0.0001283617211};
+ double q2[5] = {1.0, 0.7428795082, 0.3153932961, 0.06694219548, 0.008790609714};
- double p3[5] = {0.1788544503, 0.09359161662,0.006325387654, 0.00006611667319,-0.000002031049101};
- double q3[5] = {1.0 , 0.6097809921, 0.2560616665, 0.04746722384, 0.006957301675};
+ double p3[5] = {0.1788544503, 0.09359161662, 0.006325387654, 0.00006611667319, -0.000002031049101};
+ double q3[5] = {1.0, 0.6097809921, 0.2560616665, 0.04746722384, 0.006957301675};
- double p4[5] = {0.9874054407, 118.6723273, 849.2794360, -743.7792444, 427.0262186};
- double q4[5] = {1.0 , 106.8615961, 337.6496214, 2016.712389, 1597.063511};
+ double p4[5] = {0.9874054407, 118.6723273, 849.2794360, -743.7792444, 427.0262186};
+ double q4[5] = {1.0, 106.8615961, 337.6496214, 2016.712389, 1597.063511};
- double p5[5] = {1.003675074, 167.5702434, 4789.711289, 21217.86767, -22324.94910};
- double q5[5] = {1.0 , 156.9424537, 3745.310488, 9834.698876, 66924.28357};
+ double p5[5] = {1.003675074, 167.5702434, 4789.711289, 21217.86767, -22324.94910};
+ double q5[5] = {1.0, 156.9424537, 3745.310488, 9834.698876, 66924.28357};
- double p6[5] = {1.000827619, 664.9143136, 62972.92665, 475554.6998, -5743609.109};
- double q6[5] = {1.0 , 651.4101098, 56974.73333, 165917.4725, -2815759.939};
+ double p6[5] = {1.000827619, 664.9143136, 62972.92665, 475554.6998, -5743609.109};
+ double q6[5] = {1.0, 651.4101098, 56974.73333, 165917.4725, -2815759.939};
- double a1[3] = {0.04166666667,-0.01996527778, 0.02709538966};
+ double a1[3] = {0.04166666667, -0.01996527778, 0.02709538966};
- double a2[2] = {-1.845568670,-4.284640743};
+ double a2[2] = {-1.845568670, -4.284640743};
if (sigma <= 0) return 0.;
- double v = (x-mpv)/sigma;
+ double v = (x - mpv) / sigma;
double u, ue, us, den;
- if (v < -5.5)
- {
- u = exp(v+1.0);
- if (u < 1e-10) return 0.0;
- ue = exp(-1/u);
- us = sqrt(u);
- den = 0.3989422803*(ue/us)*(1+(a1[0]+(a1[1]+a1[2]*u)*u)*u);
- }
- else if (v < -1)
- {
- u = exp(-v-1);
- den = exp(-u)*sqrt(u)*
- (p1[0]+(p1[1]+(p1[2]+(p1[3]+p1[4]*v)*v)*v)*v)/
- (q1[0]+(q1[1]+(q1[2]+(q1[3]+q1[4]*v)*v)*v)*v);
- }
- else if (v < 1)
- {
- den = (p2[0]+(p2[1]+(p2[2]+(p2[3]+p2[4]*v)*v)*v)*v)/
- (q2[0]+(q2[1]+(q2[2]+(q2[3]+q2[4]*v)*v)*v)*v);
- }
- else if (v < 5)
- {
- den = (p3[0]+(p3[1]+(p3[2]+(p3[3]+p3[4]*v)*v)*v)*v)/
- (q3[0]+(q3[1]+(q3[2]+(q3[3]+q3[4]*v)*v)*v)*v);
- }
- else if (v < 12)
- {
- u = 1/v;
- den = u*u*(p4[0]+(p4[1]+(p4[2]+(p4[3]+p4[4]*u)*u)*u)*u)/
- (q4[0]+(q4[1]+(q4[2]+(q4[3]+q4[4]*u)*u)*u)*u);
- }
- else if (v < 50)
- {
- u = 1/v;
- den = u*u*(p5[0]+(p5[1]+(p5[2]+(p5[3]+p5[4]*u)*u)*u)*u)/
- (q5[0]+(q5[1]+(q5[2]+(q5[3]+q5[4]*u)*u)*u)*u);
- }
- else if (v < 300)
- {
- u = 1/v;
- den = u*u*(p6[0]+(p6[1]+(p6[2]+(p6[3]+p6[4]*u)*u)*u)*u)/
- (q6[0]+(q6[1]+(q6[2]+(q6[3]+q6[4]*u)*u)*u)*u);
- }
- else
- {
- u = 1/(v-v*log(v)/(v+1));
- den = u*u*(1+(a2[0]+a2[1]*u)*u);
+ if (v < -5.5) {
+ u = exp(v + 1.0);
+ if (u < 1e-10) return 0.0;
+ ue = exp(-1 / u);
+ us = sqrt(u);
+ den = 0.3989422803 * (ue / us) * (1 + (a1[0] + (a1[1] + a1[2] * u) * u) * u);
+ } else if (v < -1) {
+ u = exp(-v - 1);
+ den = exp(-u) * sqrt(u) * (p1[0] + (p1[1] + (p1[2] + (p1[3] + p1[4] * v) * v) * v) * v)
+ / (q1[0] + (q1[1] + (q1[2] + (q1[3] + q1[4] * v) * v) * v) * v);
+ } else if (v < 1) {
+ den = (p2[0] + (p2[1] + (p2[2] + (p2[3] + p2[4] * v) * v) * v) * v)
+ / (q2[0] + (q2[1] + (q2[2] + (q2[3] + q2[4] * v) * v) * v) * v);
+ } else if (v < 5) {
+ den = (p3[0] + (p3[1] + (p3[2] + (p3[3] + p3[4] * v) * v) * v) * v)
+ / (q3[0] + (q3[1] + (q3[2] + (q3[3] + q3[4] * v) * v) * v) * v);
+ } else if (v < 12) {
+ u = 1 / v;
+ den = u * u * (p4[0] + (p4[1] + (p4[2] + (p4[3] + p4[4] * u) * u) * u) * u)
+ / (q4[0] + (q4[1] + (q4[2] + (q4[3] + q4[4] * u) * u) * u) * u);
+ } else if (v < 50) {
+ u = 1 / v;
+ den = u * u * (p5[0] + (p5[1] + (p5[2] + (p5[3] + p5[4] * u) * u) * u) * u)
+ / (q5[0] + (q5[1] + (q5[2] + (q5[3] + q5[4] * u) * u) * u) * u);
+ } else if (v < 300) {
+ u = 1 / v;
+ den = u * u * (p6[0] + (p6[1] + (p6[2] + (p6[3] + p6[4] * u) * u) * u) * u)
+ / (q6[0] + (q6[1] + (q6[2] + (q6[3] + q6[4] * u) * u) * u) * u);
+ } else {
+ u = 1 / (v - v * log(v) / (v + 1));
+ den = u * u * (1 + (a2[0] + a2[1] * u) * u);
}
if (!norm) return den;
- return den/sigma;
+ return den / sigma;
}
-} // end of namespace
+} // namespace Rec