CMS 3D CMS Logo

/data/git/CMSSW_5_3_11_patch5/src/PhysicsTools/TagAndProbe/interface/ColinsSoperVariables.h

Go to the documentation of this file.
00001 #define CM_ENERGY 7000.0
00002 #include "TLorentzVector.h" 
00003 #include "TVector3.h"
00004 
00005 
00006 // calculate the Colins-Soper variables;
00007 // everything is in the lab frame
00008 void calCSVariables(TLorentzVector mu, TLorentzVector mubar, 
00009                     double *res, bool swap) {
00010 
00011   // convention. beam direction is on the positive Z direction.
00012   // beam contains quark flux.
00013   TLorentzVector Pbeam  (0, 0,  CM_ENERGY/2.0, CM_ENERGY/2.0);
00014   TLorentzVector Ptarget(0, 0, -CM_ENERGY/2.0, CM_ENERGY/2.0);
00015   
00016 
00017   TLorentzVector Q(mu+mubar);
00018   /************************************************************************
00019    *
00020    * 1) cos(theta) = 2 Q^-1 (Q^2+Qt^2)^-1 (mu^+ mubar^- - mu^- mubar^+)
00021    *
00022    *
00023    ************************************************************************/
00024   double muplus  = 1.0/sqrt(2.0) * (mu.E() + mu.Z());
00025   double muminus = 1.0/sqrt(2.0) * (mu.E() - mu.Z());
00026 
00027   double mubarplus  = 1.0/sqrt(2.0) * (mubar.E() + mubar.Z());
00028   double mubarminus = 1.0/sqrt(2.0) * (mubar.E() - mubar.Z());
00029  
00030   double costheta = 2.0 / Q.Mag() / sqrt(pow(Q.Mag(), 2) + pow(Q.Pt(), 2)) * 
00031     (muplus * mubarminus - muminus * mubarplus);
00032   if (swap) costheta = -costheta;
00033 
00034 
00035 
00036   /************************************************************************
00037    *
00038    * 2) sin2(theta) = Q^-2 Dt^2 - Q^-2 (Q^2 + Qt^2)^-1 * (Dt dot Qt)^2
00039    *
00040    ************************************************************************/
00041   TLorentzVector D(mu-mubar);
00042   double dt_qt = D.X()*Q.X() + D.Y()*Q.Y();
00043   double sin2theta = pow(D.Pt()/Q.Mag(), 2)
00044     - 1.0/pow(Q.Mag(), 2)/(pow(Q.Mag(), 2) + pow(Q.Pt(), 2))*pow(dt_qt, 2);
00045 
00046  
00047 
00048   /************************************************************************
00049    *
00050    * 3) tanphi = (Q^2 + Qt^2)^1/2 / Q (Dt dot R unit) /(Dt dot Qt unit)
00051    *
00052    ************************************************************************/
00053   // unit vector on R direction
00054   TVector3 R = Pbeam.Vect().Cross(Q.Vect());
00055   TVector3 Runit = R.Unit();
00056 
00057 
00058   // unit vector on Qt
00059   TVector3 Qt = Q.Vect(); Qt.SetZ(0);
00060   TVector3 Qtunit = Qt.Unit();
00061 
00062 
00063   TVector3 Dt = D.Vect(); Dt.SetZ(0);
00064   double tanphi = sqrt(pow(Q.Mag(), 2) + pow(Q.Pt(), 2)) / Q.Mag() * 
00065     Dt.Dot(Runit) / Dt.Dot(Qtunit);
00066   if (swap) tanphi = -tanphi;
00067 
00068   res[0] = costheta;
00069   res[1] = sin2theta;
00070   res[2] = tanphi;
00071 }
00072