Go to the documentation of this file.00001 #include "RecoVertex/KinematicFit/interface/FourMomentumKinematicConstraint.h"
00002 #include "RecoVertex/VertexPrimitives/interface/VertexException.h"
00003
00004
00005 FourMomentumKinematicConstraint::FourMomentumKinematicConstraint(const AlgebraicVector& momentum,
00006 const AlgebraicVector& deviation)
00007 {
00008 if((momentum.num_row() != 4)||(deviation.num_row() != 4))
00009 throw VertexException("FourMomentumKinematicConstraint::vector of wrong size passed as 4-Momentum or Deviations");
00010 mm = momentum;
00011 AlgebraicVector deviation_l(7,0);
00012 deviation_l(6) = momentum(3);
00013 deviation_l(5) = momentum(2);
00014 deviation_l(4) = momentum(1);
00015
00016 double mass_sq = momentum(4)*momentum(4) - momentum(3)*momentum(3)
00017 -momentum(2)*momentum(2) - momentum(1)*momentum(1);
00018
00019 if(mass_sq == 0.)throw VertexException("FourMomentumKinematicConstraint::the constraint vector passed corresponds to 0 mass");
00020
00021
00022 deviation_l(7) = momentum(4)*momentum(4)*deviation(4)/mass_sq
00023 + momentum(3)*momentum(3)*deviation(3)/mass_sq
00024 + momentum(2)*momentum(2)*deviation(2)/mass_sq
00025 + momentum(1)*momentum(1)*deviation(1)/mass_sq;
00026
00027
00028 dd = deviation_l;
00029 }
00030
00031 std::pair<AlgebraicVector,AlgebraicVector> FourMomentumKinematicConstraint::value(const AlgebraicVector& exPoint) const
00032 {
00033 if(exPoint.num_row() ==0 ) throw VertexException("FourMomentumKinematicConstraint::value requested for zero Linearization point");
00034
00035
00036 int inSize = exPoint.num_row();
00037 if((inSize%7) !=0) throw VertexException("FourMomentumKinematicConstraint::linearization point has a wrong dimension");
00038 int nStates = inSize/7;
00039 if(nStates != 1) throw VertexException("FourMomentumKinematicConstraint::Multi particle refit is not supported in this version");
00040 AlgebraicVector pr = exPoint;
00041 AlgebraicVector vl(4,0);
00042 vl(1) += (pr(4) - mm(1));
00043 vl(2) += (pr(5) - mm(2));
00044 vl(3) += (pr(6) - mm(3));
00045 vl(7) += (sqrt(pr(4)*pr(4)+pr(5)*pr(5)+pr(6)*pr(6)+pr(7)*pr(7)) - mm(4));
00046
00047 return std::pair<AlgebraicVector,AlgebraicVector>(vl,pr);
00048 }
00049
00050 std::pair<AlgebraicMatrix, AlgebraicVector> FourMomentumKinematicConstraint::derivative(const AlgebraicVector& exPoint) const
00051 {
00052 if(exPoint.num_row() ==0) throw VertexException("FourMomentumKinematicConstraint::value requested for zero Linearization point");
00053
00054
00055 int inSize = exPoint.num_row();
00056 if((inSize%7) !=0) throw VertexException("FourMomentumKinematicConstraint::linearization point has a wrong dimension");
00057 int nStates = inSize/7;
00058 if(nStates != 1) throw VertexException("FourMomentumKinematicConstraint::Multi particle refit is not supported in this version");
00059 AlgebraicVector pr = exPoint;
00060 AlgebraicMatrix dr(4,7,0);
00061
00062 dr(1,4) = 1.;
00063 dr(2,5) = 1.;
00064 dr(3,6) = 1.;
00065 dr(4,7) = pr(7)/ sqrt(pr(4)*pr(4)+pr(5)*pr(5)+pr(6)*pr(6)+pr(7)*pr(7));
00066
00067 return std::pair<AlgebraicMatrix,AlgebraicVector>(dr,pr);
00068 }
00069
00070 std::pair<AlgebraicVector, AlgebraicVector> FourMomentumKinematicConstraint::value(const std::vector<RefCountedKinematicParticle> par) const
00071 {
00072 int nStates = par.size();
00073 if(nStates == 0) throw VertexException("FourMomentumKinematicConstraint::Empty vector of particles passed");
00074 if(nStates != 1) throw VertexException("FourMomentumKinematicConstraint::Multi particle refit is not supported in this version");
00075 AlgebraicVector pr = asHepVector<7>(par.front()->currentState().kinematicParameters().vector());
00076 AlgebraicVector vl(4,0);
00077
00078 vl(1) = pr(4) - mm(1);
00079 vl(2) = pr(5) - mm(2);
00080 vl(3) = pr(6) - mm(3);
00081 vl(7) = sqrt(pr(4)*pr(4)+pr(5)*pr(5)+pr(6)*pr(6)+pr(7)*pr(7)) - mm(4);
00082
00083 return std::pair<AlgebraicVector,AlgebraicVector>(vl,pr);
00084 }
00085
00086 std::pair<AlgebraicMatrix, AlgebraicVector> FourMomentumKinematicConstraint::derivative(const std::vector<RefCountedKinematicParticle> par) const
00087 {
00088 int nStates = par.size();
00089 if(nStates == 0) throw VertexException("FourMomentumKinematicConstraint::Empty vector of particles passed");
00090 if(nStates != 1) throw VertexException("FourMomentumKinematicConstraint::Multi particle refit is not supported in this version");
00091 AlgebraicMatrix dr(4,7,0);
00092
00093 AlgebraicVector pr = asHepVector<7>(par.front()->currentState().kinematicParameters().vector());
00094 dr(1,4) = 1.;
00095 dr(2,5) = 1.;
00096 dr(3,6) = 1.;
00097 dr(4,7) = - pr(7)/sqrt(pr(4)*pr(4)+pr(5)*pr(5)+pr(6)*pr(6)+pr(7)*pr(7));
00098
00099 return std::pair<AlgebraicMatrix,AlgebraicVector>(dr,pr);
00100 }
00101
00102 AlgebraicVector FourMomentumKinematicConstraint::deviations(int nStates) const
00103 {
00104 if(nStates == 0) throw VertexException("FourMomentumKinematicConstraint::Empty vector of particles passed");
00105 if(nStates != 1) throw VertexException("FourMomentumKinematicConstraint::Multi particle refit is not supported in this version");
00106 AlgebraicVector res = dd;
00107 return res;
00108 }
00109
00110 int FourMomentumKinematicConstraint::numberOfEquations() const
00111 {return 4;}
00112