KinFitter/src/TFitParticleCart.cc

0001 // Classname: TFitParticleCart
0002 // Author: Jan E. Sundermann, Verena Klose (TU Dresden)
0003
0004 //________________________________________________________________
0005 //
0006 // TFitParticleCart::
0007 // --------------------
0008 //
0009 // Particle with cartesian 4vector parametrization and free mass
0010 // [four free parameters (px, py, pz, d) with initial values
0011 //  (px, py, pz, 1)]
0012 //
0013 // p = px*u1 + py*u2 + pz*u3
0014 // E = Sqrt( |p|^2 + d^2*m^2 )
0015 //
0016
0017 #include <iostream>
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 #include "PhysicsTools/KinFitter/interface/TFitParticleCart.h"
0020 #include "TMath.h"
0021
0022 //----------------
0023 // Constructor --
0024 //----------------
0025 TFitParticleCart::TFitParticleCart() : TAbsFitParticle() { init(nullptr, nullptr); }
0026
0027 TFitParticleCart::TFitParticleCart(const TFitParticleCart& fitParticle)
0028     : TAbsFitParticle(fitParticle.GetName(), fitParticle.GetTitle()) {
0029   _nPar = fitParticle._nPar;
0030   _u1 = fitParticle._u1;
0031   _u2 = fitParticle._u2;
0032   _u3 = fitParticle._u3;
0033   _covMatrix.ResizeTo(fitParticle._covMatrix);
0034   _covMatrix = fitParticle._covMatrix;
0035   _iniparameters.ResizeTo(fitParticle._iniparameters);
0036   _iniparameters = fitParticle._iniparameters;
0037   _parameters.ResizeTo(fitParticle._parameters);
0038   _parameters = fitParticle._parameters;
0039   _pini = fitParticle._pini;
0040   _pcurr = fitParticle._pcurr;
0041 }
0042
0043 TFitParticleCart::TFitParticleCart(TLorentzVector* pini, const TMatrixD* theCovMatrix) : TAbsFitParticle() {
0044   init(pini, theCovMatrix);
0045 }
0046
0047 TFitParticleCart::TFitParticleCart(const TString& name,
0048                                    const TString& title,
0049                                    TLorentzVector* pini,
0050                                    const TMatrixD* theCovMatrix)
0051     : TAbsFitParticle(name, title) {
0052   init(pini, theCovMatrix);
0053 }
0054
0055 TAbsFitParticle* TFitParticleCart::clone(const TString& newname) const {
0056   // Returns a copy of itself
0057
0058   TAbsFitParticle* myclone = new TFitParticleCart(*this);
0059   if (newname.Length() > 0)
0060     myclone->SetName(newname);
0061   return myclone;
0062 }
0063
0064 //--------------
0065 // Destructor --
0066 //--------------
0067 TFitParticleCart::~TFitParticleCart() {}
0068
0069 //--------------
0070 // Operations --
0071 //--------------
0072 void TFitParticleCart::init(TLorentzVector* pini, const TMatrixD* theCovMatrix) {
0073   _nPar = 4;
0074   setIni4Vec(pini);
0075   setCovMatrix(theCovMatrix);
0076 }
0077
0078 TLorentzVector* TFitParticleCart::calc4Vec(const TMatrixD* params) {
0079   // Calculates a 4vector corresponding to the given
0080   // parameter values
0081
0082   if (params == nullptr) {
0083     return nullptr;
0084   }
0085
0086   if (params->GetNcols() != 1 || params->GetNrows() != _nPar) {
0087     edm::LogError("WrongMatrixSize") << GetName() << "::calc4Vec - Parameter matrix has wrong size.";
0088     return nullptr;
0089   }
0090
0091   Double_t X = (*params)(0, 0);
0092   Double_t Y = (*params)(1, 0);
0093   Double_t Z = (*params)(2, 0);
0094   Double_t E = TMath::Sqrt(X * X + Y * Y + Z * Z + (*params)(3, 0) * (*params)(3, 0) * _pini.M2());
0095
0096   TLorentzVector* vec = new TLorentzVector(X, Y, Z, E);
0097   return vec;
0098 }
0099
0100 void TFitParticleCart::setIni4Vec(const TLorentzVector* pini) {
0101   // Set the initial 4vector. Will also set the
0102   // inital parameter values
0103
0104   if (pini == nullptr) {
0105     _u1.SetXYZ(0., 0., 0.);
0106     _u3.SetXYZ(0., 0., 0.);
0107     _u2.SetXYZ(0., 0., 0.);
0108     _pini.SetXYZT(0., 0., 0., 0.);
0109     _pcurr = _pini;
0110
0111     _iniparameters.ResizeTo(_nPar, 1);
0112     _iniparameters(0, 0) = 0;
0113     _iniparameters(1, 0) = 0;
0114     _iniparameters(2, 0) = 0;
0115     _iniparameters(3, 0) = 0.;
0116
0117     _parameters.ResizeTo(_nPar, 1);
0118     _parameters(0, 0) = 0.;
0119     _parameters(1, 0) = 0.;
0120     _parameters(2, 0) = 0.;
0121     _parameters(3, 0) = 0.;
0122
0123   } else {
0124     _pini = (*pini);
0125     _pcurr = _pini;
0126
0127     _u1.SetXYZ(1., 0., 0.);
0128     _u2.SetXYZ(0., 1., 0.);
0129     _u3.SetXYZ(0., 0., 1.);
0130
0131     _iniparameters.ResizeTo(_nPar, 1);
0132     _iniparameters(0, 0) = _pini.X();
0133     _iniparameters(1, 0) = _pini.Y();
0134     _iniparameters(2, 0) = _pini.Z();
0135     _iniparameters(3, 0) = 1.;
0136     _parameters.ResizeTo(_nPar, 1);
0137     _parameters = _iniparameters;
0138   }
0139 }
0140
0141 TMatrixD* TFitParticleCart::getDerivative() {
0142   // returns derivative dP/dy with P=(p,E) and y=(px, py, pz, d)
0143   // the free parameters of the fit. The columns of the matrix contain
0144   // (dP/dpx, dP/dpy, ...).
0145
0146   TMatrixD* DerivativeMatrix = new TMatrixD(4, 4);
0147   (*DerivativeMatrix) *= 0.;
0148
0149   //1st column: dP/dx
0150   (*DerivativeMatrix)(0, 0) = 1.;
0151   (*DerivativeMatrix)(1, 0) = 0.;
0152   (*DerivativeMatrix)(2, 0) = 0.;
0153   (*DerivativeMatrix)(3, 0) = 0.;
0154
0155   //2nd column: dP/dy
0156   (*DerivativeMatrix)(0, 1) = 0;
0157   (*DerivativeMatrix)(1, 1) = 1;
0158   (*DerivativeMatrix)(2, 1) = 0;
0159   (*DerivativeMatrix)(3, 1) = 0.;
0160
0161   //3rd column: dP/dz
0162   (*DerivativeMatrix)(0, 2) = 0.;
0163   (*DerivativeMatrix)(1, 2) = 0.;
0164   (*DerivativeMatrix)(2, 2) = 1.;
0165   (*DerivativeMatrix)(3, 2) = 0.;
0166
0167   //4th column: dP/dm
0168   (*DerivativeMatrix)(0, 3) = 0.;
0169   (*DerivativeMatrix)(1, 3) = 0.;
0170   (*DerivativeMatrix)(2, 3) = 0.;
0171   (*DerivativeMatrix)(3, 3) = _pini.M() * _pini.M() * _parameters(3, 0) / _pcurr.E();
0172
0173   return DerivativeMatrix;
0174 }
0175
0176 TMatrixD* TFitParticleCart::transform(const TLorentzVector& vec) {
0177   // Returns the parameters corresponding to the given
0178   // 4vector
0179
0180   TMatrixD* tparams = new TMatrixD(_nPar, 1);
0181   (*tparams)(0, 0) = vec.X();
0182   (*tparams)(1, 0) = vec.Y();
0183   (*tparams)(2, 0) = vec.Z();
0184   (*tparams)(3, 0) = vec.M() / _pini.M();
0185
0186   return tparams;
0187 }