- Timestamp:
- Dec 5, 2008, 9:34:25 AM (16 years ago)
- Location:
- trunk
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Resolutions.cpp
r88 r89 183 183 184 184 //read the datacard input file 185 string DetDatacard(" ");185 string DetDatacard("data/DataCardDet.dat"); 186 186 if(argc==4) DetDatacard =argv[3]; 187 187 RESOLution *DET = new RESOLution(); … … 222 222 TrackCentral.clear(); 223 223 TSimpleArray<TRootGenParticle> NFCentralQ; 224 224 225 225 input_particlesReco.clear(); 226 226 input_particlesGEN.clear(); … … 234 234 int pid = abs(particle->PID); 235 235 float eta = fabs(particle->Eta); 236 236 237 237 //input generator level particle for jet algorithm 238 238 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())); 241 241 } 242 242 243 243 //Calculate ETMIS from generated particles 244 244 if((pid == pNU1) || (pid == pNU2) || (pid == pNU3))PTmisGEN = PTmisGEN + genMomentum; 245 245 246 246 if( (particle->Status == 1) && 247 247 ((pid != pNU1) && (pid != pNU2) && (pid != pNU3)) && … … 252 252 //use of the magnetic field propagation 253 253 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 256 258 switch(pid) { 257 259 258 260 case pE: // all electrons with eta < DET->MAX_CALO_FWD 259 261 DET->SmearElectron(recoMomentum); 260 if( eta < DET->MAX_TRACKER){261 262 263 262 if(recoMomentum.E() !=0 && eta < DET->MAX_TRACKER){ 263 elementElec=(RESOLELEC*) branchelec->NewEntry(); 264 elementElec->E = genMomentum.E(); 265 elementElec->SmearedE = recoMomentum.E();} 264 266 break; // case pE 265 267 case pGAMMA: // all photons with eta < DET->MAX_CALO_FWD … … 268 270 case pMU: // all muons with eta < DET->MAX_MU 269 271 DET->SmearMu(recoMomentum); 272 if(recoMomentum.E() !=0){ 270 273 elementMuon = (RESOLMUON*) branchmuon->NewEntry(); 271 274 elementMuon->OverPT = (1/genMomentum.Pt()); 272 elementMuon->OverSmearedPT = (1/recoMomentum.Pt()); 275 elementMuon->OverSmearedPT = (1/recoMomentum.Pt());} 273 276 break; // case pMU 274 277 case pLAMBDA: // all lambdas with eta < DET->MAX_CALO_FWD … … 280 283 break; 281 284 } // switch (pid) 282 283 284 285 286 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)){ 287 290 towers.push_back(PhysicsTower(LorentzVector(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E()))); 288 291 input_particlesReco.push_back(fastjet::PseudoJet(recoMomentum.Px(),recoMomentum.Py(),recoMomentum.Pz(), recoMomentum.E())); 289 292 } 290 293 } 291 294 292 295 // all final charged particles 293 294 295 296 297 298 299 300 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 302 305 } // switch 303 306 } // 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 313 318 elementEtmis= (ETMIS*) branchetmis->NewEntry(); 314 319 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 317 326 //***************************** 318 327 -
trunk/interface/TreeClasses.h
r88 r89 73 73 }; 74 74 75 //--------------------------------------------- 76 75 77 class ETMIS : public TSortableObject 76 78 { 77 79 public: 78 80 Float_t Et; 79 Float_t SmearedEt; 81 Float_t SEt; 82 Float_t Ex; 83 Float_t EtSmeare; 84 Float_t ExSmeare; 80 85 81 86 static TCompare *fgCompare; //! … … 86 91 87 92 }; 93 88 94 89 95 -
trunk/routines/resolutions.C
r60 r89 103 103 gROOT->Reset(); 104 104 105 TFile *f1 = new TFile(" PTMIS.root","read");105 TFile *f1 = new TFile("Etmis.root","read"); 106 106 TTree *Analyze = (TTree*)f1->Get("Analysis"); 107 107 … … 132 132 char tempMin[500]; 133 133 if(i==0)binMin=5; 134 sprintf(tempMin,"ElecEResol. NonSmeareE > %d",binMin);134 sprintf(tempMin,"ElecEResol.E > %d",binMin); 135 135 string mystringMin(tempMin); 136 136 char tempMax[500]; 137 sprintf(tempMax,"ElecEResol. NonSmeareE < %d",binMax);137 sprintf(tempMax,"ElecEResol.E < %d",binMax); 138 138 string mystringMax(tempMax); 139 139 char tempName[500]; 140 sprintf(tempName,"(ElecEResol. NonSmeareE-ElecEResol.SmeareE)>>hETSoverET%d",i);140 sprintf(tempName,"(ElecEResol.E-ElecEResol.SmearedE)>>hETSoverET%d",i); 141 141 string mystringName(tempName); 142 142 … … 189 189 gROOT->Reset(); 190 190 191 TFile *f1 = new TFile(" PTMIS.root","read");191 TFile *f1 = new TFile("Etmis.root","read"); 192 192 TTree *Analyze = (TTree*)f1->Get("Analysis"); 193 193 … … 233 233 gROOT->Reset(); 234 234 235 TFile *f1 = new TFile(" PTMIS.root","read");235 TFile *f1 = new TFile("Etmis.root","read"); 236 236 TTree *Analyze = (TTree*)f1->Get("Analysis"); 237 237 … … 311 311 void General() 312 312 { 313 JetResol();313 // JetResol(); 314 314 ElecResol(); 315 315 ETmisResol();
Note:
See TracChangeset
for help on using the changeset viewer.