Fork me on GitHub

Changeset 89 in svn for trunk


Ignore:
Timestamp:
Dec 5, 2008, 9:34:25 AM (16 years ago)
Author:
severine ovyn
Message:

for resol plots

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/Resolutions.cpp

    r88 r89  
    183183 
    184184  //read the datacard input file
    185   string DetDatacard("");
     185  string DetDatacard("data/DataCardDet.dat");
    186186  if(argc==4)  DetDatacard =argv[3];
    187187  RESOLution *DET = new RESOLution();   
     
    222222      TrackCentral.clear();
    223223      TSimpleArray<TRootGenParticle> NFCentralQ;
    224    
     224      
    225225      input_particlesReco.clear();
    226226      input_particlesGEN.clear();
     
    234234          int pid = abs(particle->PID);
    235235          float eta = fabs(particle->Eta);
    236 
     236         
    237237          //input generator level particle for jet algorithm     
    238238          if(particle->Status == 1 && eta < DET->MAX_CALO_FWD)
     
    240240              input_particlesGEN.push_back(fastjet::PseudoJet(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E()));
    241241            }
    242 
     242         
    243243          //Calculate ETMIS from generated particles
    244244          if((pid == pNU1) || (pid == pNU2) || (pid == pNU3))PTmisGEN = PTmisGEN + genMomentum;
    245 
     245         
    246246          if( (particle->Status == 1)    &&
    247247              ((pid != pNU1) && (pid != pNU2) && (pid != pNU3)) &&
     
    252252              //use of the magnetic field propagation
    253253              TRACP->Propagation(particle,recoMomentum);
    254               float eta=fabs(recoMomentum.Eta());
    255 
     254              //              cout<<"eta "<<eta<<endl;
     255              eta=fabs(recoMomentum.Eta());
     256              //            cout<<"eta apres"<<eta<<endl;
     257             
    256258              switch(pid) {
    257 
     259               
    258260              case pE: // all electrons with eta < DET->MAX_CALO_FWD
    259261                DET->SmearElectron(recoMomentum);
    260                 if(eta < DET->MAX_TRACKER){
    261                 elementElec=(RESOLELEC*) branchelec->NewEntry();
    262                 elementElec->E = genMomentum.E();
    263                 elementElec->SmearedE = recoMomentum.E();}
     262                if(recoMomentum.E() !=0 && eta < DET->MAX_TRACKER){
     263                  elementElec=(RESOLELEC*) branchelec->NewEntry();
     264                  elementElec->E = genMomentum.E();
     265                  elementElec->SmearedE = recoMomentum.E();}
    264266                break; // case pE
    265267              case pGAMMA: // all photons with eta < DET->MAX_CALO_FWD
     
    268270              case pMU: // all muons with eta < DET->MAX_MU
    269271                DET->SmearMu(recoMomentum);
     272                if(recoMomentum.E() !=0){
    270273                elementMuon = (RESOLMUON*) branchmuon->NewEntry();
    271274                elementMuon->OverPT = (1/genMomentum.Pt());
    272                 elementMuon->OverSmearedPT = (1/recoMomentum.Pt());
     275                elementMuon->OverSmearedPT = (1/recoMomentum.Pt());}
    273276                break; // case pMU
    274277              case pLAMBDA: // all lambdas with eta < DET->MAX_CALO_FWD
     
    280283                break;
    281284              } // switch (pid)
    282              
    283             //information to reconstruct jets from reconstructed towers
    284             int charge=Charge(pid);
    285             if(recoMomentum.E() !=0 && pid != pMU) {
    286               if(charge == 0 || (charge !=0 && recoMomentum.Pt() >= DET->PT_TRACKS_MIN)){
     285              
     286              //information to reconstruct jets from reconstructed towers
     287              int charge=Charge(pid);
     288              if(recoMomentum.E() !=0 && pid != pMU) {
     289                if(charge == 0 || (charge !=0 && recoMomentum.Pt() >= DET->PT_TRACKS_MIN)){
    287290                  towers.push_back(PhysicsTower(LorentzVector(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E())));
    288291                  input_particlesReco.push_back(fastjet::PseudoJet(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E()));
    289292                }
    290             }
     293              }
    291294             
    292295              // all final charged particles
    293                  if(
    294                   (genMomentum.E()!=0) &&
    295                   (fabs(genMomentum.Eta()) < DET->MAX_TRACKER) &&
    296                   (genMomentum.Pt() > DET->PT_TRACKS_MIN ) &&     // pt too small to be taken into account
    297                   ((rand()%100) < DET->TRACKING_EFF)  &&
    298                   (charge!=0)
    299                   )
    300                    {TrackCentral.push_back(genMomentum);}
    301 
     296              if(
     297                (genMomentum.E()!=0) &&
     298                (fabs(genMomentum.Eta()) < DET->MAX_TRACKER) &&
     299                (genMomentum.Pt() > DET->PT_TRACKS_MIN ) &&     // pt too small to be taken into account
     300                ((rand()%100) < DET->TRACKING_EFF)  &&
     301                (charge!=0)
     302                )
     303                {TrackCentral.push_back(genMomentum);}
     304             
    302305            } // switch
    303306        } // while
    304 
    305      //compute missing transverse energy from calo towers   
    306      TLorentzVector Att(0.,0.,0.,0.);
    307      for(unsigned int i=0; i < towers.size(); i++)
    308        {
    309          Att.SetPxPyPzE(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E);
    310          PTmisReco = PTmisReco + Att;
    311        }
    312 
     307     
     308      //compute missing transverse energy from calo towers   
     309      TLorentzVector Att(0.,0.,0.,0.);
     310      float ScalarEt=0;
     311      for(unsigned int i=0; i < towers.size(); i++)
     312        {
     313          Att.SetPxPyPzE(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E);
     314          ScalarEt = ScalarEt + Att.Et();
     315          PTmisReco = PTmisReco + Att;
     316        }
     317     
    313318      elementEtmis= (ETMIS*) branchetmis->NewEntry();
    314319      elementEtmis->Et = (PTmisGEN).Pt();
    315       elementEtmis->SmearedEt = (PTmisReco).Pt()-(PTmisGEN).Pt();
    316      
     320      elementEtmis->Ex = (-PTmisGEN).Px();
     321      elementEtmis->SEt = ScalarEt;
     322
     323      elementEtmis->EtSmeare = (PTmisReco).Pt()-(PTmisGEN).Pt();
     324      elementEtmis->ExSmeare = (-PTmisReco).Px()-(PTmisGEN).Px();
     325
    317326      //*****************************
    318327     
  • trunk/interface/TreeClasses.h

    r88 r89  
    7373};
    7474
     75//---------------------------------------------
     76
    7577class ETMIS : public TSortableObject
    7678{
    7779public:
    7880Float_t Et;
    79 Float_t SmearedEt;
     81Float_t SEt;
     82Float_t Ex;
     83Float_t EtSmeare;
     84Float_t ExSmeare;
    8085
    8186 static TCompare *fgCompare; //!
     
    8691
    8792};
     93
    8894
    8995
  • trunk/routines/resolutions.C

    r60 r89  
    103103  gROOT->Reset();
    104104 
    105   TFile *f1 = new TFile("PTMIS.root","read");
     105  TFile *f1 = new TFile("Etmis.root","read");
    106106  TTree *Analyze = (TTree*)f1->Get("Analysis");
    107107 
     
    132132      char tempMin[500];
    133133      if(i==0)binMin=5;
    134       sprintf(tempMin,"ElecEResol.NonSmeareE > %d",binMin);
     134      sprintf(tempMin,"ElecEResol.E > %d",binMin);
    135135      string mystringMin(tempMin);
    136136      char tempMax[500];
    137       sprintf(tempMax,"ElecEResol.NonSmeareE < %d",binMax);
     137      sprintf(tempMax,"ElecEResol.E < %d",binMax);
    138138      string mystringMax(tempMax);
    139139      char tempName[500];
    140       sprintf(tempName,"(ElecEResol.NonSmeareE-ElecEResol.SmeareE)>>hETSoverET%d",i);
     140      sprintf(tempName,"(ElecEResol.E-ElecEResol.SmearedE)>>hETSoverET%d",i);
    141141      string mystringName(tempName);
    142142     
     
    189189  gROOT->Reset();
    190190 
    191   TFile *f1 = new TFile("PTMIS.root","read");
     191  TFile *f1 = new TFile("Etmis.root","read");
    192192  TTree *Analyze = (TTree*)f1->Get("Analysis");
    193193 
     
    233233  gROOT->Reset();
    234234 
    235   TFile *f1 = new TFile("PTMIS.root","read");
     235  TFile *f1 = new TFile("Etmis.root","read");
    236236  TTree *Analyze = (TTree*)f1->Get("Analysis");
    237237 
     
    311311void General()
    312312{
    313   JetResol();
     313//  JetResol();
    314314  ElecResol();
    315315  ETmisResol();
Note: See TracChangeset for help on using the changeset viewer.