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 {
22 declareProperty( "Z_pdg_mass", m_conMass);
23 declareProperty( "Z_pdg_width", m_conWidth);
24 declareProperty( "hasWidth", m_conHasWidth);
25 declareProperty( "ignoreInputChecks", m_ignoreInputChecks);
26 declareProperty( "EgammaCalibAndSmearingTool", m_energyRescaler = ToolHandle<CP::IEgammaCalibrationAndSmearingTool>("CP::EgammaCalibrationAndSmearingTool"));
27 declareProperty( "MuonCalibrationAndSmearingTool", m_mu_resolSFTool = ToolHandle<CP::IMuonCalibrationAndSmearingTool>("CP::MuonCalibrationAndSmearingTool"));
28 }
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 30 of file ConstraintFit.cxx.

31 {}

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 83 of file ConstraintFit.cxx.

86 {
87
88 ATH_MSG_DEBUG ("addFSRParticle: *** 4vecFsr: " << fsr4Vec.Pt() << "/" << fsr4Vec.Eta() << "/" << fsr4Vec.Phi());
89 // Get energy resolution - done differently for photons and electrons
90 float e_res = 0;
91 PATCore::ParticleType::Type phType = PATCore::ParticleType::Electron;
92 float cl_etaCalo = retrieve_eta_calo(part);
93 if (!m_energyRescaler.empty()) {
94 if (xAOD::Type::Photon == part.type()) {
95 const xAOD::Photon* fsrPf = dynamic_cast<const xAOD::Photon*>(&part);
96 phType = (xAOD::EgammaHelpers::isConvertedPhoton(fsrPf)) ?
97 PATCore::ParticleType::ConvertedPhoton :
98 PATCore::ParticleType::UnconvertedPhoton;
99 }
100 e_res = m_energyRescaler->resolution( fsr4Vec.E(), fsr4Vec.Eta(), cl_etaCalo, phType)*fsr4Vec.E();
101 }
102 else {
103 ATH_MSG_ERROR( "addFSRParticle: - cannot find EgammaCalibrationAndSmearingTool. ");
104 ATH_MSG_ERROR( "addFSRParticle: - please set property 'EgammaCalibAndSmearingToolName' with the name of your EgammaCalibrationAndSmearingTool - default is 'EgammaCalibrationAndSmearingTool'. ");
105 }
106
107 AmgMatrix(5,5) photonCovarianceMatrix;
108 photonCovarianceMatrix.setZero();
109 photonCovarianceMatrix(d0,d0) = 0.000001;
110 photonCovarianceMatrix(z0,z0) = 0.000001;
111 photonCovarianceMatrix(phi0,phi0) = 0.000001;
112 photonCovarianceMatrix(theta,theta) = 0.000001;
113 photonCovarianceMatrix(qOverP,qOverP) = e_res*e_res;
114 photonCovarianceMatrix(d0,z0) = 0;
115 photonCovarianceMatrix(d0,phi0) = 0;
116 photonCovarianceMatrix(d0,theta) = 0;
117 photonCovarianceMatrix(d0,qOverP) = 0;
118 photonCovarianceMatrix(z0,phi0) = 0;
119 photonCovarianceMatrix(z0,theta) = 0;
120 photonCovarianceMatrix(z0,qOverP) = 0;
121 photonCovarianceMatrix(phi0,theta) = 0;
122 photonCovarianceMatrix(phi0,qOverP) = 0;
123 photonCovarianceMatrix(theta,qOverP) = 0;
124
125 ATH_MSG_DEBUG ("addFSRParticle: FSR calib type,e, e_res " << phType << " " << fsr4Vec.E() << " " << e_res);
126
127 // symmetrize it
128 for (int ii = 0; ii < 5; ii++)
129 for (int jj = ii + 1; jj < 5; jj++)
130 photonCovarianceMatrix(jj, ii) = photonCovarianceMatrix(ii, jj);
131
132 bool isOK = doSanityChecksOnCovariance(fsr4Vec, photonCovarianceMatrix);
133 // Get covariance in cartesian CS
134 AmgMatrix(5,5) covarXYZ;
135 covarXYZ.setZero();
136 convertCovd0z0PhiThetaPToXYZ(fsr4Vec, photonCovarianceMatrix, covarXYZ);
137 // Save in input
138 input.addConstituent_FourVector_d0z0PhiThetaP(fsr4Vec, photonCovarianceMatrix, covarXYZ, isOK);
139
140 ATH_MSG_DEBUG( "addFSRParticle: - added fsr");
141
142 }
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:1127
ZMassConstraint::ConstraintFit::retrieve_eta_calo
float retrieve_eta_calo(const xAOD::IParticle &part) const
Definition ConstraintFit.cxx:1189
ZMassConstraint::ConstraintFit::doSanityChecksOnCovariance
bool doSanityChecksOnCovariance(const TLorentzVector &vector, const AmgMatrix(5, 5)&covar) const
Definition ConstraintFit.cxx:1090
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 310 of file ConstraintFit.cxx.

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

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

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

1161 {
1162 unsigned int matrixSize = 5 * particleList.size();
1163 Amg::MatrixX Jacobian(matrixSize, matrixSize);
1164 Jacobian.setZero();
1166 for(unsigned int ii = 0; ii < particleList.size(); ii++) {
1167 double phi = particleList.at(ii).Phi();
1168 double theta = particleList.at(ii).Theta();
1169 double P = particleList.at(ii).P();
1170 Jacobian( 5*ii, 5*ii) = 1.;
1171 Jacobian(1 + 5*ii, 1 + 5*ii) = 1.;
1172
1173 Jacobian(2 + 5*ii, 2 + 5*ii) = -P * TMath::Sin(theta) * TMath::Sin(phi);
1174 Jacobian(2 + 5*ii, 3 + 5*ii) = P * TMath::Sin(theta) * TMath::Cos(phi);
1175 Jacobian(3 + 5*ii, 2 + 5*ii) = P * TMath::Cos(theta) * TMath::Cos(phi);
1176 Jacobian(3 + 5*ii, 3 + 5*ii) = P * TMath::Cos(theta) * TMath::Sin(phi);
1177 Jacobian(3 + 5*ii, 4 + 5*ii) = -P * TMath::Sin(theta);
1178 Jacobian(4 + 5*ii, 2 + 5*ii) = TMath::Sin(theta) * TMath::Cos(phi);
1179 Jacobian(4 + 5*ii, 3 + 5*ii) = TMath::Sin(theta) * TMath::Sin(phi);
1180 Jacobian(4 + 5*ii, 4 + 5*ii) = TMath::Cos(theta);
1181 }
1182
1183 Amg::MatrixX Jacobianinverse(matrixSize, matrixSize);
1184 Jacobianinverse=Jacobian.inverse();
1185 covard0z0PhiThetaP = Jacobianinverse.transpose() * covarXYZ * Jacobianinverse;
1186 }
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 54 of file ConstraintFit.cxx.

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

◆ doSanityChecksOnCovariance()

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

Definition at line 1090 of file ConstraintFit.cxx.

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

◆ 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::TEvent 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::TEvent, 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:407

◆ 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 905 of file ConstraintFit.cxx.

907 {
908 unsigned int nPartFirst = firstInput.getNConstituents();
909 unsigned int nPartSecond = secondInput.getNConstituents();
910 unsigned int nPart = nPartFirst + nPartSecond;
911 std::vector<TLorentzVector> particles(nPart);
912 std::vector<AmgMatrix(3,3)> covariances(nPart);
913 for ( unsigned int iobj = 0; iobj < firstInput.getNConstituents(); ++iobj) {
914 covariances[iobj].setZero();
915 particles[iobj] = firstInput.getConstituentFourVector(iobj);
916 firstInput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj]);
917 }
918 for ( unsigned int iobj = 0; iobj < secondInput.getNConstituents(); ++iobj) {
919 covariances[iobj + nPartFirst].setZero();
920 particles[iobj + nPartFirst] = secondInput.getConstituentFourVector(iobj);
921 secondInput.getConstituentCovariancePhiThetaP(iobj, covariances[iobj + nPartFirst]);
922 }
923 return getMassError( particles, covariances );
924 }
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:905

◆ 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 929 of file ConstraintFit.cxx.

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

◆ 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 954 of file ConstraintFit.cxx.

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

◆ getMassError() [4/4]

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

Definition at line 976 of file ConstraintFit.cxx.

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

◆ 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::TEvent 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::TEvent, 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 34 of file ConstraintFit.cxx.

35 {
36 ATH_MSG_INFO ("initialize: retrieve tools" << m_energyRescaler.name() << " and "
37 << m_mu_resolSFTool.name());
38
39 if (!m_energyRescaler.retrieve().isSuccess()) {
40 ATH_MSG_ERROR ("initialize: unable to retrieve EgammaCalibrationAndSmearingTool");
41 return StatusCode::FAILURE;
42 }
43
44 if (!m_mu_resolSFTool.retrieve().isSuccess()) {
45 ATH_MSG_ERROR ("initialize: unable to retrieve MuonCalibrationAndSmearingTool");
46 return StatusCode::FAILURE;
47 }
48
49 // Return gracefully:
50 return StatusCode::SUCCESS;
51 }
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 434 of file ConstraintFit.cxx.

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

◆ likelihoodMass2()

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

Definition at line 502 of file ConstraintFit.cxx.

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

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

◆ massFitCalculation()

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

Definition at line 786 of file ConstraintFit.cxx.

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

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

◆ massFitRun()

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

Definition at line 584 of file ConstraintFit.cxx.

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

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 1189 of file ConstraintFit.cxx.

1190 {
1191 double eta_calo = 0;
1192 const xAOD::Egamma* eg = dynamic_cast<const xAOD::Egamma*> (&part);
1193 if (!eg) {
1194 ATH_MSG_ERROR("retrieve_eta_calo - unable to cast to Egamma");
1195 return eta_calo;
1196 }
1197 const xAOD::CaloCluster& cluster = *eg->caloCluster();
1198 static const SG::ConstAccessor<float> etaCaloAcc("etaCalo");
1199 if (cluster.retrieveMoment(xAOD::CaloCluster::ETACALOFRAME,
1200 eta_calo)) { }
1201 else if (etaCaloAcc.isAvailable(cluster)) {
1202 eta_calo = etaCaloAcc(cluster);
1203 }
1204 else {
1205 ATH_MSG_ERROR("retrieve_eta_calo - etaCalo not available as auxilliary variable");
1206 ATH_MSG_WARNING("retrieve_eta_calo - using eta as etaCalo");
1207 eta_calo = cluster.eta();
1208 }
1209 return eta_calo;
1210 }
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::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.14.0