add phi-phi moment to Jet

stadie · Jan 25, 2010 · a0c39a3 · a0c39a3
1 parent 28ce5fe
commit a0c39a3
Show file tree

Hide file tree

Showing 12 changed files with 84 additions and 59 deletions.
diff --git a/CalibData.h b/CalibData.h
@@ -1,5 +1,5 @@
 //
-// $Id: CalibData.h,v 1.77 2009/11/26 18:24:42 stadie Exp $
+// $Id: CalibData.h,v 1.78 2010/01/14 13:13:17 mschrode Exp $
 //
 #ifndef CalibData_h
 #define CalibData_h
@@ -26,19 +26,20 @@ enum DataType {Default, TrackTower, GammaJet, TrackCluster, MessMess, PtBalance,
 //!
 //!  \sa Jet, Tower, Track, JetWithTowers, JetWithTracks
 //!
-//!  \author Christian Autermann
-//!  $Id: CalibData.h,v 1.77 2009/11/26 18:24:42 stadie Exp $
+//!  \author Christian Autermann, Hartmut Stadie
+//!  $Id: CalibData.h,v 1.78 2010/01/14 13:13:17 mschrode Exp $
 class Measurement
 {
 public:
  Measurement() :
-  pt(0.),EMF(0.),HadF(0.),OutF(0.),E(0.),eta(0.),phi(0.),etaeta(0.0)
+  pt(0.),EMF(0.),HadF(0.),OutF(0.),E(0.),eta(0.),phi(0.),phiphi(0),
+   etaeta(0)
     {
     }
  Measurement(double Et,double EmEt,double HadEt,double OutEt,double E,
-	     double eta,double phi, double etaeta = 0)
+	     double eta,double phi, double phiphi = 0,double etaeta = 0)
    : pt(Et),EMF(EmEt),HadF(HadEt),OutF(OutEt),E(E),eta(eta),phi(phi),
-    etaeta(etaeta) 
+    phiphi(phiphi),etaeta(etaeta) 
   {
   }
   virtual ~Measurement() {};
@@ -50,6 +51,7 @@ class Measurement
   double E;      //!< Total energy					
   double eta;    //!< Pseudorapidity eta				
   double phi;    //!< Polar angle phi  
+  double phiphi; //!< Phi-Phi moment (width in phi) 
   double etaeta; //!< Eta-Eta moment (width in eta) 
 };
 
@@ -63,7 +65,7 @@ class Measurement
 //!  \todo Document members
 //!
 //!  \author Jan Thomsen
-//!  $Id: CalibData.h,v 1.77 2009/11/26 18:24:42 stadie Exp $
+//!  $Id: CalibData.h,v 1.78 2010/01/14 13:13:17 mschrode Exp $
 class TTrack : public Measurement
 {
 public:
@@ -120,7 +122,7 @@ class TTrack : public Measurement
 //!     The available data types are:
 //!  \author Christian Autermann
 //!  \date Wed Jul 18 13:54:50 CEST 2007
-//! $Id: CalibData.h,v 1.77 2009/11/26 18:24:42 stadie Exp $
+//! $Id: CalibData.h,v 1.78 2010/01/14 13:13:17 mschrode Exp $
 class Event
 {
 public:
@@ -235,7 +237,7 @@ class Event
 //!
 //!  \author Hartmut Stadie
 //!  \date Thu Dec 11 17:20:25 2008 UTC
-//!  $Id: CalibData.h,v 1.77 2009/11/26 18:24:42 stadie Exp $
+//!  $Id: CalibData.h,v 1.78 2010/01/14 13:13:17 mschrode Exp $
 class TAbstractData : public Event
 {
 public:

diff --git a/DiJetReader.cc b/DiJetReader.cc
@@ -1,6 +1,6 @@
 //
 //    first version: Hartmut Stadie 2008/12/12
-//    $Id: DiJetReader.cc,v 1.34 2010/01/12 19:24:48 mschrode Exp $
+//    $Id: DiJetReader.cc,v 1.35 2010/01/14 13:11:28 mschrode Exp $
 //   
 #include "DiJetReader.h"
 
@@ -300,8 +300,6 @@ int DiJetReader::createJetTruthEvents(std::vector<Event*>& data)
     Measurement tower;
     double dR        = 10;
     int closestTower = 0;
-    double seta = 0;
-    double seta2 = 0;
     double sumpt = 0;
     for(int n=0; n<nJet_->NobjTow; ++n){
       if(nJet_->Tow_jetidx[n] != calJetIdx) continue;//look for ij-jet's towers
@@ -326,8 +324,6 @@ int DiJetReader::createJetTruthEvents(std::vector<Event*>& data)
       err2     += terr[n];
       dphi      = TVector2::Phi_mpi_pi(nJet_->JetPhi[calJetIdx]-tower.phi);
       deta      = nJet_->JetEta[calJetIdx]-tower.eta;
-      seta += tower.pt * tower.eta;
-      seta2 += tower.pt * tower.eta * tower.eta;
       sumpt += tower.pt;
       double dr = sqrt( deta*deta + dphi*dphi );     
       if(dr < dR) {
@@ -351,16 +347,17 @@ int DiJetReader::createJetTruthEvents(std::vector<Event*>& data)
     tower.eta     = nJet_->JetEta[calJetIdx];
     tower.phi     = nJet_->JetPhi[calJetIdx];
     tower.E       = nJet_->JetE[calJetIdx];
-    tower.etaeta  = sqrt(seta2/sumpt - seta * seta /(sumpt * sumpt));
     double err    = jet_error_param(&tower.pt,&tower,0);
     err2         += err * err;
 
     Jet *jet;
     if(dataClass_ == 12) {
       JetWithTowers *jt = 
 	new JetWithTowers(nJet_->JetEt[calJetIdx],em * factor,had * factor,
-			  out * factor,nJet_->JetE[calJetIdx],nJet_->JetEta[calJetIdx],
-			  nJet_->JetPhi[calJetIdx],tower.etaeta, Jet::uds,
+			  out * factor,nJet_->JetE[calJetIdx],
+			  nJet_->JetEta[calJetIdx],nJet_->JetPhi[calJetIdx],
+			  nJet_->JetEtWeightedSigmaPhi[calJetIdx],
+			  nJet_->JetEtWeightedSigmaEta[calJetIdx], Jet::uds,
 			  nJet_->GenJetColEt[genJetIdx],drJetGenjet,
 			  createCorFactors(calJetIdx),
 			  par_->jet_function(nJet_->TowId_eta[closestTower],
@@ -378,8 +375,11 @@ int DiJetReader::createJetTruthEvents(std::vector<Event*>& data)
       jet = jt;
     }
     else { 
-      jet = new Jet(nJet_->JetEt[calJetIdx],em * factor,had * factor,out * factor,
-		    nJet_->JetE[calJetIdx],nJet_->JetEta[calJetIdx],nJet_->JetPhi[calJetIdx],tower.etaeta,
+      jet = new Jet(nJet_->JetCorrL2L3[calJetIdx]*nJet_->JetEt[calJetIdx],em * factor,had * factor,out * factor,
+		    nJet_->JetCorrL2L3[calJetIdx]*nJet_->JetE[calJetIdx],
+		    nJet_->JetEta[calJetIdx],nJet_->JetPhi[calJetIdx],
+		    nJet_->JetEtWeightedSigmaPhi[calJetIdx],
+		    nJet_->JetEtWeightedSigmaEta[calJetIdx],
 		    Jet::uds,nJet_->GenJetColEt[genJetIdx],drJetGenjet,
 		    createCorFactors(calJetIdx),
 		    par_->jet_function(nJet_->TowId_eta[closestTower],
@@ -721,7 +721,8 @@ std::vector<Jet*> DiJetReader::readCaloJets(int nJets) const {
 			  nJet_->JetE[j],
 			  nJet_->JetEta[j],
 			  nJet_->JetPhi[j],
-			  0.,
+			  nJet_->JetEtWeightedSigmaPhi[j],
+			  nJet_->JetEtWeightedSigmaEta[j],
 			  Jet::uds,
 			  nJet_->GenJetPt[j],
 			  drJetGenjet,

diff --git a/Jet.cc b/Jet.cc
@@ -2,7 +2,7 @@
 //    Class for basic jets 
 //
 //    first version: Hartmut Stadie 2008/12/14
-//    $Id: Jet.cc,v 1.35 2009/11/26 18:24:41 stadie Exp $
+//    $Id: Jet.cc,v 1.36 2009/11/27 15:28:12 stadie Exp $
 //   
 #include "Jet.h"  
 
@@ -12,13 +12,13 @@
 #include <iomanip>
 
 Jet::Jet(double Et, double EmEt, double HadEt ,double OutEt, double E,
-         double eta,double phi, double etaeta, Flavor flavor, double genPt, 
-	 double dR, CorFactors* corFactors, const Function& f, 
+         double eta,double phi, double phiphi, double etaeta, Flavor flavor, 
+	 double genPt, double dR, CorFactors* corFactors, const Function& f, 
 	 double (*errfunc)(const double *x, const Measurement *xorig, double err), 
 	 const Function& gf, double Etmin) 
-  : Measurement(Et,EmEt,HadEt,OutEt,E,eta,phi,etaeta),flavor_(flavor), genPt_(genPt), 
-    dR_(dR), ptHat_(0.), corFactors_(corFactors),f(f),gf(gf),errf(errfunc),
-    etmin(Etmin),EoverPt(E/Et),gsl_impl(this)
+  : Measurement(Et,EmEt,HadEt,OutEt,E,eta,phi,phiphi,etaeta),flavor_(flavor), 
+    genPt_(genPt),dR_(dR), ptHat_(0.), corFactors_(corFactors),f(f),gf(gf),
+    errf(errfunc),etmin(Etmin),EoverPt(E/Et),gsl_impl(this)
 {
   temp = *this;
   varcoll.resize(f.nPars() + gf.nPars());

diff --git a/Jet.h b/Jet.h
@@ -4,7 +4,7 @@
 //!
 //!    \date 2008/12/14
 //!
-//!    $Id: Jet.h,v 1.28 2009/11/27 15:28:12 stadie Exp $
+//!    $Id: Jet.h,v 1.29 2010/01/12 19:24:48 mschrode Exp $
 #ifndef JET_H
 #define JET_H
 
@@ -45,8 +45,8 @@ class Jet : public Measurement
 
  public:
   Jet(double Et, double EmEt, double HadEt ,double OutEt, double E,
-      double eta,double phi, double etaeta, Flavor flavor, double genPt, 
-      double dR, CorFactors* corFactors, const Function& f,
+      double eta,double phi, double phiphi, double etaeta, Flavor flavor, 
+      double genPt, double dR, CorFactors* corFactors, const Function& f,
       double (*errfunc)(const double *x, const Measurement *xorig, double err), 
       const Function& gf, double Etmin = 0); 
   virtual ~Jet();
@@ -59,6 +59,7 @@ class Jet : public Measurement
   double E()      const {return Measurement::E;}    //!< Return energy
   double eta()    const {return Measurement::eta;}  //!< Return pseudorapidity
   double phi()    const {return Measurement::phi;}  //!< Return azimuthal angle
+  double momentPhiPhi() const {return Measurement::phiphi;}  //!< Return phi-phi moment (width of jet in phi)
   double momentEtaEta() const {return Measurement::etaeta;}  //!< Return eta-eta moment (width of jet in eta)
   Flavor flavor() const {return flavor_;}       //!< Return jet flavor
   double genPt()  const {return genPt_;}        //!< Return Pt for corresponding GenJet 

diff --git a/JetConstraintsReader.cc b/JetConstraintsReader.cc
@@ -4,7 +4,7 @@
 //    This class add user defined jet constraints
 //
 //    first version: Hartmut Stadie 2008/12/12
-//    $Id: JetConstraintsReader.cc,v 1.7 2009/11/25 13:07:45 stadie Exp $
+//    $Id: JetConstraintsReader.cc,v 1.8 2009/11/26 18:24:41 stadie Exp $
 //   
 
 #include "JetConstraintsReader.h"
@@ -49,7 +49,7 @@ int JetConstraintsReader::readEvents(std::vector<Event*>& data)
     JetConstraintEvent* jce = new JetConstraintEvent(ic->Et,ic->weight);
     for(int ideta = ic->mineta ; ideta <= ic->maxeta  ; ++ideta) {
       if(ideta == 0) ideta = 1;
-      jce->addJet(new  Jet(ic->Et,0,ic->Et,0,ic->Et,0,0,0.13,Jet::uds,ic->Et,0,
+      jce->addJet(new  Jet(ic->Et,0,ic->Et,0,ic->Et,0,0,0.13,0.13,Jet::uds,ic->Et,0,
 			   new CorFactors(1,1,1,1,1,1,1),par_->jet_function(ideta,1),
 			   jet_error_param,Function(&Parametrization::correctedJetEt,0,0,0,0,cp),0));
     }

diff --git a/JetWithTowers.cc b/JetWithTowers.cc
@@ -2,18 +2,19 @@
 //    Class for jets with towers 
 //
 //    first version: Hartmut Stadie 2008/12/25
-//    $Id: JetWithTowers.cc,v 1.20 2009/11/25 13:07:45 stadie Exp $
+//    $Id: JetWithTowers.cc,v 1.21 2009/11/26 18:24:41 stadie Exp $
 //   
 #include"JetWithTowers.h"
 
 #include "TLorentzVector.h"
 
 JetWithTowers::JetWithTowers(double Et, double EmEt, double HadEt ,double OutEt, double E,
-			     double eta,double phi, double etaeta, Flavor flavor, double genPt, double dR,
+			     double eta,double phi, double phiphi, double etaeta, 
+			     Flavor flavor, double genPt, double dR,
 			     CorFactors* corFactors, const Function& func, 
 			     double (*errfunc)(const double *x, const Measurement *xorig, double err), 
 			     const Function& gfunc, double Etmin) 
-  :  Jet(Et,EmEt,HadEt,OutEt,E,eta,phi,etaeta,flavor,genPt,dR,corFactors,
+  :  Jet(Et,EmEt,HadEt,OutEt,E,eta,phi,phiphi,etaeta,flavor,genPt,dR,corFactors,
 	 func,errfunc,gfunc,Etmin),
      ntowerpars(0)
 {

diff --git a/JetWithTowers.h b/JetWithTowers.h
@@ -12,14 +12,15 @@
 //!
 //!    \author Hartmut Stadie
 //!    \date 2008/12/25
-//!    $Id: JetWithTowers.h,v 1.17 2009/11/25 13:07:45 stadie Exp $
+//!    $Id: JetWithTowers.h,v 1.18 2009/11/26 18:24:42 stadie Exp $
 // ----------------------------------------------------------------   
 class JetWithTowers : public Jet
 {
  public:
   JetWithTowers(double Et, double EmEt, double HadEt ,double OutEt, double E,
-		double eta,double phi, double etaeta, Flavor flavor,
-		double genPt, double dR, CorFactors* corFactors, const Function& f,
+		double eta,double phi, double phiphi, double etaeta, 
+		Flavor flavor, double genPt, double dR, CorFactors* corFactors,
+		const Function& f,
 		double (*errfunc)(const double *x, const Measurement *xorig, double err), 
 		const Function& gf, double Etmin = 0); 
   virtual ~JetWithTowers(); 

diff --git a/JetWithTracks.cc b/JetWithTracks.cc
@@ -2,18 +2,18 @@
 //    Class for jets with tracks 
 //
 //    first version: Hartmut Stadie 2009/04/08
-//    $Id: JetWithTracks.cc,v 1.7 2009/11/25 13:07:45 stadie Exp $
+//    $Id: JetWithTracks.cc,v 1.8 2009/11/26 18:24:41 stadie Exp $
 //   
 #include"JetWithTracks.h"
 
 #include "TLorentzVector.h"
 
 JetWithTracks::JetWithTracks(double Et, double EmEt, double HadEt ,double OutEt, double E,
-			     double eta,double phi, double etaeta, Flavor flavor, double genPt, double dR,
+			     double eta,double phi, double phiphi, double etaeta, Flavor flavor, double genPt, double dR,
 			     CorFactors* corFactors, const Function& func, 
 			     double (*errfunc)(const double *x, const Measurement *xorig, double err), 
 			     const Function& gfunc, double Etmin) 
-  :  Jet(Et,EmEt,HadEt,OutEt,E,eta,phi,etaeta,flavor,genPt,dR,corFactors,
+  :  Jet(Et,EmEt,HadEt,OutEt,E,eta,phi,phiphi,etaeta,flavor,genPt,dR,corFactors,
 	 func,errfunc,gfunc,Etmin),
      ntrackpars(0), expectedCaloEt(0), trackPt(0)
 {

diff --git a/JetWithTracks.h b/JetWithTracks.h
@@ -12,14 +12,15 @@
 //!
 //!    \author Hartmut Stadie
 //!    \date 2008/12/25
-//!    $Id: JetWithTracks.h,v 1.6 2009/11/25 13:07:45 stadie Exp $
+//!    $Id: JetWithTracks.h,v 1.7 2009/11/26 18:24:42 stadie Exp $
 // ---------------------------------------------------------------   
 class JetWithTracks : public Jet
 {
  public:
   JetWithTracks(double Et, double EmEt, double HadEt ,double OutEt, double E,
-		double eta,double phi, double etaeta, Flavor flavor,double genPt, double dR,
-		CorFactors* corFactors, const Function& f,
+		double eta,double phi, double phiphi, double etaeta, 
+		Flavor flavor,double genPt, double dR, CorFactors* corFactors,
+		const Function& f,
 		double (*errfunc)(const double *x, const Measurement *xorig, double err), 
 		const Function& gf, double Etmin = 0); 
   virtual ~JetWithTracks(); 

diff --git a/PhotonJetReader.cc b/PhotonJetReader.cc
@@ -1,5 +1,5 @@
 //
-//  $Id: PhotonJetReader.cc,v 1.26 2009/11/26 18:24:41 stadie Exp $
+//  $Id: PhotonJetReader.cc,v 1.27 2009/11/27 15:28:12 stadie Exp $
 //
 #include "PhotonJetReader.h"
 
@@ -214,7 +214,9 @@ Event* PhotonJetReader::createJetTruthEvent()
   // Loop over towers, find closest tower to jet axis,
   // and sum up emf, hadf, outf
   double seta = 0;
-  double seta2 = 0;
+  double seta2 = 0; 
+  double sphi = 0;
+  double sphi2 = 0;
   double sumpt = 0;
   for(int n = 0; n < gammaJet_->NobjTowCal; ++n) {
     em          += gammaJet_->TowEm[n];
@@ -224,7 +226,9 @@ Event* PhotonJetReader::createJetTruthEvent()
     double phi =  gammaJet_->TowPhi[n];
     double eta =  gammaJet_->TowEta[n];
     seta += gammaJet_->TowEt[n] * eta;
-    seta2 += gammaJet_->TowEt[n] * eta * eta;
+    seta2 += gammaJet_->TowEt[n] * eta * eta; 
+    sphi += gammaJet_->TowEt[n] * phi;
+    sphi2 += gammaJet_->TowEt[n] * phi * phi;
     sumpt += gammaJet_->TowEt[n];
     double dphi  = TVector2::Phi_mpi_pi(gammaJet_->JetCalPhi-phi);
     double dr    = sqrt((gammaJet_->JetCalEta-eta)*(gammaJet_->JetCalEta-eta)+
@@ -246,12 +250,13 @@ Event* PhotonJetReader::createJetTruthEvent()
 
   double factor = gammaJet_->JetCalEt /  gammaJet_->JetCalE;
   double etaeta = sqrt(seta2/sumpt - seta * seta /(sumpt * sumpt));
+  double phiphi = sqrt(sphi2/sumpt - sphi * sphi /(sumpt * sumpt));
   Jet *j;
   if(dataClass_ == 2) {
     JetWithTowers *jt = 
       new JetWithTowers(gammaJet_->JetCalEt,em * factor,had * factor,
 			out * factor,gammaJet_->JetCalE,gammaJet_->JetCalEta,
-			gammaJet_->JetCalPhi,etaeta,Jet::uds,gammaJet_->JetGenEt,
+			gammaJet_->JetCalPhi,phiphi,etaeta,Jet::uds,gammaJet_->JetGenEt,
 			LJet.DeltaR(LGenJet),createCorFactors(0),
 			par_->jet_function(gammaJet_->TowId_eta[closestTower],
 					   gammaJet_->TowId_phi[closestTower]),
@@ -269,7 +274,7 @@ Event* PhotonJetReader::createJetTruthEvent()
   else { 
     j = new Jet(gammaJet_->JetCalEt,em * factor,had * factor,out * factor,
 		gammaJet_->JetCalE,gammaJet_->JetCalEta,gammaJet_->JetCalPhi,
-		etaeta,Jet::uds,gammaJet_->JetGenEt,LJet.DeltaR(LGenJet),
+		phiphi,etaeta,Jet::uds,gammaJet_->JetGenEt,LJet.DeltaR(LGenJet),
 		createCorFactors(0),
 		par_->jet_function(gammaJet_->TowId_eta[closestTower],
 				   gammaJet_->TowId_phi[closestTower]),