Changeset 152 in svn

Timestamp:

Jan 7, 2009, 12:23:20 PM (16 years ago)

Author:

severine ovyn

Message:

used for the publi

File:

• trunk/Resolutions.cpp (modified) (9 diffs)

trunk/Resolutions.cpp

@@ -181 +181 @@
   ETMIS *elementEtmis;
   
+int numTau=0;
+int numTauRec=0;
   
   //read the datacard input file
@@ -195 +197 @@
    
   TLorentzVector genMomentum(0,0,0,0);//TLorentzVector containing generator level information
+  TLorentzVector recoMomentumCalo(0,0,0,0);
   TLorentzVector recoMomentum(0,0,0,0);//TLorentzVector containing Reco level information
   LorentzVector jetMomentum;
@@ -206 +209 @@
 
   vector<PhysicsTower> towers;
-  
+ 
+ float iPhi=0,iEta=0;
+ 
   // Loop over all events
   Long64_t entry, allEntries = treeReader->GetEntries();
@@ -251 +256 @@
               recoMomentum = genMomentum;          
               //use of the magnetic field propagation
-              TRACP->Propagation(particle,recoMomentum);
+              //TRACP->Propagation(particle,recoMomentum);
               //              cout<<"eta "<<eta<<endl;
               eta=fabs(recoMomentum.Eta());
@@ -285 +290 @@
               
               //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->TRACK_ptmin)){
-                  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()));
-                }
-              }
+               int charge=Charge(pid);
+                if(recoMomentum.E() !=0 && pid != pMU) {
+                if(charge == 0 || (charge !=0 && recoMomentum.Pt() >= DET->TRACK_ptmin)){
+                  DET->BinEtaPhi(recoMomentum.Phi(), recoMomentum.Eta(), iPhi, iEta);
+                  if(iEta != -100 && iPhi != -100)
+                    {
+                      recoMomentumCalo.SetPtEtaPhiE(recoMomentum.Pt(),iEta,iPhi,recoMomentum.E());
+
+                      PhysicsTower Tower(LorentzVector(recoMomentumCalo.Px(),recoMomentumCalo.Py(),recoMomentumCalo.Pz(), recoMomentumCalo.E()));
+                      towers.push_back(Tower);
+                    }
+                }
+              }
               
               // all final charged particles
               if(
-                 (genMomentum.E()!=0) &&
-                 (fabs(genMomentum.Eta()) < DET->CEN_max_tracker) &&
-                 (genMomentum.Pt() > DET->TRACK_ptmin ) &&     // pt too small to be taken into account
+                 (recoMomentum.E()!=0) &&
+                 (fabs(recoMomentum.Eta()) < DET->CEN_max_tracker) &&
+                 (recoMomentum.Pt() > DET->TRACK_ptmin ) &&     // pt too small to be taken into account
                  ((rand()%100) < DET->TRACK_eff)  &&
                  (charge!=0)
                  )
-                {TrackCentral.push_back(genMomentum);}
+                {TrackCentral.push_back(recoMomentum);}
               
             } // switch
@@ -309 +320 @@
       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;
-        }
+      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);
+          if(fabs(Att.Eta()) < DET->CEN_max_calo_fwd)
+            {
+               ScalarEt = ScalarEt + Att.Et();
+               PTmisReco = PTmisReco + Att;
+               // create a fastjet::PseudoJet with these components and put it onto
+               // back of the input_particles vector
+               input_particlesReco.push_back(fastjet::PseudoJet(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E));
+            }
+        }
       
       elementEtmis= (ETMIS*) branchetmis->NewEntry();
@@ -343 +360 @@
           }
       }
+      numTau = numTau+TausHadr.GetEntries();
       for (unsigned int i = 0; i < sorted_jetsReco.size(); i++) {
         TLorentzVector JETT(0,0,0,0);
@@ -355 +373 @@
                     elementTaujet->EnergieCen = (DET->EnergySmallCone(towers,JETT.Eta(),JETT.Phi())/JETT.E());
                     elementTaujet->NumTrack = DET->NumTracks(TrackCentral,DET->TAU_track_pt,JETT.Eta(),JETT.Phi());
+                    if( (DET->EnergySmallCone(towers,JETT.Eta(),JETT.Phi())/JETT.E()) > 0.95 
+                     &&  (DET->NumTracks(TrackCentral,DET->TAU_track_pt,JETT.Eta(),JETT.Phi()))==1)numTauRec++;
+                    
                  }
               }
@@ -365 +386 @@
     } // Loop over all events
   treeWriter->Write();
-  
+float frac = numTauRec/numTau;  
+cout<<numTauRec<<endl;
+cout<<numTau<<endl;
+
   cout << "** Exiting..." << endl;
+  cout<<frac<<endl;
+
   
   delete treeWriter;