KinFitter/src/TFitConstraintEp.cc

0001 // Classname: TFitConstraintEp
0002 // Author: Jan E. Sundermann, Verena Klose (TU Dresden)
0003
0004 //________________________________________________________________
0005 //
0006 // TFitConstraintEp::
0007 // --------------------
0008 //
0009 // Fit constraint: energy and momentum conservation
0010 //
0011
0012 #include "PhysicsTools/KinFitter/interface/TFitConstraintEp.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include <iostream>
0015 #include <iomanip>
0016 #include "TClass.h"
0017
0018 //----------------
0019 // Constructor --
0020 //----------------
0021
0022 TFitConstraintEp::TFitConstraintEp()
0023     : TAbsFitConstraint(), _particles(0), _constraint(0.), _component(TFitConstraintEp::pX) {}
0024
0025 TFitConstraintEp::TFitConstraintEp(const TString& name,
0026                                    const TString& title,
0027                                    TFitConstraintEp::component thecomponent,
0028                                    Double_t constraint)
0029     : TAbsFitConstraint(name, title), _constraint(constraint), _component(thecomponent) {}
0030
0031 TFitConstraintEp::TFitConstraintEp(std::vector<TAbsFitParticle*>* particles,
0032                                    TFitConstraintEp::component thecomponent,
0033                                    Double_t constraint)
0034     : TAbsFitConstraint(), _particles(0), _constraint(constraint), _component(thecomponent) {
0035   // particles: vector containing pointer to TAbsFitParticle objects.
0036   //            Energy or momentum conservation will be calculated for
0037   //            those particles.
0038   // thecomponent: conserved 4vector component ( pX, pY, pZ, E ). For
0039   //            full 4vector conservation four objects of type TFitConstraintEp
0040   //            are needed (four constraints)
0041   // constraint:value of energy or momentum constraint ( e.g. sum[ pX_i ] = constraint )
0042
0043   if (particles) {
0044     _particles = (*particles);
0045   }
0046 }
0047
0048 TFitConstraintEp::TFitConstraintEp(const TString& name,
0049                                    const TString& title,
0050                                    std::vector<TAbsFitParticle*>* particles,
0051                                    TFitConstraintEp::component thecomponent,
0052                                    Double_t constraint)
0053     : TAbsFitConstraint(name, title), _particles(0), _constraint(constraint), _component(thecomponent) {
0054   // particles: vector containing pointer to TAbsFitParticle objects.
0055   //            Energy or momentum conservation will be calculated for
0056   //            those particles.
0057   // thecomponent: conserved 4vector component ( pX, pY, pZ, E ). For
0058   //            full 4vector conservation four objects of type TFitConstraintEp
0059   //            are needed (four constraints)
0060   // constraint:value of energy or momentum constraint ( e.g. sum[ pX_i ] = constraint )
0061
0062   if (particles) {
0063     _particles = (*particles);
0064   }
0065 }
0066
0067 //--------------
0068 // Destructor --
0069 //--------------
0070 TFitConstraintEp::~TFitConstraintEp() {}
0071
0072 void TFitConstraintEp::addParticle(TAbsFitParticle* particle) {
0073   // Add one particles to list of constrained particles
0074
0075   _particles.push_back(particle);
0076 }
0077
0078 void TFitConstraintEp::addParticles(TAbsFitParticle* p1,
0079                                     TAbsFitParticle* p2,
0080                                     TAbsFitParticle* p3,
0081                                     TAbsFitParticle* p4,
0082                                     TAbsFitParticle* p5,
0083                                     TAbsFitParticle* p6,
0084                                     TAbsFitParticle* p7,
0085                                     TAbsFitParticle* p8,
0086                                     TAbsFitParticle* p9,
0087                                     TAbsFitParticle* p10) {
0088   // Add many particles to list of constrained particles
0089
0090   if (p1)
0091     addParticle(p1);
0092   if (p2)
0093     addParticle(p2);
0094   if (p3)
0095     addParticle(p3);
0096   if (p4)
0097     addParticle(p4);
0098   if (p5)
0099     addParticle(p5);
0100   if (p6)
0101     addParticle(p6);
0102   if (p7)
0103     addParticle(p7);
0104   if (p8)
0105     addParticle(p8);
0106   if (p9)
0107     addParticle(p9);
0108   if (p10)
0109     addParticle(p10);
0110 }
0111
0112 //--------------
0113 // Operations --
0114 //--------------
0115 TMatrixD* TFitConstraintEp::getDerivative(TAbsFitParticle* particle) {
0116   // returns derivative df/dP with P=(p,E) and f the constraint (f=0).
0117   // The matrix contains one row (df/dp, df/dE).
0118
0119   TMatrixD* DerivativeMatrix = new TMatrixD(1, 4);
0120   (*DerivativeMatrix) *= 0.;
0121   (*DerivativeMatrix)(0, (int)_component) = 1.;
0122   return DerivativeMatrix;
0123 }
0124
0125 Double_t TFitConstraintEp::getInitValue() {
0126   // Get initial value of constraint (before the fit)
0127
0128   Double_t InitValue(0);
0129   UInt_t Npart = _particles.size();
0130   for (unsigned int i = 0; i < Npart; i++) {
0131     const TLorentzVector* FourVec = _particles[i]->getIni4Vec();
0132     InitValue += (*FourVec)[(int)_component];
0133   }
0134   InitValue -= _constraint;
0135   return InitValue;
0136 }
0137
0138 Double_t TFitConstraintEp::getCurrentValue() {
0139   // Get value of constraint after the fit
0140
0141   Double_t CurrentValue(0);
0142   UInt_t Npart = _particles.size();
0143   for (unsigned int i = 0; i < Npart; i++) {
0144     const TLorentzVector* FourVec = _particles[i]->getCurr4Vec();
0145     CurrentValue += (*FourVec)[(int)_component];
0146   }
0147   CurrentValue -= _constraint;
0148   return CurrentValue;
0149 }
0150
0151 TString TFitConstraintEp::getInfoString() {
0152   // Collect information to be used for printout
0153
0154   std::stringstream info;
0155   info << std::scientific << std::setprecision(6);
0156
0157   info << "__________________________" << std::endl << std::endl;
0158   info << "OBJ: " << IsA()->GetName() << "\t" << GetName() << "\t" << GetTitle() << std::endl;
0159
0160   info << "initial value: " << getInitValue() << std::endl;
0161   info << "current value: " << getCurrentValue() << std::endl;
0162   info << "component: " << _component << std::endl;
0163   info << "constraint: " << _constraint << std::endl;
0164
0165   return info.str();
0166 }
0167
0168 void TFitConstraintEp::print() {
0169   // Print constraint contents
0170
0171   edm::LogVerbatim("KinFitter") << this->getInfoString();
0172 }