ATLAS Offline Software
Loading...
Searching...
No Matches
Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
ZMassConstraint::ConstraintFit Class Reference

#include <ConstraintFit.h>

Inheritance diagram for ZMassConstraint::ConstraintFit:
Collaboration diagram for ZMassConstraint::ConstraintFit:

Public Member Functions

 ConstraintFit (const std::string &name)
 Create a proper constructor for Athena.
 ~ConstraintFit (void)
StatusCode initialize ()
 Initialize constraint fit.
StatusCode doMassFit (const ConstraintFitInput &input, ConstraintFitOutput &output)
 Perform the constrained mass fit.
double getMassError (const ConstraintFitInput &firstInput, const ConstraintFitInput &secondInput=ConstraintFitInput())
 Calculate the mass error without fit - use just the inputs.
double getMassError (const ConstraintFitOutput &fitOutput, const ConstraintFitInput &extraInput=ConstraintFitInput())
 Calculate the mass error with the result of fit and possible additional input (4-vecs/cov mat) that has NOT been fit.
double getMassError (const ConstraintFitOutput &fitOutput, const ConstraintFitOutput &secondFitOutput)
 Calculate the mass error with the result of two mass fits.
void addParticle (const xAOD::Muon &part, ConstraintFitInput &input, MassConstraintMuonType muonType=isCombMCMT)
 Add muon to input, must provide the resolution Scale Factor.
void addParticle (const xAOD::Electron &part, float elEnergyRes, ConstraintFitInput &input)
 Add electron to input, must provide the electron energy resolution.
void addFSRParticle (const xAOD::IParticle &part, const TLorentzVector &fsr4Vec, ConstraintFitInput &input)
 Add in FSR photon to input, (energy resolution is obtain in method).
virtual void print () const
 Print the state of the tool.
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysInitialize () override
 Perform system initialization for an algorithm.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const
Additional helper functions, not directly mimicking Athena
template<class T>
const T * getProperty (const std::string &name) const
 Get one of the tool's properties.
const std::string & msg_level_name () const __attribute__((deprecated))
 A deprecated function for getting the message level's name.
const std::string & getName (const void *ptr) const
 Get the name of an object that is / should be in the event store.
SG::sgkey_t getKey (const void *ptr) const
 Get the (hashed) key of an object that is in the event store.

Protected Member Functions

void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.

Private Types

enum  TrackParameters {
  d0 = 0 , z0 = 1 , phi0 = 2 , theta = 3 ,
  qOverP = 4 , qP = 4 , x = 0 , y = 1 ,
  z = 2
}
typedef ServiceHandle< StoreGateSvc > StoreGateSvc_t

Private Member Functions

int massFitInterface (const ConstraintFitInput &theInput)
void massFitRun (ConstraintFitOutput &output, double zresol=-1.)
double getMassError (const std::vector< TLorentzVector > &particles, const std::vector< AmgMatrix(3, 3)> &covariances)
double likelihoodMass2 (void)
double likelihoodMass (double)
void setIgnore (bool val)
bool doSanityChecksOnCovariance (const TLorentzVector &vector, const AmgMatrix(5, 5)&covar) const
void convertCovd0z0PhiThetaPToXYZ (const TLorentzVector &fourVec, const AmgMatrix(5, 5)&covard0z0PhiThetaP, AmgMatrix(5, 5)&covarXYZ) const
void convertCovXYZTod0z0PhiThetaP (const std::vector< TLorentzVector > &particleList, const Amg::MatrixX &covarXYZ, Amg::MatrixX &covard0z0PhiThetaP) const
double massFitCalculation (const Amg::MatrixX &var, double mass, Amg::MatrixX &p0)
double massFit (const Amg::MatrixX &, const Amg::MatrixX &var, double mass, Amg::MatrixX &pOut, Amg::MatrixX &)
float retrieve_eta_calo (const xAOD::IParticle &part) const
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

bool m_conHasWidth
double m_conMass
double m_conWidth
double m_resolution
bool m_ignoreInputChecks
ToolHandle< CP::IEgammaCalibrationAndSmearingTool > m_energyRescaler
ToolHandle< CP::IMuonCalibrationAndSmearingTool > m_mu_resolSFTool
ConstraintFitInput m_theInput
unsigned int m_parameters
unsigned int m_nobj
std::vector< double > m_objmass
Amg::MatrixX m_parametersInit
Amg::MatrixX m_covarianceInit
Amg::MatrixX m_parametersFinal
Amg::MatrixX m_covarianceFinal
StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default).
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default).
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared

Detailed Description

Definition at line 30 of file ConstraintFit.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Member Enumeration Documentation

◆ TrackParameters

enum ZMassConstraint::ConstraintFit::TrackParameters
private
Enumerator
d0 
z0 
phi0 
theta 
qOverP 
qP 
x 
y 
z 

Definition at line 82 of file ConstraintFit.h.

82 {
83 d0 = 0,
84 z0 = 1,
85 phi0 = 2,
86 theta = 3,
87 qOverP = 4,
88 qP = 4,
89 x = 0,
90 y = 1,
91 z = 2
92 };

Constructor & Destructor Documentation

◆ ConstraintFit()

ZMassConstraint::ConstraintFit::ConstraintFit ( const std::string & name)

Create a proper constructor for Athena.

Create a constructor for standalone usage

Definition at line 13 of file ConstraintFit.cxx.

13 :
14 asg::AsgTool( name ),
15 m_conHasWidth(true),
16 m_conMass(91187.6),
17 m_conWidth(2495.2),
18 m_resolution(0.01),
19 m_ignoreInputChecks(false),
20 m_parameters(3),
21 m_nobj(0)
22 {
23 declareProperty( "Z_pdg_mass", m_conMass);
24 declareProperty( "Z_pdg_width", m_conWidth);
25 declareProperty( "hasWidth", m_conHasWidth);
26 declareProperty( "ignoreInputChecks", m_ignoreInputChecks);
27 declareProperty( "EgammaCalibAndSmearingTool", m_energyRescaler = ToolHandle<CP::IEgammaCalibrationAndSmearingTool>("CP::EgammaCalibrationAndSmearingTool"));
28 declareProperty( "MuonCalibrationAndSmearingTool", m_mu_resolSFTool = ToolHandle<CP::IMuonCalibrationAndSmearingTool>("CP::MuonCalibrationAndSmearingTool"));
29 }
AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
Definition AthCommonDataStore.h:145
ZMassConstraint::ConstraintFit::m_energyRescaler
ToolHandle< CP::IEgammaCalibrationAndSmearingTool > m_energyRescaler
Definition ConstraintFit.h:116
ZMassConstraint::ConstraintFit::m_mu_resolSFTool
ToolHandle< CP::IMuonCalibrationAndSmearingTool > m_mu_resolSFTool
Definition ConstraintFit.h:117

◆ ~ConstraintFit()

ZMassConstraint::ConstraintFit::~ConstraintFit ( void )

Definition at line 31 of file ConstraintFit.cxx.

32 {}

Member Function Documentation

◆ addFSRParticle()

void ZMassConstraint::ConstraintFit::addFSRParticle ( const xAOD::IParticle & part,
const TLorentzVector & fsr4Vec,
ConstraintFitInput & input )
virtual

Add in FSR photon to input, (energy resolution is obtain in method).

Implements ZMassConstraint::IConstraintFit.

Definition at line 84 of file ConstraintFit.cxx.

87 {
88
89 ATH_MSG_DEBUG ("addFSRParticle: *** 4vecFsr: " << fsr4Vec.Pt() << "/" << fsr4Vec.Eta() << "/" << fsr4Vec.Phi());
90 // Get energy resolution - done differently for photons and electrons
91 float e_res = 0;
92 PATCore::ParticleType::Type phType = PATCore::ParticleType::Electron;
93 float cl_etaCalo = retrieve_eta_calo(part);
94 if (!m_energyRescaler.empty()) {
95 if (xAOD::Type::Photon == part.type()) {
96 const xAOD::Photon* fsrPf = dynamic_cast<const xAOD::Photon*>(&part);
97 phType = (xAOD::EgammaHelpers::isConvertedPhoton(fsrPf)) ?
98 PATCore::ParticleType::ConvertedPhoton :
99 PATCore::ParticleType::UnconvertedPhoton;
100 }
101 e_res = m_energyRescaler->resolution( fsr4Vec.E(), fsr4Vec.Eta(), cl_etaCalo, phType)*fsr4Vec.E();
102 }
103 else {
104 ATH_MSG_ERROR( "addFSRParticle: - cannot find EgammaCalibrationAndSmearingTool. ");
105 ATH_MSG_ERROR( "addFSRParticle: - please set property 'EgammaCalibAndSmearingToolName' with the name of your EgammaCalibrationAndSmearingTool - default is 'EgammaCalibrationAndSmearingTool'. ");
106 }
107
108 AmgMatrix(5,5) photonCovarianceMatrix;
109 photonCovarianceMatrix.setZero();
110 photonCovarianceMatrix(d0,d0) = 0.000001;
111 photonCovarianceMatrix(z0,z0) = 0.000001;
112 photonCovarianceMatrix(phi0,phi0) = 0.000001;
113 photonCovarianceMatrix(theta,theta) = 0.000001;
114 photonCovarianceMatrix(qOverP,qOverP) = e_res*e_res;
115 photonCovarianceMatrix(d0,z0) = 0;
116 photonCovarianceMatrix(d0,phi0) = 0;
117 photonCovarianceMatrix(d0,theta) = 0;
118 photonCovarianceMatrix(d0,qOverP) = 0;
119 photonCovarianceMatrix(z0,phi0) = 0;
120 photonCovarianceMatrix(z0,theta) = 0;
121 photonCovarianceMatrix(z0,qOverP) = 0;
122 photonCovarianceMatrix(phi0,theta) = 0;
123 photonCovarianceMatrix(phi0,qOverP) = 0;
124 photonCovarianceMatrix(theta,qOverP) = 0;
125
126 ATH_MSG_DEBUG ("addFSRParticle: FSR calib type,e, e_res " << phType << " " << fsr4Vec.E() << " " << e_res);
127
128 // symmetrize it
129 for (int ii = 0; ii < 5; ii++)
130 for (int jj = ii + 1; jj < 5; jj++)
131 photonCovarianceMatrix(jj, ii) = photonCovarianceMatrix(ii, jj);
132
133 bool isOK = doSanityChecksOnCovariance(fsr4Vec, photonCovarianceMatrix);
134 // Get covariance in cartesian CS
135 AmgMatrix(5,5) covarXYZ;
136 covarXYZ.setZero();
137 convertCovd0z0PhiThetaPToXYZ(fsr4Vec, photonCovarianceMatrix, covarXYZ);
138 // Save in input
139 input.addConstituent_FourVector_d0z0PhiThetaP(fsr4Vec, photonCovarianceMatrix, covarXYZ, isOK);
140
141 ATH_MSG_DEBUG( "addFSRParticle: - added fsr");
142
143 }
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition AthMsgStreamMacros.h:33
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition AthMsgStreamMacros.h:29
AmgMatrix
#define AmgMatrix(rows, cols)
Definition EventPrimitives.h:49
ParticleType
ParticleType
Definition TruthClasses.h:8
ZMassConstraint::ConstraintFit::convertCovd0z0PhiThetaPToXYZ
void convertCovd0z0PhiThetaPToXYZ(const TLorentzVector &fourVec, const AmgMatrix(5, 5)&covard0z0PhiThetaP, AmgMatrix(5, 5)&covarXYZ) const
Definition ConstraintFit.cxx:1128
ZMassConstraint::ConstraintFit::retrieve_eta_calo
float retrieve_eta_calo(const xAOD::IParticle &part) const
Definition ConstraintFit.cxx:1190
ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance
bool doSanityChecksOnCovariance(const TLorentzVector &vector, const AmgMatrix(5, 5)&covar) const
Definition ConstraintFit.cxx:1091
xAODType::Photon
@ Photon
The object is a photon.
Definition ObjectType.h:47
xAOD::EgammaHelpers::isConvertedPhoton
bool isConvertedPhoton(const xAOD::Egamma *eg, bool excludeTRT=false)
is the object a converted photon
Definition EgammaxAODHelpers.cxx:25
xAOD::Photon
Photon_v1 Photon
Definition of the current "egamma version".
Definition Event/xAOD/xAODEgamma/xAODEgamma/Photon.h:17

◆ addParticle() [1/2]

void ZMassConstraint::ConstraintFit::addParticle ( const xAOD::Electron & part,
float elEnergyRes,
ConstraintFitInput & input )
virtual

Add electron to input, must provide the electron energy resolution.

Implements ZMassConstraint::IConstraintFit.

Definition at line 311 of file ConstraintFit.cxx.

314 {
315 float el_E_res = elEnergyRes;
316 // if (!m_energyRescaler.empty()) {
317 // // Use el.e() which is assumed to be equal to the corrected cluster energy
318 // el_E_res = m_energyRescaler->resolution(el.caloCluster()->e(), el.caloCluster()->eta())*
319 // el.caloCluster()->e();
320 // }
321 // else ATH_MSG_ERROR( "addParticle: - energyRescaler is NOT AVAILABLE. ");
322
323 const xAOD::TrackParticle* track = el.trackParticle();
324
325 if (!track) {
326 ATH_MSG_ERROR( "addParticle: - ERROR could not get track particle for electron");
327 return;
328 }
329
330 ATH_MSG_DEBUG ("addParticle: *** 4vecEl: " << el.pt() << "/" << el.eta() << "/" << el.phi());
331
332 AmgMatrix(5,5) covmatrix;
333 covmatrix.setZero();
334 covmatrix(d0,d0) = track->definingParametersCovMatrix()(d0,d0);
335 covmatrix(z0,z0) = track->definingParametersCovMatrix()(z0,z0);
336 covmatrix(phi0,phi0) = track->definingParametersCovMatrix()(phi0,phi0);
337 covmatrix(theta,theta) = track->definingParametersCovMatrix()(theta,theta);
338 covmatrix(qOverP,qOverP) = el_E_res*el_E_res;
339 covmatrix(d0,z0) = track->definingParametersCovMatrix()(d0,z0);
340 covmatrix(d0,phi0) = track->definingParametersCovMatrix()(d0,phi0);
341 covmatrix(d0,theta) = track->definingParametersCovMatrix()(d0,theta);
342 covmatrix(d0,qOverP) = 0;
343 covmatrix(z0,phi0) = track->definingParametersCovMatrix()(z0,phi0);
344 covmatrix(z0,theta) = track->definingParametersCovMatrix()(z0,theta);
345 covmatrix(z0,qOverP) = 0;
346 covmatrix(phi0,theta) = track->definingParametersCovMatrix()(phi0,theta);
347 covmatrix(phi0,qOverP) = 0;
348 covmatrix(theta,qOverP) = 0;
349
350 ATH_MSG_DEBUG( "addParticle: el p,cov"
351 << " " << el_E_res << " " << el.e()
352 << " " << track->definingParametersCovMatrix()(d0, d0)
353 << " " << track->definingParametersCovMatrix()(z0, z0)
354 << " " << track->definingParametersCovMatrix()(phi0, phi0)
355 << " " << track->definingParametersCovMatrix()(theta, theta)
356 << " " << track->definingParametersCovMatrix()(d0, z0)
357 << " " << track->definingParametersCovMatrix()(d0, phi0)
358 << " " << track->definingParametersCovMatrix()(d0, theta)
359 << " " << track->definingParametersCovMatrix()(z0, phi0)
360 << " " << track->definingParametersCovMatrix()(z0, theta)
361 << " " << track->definingParametersCovMatrix()(phi0, theta));
362
363
364
365 for(int ii=0;ii<5;ii++)
366 for(int jj=ii+1;jj<5;jj++)
367 covmatrix(jj,ii) = covmatrix(ii,jj);
368
369 bool isOK = doSanityChecksOnCovariance(el.p4(), covmatrix);
370 // Get covariance in cartesian CS
371 AmgMatrix(5,5) covarXYZ;
372 covarXYZ.setZero();
373 convertCovd0z0PhiThetaPToXYZ(el.p4(), covmatrix, covarXYZ);
374 // Save in input
375 input.addConstituent_FourVector_d0z0PhiThetaP(el.p4(), covmatrix, covarXYZ, isOK);
376
377 ATH_MSG_DEBUG( "addParticle: - el \n" << covmatrix);
378
379 ATH_MSG_VERBOSE( "addParticle: - el\n" << input.getCovarianceCartesian(input.getNConstituents()-1));
380
381 ATH_MSG_VERBOSE( "addParticle: - added electron");
382 }
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition AthMsgStreamMacros.h:28
xAOD::TrackParticle
TrackParticle_v1 TrackParticle
Reference the current persistent version:
Definition Event/xAOD/xAODTracking/xAODTracking/TrackParticle.h:13

◆ addParticle() [2/2]

void ZMassConstraint::ConstraintFit::addParticle ( const xAOD::Muon & part,
ConstraintFitInput & input,
MassConstraintMuonType muonType = isCombMCMT )
virtual

Add muon to input, must provide the resolution Scale Factor.

Implements ZMassConstraint::IConstraintFit.

Definition at line 146 of file ConstraintFit.cxx.

149 {
150
151 // For muons, we allow the mass constraint to be done in one of three ways:
152 // 1) standard way - use muon momentum and primary track for cov matrix, or
153 // 2,3) ID or MS - use the ID or MS tracks for both the momentum and cov matrix
154
155 TLorentzVector mu4vec; // to be filled according to muon type
156
157 // setup accessors for muon decorations - set in MuonCalibrationAndSmearingTool
158 static const SG::AuxElement::ConstAccessor<float> muonSpectrometerPt ("MuonSpectrometerPt");
159 static const SG::AuxElement::ConstAccessor<float> innerDetectorPt ("InnerDetectorPt");
160
161 // For the momentum,
162 // Get the track particle according to the requested muon type for the covariance matrix
163 const xAOD::TrackParticle* track = mu.primaryTrackParticle();
164 bool set4vec = false;
165 if (((isMS_MCMT == muonType) || (isID_MCMT == muonType)) && xAOD::Muon::Combined == mu.muonType()) {
166 if (isMS_MCMT == muonType) {
167 track = mu.trackParticle(xAOD::Muon::MuonSpectrometerTrackParticle);
168 if (!muonSpectrometerPt.isAvailable(mu))
169 ATH_MSG_ERROR( "addParticle: - could not get muonSpectrometerPt from muon. Please applyCorrection with the MuonCalibAndSmearTool");
170 if (!track) {
171 ATH_MSG_ERROR( "addParticle: - Combined muon is missing MS track particle. Using combined track particle");
172 track = mu.primaryTrackParticle();
173 }
174 mu4vec.SetPtEtaPhiM(muonSpectrometerPt(mu), track->eta(), track->phi(), mu.m());
175 set4vec = true;
176 ATH_MSG_DEBUG ("addParticle: set muon track to MS ");
177 }
178 else {
179 track = mu.trackParticle(xAOD::Muon::InnerDetectorTrackParticle);
180 if (!innerDetectorPt.isAvailable(mu))
181 ATH_MSG_ERROR( "addParticle: - could not get innerDetectorPt from muon. Please applyCorrection with the MuonCalibAndSmearTool");
182 mu4vec.SetPtEtaPhiM(innerDetectorPt(mu), track->eta(), track->phi(), mu.m());
183 set4vec = true;
184 ATH_MSG_DEBUG ("addParticle: set muon track to ID ");
185 }
186 }
187
188 if (!set4vec) mu4vec.SetPtEtaPhiM(mu.pt(), mu.eta(), mu.phi(), mu.m());
189
190 if (!track) {
191 ATH_MSG_ERROR( "addParticle: - could not get track particle for muon");
192 return;
193 }
194
195 ATH_MSG_DEBUG ("addParticle: *** 4vecMu: " << mu4vec.Pt() << "/" << mu4vec.Eta() << "/" << mu4vec.Phi());
196
197
198 float muSF = 1;
199 int type = -1;
200
201 // Comment out use of MuonCalibrationAndSmearingTool until the SF is implemented by Giacomo
202
203 // if (m_mu_resolSFTool) {
204 // // Get scale factor for the muon covariance qoverp resolution.
205 // // For muonType == isCombMCMT, we get a SF differently according to the type of muon: Comb/ST/CT/SA
206 // // muonType == isID_MCMT, we use the type 2 for ID tracks
207 // // muonType == isMS_MCMT, we use the type 3 for MS tracks
208
209 // if (muonType == isCombMCMT) {
210 // if (xAOD::Muon::Combined == mu.muonType()) type = 1;
211 // else if (xAOD::Muon::SegmentTagged == mu.muonType()) type = 2;
212 // else if (xAOD::Muon::CaloTagged == mu.muonType()) type = 2;
213 // else if (xAOD::Muon::MuonStandAlone == mu.muonType()) type = 3;
214 // else ATH_MSG_ERROR("addParticle: - **** MUON type is not CB/ST/SA/Calo **** ");
215 // }
216 // else if (muonType == isID_MCMT) {
217 // type = 2;
218 // if (xAOD::Muon::MuonStandAlone == mu.muonType()) type = 3;
219 // }
220 // else if (muonType == isMS_MCMT) {
221 // type = 3;
222 // if (xAOD::Muon::SegmentTagged == mu.muonType() ||
223 // xAOD::Muon::CaloTagged == mu.muonType()) type = 2;
224 // }
225
226 // if (muonType < isCombMCMT) {
227 // ATH_MSG_ERROR( "addParticle: - Invalid muon type requested, please use .");
228 // }
229 // else {
230 // muSF = m_mu_resolSF->getResolutionScaleFactor(mu4vec, type);
231 // }
232
233 // ATH_MSG_DEBUG ("addParticle: 4v,type,errSF " << mu4vec.Eta() << " " << type << " " << muSF);
234
235
236 // }
237 // else {
238 // ATH_MSG_ERROR( "addParticle: - MuonResSFTool is NOT set. Please call setMuonResSFTool with MuonResolutionAndMomentumScaleFactors tool.");
239 // }
240
241
242 AmgMatrix(5,5) covmatrix;
243 covmatrix.setZero();
244 covmatrix(d0,d0) = track->definingParametersCovMatrix()(d0,d0);
245 covmatrix(z0,z0) = track->definingParametersCovMatrix()(z0,z0);
246 covmatrix(phi0,phi0) = track->definingParametersCovMatrix()(phi0,phi0);
247 covmatrix(theta,theta) = track->definingParametersCovMatrix()(theta,theta);
248 covmatrix(qOverP,qOverP) = track->definingParametersCovMatrix()(qOverP,qOverP)*muSF*muSF;
249 covmatrix(d0,z0) = track->definingParametersCovMatrix()(d0,z0);
250 covmatrix(d0,phi0) = track->definingParametersCovMatrix()(d0,phi0);
251 covmatrix(d0,theta) = track->definingParametersCovMatrix()(d0,theta);
252 covmatrix(d0,qOverP) = track->definingParametersCovMatrix()(d0,qOverP)*muSF;
253 covmatrix(z0,phi0) = track->definingParametersCovMatrix()(z0,phi0);
254 covmatrix(z0,theta) = track->definingParametersCovMatrix()(z0,theta);
255 covmatrix(z0,qOverP) = track->definingParametersCovMatrix()(z0,qOverP)*muSF;
256 covmatrix(phi0,theta) = track->definingParametersCovMatrix()(phi0,theta);
257 covmatrix(phi0,qOverP) = track->definingParametersCovMatrix()(phi0,qOverP)*muSF;
258 covmatrix(theta,qOverP) = track->definingParametersCovMatrix()(theta,qOverP)*muSF;
259
260 ATH_MSG_DEBUG ("addParticle: type,errSF,mu p,cov " << type << " " << muSF
261 << " " << mu4vec.E()
262 << " " << track->definingParametersCovMatrix()(d0, d0)
263 << " " << track->definingParametersCovMatrix()(z0, z0)
264 << " " << track->definingParametersCovMatrix()(phi0, phi0)
265 << " " << track->definingParametersCovMatrix()(theta, theta)
266 << " " << track->definingParametersCovMatrix()(qOverP,qOverP)*muSF*muSF
267 << " " << track->definingParametersCovMatrix()(d0, z0)
268 << " " << track->definingParametersCovMatrix()(d0, phi0)
269 << " " << track->definingParametersCovMatrix()(d0, theta)
270 << " " << track->definingParametersCovMatrix()(d0,qOverP)*muSF
271 << " " << track->definingParametersCovMatrix()(z0, phi0)
272 << " " << track->definingParametersCovMatrix()(z0, theta)
273 << " " << track->definingParametersCovMatrix()(z0,qOverP)*muSF
274 << " " << track->definingParametersCovMatrix()(phi0, theta)
275 << " " << track->definingParametersCovMatrix()(phi0,qOverP)*muSF
276 << " " << track->definingParametersCovMatrix()(theta,qOverP)*muSF);
277
278 for (int ii = 0; ii < 5; ii++)
279 for (int jj = ii+1; jj < 5; jj++)
280 covmatrix(jj,ii) = covmatrix(ii,jj);
281
282 ATH_MSG_VERBOSE( "addParticle: cov matrix \n" << covmatrix );
283
284 // Transform cov of q/p to cov p in covariance matrix
285 double P = mu4vec.P();
286 AmgMatrix(5,5) jacobian0;
287 jacobian0.setZero();
288 jacobian0(0,0) = 1.;
289 jacobian0(1,1) = 1.;
290 jacobian0(2,2) = 1.;
291 jacobian0(3,3) = 1.;
292 jacobian0(4,4) = -P*P;
293
294 AmgMatrix(5,5) newmatrix = jacobian0.transpose() * covmatrix * jacobian0;
295
296 ATH_MSG_DEBUG( "addParticle: new matrix \n" << newmatrix );
297
298 bool isOK = doSanityChecksOnCovariance(mu4vec, newmatrix);
299 // Get covariance in cartesian CS
300 AmgMatrix(5,5) covarXYZ;
301 covarXYZ.setZero();
302 convertCovd0z0PhiThetaPToXYZ(mu4vec, newmatrix, covarXYZ);
303 // Save in input
304 input.addConstituent_FourVector_d0z0PhiThetaP(mu4vec, newmatrix, covarXYZ, isOK);
305
306 ATH_MSG_VERBOSE( "addParticle: - added muon");
307
308 }
P
static Double_t P(Double_t *tt, Double_t *par)
Definition LArPhysWaveHECTool.cxx:280
CaloNoise_fillDB.mu
mu
Definition CaloNoise_fillDB.py:51

◆ convertCovd0z0PhiThetaPToXYZ()

void ZMassConstraint::ConstraintFit::convertCovd0z0PhiThetaPToXYZ ( const TLorentzVector & fourVec,
const AmgMatrix(5, 5)& covard0z0PhiThetaP,
AmgMatrix(5, 5)& covarXYZ ) const
private

Definition at line 1128 of file ConstraintFit.cxx.

1131 {
1132 //going from d0,z0,phi,theta,P --> d0,z0,px,py,pz
1133 double phi = fourVec.Phi();
1134 double theta = fourVec.Theta();
1135 double P = fourVec.P();
1136 //std::cout << ":::::L \t"<<P<<"\t"<< phi<<"\t"<<theta<<std::endl;
1137 AmgMatrix(5,5) jacobian;
1138 jacobian.setZero();
1139 jacobian(0,0) = 1.;
1140 jacobian(1,1) = 1.;
1141 jacobian(2,2) = -P * TMath::Sin(theta) * TMath::Sin(phi);
1142 jacobian(2,3) = P * TMath::Sin(theta) * TMath::Cos(phi);
1143 jacobian(3,2) = P * TMath::Cos(theta) * TMath::Cos(phi);
1144 jacobian(3,3) = P * TMath::Cos(theta) * TMath::Sin(phi);
1145 jacobian(3,4) = -P * TMath::Sin(theta);
1146 jacobian(4,2) = TMath::Sin(theta) * TMath::Cos(phi);
1147 jacobian(4,3) = TMath::Sin(theta) * TMath::Sin(phi);
1148 jacobian(4,4) = TMath::Cos(theta);
1149
1150 covarXYZ = jacobian.transpose() * covard0z0PhiThetaP * jacobian;
1151
1152 // std::cout << "initial\n"<< covar << "final \n" << newcovariance << "jac\n"
1153 // << jacobian << std::endl;
1154
1155 }
phi
Scalar phi() const
phi method
Definition AmgMatrixBasePlugin.h:67

◆ convertCovXYZTod0z0PhiThetaP()

void ZMassConstraint::ConstraintFit::convertCovXYZTod0z0PhiThetaP ( const std::vector< TLorentzVector > & particleList,
const Amg::MatrixX & covarXYZ,
Amg::MatrixX & covard0z0PhiThetaP ) const
private

build the Jacobian of the phi,theta,P --> px,py,pz

Definition at line 1159 of file ConstraintFit.cxx.

1162 {
1163 unsigned int matrixSize = 5 * particleList.size();
1164 Amg::MatrixX Jacobian(matrixSize, matrixSize);
1165 Jacobian.setZero();
1167 for(unsigned int ii = 0; ii < particleList.size(); ii++) {
1168 double phi = particleList.at(ii).Phi();
1169 double theta = particleList.at(ii).Theta();
1170 double P = particleList.at(ii).P();
1171 Jacobian( 5*ii, 5*ii) = 1.;
1172 Jacobian(1 + 5*ii, 1 + 5*ii) = 1.;
1173
1174 Jacobian(2 + 5*ii, 2 + 5*ii) = -P * TMath::Sin(theta) * TMath::Sin(phi);
1175 Jacobian(2 + 5*ii, 3 + 5*ii) = P * TMath::Sin(theta) * TMath::Cos(phi);
1176 Jacobian(3 + 5*ii, 2 + 5*ii) = P * TMath::Cos(theta) * TMath::Cos(phi);
1177 Jacobian(3 + 5*ii, 3 + 5*ii) = P * TMath::Cos(theta) * TMath::Sin(phi);
1178 Jacobian(3 + 5*ii, 4 + 5*ii) = -P * TMath::Sin(theta);
1179 Jacobian(4 + 5*ii, 2 + 5*ii) = TMath::Sin(theta) * TMath::Cos(phi);
1180 Jacobian(4 + 5*ii, 3 + 5*ii) = TMath::Sin(theta) * TMath::Sin(phi);
1181 Jacobian(4 + 5*ii, 4 + 5*ii) = TMath::Cos(theta);
1182 }
1183
1184 Amg::MatrixX Jacobianinverse(matrixSize, matrixSize);
1185 Jacobianinverse=Jacobian.inverse();
1186 covard0z0PhiThetaP = Jacobianinverse.transpose() * covarXYZ * Jacobianinverse;
1187 }
Amg::MatrixX
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Definition EventPrimitives.h:27

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
159 return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
160 hndl.value(),
161 hndl.documentation());
162
163 }

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
147 return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
148 }
AthCommonDataStore::declareGaudiProperty
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>
Definition AthCommonDataStore.h:156

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

95{ return m_detStore; }

◆ doMassFit()

StatusCode ZMassConstraint::ConstraintFit::doMassFit ( const ConstraintFitInput & input,
ConstraintFitOutput & output )
virtual

Perform the constrained mass fit.

Implements ZMassConstraint::IConstraintFit.

Definition at line 55 of file ConstraintFit.cxx.

56 {
57
58 // For the moment, call the 'old' methods
59 if (massFitInterface(input)) {
60 ATH_MSG_ERROR( "Could not add in input" );
61 return StatusCode::FAILURE;
62 }
63
64 massFitRun(output);
65
66 // Return gracefully:
67 return StatusCode::SUCCESS;
68 }
ZMassConstraint::ConstraintFit::massFitInterface
int massFitInterface(const ConstraintFitInput &theInput)
Definition ConstraintFit.cxx:385
ZMassConstraint::ConstraintFit::massFitRun
void massFitRun(ConstraintFitOutput &output, double zresol=-1.)
Definition ConstraintFit.cxx:585

◆ doSanityChecksOnCovariance()

bool ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance ( const TLorentzVector & vector,
const AmgMatrix(5, 5)& covar ) const
private

Definition at line 1091 of file ConstraintFit.cxx.

1093 {
1094 bool isOK = true;
1095 if ( sqrt(covar(qP,qP))/vector.P()>1.0 ) {
1096 isOK = false;
1097 ATH_MSG_WARNING("ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance:: Fractional uncertainty on P = "
1098 << sqrt(covar(qP,qP))/vector.P() << " is > 1 ...this is not ok...");
1099
1100 }
1101
1102 if( covar(theta,theta) > 1.e-2) {
1103 isOK = false;
1104 ATH_MSG_WARNING("ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance:: Uncertainty on Theta = "
1105 << covar(theta,theta) << " is > 10^-2 ...this is not ok");
1106 }
1107 if( covar(phi0,phi0) > 1.e-2) {
1108 isOK = false;
1109 ATH_MSG_WARNING("ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance:: Uncertainty on Phi = "
1110 << covar(phi0,phi0) << " is > 10^-2 ...this is not ok");
1111 }
1112
1113 for(int i=0;i<5;i++) {
1114 for(int j=0;j<5;j++) {
1115 if(i==j) continue;
1116 if(fabs(covar(i,j))/sqrt(covar(i,i))/sqrt(covar(j,j)) > 1.) {
1117 ATH_MSG_WARNING( "ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance:: Off-diagonal term "
1118 << i << " " << j << " is > 1 - doesn't look ok");
1119 isOK = false;
1120 }
1121 }
1122 }
1123
1124 return isOK;
1125 }
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition AthMsgStreamMacros.h:32
TauHitVars::j
float j(const xAOD::IParticle &, const xAOD::TrackMeasurementValidation &hit, const Eigen::Matrix3d &jab_inv)
Definition ConstituentLoaderTauHit.cxx:102

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

85{ return m_evtStore; }

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ getKey()

SG::sgkey_t asg::AsgTool::getKey ( const void * ptr) const
inherited

Get the (hashed) key of an object that is in the event store.

This is a bit of a special one. StoreGateSvc and xAOD::Event both provide ways for getting the SG::sgkey_t key for an object that is in the store, based on a bare pointer. But they provide different interfaces for doing so.

In order to allow tools to efficiently perform this operation, they can use this helper function.

See also
asg::AsgTool::getName
Parameters
ptrThe bare pointer to the object that the event store should know about
Returns
The hashed key of the object in the store. If not found, an invalid (zero) key.

Definition at line 119 of file AsgTool.cxx.

119 {
120
121#ifdef XAOD_STANDALONE
122 // In case we use @c xAOD::Event, we have a direct function call
123 // for this.
124 return evtStore()->event()->getKey( ptr );
125#else
126 const SG::DataProxy* proxy = evtStore()->proxy( ptr );
127 return ( proxy == nullptr ? 0 : proxy->sgkey() );
128#endif // XAOD_STANDALONE
129 }
AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore
ServiceHandle< StoreGateSvc > & evtStore()
Definition AthCommonDataStore.h:85
StateLessPT_NewConfig.proxy
proxy
Definition StateLessPT_NewConfig.py:415

◆ getMassError() [1/4]

double ZMassConstraint::ConstraintFit::getMassError ( const ConstraintFitInput & firstInput,
const ConstraintFitInput & secondInput = ConstraintFitInput() )
virtual

Calculate the mass error without fit - use just the inputs.

Implements ZMassConstraint::IConstraintFit.

Definition at line 906 of file ConstraintFit.cxx.

908 {
909 unsigned int nPartFirst = firstInput.getNConstituents();
910 unsigned int nPartSecond = secondInput.getNConstituents();
911 unsigned int nPart = nPartFirst + nPartSecond;
912 std::vector<TLorentzVector> particles(nPart);
913 std::vector<AmgMatrix(3,3)> covariances(nPart);
914 for ( unsigned int iobj = 0; iobj < firstInput.getNConstituents(); ++iobj) {
915 covariances[iobj].setZero();
916 particles[iobj] = firstInput.getConstituentFourVector(iobj);
917 firstInput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj]);
918 }
919 for ( unsigned int iobj = 0; iobj < secondInput.getNConstituents(); ++iobj) {
920 covariances[iobj + nPartFirst].setZero();
921 particles[iobj + nPartFirst] = secondInput.getConstituentFourVector(iobj);
922 secondInput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj + nPartFirst]);
923 }
924 return getMassError( particles, covariances );
925 }
ZMassConstraint::ConstraintFit::getMassError
double getMassError(const ConstraintFitInput &firstInput, const ConstraintFitInput &secondInput=ConstraintFitInput())
Calculate the mass error without fit - use just the inputs.
Definition ConstraintFit.cxx:906

◆ getMassError() [2/4]

double ZMassConstraint::ConstraintFit::getMassError ( const ConstraintFitOutput & fitOutput,
const ConstraintFitInput & extraInput = ConstraintFitInput() )
virtual

Calculate the mass error with the result of fit and possible additional input (4-vecs/cov mat) that has NOT been fit.

Implements ZMassConstraint::IConstraintFit.

Definition at line 930 of file ConstraintFit.cxx.

932 {
933 unsigned int nPartFirst = fitOutput.getNConstituents();
934 unsigned int nPartSecond = extraInput.getNConstituents();
935 unsigned int nPart = nPartFirst + nPartSecond;
936 std::vector<TLorentzVector> particles(nPart);
937 std::vector<AmgMatrix(3,3)> covariances(nPart);
938 for ( unsigned int iobj = 0; iobj < fitOutput.getNConstituents(); ++iobj) {
939 covariances[iobj].setZero();
940 particles[iobj] = fitOutput.getConstituentFourVector(iobj);
941 fitOutput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj]);
942 }
943 for ( unsigned int iobj = 0; iobj < extraInput.getNConstituents(); ++iobj) {
944 covariances[iobj + nPartFirst].setZero();
945 particles[iobj + nPartFirst] = extraInput.getConstituentFourVector(iobj);
946 extraInput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj + nPartFirst]);
947 }
948 return getMassError( particles, covariances );
949 }

◆ getMassError() [3/4]

double ZMassConstraint::ConstraintFit::getMassError ( const ConstraintFitOutput & fitOutput,
const ConstraintFitOutput & secondFitOutput )
virtual

Calculate the mass error with the result of two mass fits.

Implements ZMassConstraint::IConstraintFit.

Definition at line 955 of file ConstraintFit.cxx.

957 {
958 unsigned int nPartFirst = fitOutput.getNConstituents();
959 unsigned int nPartSecond = secondFitOutput.getNConstituents();
960 unsigned int nPart = nPartFirst + nPartSecond;
961 std::vector<TLorentzVector> particles(nPart);
962 std::vector<AmgMatrix(3,3)> covariances(nPart);
963 for ( unsigned int iobj = 0; iobj < fitOutput.getNConstituents(); ++iobj) {
964 covariances[iobj].setZero();
965 particles[iobj] = fitOutput.getConstituentFourVector(iobj);
966 fitOutput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj]);
967 }
968 for ( unsigned int iobj = 0; iobj < secondFitOutput.getNConstituents(); ++iobj) {
969 covariances[iobj + nPartFirst].setZero();
970 particles[iobj + nPartFirst] = secondFitOutput.getConstituentFourVector(iobj);
971 secondFitOutput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj + nPartFirst]);
972 }
973 return getMassError( particles, covariances );
974 }

◆ getMassError() [4/4]

double ZMassConstraint::ConstraintFit::getMassError ( const std::vector< TLorentzVector > & particles,
const std::vector< AmgMatrix(3, 3)> & covariances )
private

Definition at line 977 of file ConstraintFit.cxx.

979 {
980 double masserror=0;
981
982 if (particles.size() != covariances.size()) {
983 ATH_MSG_ERROR( "getMassError: Number of particles is not equal to the number of covariance matrices");
984 return 1E11;
985 }
986
987 ATH_MSG_DEBUG ("getMassError: 1 " << particles.size() << "/" << covariances.size());
988
989 //composite lorentz vector and its invariant mass
990 TLorentzVector combv;
991 for ( auto lv : particles ) combv += lv;
992 double invmass = combv.M();
993
994 ATH_MSG_VERBOSE ("getMassError: 1.1 invmass " << invmass);
995
996 //Calculation of 1x4 jacobian - used in derivation of var(M)
997 Amg::MatrixX jacobianNPtoM(1,4);
998 jacobianNPtoM.setZero();
999 jacobianNPtoM(0,0) = -1.*combv.Px()/invmass; // dM/dpxN
1000 jacobianNPtoM(0,1) = -1.*combv.Py()/invmass; // dM/dpyN
1001 jacobianNPtoM(0,2) = -1.*combv.Pz()/invmass; // dM/dpzN
1002 jacobianNPtoM(0,3) = combv.E()/invmass; // dM/dEN
1003
1004 ATH_MSG_VERBOSE ("getMassError: 2\n" << jacobianNPtoM);
1005
1006 int iobj = -1;
1007 for ( auto lv : particles ) {
1008 ++iobj;
1009
1010 ATH_MSG_VERBOSE ("getMassError: 1.2 iobj " << iobj);
1011
1012 // Constituent 3*3 covariance matrix (Phi, Theta, P)
1013 const AmgMatrix(3,3)& covPhiThetaP = covariances[iobj];
1014
1015 ATH_MSG_VERBOSE ("getMassError: 1.3 covPhiThetaP\n " << covPhiThetaP);
1016
1017 double theta = lv.Theta();
1018 double phi = lv.Phi();
1019 double p = lv.P();
1020 double e = lv.E();
1021 double m = lv.M();
1022
1023 ATH_MSG_VERBOSE ("getMassError: 2.1 iobj " << iobj << "covPhiThetaP\n" << covPhiThetaP);
1024
1025 // Convert (Phi, Theta, P) to (Phi, eta, P), as covPhiEtaP
1026 AmgMatrix(3,3) jacobian1;
1027 jacobian1.setZero();
1028 jacobian1(0,0) = 1.; //
1029 jacobian1(1,1) = -1*(1./sin(theta)); // deta/dtheta
1030 jacobian1(2,2) = 1.; //
1031 AmgMatrix(3,3) covPhiEtaP = jacobian1 * covPhiThetaP * jacobian1.transpose();
1032
1033 ATH_MSG_VERBOSE ("getMassError: 2.2 jacobian1\n" << jacobian1);
1034 ATH_MSG_VERBOSE ("getMassError: 2.4 covPhiEtaP\n" << covPhiEtaP);
1035
1036
1037 // Rearrange things to get (E, eta, phi, M), as covEEtaPhiM
1038 AmgMatrix(4,4) covEEtaPhiM;
1039 covEEtaPhiM.setZero();
1040 for(int i = 0; i < 3; i++) {
1041 for(int j = 0; j < 3; j++) {
1042 covEEtaPhiM(i,j) = covEEtaPhiM(j,i) = covPhiEtaP(2-i,2-j);
1043 }
1044 }
1045
1046 ATH_MSG_VERBOSE ("getMassError: 3 covEEtaPhiM\n" << covEEtaPhiM);
1047
1048 //convert (E, eta, phi, M) to (px,py,pz,E) representation, as covPxPyPzE
1049 AmgMatrix(4,4) jacobian2;
1050 jacobian2.setZero();
1051 jacobian2(0,0) = e/p * sin(theta)*cos(phi); // dpx/dE
1052 jacobian2(1,0) = e/p * sin(theta)*sin(phi); // dpy/dE
1053 jacobian2(2,0) = e/p * cos(theta); // dpy/dE
1054 jacobian2(3,0) = 1; // dE/dE
1055
1056 jacobian2(0,1) = -p * cos(phi) * cos(theta)*sin(theta); // dpx/deta
1057 jacobian2(1,1) = -p * sin(phi) * cos(theta)*sin(theta); // dpy/deta
1058 jacobian2(2,1) = p * sin(theta)*sin(theta); // dpz/deta;
1059
1060 jacobian2(0,2) = -p * sin(theta) * sin(phi); // dpx/dphi
1061 jacobian2(1,2) = p * sin(theta) * cos(phi); // dpy/dphi
1062
1063 jacobian2(0,3) = -m/p * sin(theta)*cos(phi); // dpx/dM
1064 jacobian2(1,3) = -m/p * sin(theta)*sin(phi); // dpy/dM
1065 jacobian2(2,3) = -m/p * cos(theta); // dpz/dM
1066
1067 AmgMatrix(4,4) covPxPyPzE = jacobian2 * covEEtaPhiM * jacobian2.transpose();
1068
1069 ATH_MSG_VERBOSE ("getMassError: 4 jacobian2\n" << jacobian2);
1070 ATH_MSG_VERBOSE ("getMassError: 5 covPxPyPzE\n" << covPxPyPzE);
1071
1072 //Get mass variance
1073 // (Similarity transform: C = A*B*A^T (from egammaFourMomentumError/GeneralUtils.cxx))
1074 Amg::MatrixX em = jacobianNPtoM * covPxPyPzE * jacobianNPtoM.transpose();
1075 masserror += em(0,0);
1076
1077 ATH_MSG_VERBOSE ("getMassError: 6 em\n" << em);
1078
1079 }
1080
1081 //Return square root of variance.
1082 if (masserror < 0.) {
1083 ATH_MSG_WARNING( "getMassError: Mass covariance element less than zero! Returning 1E11 ...");
1084 masserror = 1E11;
1085 }
1086 return sqrt(masserror);
1087
1088 }
if
if(pathvar)
Definition SealSharedLib.cxx:167

◆ getName()

const std::string & asg::AsgTool::getName ( const void * ptr) const
inherited

Get the name of an object that is / should be in the event store.

This is a bit of a special one. StoreGateSvc and xAOD::Event both provide ways for getting the std::string name for an object that is in the store, based on a bare pointer. But they provide different interfaces for doing so.

In order to allow tools to efficiently perform this operation, they can use this helper function.

See also
asg::AsgTool::getKey
Parameters
ptrThe bare pointer to the object that the event store should know about
Returns
The string name of the object in the store. If not found, an empty string.

Definition at line 106 of file AsgTool.cxx.

106 {
107
108#ifdef XAOD_STANDALONE
109 // In case we use @c xAOD::Event, we have a direct function call
110 // for this.
111 return evtStore()->event()->getName( ptr );
112#else
113 const SG::DataProxy* proxy = evtStore()->proxy( ptr );
114 static const std::string dummy = "";
115 return ( proxy == nullptr ? dummy : proxy->name() );
116#endif // XAOD_STANDALONE
117 }

◆ getProperty()

template<class T>
const T * asg::AsgTool::getProperty ( const std::string & name) const
inherited

Get one of the tool's properties.

◆ initialize()

StatusCode ZMassConstraint::ConstraintFit::initialize ( void )
virtual

Initialize constraint fit.

Reimplemented from asg::AsgTool.

Definition at line 35 of file ConstraintFit.cxx.

36 {
37 ATH_MSG_INFO ("initialize: retrieve tools" << m_energyRescaler.name() << " and "
38 << m_mu_resolSFTool.name());
39
40 if (!m_energyRescaler.retrieve().isSuccess()) {
41 ATH_MSG_ERROR ("initialize: unable to retrieve EgammaCalibrationAndSmearingTool");
42 return StatusCode::FAILURE;
43 }
44
45 if (!m_mu_resolSFTool.retrieve().isSuccess()) {
46 ATH_MSG_ERROR ("initialize: unable to retrieve MuonCalibrationAndSmearingTool");
47 return StatusCode::FAILURE;
48 }
49
50 // Return gracefully:
51 return StatusCode::SUCCESS;
52 }
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition AthMsgStreamMacros.h:31

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ likelihoodMass()

double ZMassConstraint::ConstraintFit::likelihoodMass ( double MassResol)
private

Definition at line 435 of file ConstraintFit.cxx.

436 {
437 const int rows =4;
438 const int columns = m_nobj;
439 double **p = 0;
440 p = new double *[rows] ;
441 for( int i = 0 ; i < rows ; i++ )
442 p[i] = new double[columns];
443
444 for(unsigned int iobj=0; iobj<m_nobj; iobj++)
445 {
446 p[3][iobj] = m_objmass[iobj]*m_objmass[iobj];
447 for(int i=0;i<3;i++)
448 {
449 p[i][iobj] = m_parametersInit(i + 3*iobj, 0);
450 p[3][iobj] += p[i][iobj]*p[i][iobj];
451 }
452 p[3][iobj] = sqrt(p[3][iobj]);
453 }
454 double Etot = 0.;
455 double Xtot = 0.;
456 double Ytot = 0.;
457 double Ztot = 0.;
458 for(unsigned int ij=0;ij<m_nobj;ij++)
459 {
460 Etot += p[3][ij];
461 Xtot += p[0][ij];
462 Ytot += p[1][ij];
463 Ztot += p[2][ij];
464 }
465 for( int i = 0 ; i < rows ; i++ )
466 delete [] p[i];
467 delete [] p;
468
469 double initMass = Etot * Etot - Xtot * Xtot - Ytot * Ytot - Ztot * Ztot;
470 initMass = sqrt(initMass);
471
472 double sig = MassResol;
473 double xLeft = initMass;
474 double xRight = m_conMass;
475 if (m_conMass < initMass) {
476 xLeft = m_conMass;
477 xRight = initMass;
478 }
479 double dLinitL = (initMass - xLeft) / sig / sig - 2.*(xLeft - m_conMass) / ((xLeft - m_conMass) * (xLeft - m_conMass) + m_conWidth * m_conWidth);
480 double dLinitR = (initMass - xRight) / sig / sig - 2.*(xRight - m_conMass) / ((xRight - m_conMass) * (xRight - m_conMass) + m_conWidth * m_conWidth);
481 if (dLinitL * dLinitR < 0.) {
482 while (xRight - xLeft > 1.) { //1 MeV
483 double xM = (xRight + xLeft) / 2.;
484 double dL = (initMass - xM) / sig / sig - 2.*(xM - m_conMass) / ((xM - m_conMass) * (xM - m_conMass) + m_conWidth * m_conWidth);
485 if (dL * dLinitL < 0.) {
486 xRight = xM;
487 dLinitR = dL;
488 } else {
489 xLeft = xM;
490 dLinitL = dL;
491 }
492 }
493 return (xLeft + xRight) / 2.;
494 } else {
495 if (dLinitL > dLinitR)
496 return xLeft;
497 else
498 return xRight;
499 }
500 }
ZMassConstraint::ConstraintFit::m_objmass
std::vector< double > m_objmass
Definition ConstraintFit.h:122

◆ likelihoodMass2()

double ZMassConstraint::ConstraintFit::likelihoodMass2 ( void )
private

Definition at line 503 of file ConstraintFit.cxx.

504 {
505 ATH_MSG_VERBOSE( "likelihoodMass2: entering" );
506 const int rows =4;
507 const int columns = m_nobj;
508 double **p = 0;
509 p = new double *[rows] ;
510 for( int i = 0 ; i < rows ; i++ )
511 p[i] = new double[columns];
512
513 for(unsigned int iobj=0; iobj<m_nobj; iobj++) {
514 p[3][iobj] = m_objmass[iobj]*m_objmass[iobj];
515 for(int i=0;i<3;i++) {
516 p[i][iobj] = m_parametersInit(i + 3*iobj, 0);
517 p[3][iobj] += p[i][iobj]*p[i][iobj];
518 }
519 p[3][iobj] = sqrt(p[3][iobj]);
520 }
521 double etot = 0.;
522 double xtot = 0.;
523 double ytot = 0.;
524 double ztot = 0.;
525 for(unsigned int ij=0;ij<m_nobj;ij++) {
526 etot += p[3][ij];
527 xtot += p[0][ij];
528 ytot += p[1][ij];
529 ztot += p[2][ij];
530 }
531
532 double initMass = etot * etot - xtot * xtot - ytot * ytot - ztot * ztot;
533 initMass = sqrt(initMass);
534 //what is the uncertainty on initMass?
535 Amg::MatrixX jacobianMass(3*m_nobj, 1);
536 jacobianMass.setZero();
537 for(unsigned int ii=0;ii<m_nobj;ii++) {
538 jacobianMass(3*ii + 0, 0) = (2.*(etot) * p[0][ii] / p[3][ii] - 2.* xtot);
539 jacobianMass(3*ii + 1, 0) = (2.*(etot) * p[1][ii] / p[3][ii] - 2.* ytot);
540 jacobianMass(3*ii + 2, 0) = (2.*(etot) * p[2][ii] / p[3][ii] - 2.* ztot);
541 }
542 ATH_MSG_VERBOSE( "likelihoodMass2: the jacobian mass/matrix \n" << jacobianMass << "\n" << m_covarianceInit );
543
544
545 Amg::MatrixX sigm = jacobianMass.transpose() * m_covarianceInit * jacobianMass;
546
547 ATH_MSG_VERBOSE( "likelihoodMass2: sigm \n" << sigm );
548
549 double sig = sigm(0, 0);
550 sig = sqrt(sig);
551 m_resolution = sig;
552 sig /= 2 * initMass;
553 double MassResol = sig;
554 double maxmass = initMass;
555
556 ATH_MSG_VERBOSE( "likelihoodMass2: 2 " );
557
558 double max = -(maxmass - initMass) * (maxmass - initMass) / 2. / MassResol / MassResol
559 - log((maxmass * maxmass - m_conMass * m_conMass) * (maxmass * maxmass - m_conMass * m_conMass) + m_conMass * m_conMass * m_conWidth * m_conWidth);
560
561 for (int i = 1; i < 401; i++) {
562 double ytest = initMass + (m_conMass - initMass) / 400 * i;
563
564 double val = -(ytest - initMass) * (ytest - initMass) / 2. / MassResol / MassResol
565 - log((ytest * ytest - m_conMass * m_conMass) * (ytest * ytest - m_conMass * m_conMass) + m_conMass * m_conMass * m_conWidth * m_conWidth);
566
567 if (val > max) {
568 max = val;
569 maxmass = ytest;
570 }
571 }
572
573 ATH_MSG_VERBOSE( "likelihoodMass2: 3 " );
574
575 for( int i = 0 ; i < rows ; i++ )
576 delete [] p[i];
577 delete [] p;
578
579 ATH_MSG_VERBOSE( "likelihoodMass2: maxmass " << maxmass );
580
581 return maxmass;
582 }
max
#define max(a, b)
Definition cfImp.cxx:41
Pythia8_RapidityOrderMPI.val
val
Definition Pythia8_RapidityOrderMPI.py:14

◆ massFit()

double ZMassConstraint::ConstraintFit::massFit ( const Amg::MatrixX & ,
const Amg::MatrixX & var,
double mass,
Amg::MatrixX & pOut,
Amg::MatrixX &  )
private

Definition at line 857 of file ConstraintFit.cxx.

859 {
860 Amg::MatrixX ivar = var.inverse();
861
862 bool doIter = true;
863 int maxIterations = 20;
864 int iIter = 0;
865 double constraintValue = -1.e10;
866 ATH_MSG_VERBOSE( "massFit: transpose matrix \n" << pOut.transpose()
867 << " Iterations " << iIter
868 << " parameters \n" << m_parametersFinal
869 << " cov matrix \n" << m_covarianceFinal
870 << " Iteration " << iIter);
871 while (doIter) {
872 Amg::MatrixX p_old(pOut);
873 constraintValue = massFitCalculation(var, mass, pOut); //call the fit
874
875 // Convergence criteria
876 // 1. the parameters should not change too much (<1.e-6 relative change)
877 // 2. the constraint should be satisfied very well (<1.e-6 absolute)
878 double maxDiff = 5.e15;
879 Amg::MatrixX diff = (pOut - p_old);
880 for (unsigned int i = 0; i < m_parameters * m_nobj; i++) {
881 diff(i, 0) = diff(i, 0) / p_old(i, 0);
882 if (maxDiff > diff(i, 0))
883 maxDiff = diff(i, 0);
884 }
885 if ((maxDiff < 1.e-4 && fabs(constraintValue) < 5.e-5) || maxIterations <= iIter)
886 doIter = false;
887 ATH_MSG_VERBOSE( "massFit: transpose matrix \n" << pOut.transpose()
888 << " Iterations " << iIter
889 << " parameters \n" << m_parametersFinal
890 << " cov matrix \n" << m_covarianceFinal
891 << " \n Iteration " << iIter << " doIter " << doIter);
892 iIter++;
893
894 }
895 ATH_MSG_VERBOSE( "massFit: transpose matrix \n" << pOut.transpose()
896 << " Iterations " << iIter
897 << " parameters \n" << m_parametersFinal
898 << " cov matrix \n" << m_covarianceFinal
899 << " Iteration " << iIter);
900 // double chi2 = 0.;//calculateChi2(p0, var, mass);
901 return 0;
902 }
diff
void diff(const Jet &rJet1, const Jet &rJet2, std::map< std::string, double > varDiff)
Difference between jets - Non-Class function required by trigger.
Definition Jet.cxx:631
ZMassConstraint::ConstraintFit::massFitCalculation
double massFitCalculation(const Amg::MatrixX &var, double mass, Amg::MatrixX &p0)
Definition ConstraintFit.cxx:787

◆ massFitCalculation()

double ZMassConstraint::ConstraintFit::massFitCalculation ( const Amg::MatrixX & var,
double mass,
Amg::MatrixX & p0 )
private

Definition at line 787 of file ConstraintFit.cxx.

788 {
789 const int rows =4;
790 const int columns = m_nobj;
791
792 double **p = 0;
793 p = new double *[rows] ;
794 for( int i = 0 ; i < rows ; i++ )
795 p[i] = new double[columns];
796
797 for(unsigned int iobj=0; iobj<m_nobj; iobj++) {
798 p[3][iobj] = m_objmass[iobj]*m_objmass[iobj];
799 for(int i=0;i<3;i++)
800 {
801 p[i][iobj] = p0(i+3*iobj, 0);
802 p[3][iobj] += p[i][iobj]*p[i][iobj];
803 }
804 p[3][iobj] = sqrt(p[3][iobj]);
805 }
806
807 double Etot = 0.;
808 double Xtot = 0.;
809 double Ytot = 0.;
810 double Ztot = 0.;
811 for(unsigned int ij=0;ij<m_nobj;ij++)
812 {
813 Etot += p[3][ij];
814 Xtot += p[0][ij];
815 Ytot += p[1][ij];
816 Ztot += p[2][ij];
817 }
818 double XYZtot[3] = {Xtot, Ytot, Ztot};
819 std::vector<double> constraintD;
820 for(unsigned int ij=0;ij<m_nobj;ij++)
821 for(int ik=0;ik<3;ik++)
822 constraintD.push_back((2.*(Etot) * p[ik][ij] / p[3][ij] - 2.* XYZtot[ik]));
823
824 double constraintd = Etot * Etot - Xtot * Xtot - Ytot * Ytot - Ztot * Ztot;
825 constraintd = constraintd - mass * mass;
826 //std::cout << "constraint = " <<constraintd << " mass value " <<mass<<std::endl;
827 Amg::MatrixX D(1,3*m_nobj);
828 for(unsigned int ij=0;ij<3*m_nobj;ij++)
829 D(0, ij) = constraintD.at(ij);
830 AmgMatrix(1,1) d;
831 d(0, 0) = constraintd;
832
833 Amg::MatrixX DVD(D * var * D.transpose());
834
835 Amg::MatrixX DVDinverse(1,1);
836 DVDinverse=DVD.inverse();
837
838 Amg::MatrixX VD(DVDinverse);
839 Amg::MatrixX lambda(VD * d);
840 Amg::MatrixX test(p0 - var * D.transpose() * lambda);
841 Amg::MatrixX Vp(var - var * D.transpose() * VD * D * var);
842 for(unsigned int ij=0;ij<3*m_nobj;ij++) {
843 p0(ij, 0) = test(ij, 0);
844 m_parametersFinal(ij,0) = p0(ij,0);
845 for(unsigned int jk = 0; jk < 3 * m_nobj; jk++)
846 m_covarianceFinal(ij,jk) = Vp(ij,jk);
847 }
848 double constraintValue = d(0, 0);
849
850 for( int i = 0 ; i < rows ; i++ )
851 delete [] p[i];
852 delete [] p;
853 return constraintValue;
854 }
Base_Fragment.mass
mass
Definition Sherpa_i/share/common/Base_Fragment.py:59
TRTCalib_cfilter.p0
p0
Definition TRTCalib_cfilter.py:129
TrigInDetValidation_Base.test
test
Definition TrigInDetValidation_Base.py:142
hist_file_dump.d
d
Definition hist_file_dump.py:142

◆ massFitInterface()

int ZMassConstraint::ConstraintFit::massFitInterface ( const ConstraintFitInput & theInput)
private

Definition at line 385 of file ConstraintFit.cxx.

386 {
387 m_objmass.clear();
388 m_theInput = theInput;
389
390 // the m_parameters first parameters of pin are the fitted parameters
391 // the next one is the mass of the particle
392 if(theInput.getNConstituents() <2 ) {
393 ATH_MSG_ERROR( "massFitInterface: Found number of input particles less than 2");
394 return 1;
395 }
396 else {
397 ATH_MSG_DEBUG( "massFitInterface: nparticles " << theInput.getNConstituents());
398 }
399
400 m_nobj = theInput.getNConstituents();
401 unsigned int dimension = m_parameters * m_nobj;
402 m_parametersInit = Amg::MatrixX(dimension, 1);
403 m_parametersInit.setZero();
404 m_covarianceInit = Amg::MatrixX(dimension, dimension);
405 m_covarianceInit.setZero();
406 m_parametersFinal= Amg::MatrixX(dimension, 1);
407 m_parametersFinal.setZero();
408 m_covarianceFinal= Amg::MatrixX(dimension, dimension);
409 m_covarianceFinal.setZero();
410
411 for (unsigned int iobj = 0; iobj < m_nobj; iobj++)
412 {
413 m_objmass.push_back(theInput.getConstituentFourVector(iobj).M());
414 m_parametersInit(0 + iobj*3, 0) = theInput.getConstituentFourVector(iobj).Px();
415 m_parametersInit(1 + iobj*3, 0) = theInput.getConstituentFourVector(iobj).Py();
416 m_parametersInit(2 + iobj*3, 0) = theInput.getConstituentFourVector(iobj).Pz();
417
418 for (unsigned int i = 0; i < m_parameters; i++)
419 for (unsigned int j = 0; j < m_parameters; j++)
420 m_covarianceInit(i + 3*iobj, j + 3*iobj) = (theInput.getCovarianceCartesian(iobj))(i + 2, j + 2); //keep only the lower right 3x3 part
421 }
422
423 for(unsigned int i=0;i<dimension;i++)
424 {
425 m_parametersFinal(i,0) = m_parametersInit(i,0);
426 for(unsigned int j=0;j<dimension;j++)
427 m_covarianceFinal(i,j) = m_covarianceInit(i,j);
428 }
429 ATH_MSG_VERBOSE( "massFitInterface: parameters \n" << m_parametersInit );
430 ATH_MSG_VERBOSE( "massFitInterface: cov matrix \n" << m_covarianceInit );
431 return 0;
432 }

◆ massFitRun()

void ZMassConstraint::ConstraintFit::massFitRun ( ConstraintFitOutput & output,
double zresol = -1. )
private

Definition at line 585 of file ConstraintFit.cxx.

586 {
587 ATH_MSG_VERBOSE( "massFitRun: entering ");
588
589 if(!m_ignoreInputChecks && !m_theInput.isOK()) {
590 ATH_MSG_WARNING("ConstraintFit::massFitRun() the input is not OK.");
591 ATH_MSG_WARNING("ConstraintFit::massFitRun()(cont) Following H->4l group decision (Feb 1, 2013), this event will not be fit");
592 ATH_MSG_WARNING("ConstraintFit::massFitRun()(cont) To bypass this behaviour set the ingore flag");
593
594 std::vector<TLorentzVector> particleList;
595 for(unsigned int i=0;i<m_nobj;i++)
596 {
597 TLorentzVector particle;
598 particle.SetXYZM(m_parametersInit(0 + i*m_parameters,0),
599 m_parametersInit(1 + i*m_parameters,0),
600 m_parametersInit(2 + i*m_parameters,0),
601 m_objmass[i]);
602 particleList.push_back(particle);
603
604 }
605 Amg::MatrixX covOut(5*m_nobj,5*m_nobj);
606 covOut.setZero();
607 for(unsigned int iobj=0;iobj<m_nobj;iobj++)
608 for(int ii=0;ii<5;ii++)
609 for(int jj=0;jj<5;jj++)
610 covOut(ii+5*iobj,jj+5*iobj) = m_theInput.getCovarianceCartesian(iobj)(ii,jj);
611
612 ATH_MSG_VERBOSE( "massFitRun: cov 1 ");
613
614 // convert covariance matrix from cartesian to d0z0PhiThetaP
615 Amg::MatrixX covard0z0PhiThetaP;
616 covard0z0PhiThetaP.setZero();
617 convertCovXYZTod0z0PhiThetaP(particleList, covOut, covard0z0PhiThetaP);
618 // save output
619 output.setFitOutput(particleList, covOut, covard0z0PhiThetaP);
620
621 // ATH_MSG_VERBOSE( "massFitRun: cov 1 " << covOut);
622
623 return;
624 }
625
626 ATH_MSG_VERBOSE( "massFitRun: 1 ");
627
628 Amg::MatrixX parametersFit = m_parametersInit;
629 double mass = -1.;
630 if (!m_conHasWidth)
631 mass = m_conMass;
632 else if(zresol<0.)
633 mass = likelihoodMass2();
634 else
635 mass = likelihoodMass(zresol);
636
637 ATH_MSG_VERBOSE( "massFitRun: resolution/mass " << m_resolution << "/" << mass);
638
639 // do not use chi2:
640 // double chi2 = massFit(&m_parametersInit, &m_covarianceInit, mass, &parametersFit, &m_covarianceFinal);
641 massFit(m_parametersInit, m_covarianceInit, mass, parametersFit, m_covarianceFinal);
642
643 ATH_MSG_VERBOSE( "massFitRun: 2");
644
645 std::vector<TLorentzVector> particleList;
646 for(unsigned int i=0;i<m_nobj;i++) {
647 TLorentzVector particle;
648 particle.SetXYZM(m_parametersFinal(0 + i*m_parameters,0),
649 m_parametersFinal(1 + i*m_parameters,0),
650 m_parametersFinal(2 + i*m_parameters,0),
651 m_objmass[i]);
652 particleList.push_back(particle);
653 }
654
655 ATH_MSG_VERBOSE( "massFitRun: 3 ");
656
657 Amg::MatrixX covOut(5*m_nobj,5*m_nobj);
658 covOut.setZero();
659 for(unsigned int iobj=0;iobj<m_nobj;iobj++)
660 for(int ii=0;ii<2;ii++)
661 for(int jj=0;jj<2;jj++)
662 covOut(ii+5*iobj,jj+5*iobj) = m_theInput.getCovarianceCartesian(iobj)(ii,jj);
663 for(unsigned int iobj=0;iobj<m_nobj;iobj++)
664 for(int ii=2;ii<5;ii++)
665 for(int jj=2;jj<5;jj++)
666 covOut(ii+5*iobj,jj+5*iobj) = m_covarianceFinal(ii-2+iobj*3,jj-2+iobj*3);
667
668 ATH_MSG_VERBOSE( "massFitRun: 4 ");
669
670 // empirical corrections for constraint width
671 if(!m_conHasWidth) {
672 // correlation terms
673 for(unsigned int iobj1=0;iobj1<m_nobj;iobj1++) {
674 for(unsigned int iobj2=0;iobj2<m_nobj;iobj2++) {
675 if(iobj1 == iobj2) continue;
676 for(int ii=2;ii<5;ii++)
677 for(int jj=2;jj<5;jj++)
678 covOut(ii + 5*iobj1, jj+5*iobj2) = m_covarianceFinal(ii-2+iobj1*3,jj-2+iobj2*3);
679 }
680 }
681 }
682
683 // ATH_MSG_DEBUG( "massFitRun: cov 2 ");
684
685 // convert covariance matrix from cartesian to d0z0PhiThetaP
686 Amg::MatrixX covard0z0PhiThetaP;
687 covard0z0PhiThetaP.setZero();
688 convertCovXYZTod0z0PhiThetaP(particleList, covOut, covard0z0PhiThetaP);
689 // save output
690 output.setFitOutput(particleList, covOut, covard0z0PhiThetaP);
691
692 // ATH_MSG_DEBUG( "massFitRun: cov 2 " << covOut);
693
694
695 ATH_MSG_VERBOSE( "massFitRun: 5 conHasWidth " << m_conHasWidth);
696
697 // ATH_MSG_DEBUG( "massFitRun: 5.5 conHasWidth " << m_conHasWidth);
698
699 if(m_conHasWidth) {
700 TLorentzVector aa;
701 for(unsigned int iobj = 0; iobj < m_nobj; iobj++) {
702 aa += m_theInput.getConstituentFourVector(iobj);
703 }
704 double massInit = aa.M();
705
706 // ATH_MSG_DEBUG( "massFitRun: 5.6 massInit " << massInit);
707
708 if(massInit<83.e3) {
709
710 ATH_MSG_VERBOSE( "massFitRun: 6 \n" << covOut);
711
712 for(unsigned int iobj = 0; iobj < m_nobj; iobj++) {
713 for(int ii = 0; ii < 5; ii++) {
714 for(int jj = 0; jj < 5; jj++) {
715 covOut(ii+5*iobj, jj+5*iobj) = m_theInput.getCovarianceCartesian(iobj)(ii, jj);
716 }
717 }
718 }
719
720 ATH_MSG_VERBOSE( "massFitRun: 6.1 \n" << covOut);
721 }
722 else {
723
724 // ATH_MSG_DEBUG( "massFitRun: 7 ");
725
726 std::vector<double> resolInit;
727 std::vector<double> resolFinal;
728 std::vector<double> resolFinalnew;
729 resolInit.reserve(m_nobj);
730 resolFinal.reserve(m_nobj);
731 resolFinalnew.reserve(m_nobj);
732 double sumResolInit=0.;
733 const double constraintWidthSquare = /*m_resolution*m_resolution + */m_conWidth*m_conWidth/2.35/2.35;
734 for(unsigned int iobj = 0; iobj < m_nobj; iobj++) {
735
736 // ATH_MSG_DEBUG( "massFitRun: 7 iobj " << iobj << " " << m_theInput.getNConstituents());
737
738 resolInit.push_back((m_theInput.getConstituentCovariance(iobj))(4,4));
739
740 // ATH_MSG_DEBUG( "massFitRun: 7 1 ");
741
742 AmgMatrix(5,5) tmpMtx;
743 tmpMtx.setZero();
744 output.getConstituentCovarianced0z0PhiThetaP(iobj, tmpMtx);
745 resolFinal.push_back(tmpMtx(4,4));
746 // resolFinal.push_back((output.getConstituentCovarianced0z0PhiThetaP(iobj))(4,4));
747
748 // ATH_MSG_DEBUG( "massFitRun: 7 2 ");
749
750 sumResolInit += resolInit[iobj];
751 //std::cout << " ### "<< iobj <<"\t"<< resolInit[iobj] <<"\t"<<resolFinal[iobj]<<std::endl;
752 }
753
754
755 ATH_MSG_VERBOSE( "massFitRun: 8 \n" << covOut);
756
757 for(unsigned int iobj=0;iobj<m_nobj;iobj++) {
758 resolFinalnew.push_back( resolFinal[iobj] + constraintWidthSquare *
759 resolInit[iobj]*resolInit[iobj]
760 /sumResolInit/sumResolInit );
761 covOut(4+5*iobj,4+5*iobj) = covOut(4+5*iobj,4+5*iobj)*resolFinalnew[iobj]/resolFinal[iobj];
762 covOut(3+5*iobj,3+5*iobj) = covOut(3+5*iobj,3+5*iobj)*resolFinalnew[iobj]/resolFinal[iobj];
763 covOut(2+5*iobj,2+5*iobj) = covOut(2+5*iobj,2+5*iobj)*resolFinalnew[iobj]/resolFinal[iobj];
764 //std::cout << " ### "<< iobj <<"\t"<< resolInit[iobj] <<"\t"<<resolFinal[iobj]<<"\t"<< resolFinalnew[iobj]<<std::endl;
765 }
766 ATH_MSG_VERBOSE( "massFitRun: 8.1 \n" << covOut);
767 }
768
769 ATH_MSG_VERBOSE( "massFitRun: 9 ");
770
771 // ATH_MSG_DEBUG( "massFitRun: cov 3 ");
772
773 // convert covariance matrix from cartesian to d0z0PhiThetaP
774 Amg::MatrixX covard0z0PhiThetaP;
775 convertCovXYZTod0z0PhiThetaP(particleList, covOut, covard0z0PhiThetaP);
776 // save output
777 output.setFitOutput(particleList, covOut, covard0z0PhiThetaP);
778
779 ATH_MSG_VERBOSE( "massFitRun: cov 3 \n" << covOut);
780
781 }
782
783 ATH_MSG_VERBOSE( "massFitRun: end ");
784 }
ZMassConstraint::ConstraintFit::convertCovXYZTod0z0PhiThetaP
void convertCovXYZTod0z0PhiThetaP(const std::vector< TLorentzVector > &particleList, const Amg::MatrixX &covarXYZ, Amg::MatrixX &covard0z0PhiThetaP) const
Definition ConstraintFit.cxx:1159
ZMassConstraint::ConstraintFit::massFit
double massFit(const Amg::MatrixX &, const Amg::MatrixX &var, double mass, Amg::MatrixX &pOut, Amg::MatrixX &)
Definition ConstraintFit.cxx:857
Trk::ParticleSwitcher::particle
constexpr ParticleHypothesis particle[PARTICLEHYPOTHESES]
the array of masses
Definition ParticleHypothesis.h:79
merge.output
output
Definition merge.py:16
test_AnalysisBaseEventLoopJob.aa
aa
Definition test_AnalysisBaseEventLoopJob.py:37

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msg_level_name()

const std::string & asg::AsgTool::msg_level_name ( ) const
inherited

A deprecated function for getting the message level's name.

Instead of using this, weirdly named function, user code should get the string name of the current minimum message level (in case they really need it...), with:

MSG::name( msg().level() )

This function's name doesn't follow the ATLAS coding rules, and as such will be removed in the not too distant future.

Returns
The string name of the current minimum message level that's printed

Definition at line 101 of file AsgTool.cxx.

101 {
102
103 return MSG::name( msg().level() );
104 }
AthCommonMsg< AlgTool >::msg
MsgStream & msg() const
Definition AthCommonMsg.h:24
MSG::name
const std::string & name(Level lvl)
Convenience function for translating message levels to strings.
Definition MsgLevel.cxx:19

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ print()

void asg::AsgTool::print ( ) const
virtualinherited

Print the state of the tool.

Implements asg::IAsgTool.

Reimplemented in AsgHelloTool, HI::HIPileupTool, JetBottomUpSoftDrop, JetConstituentsRetriever, JetDumper, JetFinder, JetFromPseudojet, JetModifiedMassDrop, JetPileupLabelingTool, JetPruner, JetPseudojetRetriever, JetReclusterer, JetReclusteringTool, JetRecTool, JetRecursiveSoftDrop, JetSoftDrop, JetSplitter, JetSubStructureMomentToolsBase, JetToolRunner, JetTrimmer, JetTruthLabelingTool, and KtDeltaRTool.

Definition at line 131 of file AsgTool.cxx.

131 {
132
133 ATH_MSG_INFO( "AsgTool " << name() << " @ " << this );
134 return;
135 }

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
383 PBASE::renounce (h);
384 }
AthCommonDataStore::renounce
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)
Definition AthCommonDataStore.h:380

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
365 handlesArray.renounce();
366 }

◆ retrieve_eta_calo()

float ZMassConstraint::ConstraintFit::retrieve_eta_calo ( const xAOD::IParticle & part) const
private

Definition at line 1190 of file ConstraintFit.cxx.

1191 {
1192 double eta_calo = 0;
1193 const xAOD::Egamma* eg = dynamic_cast<const xAOD::Egamma*> (&part);
1194 if (!eg) {
1195 ATH_MSG_ERROR("retrieve_eta_calo - unable to cast to Egamma");
1196 return eta_calo;
1197 }
1198 const xAOD::CaloCluster& cluster = *eg->caloCluster();
1199 static const SG::ConstAccessor<float> etaCaloAcc("etaCalo");
1200 if (cluster.retrieveMoment(xAOD::CaloCluster::ETACALOFRAME,
1201 eta_calo)) { }
1202 else if (etaCaloAcc.isAvailable(cluster)) {
1203 eta_calo = etaCaloAcc(cluster);
1204 }
1205 else {
1206 ATH_MSG_ERROR("retrieve_eta_calo - etaCalo not available as auxilliary variable");
1207 ATH_MSG_WARNING("retrieve_eta_calo - using eta as etaCalo");
1208 eta_calo = cluster.eta();
1209 }
1210 return eta_calo;
1211 }
xAOD::CaloCluster_v1::retrieveMoment
bool retrieveMoment(MomentType type, double &value) const
Retrieve individual moment.
Definition CaloCluster_v1.cxx:662
xAOD::CaloCluster_v1::eta
virtual double eta() const
The pseudorapidity ( ) of the particle.
Definition CaloCluster_v1.cxx:251
xAOD::CaloCluster_v1::ETACALOFRAME
@ ETACALOFRAME
Eta in the calo frame (for egamma).
Definition CaloCluster_v1.h:190
xAOD::Egamma_v1::caloCluster
const xAOD::CaloCluster * caloCluster(size_t index=0) const
Pointer to the xAOD::CaloCluster/s that define the electron candidate.
Definition Egamma_v1.cxx:377
xAOD::CaloCluster
CaloCluster_v1 CaloCluster
Define the latest version of the calorimeter cluster class.
Definition Event/xAOD/xAODCaloEvent/xAODCaloEvent/CaloCluster.h:19
xAOD::Egamma
Egamma_v1 Egamma
Definition of the current "egamma version".
Definition Egamma.h:17

◆ setIgnore()

void ZMassConstraint::ConstraintFit::setIgnore ( bool val)
inlineprivate

Definition at line 100 of file ConstraintFit.h.

100{ m_ignoreInputChecks=val; }

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )
overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in asg::AsgMetadataTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and DerivationFramework::CfAthAlgTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
312 std::vector<SG::VarHandleKey*> keys = a->keys();
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }
AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
std::vector< SG::VarHandleKeyArray * > m_vhka
Definition AthCommonDataStore.h:398

Member Data Documentation

◆ m_conHasWidth

bool ZMassConstraint::ConstraintFit::m_conHasWidth
private

Definition at line 109 of file ConstraintFit.h.

◆ m_conMass

double ZMassConstraint::ConstraintFit::m_conMass
private

Definition at line 110 of file ConstraintFit.h.

◆ m_conWidth

double ZMassConstraint::ConstraintFit::m_conWidth
private

Definition at line 111 of file ConstraintFit.h.

◆ m_covarianceFinal

Amg::MatrixX ZMassConstraint::ConstraintFit::m_covarianceFinal
private

Definition at line 126 of file ConstraintFit.h.

◆ m_covarianceInit

Amg::MatrixX ZMassConstraint::ConstraintFit::m_covarianceInit
private

Definition at line 124 of file ConstraintFit.h.

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default).

Definition at line 393 of file AthCommonDataStore.h.

◆ m_energyRescaler

ToolHandle<CP::IEgammaCalibrationAndSmearingTool> ZMassConstraint::ConstraintFit::m_energyRescaler
private

Definition at line 116 of file ConstraintFit.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default).

Definition at line 390 of file AthCommonDataStore.h.

◆ m_ignoreInputChecks

bool ZMassConstraint::ConstraintFit::m_ignoreInputChecks
private

Definition at line 113 of file ConstraintFit.h.

◆ m_mu_resolSFTool

ToolHandle<CP::IMuonCalibrationAndSmearingTool> ZMassConstraint::ConstraintFit::m_mu_resolSFTool
private

Definition at line 117 of file ConstraintFit.h.

◆ m_nobj

unsigned int ZMassConstraint::ConstraintFit::m_nobj
private

Definition at line 121 of file ConstraintFit.h.

◆ m_objmass

std::vector<double> ZMassConstraint::ConstraintFit::m_objmass
private

Definition at line 122 of file ConstraintFit.h.

◆ m_parameters

unsigned int ZMassConstraint::ConstraintFit::m_parameters
private

Definition at line 120 of file ConstraintFit.h.

◆ m_parametersFinal

Amg::MatrixX ZMassConstraint::ConstraintFit::m_parametersFinal
private

Definition at line 125 of file ConstraintFit.h.

◆ m_parametersInit

Amg::MatrixX ZMassConstraint::ConstraintFit::m_parametersInit
private

Definition at line 123 of file ConstraintFit.h.

◆ m_resolution

double ZMassConstraint::ConstraintFit::m_resolution
private

Definition at line 112 of file ConstraintFit.h.

◆ m_theInput

ConstraintFitInput ZMassConstraint::ConstraintFit::m_theInput
private

Definition at line 119 of file ConstraintFit.h.

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

The documentation for this class was generated from the following files:
Generated on for ATLAS Offline Software by doxygen 1.16.1