Changeset 89 in svn

Timestamp:

Dec 5, 2008, 9:34:25 AM (16 years ago)

Author:

severine ovyn

Message:

for resol plots

Location:

trunk

Files:

• Resolutions.cpp (modified) (7 diffs)
• interface/TreeClasses.h (modified) (2 diffs)
• routines/resolutions.C (modified) (5 diffs)

trunk/Resolutions.cpp

@@ -183 +183 @@
   
   //read the datacard input file
-  string DetDatacard("");
+  string DetDatacard("data/DataCardDet.dat");
   if(argc==4)  DetDatacard =argv[3];
   RESOLution *DET = new RESOLution();   
@@ -222 +222 @@
       TrackCentral.clear();
       TSimpleArray<TRootGenParticle> NFCentralQ;
-   
+      
       input_particlesReco.clear();
       input_particlesGEN.clear();
@@ -234 +234 @@
           int pid = abs(particle->PID);
           float eta = fabs(particle->Eta);
-
+          
           //input generator level particle for jet algorithm      
           if(particle->Status == 1 && eta < DET->MAX_CALO_FWD)
@@ -240 +240 @@
               input_particlesGEN.push_back(fastjet::PseudoJet(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E()));
             }
-
+          
           //Calculate ETMIS from generated particles
           if((pid == pNU1) || (pid == pNU2) || (pid == pNU3))PTmisGEN = PTmisGEN + genMomentum;
-
+          
           if( (particle->Status == 1)    &&
               ((pid != pNU1) && (pid != pNU2) && (pid != pNU3)) &&
@@ -252 +252 @@
               //use of the magnetic field propagation
               TRACP->Propagation(particle,recoMomentum);
-              float eta=fabs(recoMomentum.Eta());
-
+              //              cout<<"eta "<<eta<<endl;
+              eta=fabs(recoMomentum.Eta());
+              //            cout<<"eta apres"<<eta<<endl;
+              
               switch(pid) {
-
+                
               case pE: // all electrons with eta < DET->MAX_CALO_FWD
                 DET->SmearElectron(recoMomentum);
-                if(eta < DET->MAX_TRACKER){
-                elementElec=(RESOLELEC*) branchelec->NewEntry();
-                elementElec->E = genMomentum.E();
-                elementElec->SmearedE = recoMomentum.E();}
+                if(recoMomentum.E() !=0 && eta < DET->MAX_TRACKER){
+                  elementElec=(RESOLELEC*) branchelec->NewEntry();
+                  elementElec->E = genMomentum.E();
+                  elementElec->SmearedE = recoMomentum.E();}
                 break; // case pE
               case pGAMMA: // all photons with eta < DET->MAX_CALO_FWD
@@ -268 +270 @@
               case pMU: // all muons with eta < DET->MAX_MU
                 DET->SmearMu(recoMomentum);
+                if(recoMomentum.E() !=0){
                 elementMuon = (RESOLMUON*) branchmuon->NewEntry();
                 elementMuon->OverPT = (1/genMomentum.Pt());
-                elementMuon->OverSmearedPT = (1/recoMomentum.Pt());
+                elementMuon->OverSmearedPT = (1/recoMomentum.Pt());}
                 break; // case pMU
               case pLAMBDA: // all lambdas with eta < DET->MAX_CALO_FWD
@@ -280 +283 @@
                 break;
               } // switch (pid)
-             
-            //information to reconstruct jets from reconstructed towers
-            int charge=Charge(pid);
-            if(recoMomentum.E() !=0 && pid != pMU) {
-              if(charge == 0 || (charge !=0 && recoMomentum.Pt() >= DET->PT_TRACKS_MIN)){
+              
+              //information to reconstruct jets from reconstructed towers
+              int charge=Charge(pid);
+              if(recoMomentum.E() !=0 && pid != pMU) {
+                if(charge == 0 || (charge !=0 && recoMomentum.Pt() >= DET->PT_TRACKS_MIN)){
                   towers.push_back(PhysicsTower(LorentzVector(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E())));
                   input_particlesReco.push_back(fastjet::PseudoJet(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E()));
                 }
-            }
+              }
               
               // all final charged particles
-                 if(
-                   (genMomentum.E()!=0) &&
-                   (fabs(genMomentum.Eta()) < DET->MAX_TRACKER) &&
-                   (genMomentum.Pt() > DET->PT_TRACKS_MIN ) &&     // pt too small to be taken into account
-                   ((rand()%100) < DET->TRACKING_EFF)  &&
-                   (charge!=0)
-                   )
-                   {TrackCentral.push_back(genMomentum);}
-
+              if(
+                 (genMomentum.E()!=0) &&
+                 (fabs(genMomentum.Eta()) < DET->MAX_TRACKER) &&
+                 (genMomentum.Pt() > DET->PT_TRACKS_MIN ) &&     // pt too small to be taken into account
+                 ((rand()%100) < DET->TRACKING_EFF)  &&
+                 (charge!=0)
+                 )
+                {TrackCentral.push_back(genMomentum);}
+              
             } // switch
         } // while
-
-     //compute missing transverse energy from calo towers    
-     TLorentzVector Att(0.,0.,0.,0.);
-     for(unsigned int i=0; i < towers.size(); i++) 
-       {
-         Att.SetPxPyPzE(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E);
-         PTmisReco = PTmisReco + Att;
-       }
-
+      
+      //compute missing transverse energy from calo towers    
+      TLorentzVector Att(0.,0.,0.,0.);
+      float ScalarEt=0;
+      for(unsigned int i=0; i < towers.size(); i++) 
+        {
+          Att.SetPxPyPzE(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E);
+          ScalarEt = ScalarEt + Att.Et();
+          PTmisReco = PTmisReco + Att;
+        }
+      
       elementEtmis= (ETMIS*) branchetmis->NewEntry();
       elementEtmis->Et = (PTmisGEN).Pt();
-      elementEtmis->SmearedEt = (PTmisReco).Pt()-(PTmisGEN).Pt();
-      
+      elementEtmis->Ex = (-PTmisGEN).Px();
+      elementEtmis->SEt = ScalarEt;
+
+      elementEtmis->EtSmeare = (PTmisReco).Pt()-(PTmisGEN).Pt();
+      elementEtmis->ExSmeare = (-PTmisReco).Px()-(PTmisGEN).Px();
+
       //*****************************
       

trunk/interface/TreeClasses.h

@@ -73 +73 @@
 };
 
+//---------------------------------------------
+
 class ETMIS : public TSortableObject
 {
 public:
 Float_t Et;
-Float_t SmearedEt;
+Float_t SEt;
+Float_t Ex;
+Float_t EtSmeare;
+Float_t ExSmeare;
 
  static TCompare *fgCompare; //!
@@ -86 +91 @@
 
 };
+
 
 

trunk/routines/resolutions.C

@@ -103 +103 @@
   gROOT->Reset();
   
-  TFile *f1 = new TFile("PTMIS.root","read");
+  TFile *f1 = new TFile("Etmis.root","read");
   TTree *Analyze = (TTree*)f1->Get("Analysis");
   
@@ -132 +132 @@
       char tempMin[500];
       if(i==0)binMin=5;
-      sprintf(tempMin,"ElecEResol.NonSmeareE > %d",binMin);
+      sprintf(tempMin,"ElecEResol.E > %d",binMin);
       string mystringMin(tempMin);
       char tempMax[500];
-      sprintf(tempMax,"ElecEResol.NonSmeareE < %d",binMax);
+      sprintf(tempMax,"ElecEResol.E < %d",binMax);
       string mystringMax(tempMax);
       char tempName[500];
-      sprintf(tempName,"(ElecEResol.NonSmeareE-ElecEResol.SmeareE)>>hETSoverET%d",i);
+      sprintf(tempName,"(ElecEResol.E-ElecEResol.SmearedE)>>hETSoverET%d",i);
       string mystringName(tempName);
       
@@ -189 +189 @@
   gROOT->Reset();
   
-  TFile *f1 = new TFile("PTMIS.root","read");
+  TFile *f1 = new TFile("Etmis.root","read");
   TTree *Analyze = (TTree*)f1->Get("Analysis");
   
@@ -233 +233 @@
   gROOT->Reset();
   
-  TFile *f1 = new TFile("PTMIS.root","read");
+  TFile *f1 = new TFile("Etmis.root","read");
   TTree *Analyze = (TTree*)f1->Get("Analysis");
   
@@ -311 +311 @@
 void General()
 {
-  JetResol();
+//  JetResol();
   ElecResol();
   ETmisResol();