KinFitter/src/TFitParticleSpher.cc

0001 // Classname: TFitParticleSpher
0002 // Author: Jan E. Sundermann, Verena Klose (TU Dresden)
0003
0004 //________________________________________________________________
0005 //
0006 // TFitParticleSpher::
0007 // --------------------
0008 //
0009 // Particle with spherical  parametrization of the momentum 4vector and
0010 // free mass (4 free parameters). The parametrization is chosen as
0011 // follows:
0012 //
0013 // p = (r, theta, phi)
0014 // E(fit) =  Sqrt( |p|^2 + d^2*m^2 )
0015 //
0016
0017 #include <iostream>
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 #include "PhysicsTools/KinFitter/interface/TFitParticleSpher.h"
0020 #include "TMath.h"
0021
0022 //----------------
0023 // Constructor --
0024 //----------------
0025 TFitParticleSpher::TFitParticleSpher() : TAbsFitParticle() { init(nullptr, nullptr); }
0026
0027 TFitParticleSpher::TFitParticleSpher(const TFitParticleSpher& 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 TFitParticleSpher::TFitParticleSpher(TLorentzVector* pini, const TMatrixD* theCovMatrix) : TAbsFitParticle() {
0044   init(pini, theCovMatrix);
0045 }
0046
0047 TFitParticleSpher::TFitParticleSpher(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* TFitParticleSpher::clone(const TString& newname) const {
0056   // Returns a copy of itself
0057
0058   TAbsFitParticle* myclone = new TFitParticleSpher(*this);
0059   if (newname.Length() > 0)
0060     myclone->SetName(newname);
0061   return myclone;
0062 }
0063
0064 //--------------
0065 // Destructor --
0066 //--------------
0067 TFitParticleSpher::~TFitParticleSpher() {}
0068
0069 //--------------
0070 // Operations --
0071 //--------------
0072 void TFitParticleSpher::init(TLorentzVector* pini, const TMatrixD* theCovMatrix) {
0073   _nPar = 4;
0074   setIni4Vec(pini);
0075   setCovMatrix(theCovMatrix);
0076 }
0077
0078 TLorentzVector* TFitParticleSpher::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 r = (*params)(0, 0);
0092   Double_t theta = (*params)(1, 0);
0093   Double_t phi = (*params)(2, 0);
0094   Double_t d = (*params)(3, 0);
0095
0096   Double_t X = r * TMath::Cos(phi) * TMath::Sin(theta);
0097   Double_t Y = r * TMath::Sin(phi) * TMath::Sin(theta);
0098   Double_t Z = r * TMath::Cos(theta);
0099   Double_t E = TMath::Sqrt(X * X + Y * Y + Z * Z + d * d * _pini.M2());
0100
0101   TLorentzVector* vec = new TLorentzVector(X, Y, Z, E);
0102   return vec;
0103 }
0104
0105 void TFitParticleSpher::setIni4Vec(const TLorentzVector* pini) {
0106   // Set the initial 4vector. Will also set the
0107   // inital parameter values
0108
0109   if (pini == nullptr) {
0110     _u1.SetXYZ(0., 0., 0.);
0111     _u3.SetXYZ(0., 0., 0.);
0112     _u2.SetXYZ(0., 0., 0.);
0113     _pini.SetXYZT(0., 0., 0., 0.);
0114     _pcurr = _pini;
0115
0116     _iniparameters.ResizeTo(_nPar, 1);
0117     _iniparameters(0, 0) = 0.;
0118     _iniparameters(1, 0) = 0.;
0119     _iniparameters(2, 0) = 0.;
0120     _iniparameters(3, 0) = 1.;
0121
0122     _parameters.ResizeTo(_nPar, 1);
0123     _parameters(0, 0) = 0.;
0124     _parameters(1, 0) = 0.;
0125     _parameters(2, 0) = 0.;
0126     _parameters(3, 0) = 1.;
0127
0128   } else {
0129     Double_t r = pini->P();
0130     Double_t theta = pini->Theta();
0131     Double_t phi = pini->Phi();
0132
0133     _pini = (*pini);
0134     _pcurr = _pini;
0135
0136     _iniparameters.ResizeTo(_nPar, 1);
0137     _iniparameters(0, 0) = r;
0138     _iniparameters(1, 0) = theta;
0139     _iniparameters(2, 0) = phi;
0140     _iniparameters(3, 0) = 1.;
0141
0142     _parameters.ResizeTo(_nPar, 1);
0143     _parameters = _iniparameters;
0144
0145     _u1.SetXYZ(TMath::Cos(phi) * TMath::Sin(theta), TMath::Sin(phi) * TMath::Sin(theta), TMath::Cos(theta));
0146     _u2.SetXYZ(TMath::Cos(phi) * TMath::Cos(theta), TMath::Sin(phi) * TMath::Cos(theta), -1. * TMath::Sin(theta));
0147     _u3.SetXYZ(-1. * TMath::Sin(phi), TMath::Cos(phi), 0.);
0148   }
0149 }
0150
0151 TMatrixD* TFitParticleSpher::getDerivative() {
0152   // returns derivative dP/dy with P=(p,E) and y=(r, theta, phi, d)
0153   // the free parameters of the fit. The columns of the matrix contain
0154   // (dP/dr, dP/dtheta, ...).
0155
0156   TMatrixD* DerivativeMatrix = new TMatrixD(4, 4);
0157   (*DerivativeMatrix) *= 0.;
0158
0159   Double_t r = _parameters(0, 0);
0160   Double_t theta = _parameters(1, 0);
0161   Double_t phi = _parameters(2, 0);
0162   Double_t d = _parameters(3, 0);
0163
0164   //1st column: dP/dr
0165   (*DerivativeMatrix)(0, 0) = TMath::Cos(phi) * TMath::Sin(theta);
0166   (*DerivativeMatrix)(1, 0) = TMath::Sin(phi) * TMath::Sin(theta);
0167   (*DerivativeMatrix)(2, 0) = TMath::Cos(theta);
0168   (*DerivativeMatrix)(3, 0) = 0.;
0169
0170   //2nd column: dP/dtheta
0171   (*DerivativeMatrix)(0, 1) = r * TMath::Cos(phi) * TMath::Cos(theta);
0172   (*DerivativeMatrix)(1, 1) = r * TMath::Sin(phi) * TMath::Cos(theta);
0173   (*DerivativeMatrix)(2, 1) = -1. * r * TMath::Sin(theta);
0174   (*DerivativeMatrix)(3, 1) = 0.;
0175
0176   //3rd column: dP/dphi
0177   (*DerivativeMatrix)(0, 2) = -1. * r * TMath::Sin(phi) * TMath::Sin(theta);
0178   (*DerivativeMatrix)(1, 2) = r * TMath::Cos(phi) * TMath::Sin(theta);
0179   ;
0180   (*DerivativeMatrix)(2, 2) = 0.;
0181   (*DerivativeMatrix)(3, 2) = 0.;
0182
0183   //4th column: dP/dm
0184   (*DerivativeMatrix)(0, 3) = 0.;
0185   (*DerivativeMatrix)(1, 3) = 0.;
0186   (*DerivativeMatrix)(2, 3) = 0.;
0187   (*DerivativeMatrix)(3, 3) = _pini.M() * _pini.M() * d / _pcurr.E();
0188
0189   return DerivativeMatrix;
0190 }
0191
0192 TMatrixD* TFitParticleSpher::transform(const TLorentzVector& vec) {
0193   // Returns the parameters corresponding to the given
0194   // 4vector
0195
0196   // retrieve parameters
0197   TMatrixD* tparams = new TMatrixD(_nPar, 1);
0198   (*tparams)(0, 0) = vec.P();
0199   (*tparams)(1, 0) = vec.Theta();
0200   (*tparams)(2, 0) = vec.Phi();
0201   (*tparams)(3, 0) = vec.M() / _pini.M();
0202
0203   return tparams;
0204 }