|
|
|
#include <Calibration/HcalCalibAlgos/src/HcalIsoTrkAnalyzer.cc>
Implementation: <Notes on="" implementation>="">
Definition at line 101 of file HcalIsoTrkAnalyzer.cc.
| HcalIsoTrkAnalyzer::HcalIsoTrkAnalyzer | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 193 of file HcalIsoTrkAnalyzer.cc.
References allowMissingInputs_, associationConeSize_, AxB_, eLabel_, energyECALmip, energyMaxIso, energyMinIso, eventWeight, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), hbheLabel_, hoLabel_, TrackAssociatorParameters::loadParameters(), m_ebInstance, m_ecalLabel, m_eeInstance, m_hcalLabel, m_histoFlag, maxPNear, MinNTECHitsEndcap, MinNTrackHitsBarrel, nIterations, outputFileName_, python::trackProbabilityAnalysis_cff::parameters, parameters_, trackAssociator_, trackLabel1_, trackLabel_, and TrackDetectorAssociator::useDefaultPropagator().
00194 { 00195 //now do what ever initialization is needed 00196 00197 m_ecalLabel = iConfig.getUntrackedParameter<std::string> ("ecalRecHitsLabel","ecalRecHit"); 00198 m_ebInstance = iConfig.getUntrackedParameter<std::string> ("ebRecHitsInstance","EcalRecHitsEB"); 00199 m_eeInstance = iConfig.getUntrackedParameter<std::string> ("eeRecHitsInstance","EcalRecHitsEE"); 00200 m_hcalLabel = iConfig.getUntrackedParameter<std::string> ("hcalRecHitsLabel","hbhereco"); 00201 m_histoFlag = iConfig.getUntrackedParameter<int>("histoFlag",0); 00202 00203 hbheLabel_= iConfig.getParameter<edm::InputTag>("hbheInput"); 00204 hoLabel_=iConfig.getParameter<edm::InputTag>("hoInput"); 00205 eLabel_=iConfig.getParameter<edm::InputTag>("eInput"); 00206 trackLabel_ = iConfig.getParameter<edm::InputTag>("HcalIsolTrackInput"); 00207 trackLabel1_ = iConfig.getParameter<edm::InputTag>("trackInput"); 00208 associationConeSize_=iConfig.getParameter<double>("associationConeSize"); 00209 allowMissingInputs_=iConfig.getUntrackedParameter<bool>("allowMissingInputs",true); 00210 outputFileName_=iConfig.getParameter<std::string>("outputFileName"); 00211 00212 AxB_=iConfig.getParameter<std::string>("AxB"); 00213 00214 nIterations = iConfig.getParameter<int>("noOfIterations"); 00215 eventWeight = iConfig.getParameter<double>("eventWeight"); 00216 energyMinIso = iConfig.getParameter<double>("energyMinIso"); 00217 energyMaxIso = iConfig.getUntrackedParameter<double>("energyMaxIso",1000.); 00218 energyECALmip = iConfig.getParameter<double>("energyECALmip"); 00219 maxPNear = iConfig.getParameter<double>("maxPNear"); 00220 MinNTrackHitsBarrel = iConfig.getParameter<int>("MinNTrackHitsBarrel"); 00221 MinNTECHitsEndcap = iConfig.getParameter<int>("MinNTECHitsEndcap"); 00222 00223 edm::ParameterSet parameters = iConfig.getParameter<edm::ParameterSet>("TrackAssociatorParameters"); 00224 parameters_.loadParameters( parameters ); 00225 trackAssociator_.useDefaultPropagator(); 00226 00227 }
| HcalIsoTrkAnalyzer::~HcalIsoTrkAnalyzer | ( | ) |
| void HcalIsoTrkAnalyzer::analyze | ( | const edm::Event & | iEvent, | |
| const edm::EventSetup & | iSetup | |||
| ) | [private, virtual] |
Implements edm::EDAnalyzer.
Definition at line 237 of file HcalIsoTrkAnalyzer.cc.
References funct::abs(), allowMissingInputs_, TrackDetectorAssociator::associate(), associationConeSize_, AxB_, edm::SortedCollection< T, SORT >::begin(), cell, cellEnergy, cells, cells3x3, cellsPF, TrackAssociatorParameters::dREcal, TrackAssociatorParameters::dRHcal, e, ReconstructionGR_cff::ecal, TrackDetMatchInfo::EcalRecHits, ecRHen, ecRHeta, ecRHphi, eecal, ehcal, eLabel_, emEnergy, edm::SortedCollection< T, SORT >::end(), lat::endl(), energyECALmip, energyMaxIso, energyMinIso, EnergyVector, PV3DBase< T, PVType, FrameType >::eta(), edm::EventID::event(), EventIds, EventMatrix, eventNumber, exampleP4, MergeTrackCollections_cff::generalTracks, geo, edm::EventSetup::get(), edm::Event::getByLabel(), TrackDetectorAssociator::getFreeTrajectoryState(), CaloGeometry::getPosition(), hbheLabel_, HCAL3x3, HCAL5x5, HcalAxBxDepthEnergy, TrackDetMatchInfo::HcalRecHits, hcRHdepth, hcRHen, hcRHeta, hcRHieta, hcRHiphi, hcRHphi, hoLabel_, hoRHdepth, hoRHen, hoRHeta, hoRHieta, hoRHiphi, hoRHphi, i, edm::Event::id(), iEtaHit, info, input_to_L3, iPhiHit, it, k, m_histoFlag, MaxHhitEnergy, maxPNear, MinNTECHitsEndcap, MinNTrackHitsBarrel, nECRecHits, nHCRecHits, nHORecHits, numberOfCells, numbers, numbers2, parameters_, PV3DBase< T, PVType, FrameType >::phi(), edm::ESHandle< T >::product(), edm::Handle< T >::product(), ptrack, edm::EventID::run(), runNumber, rvert, funct::sqrt(), std, targetE, trackAssociator_, trackEta, trackLabel1_, trackLabel_, trackPhi, trackPt, tree, TrackDetMatchInfo::trkGlobPosAtEcal, TrackAssociatorParameters::useCalo, TrackAssociatorParameters::useEcal, TrackAssociatorParameters::useHcal, and TrackAssociatorParameters::useMuon.
00238 { 00239 using namespace edm; 00240 using namespace std; 00241 00242 vector<float> rawEnergyVec; 00243 vector<int> detiphi; 00244 vector<int> detieta; 00245 vector<int> i3i5; 00246 vector<HcalDetId> detidvec; 00247 float calEnergy; 00248 00249 edm::Handle<reco::TrackCollection> generalTracks; 00250 iEvent.getByLabel(trackLabel1_,generalTracks); 00251 00252 edm::Handle<reco::IsolatedPixelTrackCandidateCollection> trackCollection; 00253 iEvent.getByLabel(trackLabel_,trackCollection); 00254 //const reco::IsolatedPixelTrackCandidateCollection isoTrack = *(trackCollection.product()); 00255 // LogInfo("IsoTracks: ")<<" IsoTracks size "<<(isoTrack).size(); 00256 //cout << " IsoTracks size "<<(isoTrack).size() << endl; 00257 00258 edm::Handle<EcalRecHitCollection> ecal; 00259 iEvent.getByLabel(eLabel_,ecal); 00260 const EcalRecHitCollection Hitecal = *(ecal.product()); 00261 // LogInfo("ECAL: ")<<" Size of ECAl "<<(Hitecal).size(); 00262 // cout << " Size of ECAl "<<(Hitecal).size() << endl; 00263 00264 edm::Handle<HBHERecHitCollection> hbhe; 00265 iEvent.getByLabel(hbheLabel_,hbhe); 00266 const HBHERecHitCollection Hithbhe = *(hbhe.product()); 00267 // LogInfo("HBHE: ")<<" Size of HBHE "<<(Hithbhe).size(); 00268 00269 edm::ESHandle<CaloGeometry> pG; 00270 iSetup.get<CaloGeometryRecord>().get(pG); 00271 geo = pG.product(); 00272 00273 // rof 16.05.2008 start: include the possibility for recalibration (use "recalibrate" label for safety) 00274 /* 00275 edm::ESHandle <HcalRespCorrs> recalibCorrs; 00276 iSetup.get<HcalRespCorrsRcd>().get("recalibrate",recalibCorrs); 00277 const HcalRespCorrs* myRecalib = recalibCorrs.product(); 00278 */ 00279 // rof end 00280 00281 parameters_.useEcal = true; 00282 parameters_.useHcal = true; 00283 parameters_.useCalo = false; 00284 parameters_.useMuon = false; 00285 parameters_.dREcal = 0.5; 00286 parameters_.dRHcal = 0.6; 00287 00288 if (trackCollection->size()==0) return; 00289 00290 for (reco::TrackCollection::const_iterator trit=generalTracks->begin(); trit!=generalTracks->end(); trit++) 00291 { 00292 reco::IsolatedPixelTrackCandidateCollection::const_iterator isoMatched=trackCollection->begin(); 00293 bool matched=false; 00294 for (reco::IsolatedPixelTrackCandidateCollection::const_iterator it = trackCollection->begin(); it!=trackCollection->end(); it++) 00295 { 00296 if (floor(100000*trit->pt())==floor(100000*it->pt())) 00297 { 00298 isoMatched=it; 00299 matched=true; 00300 } 00301 } 00302 if (!matched) continue; 00303 00304 if (trit->hitPattern().numberOfValidHits()<MinNTrackHitsBarrel) continue; 00305 if (fabs(trit->eta())>1.47&&trit->hitPattern().numberOfValidStripTECHits()<MinNTECHitsEndcap) 00306 00307 ptrack = sqrt(trit->px()*trit->px()+trit->py()*trit->py()+trit->pz()*trit->pz()); 00308 calEnergy = sqrt(trit->px()*trit->px()+trit->py()*trit->py()+trit->pz()*trit->pz()+0.14*0.14); 00309 00310 trackPt = trit->pt(); 00311 trackEta = trit->eta(); 00312 trackPhi = trit->phi(); 00313 00314 00315 double corrHCAL = 1.; 00316 /* 00317 if (fabs(trackEta)<1.47) { 00318 00319 if (calEnergy < 5.) corrHCAL = 1.55; 00320 if (calEnergy > 5. && calEnergy < 9.) corrHCAL = 1.55 - 0.18*(calEnergy-5.)/4.; 00321 if (calEnergy > 9. && calEnergy < 20.) corrHCAL = 1.37 - 0.18*(calEnergy-9.)/11.; 00322 if (calEnergy > 20. && calEnergy < 30.) corrHCAL = 1.19 - 0.06*(calEnergy-20.)/10.; 00323 if (calEnergy > 30. && calEnergy < 50.) corrHCAL = 1.13 - 0.02*(calEnergy-30.)/20.; 00324 if (calEnergy > 50. && calEnergy < 100.) corrHCAL = 1.11 - 0.02*(calEnergy-50.)/50.; 00325 if (calEnergy > 100. && calEnergy < 1000.) corrHCAL = 1.09 - 0.09*(calEnergy-100.)/900.; 00326 00327 } 00328 00329 if (fabs(trackEta)>1.47) { 00330 00331 if (calEnergy < 5.) corrHCAL = 1.49; 00332 if (calEnergy > 5. && calEnergy < 9.) corrHCAL = 1.49 - 0.08*(calEnergy-5.)/4.; 00333 if (calEnergy > 9. && calEnergy < 20.) corrHCAL = 1.41- 0.15*(calEnergy-9.)/11.; 00334 if (calEnergy > 20. && calEnergy < 30.) corrHCAL = 1.26 - 0.07*(calEnergy-20.)/10.; 00335 if (calEnergy > 30. && calEnergy < 50.) corrHCAL = 1.19 - 0.04*(calEnergy-30.)/20.; 00336 if (calEnergy > 50. && calEnergy < 100.) corrHCAL = 1.15 - 0.03*(calEnergy-50.)/50.; 00337 if (calEnergy > 100. && calEnergy < 1000.) corrHCAL = 1.12 - 0.12*(calEnergy-100.)/900.; 00338 00339 } 00340 */ 00341 00342 // cout << endl << " ISO TRACK E = "<< calEnergy << " ETA = " << trackEta<< " PHI = " << trackPhi << " Correction " << corrHCAL<< endl; 00343 00344 rvert = sqrt(trit->vx()*trit->vx()+trit->vy()*trit->vy()+trit->vz()*trit->vz()); 00345 00346 //Associate track with a calorimeter 00347 TrackDetMatchInfo info = trackAssociator_.associate(iEvent, iSetup,trackAssociator_.getFreeTrajectoryState(iSetup, *trit),parameters_); 00348 //*(it->track().get()) 00349 double etaecal=info.trkGlobPosAtEcal.eta(); 00350 double phiecal=info.trkGlobPosAtEcal.phi(); 00351 00352 eecal=info.coneEnergy(parameters_.dREcal, TrackDetMatchInfo::EcalRecHits); 00353 ehcal=info.coneEnergy(parameters_.dRHcal, TrackDetMatchInfo::HcalRecHits); 00354 00355 double rmin = 0.07; 00356 if( fabs(etaecal) > 1.47 ) rmin = 0.07*(fabs(etaecal)-0.47)*1.2; 00357 if( fabs(etaecal) > 2.2 ) rmin = 0.07*(fabs(etaecal)-0.47)*1.4; 00358 00359 struct 00360 { 00361 int iphihitm; 00362 int ietahitm; 00363 int depthhit; 00364 double hitenergy; 00365 } MaxHit; 00366 00367 // Find Ecal RecHit with maximum energy and collect other information 00368 MaxHit.hitenergy=-100; 00369 nECRecHits=0; 00370 00371 double econus = 0.; 00372 float ecal_cluster = 0.; 00373 00374 for (std::vector<EcalRecHit>::const_iterator ehit=Hitecal.begin(); ehit!=Hitecal.end(); ehit++) 00375 { 00376 00377 if((*ehit).energy() > MaxHit.hitenergy) 00378 { 00379 MaxHit.hitenergy = (*ehit).energy(); 00380 } 00381 00382 GlobalPoint pos = geo->getPosition((*ehit).detid()); 00383 double phihit = pos.phi(); 00384 double etahit = pos.eta(); 00385 00386 double dphi = fabs(phiecal - phihit); 00387 if(dphi > 4.*atan(1.)) dphi = 8.*atan(1.) - dphi; 00388 double deta = fabs(etaecal - etahit); 00389 double dr = sqrt(dphi*dphi + deta*deta); 00390 00391 // cout << " eta "<< etahit << " phi "<< phihit << " en " << (*ehit).energy() << " dr "<< dr << endl; 00392 00393 if (dr < rmin) { 00394 econus = econus + (*ehit).energy(); 00395 } 00396 if (dr < 0.13) ecal_cluster += (*ehit).energy(); 00397 00398 ecRHen [nECRecHits] = (*ehit).energy(); 00399 ecRHeta[nECRecHits] = etahit; 00400 ecRHphi[nECRecHits] = phihit; 00401 nECRecHits++; 00402 } 00403 00404 MaxHit.hitenergy=-100; 00405 nHCRecHits=0; 00406 00407 float dddeta = 1000.; 00408 float dddphi = 1000.; 00409 int iphitrue = 1234; 00410 int ietatrue = 1234; 00411 00412 for (HBHERecHitCollection::const_iterator hhit=Hithbhe.begin(); hhit!=Hithbhe.end(); hhit++) 00413 { 00414 00415 // rof 16.05.2008 start: include the possibility for recalibration 00416 float recal = 1; 00417 // rof end 00418 00419 GlobalPoint pos = geo->getPosition(hhit->detid()); 00420 double phihit = pos.phi(); 00421 double etahit = pos.eta(); 00422 00423 int iphihitm = (hhit->id()).iphi(); 00424 int ietahitm = (hhit->id()).ieta(); 00425 int depthhit = (hhit->id()).depth(); 00426 00427 double dphi = fabs(info.trkGlobPosAtHcal.phi() - phihit); 00428 if(dphi > 4.*atan(1.)) dphi = 8.*atan(1.) - dphi; 00429 double deta = fabs(info.trkGlobPosAtHcal.eta() - etahit); 00430 double dr = sqrt(dphi*dphi + deta*deta); 00431 00432 if (deta<dddeta) { 00433 ietatrue = ietahitm; 00434 dddeta=deta; 00435 } 00436 00437 if (dphi<dddphi) { 00438 iphitrue = iphihitm; 00439 dddphi=dphi; 00440 } 00441 00442 if(dr<associationConeSize_) 00443 { 00444 hcRHen[nHCRecHits] = hhit->energy() * recal; 00445 hcRHeta[nHCRecHits] = etahit; 00446 hcRHphi[nHCRecHits] = phihit; 00447 hcRHieta[nHCRecHits] = ietahitm; 00448 hcRHiphi[nHCRecHits] = iphihitm; 00449 hcRHdepth[nHCRecHits] = depthhit; 00450 nHCRecHits++; 00451 00452 if(hhit->energy() * recal > MaxHit.hitenergy) 00453 { 00454 MaxHit.hitenergy = hhit->energy() * recal; 00455 MaxHit.ietahitm = (hhit->id()).ieta(); 00456 MaxHit.iphihitm = (hhit->id()).iphi(); 00457 MaxHit.depthhit = (hhit->id()).depth(); 00458 } 00459 } 00460 } 00461 MaxHhitEnergy = MaxHit.hitenergy; 00462 00463 if(m_histoFlag==1) 00464 { 00465 int MinIETA= 999; 00466 int MaxIETA= -999; 00467 int MinIPHI= 999; 00468 int MaxIPHI= -999; 00469 for (int k=0; k<nHCRecHits; k++) 00470 { 00471 00472 MinIETA = MinIETA > hcRHieta[k] ? hcRHieta[k] : MinIETA; 00473 MaxIETA = MaxIETA < hcRHieta[k] ? hcRHieta[k] : MaxIETA; 00474 MinIPHI = MinIPHI > hcRHiphi[k] ? hcRHiphi[k] : MinIPHI; 00475 MaxIPHI = MaxIPHI < hcRHiphi[k] ? hcRHiphi[k] : MaxIPHI; 00476 } 00477 } 00478 00479 if(AxB_=="5x5" || AxB_=="3x3") 00480 { 00481 HcalAxBxDepthEnergy=0.; 00482 for(int ih=0; ih<9; ih++){HCAL3x3[ih]=0.;} 00483 00484 for (HBHERecHitCollection::const_iterator hhit=Hithbhe.begin(); hhit!=Hithbhe.end(); hhit++) 00485 { 00486 00487 int DIETA = 100; 00488 if(MaxHit.ietahitm*(hhit->id()).ieta()>0) 00489 { 00490 DIETA = MaxHit.ietahitm - (hhit->id()).ieta(); 00491 } 00492 if(MaxHit.ietahitm*(hhit->id()).ieta()<0) 00493 { 00494 DIETA = MaxHit.ietahitm - (hhit->id()).ieta(); 00495 DIETA = DIETA>0 ? DIETA-1 : DIETA+1; 00496 } 00497 00498 int DIPHI = MaxHit.iphihitm - (hhit->id()).iphi(); 00499 DIPHI = DIPHI>36 ? 72-DIPHI : DIPHI; 00500 00501 int numbercell=0; 00502 if(AxB_=="3x3") numbercell = 1; 00503 if(AxB_=="5x5") numbercell = 2; 00504 00505 if( abs(DIETA)<=numbercell && abs(DIPHI)<=numbercell) { 00506 00507 // rof 16.05.2008 start: include the possibility for recalibration 00508 float recal = 1; 00509 00510 int iii3i5 = 0; 00511 00512 if(DIPHI==-1 && DIETA== 1) {HCAL3x3[0] += hhit->energy() * recal; iii3i5 = 1;} 00513 if(DIPHI==-1 && DIETA== 0) {HCAL3x3[1] += hhit->energy() * recal; iii3i5 = 1;} 00514 if(DIPHI==-1 && DIETA==-1) {HCAL3x3[2] += hhit->energy() * recal; iii3i5 = 1;} 00515 if(DIPHI== 0 && DIETA== 1) {HCAL3x3[3] += hhit->energy() * recal; iii3i5 = 1;} 00516 if(DIPHI== 0 && DIETA== 0) {HCAL3x3[4] += hhit->energy() * recal; iii3i5 = 1;} 00517 if(DIPHI== 0 && DIETA==-1) {HCAL3x3[5] += hhit->energy() * recal; iii3i5 = 1;} 00518 if(DIPHI== 1 && DIETA== 1) {HCAL3x3[6] += hhit->energy() * recal; iii3i5 = 1;} 00519 if(DIPHI== 1 && DIETA== 0) {HCAL3x3[7] += hhit->energy() * recal; iii3i5 = 1;} 00520 if(DIPHI== 1 && DIETA==-1) {HCAL3x3[8] += hhit->energy() * recal; iii3i5 = 1;} 00521 00522 if(DIPHI==-2 && DIETA== 2) {HCAL5x5[0] += hhit->energy() * recal;} 00523 if(DIPHI==-2 && DIETA== 1) {HCAL5x5[1] += hhit->energy() * recal;} 00524 if(DIPHI==-2 && DIETA== 0) {HCAL5x5[2] += hhit->energy() * recal;} 00525 if(DIPHI==-2 && DIETA==-1) {HCAL5x5[3] += hhit->energy() * recal;} 00526 if(DIPHI==-2 && DIETA==-2) {HCAL5x5[4] += hhit->energy() * recal;} 00527 00528 if(DIPHI==-1 && DIETA== 2) {HCAL5x5[5] += hhit->energy() * recal;} 00529 if(DIPHI==-1 && DIETA== 1) {HCAL5x5[6] += hhit->energy() * recal;} 00530 if(DIPHI==-1 && DIETA== 0) {HCAL5x5[7] += hhit->energy() * recal;} 00531 if(DIPHI==-1 && DIETA==-1) {HCAL5x5[8] += hhit->energy() * recal;} 00532 if(DIPHI==-1 && DIETA==-2) {HCAL5x5[9] += hhit->energy() * recal;} 00533 00534 if(DIPHI== 0 && DIETA== 2) {HCAL5x5[10] += hhit->energy() * recal;} 00535 if(DIPHI== 0 && DIETA== 1) {HCAL5x5[11] += hhit->energy() * recal;} 00536 if(DIPHI== 0 && DIETA== 0) {HCAL5x5[12] += hhit->energy() * recal;} 00537 if(DIPHI== 0 && DIETA==-1) {HCAL5x5[13] += hhit->energy() * recal;} 00538 if(DIPHI== 0 && DIETA==-2) {HCAL5x5[14] += hhit->energy() * recal;} 00539 00540 if(DIPHI== 1 && DIETA== 2) {HCAL5x5[15] += hhit->energy() * recal;} 00541 if(DIPHI== 1 && DIETA== 1) {HCAL5x5[16] += hhit->energy() * recal;} 00542 if(DIPHI== 1 && DIETA== 0) {HCAL5x5[17] += hhit->energy() * recal;} 00543 if(DIPHI== 1 && DIETA==-1) {HCAL5x5[18] += hhit->energy() * recal;} 00544 if(DIPHI== 1 && DIETA==-2) {HCAL5x5[19] += hhit->energy() * recal;} 00545 00546 if(DIPHI== 2 && DIETA== 2) {HCAL5x5[20] += hhit->energy() * recal;} 00547 if(DIPHI== 2 && DIETA== 1) {HCAL5x5[21] += hhit->energy() * recal;} 00548 if(DIPHI== 2 && DIETA== 0) {HCAL5x5[22] += hhit->energy() * recal;} 00549 if(DIPHI== 2 && DIETA==-1) {HCAL5x5[23] += hhit->energy() * recal;} 00550 if(DIPHI== 2 && DIETA==-2) {HCAL5x5[24] += hhit->energy() * recal;} 00551 00552 if(calEnergy > energyMinIso && calEnergy < energyMaxIso 00553 && isoMatched->energyIn() < energyECALmip && isoMatched->maxPtPxl() < maxPNear 00554 && MaxHit.hitenergy > 0.){ 00555 rawEnergyVec.push_back(hhit->energy() * recal * corrHCAL); 00556 detidvec.push_back(hhit->id()); 00557 detiphi.push_back((hhit->id()).iphi()); 00558 detieta.push_back((hhit->id()).ieta()); 00559 i3i5.push_back(iii3i5); 00560 numbers2[(hhit->id()).ieta()+29][(hhit->id()).iphi()] = numbers2[(hhit->id()).ieta()+29][(hhit->id()).iphi()] + 1; 00561 00562 } 00563 } 00564 } 00565 } 00566 00567 if(AxB_!="3x3" && AxB_!="5x5") LogWarning(" AxB ")<<" Not supported: "<< AxB_; 00568 00569 if(calEnergy > energyMinIso && calEnergy < energyMaxIso && 00570 isoMatched->energyIn() < energyECALmip && isoMatched->maxPtPxl() < maxPNear && 00571 (MaxHit.ietahitm+29) > -1 && (MaxHit.ietahitm+29) < 60 && MaxHit.hitenergy > 0.) 00572 { 00573 EventMatrix.push_back(rawEnergyVec); 00574 EventIds.push_back(detidvec); 00575 EnergyVector.push_back(calEnergy); 00576 00577 numbers[MaxHit.ietahitm+29][MaxHit.iphihitm] = numbers[MaxHit.ietahitm+29][MaxHit.iphihitm] + 1; 00578 } 00579 00580 if (calEnergy > 10. && isoMatched->energyIn() < energyECALmip && isoMatched->maxPtPxl() < maxPNear && 00581 (MaxHit.ietahitm+29) > -1 && (MaxHit.ietahitm+29) < 60 && MaxHit.hitenergy > 0.){ 00582 // input_to_L3 << rawEnergyVec.size() << " " << calEnergy << " " << ietatrue << " " << iphitrue; 00583 input_to_L3 << rawEnergyVec.size() << " " << calEnergy; 00584 00585 for (unsigned int i=0; i<rawEnergyVec.size(); i++) 00586 { 00587 // input_to_L3 << " " << rawEnergyVec.at(i) << " " << detidvec.at(i).rawId() << " " << detiphi.at(i) << " " << detieta.at(i); 00588 input_to_L3 << " " << rawEnergyVec.at(i) << " " << detidvec.at(i).rawId() ; 00589 } 00590 input_to_L3 <<endl; 00591 00592 // set dummy values for - not in ascii file 00593 eventNumber = iEvent.id().event(); 00594 runNumber = iEvent.id().run(); 00595 iEtaHit = ietatrue; 00596 iPhiHit = iphitrue; 00597 emEnergy = econus; 00598 exampleP4->SetPxPyPzE(2, 1, 1, 10); 00599 00600 numberOfCells=rawEnergyVec.size(); 00601 targetE = calEnergy; 00602 00603 for (unsigned int ia=0; ia<numberOfCells; ++ia) { 00604 cellEnergy = rawEnergyVec.at(ia); 00605 cell = detidvec.at(ia).rawId(); 00606 00607 new((*cells)[ia]) TCell(cell, cellEnergy); 00608 00609 if (i3i5.at(ia)==1) { 00610 cells3x3->Add((*cells)[ia]); 00611 // input_to_L3 << " 3x3 " << detiphi.at(ia) << " " <<detieta.at(ia) << endl; 00612 } 00613 if (i3i5.at(ia)==1) { 00614 cellsPF->Add((*cells)[ia]); 00615 } 00616 00617 } 00618 00619 tree->Fill(); 00620 00621 cells3x3->Clear(); 00622 cellsPF->Clear(); 00623 cells->Clear(); 00624 00625 } 00626 00627 rawEnergyVec.clear(); 00628 detidvec.clear(); 00629 detiphi.clear(); 00630 detieta.clear(); 00631 i3i5.clear(); 00632 00633 nHORecHits=0; 00634 try { 00635 Handle<HORecHitCollection> ho; 00636 iEvent.getByLabel(hoLabel_,ho); 00637 const HORecHitCollection Hitho = *(ho.product()); 00638 00639 for(HORecHitCollection::const_iterator hoItr=Hitho.begin(); hoItr!=Hitho.end(); hoItr++) 00640 { 00641 00642 // rof 16.05.2008 start: include the possibility for recalibration 00643 float recal = 1; 00644 // rof end 00645 00646 GlobalPoint pos = geo->getPosition(hoItr->detid()); 00647 double phihit = pos.phi(); 00648 double etahit = pos.eta(); 00649 00650 int iphihitm = (hoItr->id()).iphi(); 00651 int ietahitm = (hoItr->id()).ieta(); 00652 int depthhit = (hoItr->id()).depth(); 00653 00654 double dphi = fabs(trackPhi - phihit); 00655 if(dphi > 4.*atan(1.)) dphi = 8.*atan(1.) - dphi; 00656 double deta = fabs(trackEta - etahit); 00657 double dr = sqrt(dphi*dphi + deta*deta); 00658 00659 if(dr<associationConeSize_) { 00660 hoRHen[nHORecHits] = hoItr->energy() * recal; 00661 hoRHeta[nHORecHits] = etahit; 00662 hoRHphi[nHORecHits] = phihit; 00663 hoRHieta[nHORecHits] = ietahitm; 00664 hoRHiphi[nHORecHits] = iphihitm; 00665 hoRHdepth[nHORecHits] = depthhit; 00666 nHORecHits++; 00667 00668 } 00669 } 00670 } catch (cms::Exception& e) { // can't find it! 00671 if (!allowMissingInputs_) throw e; 00672 } 00673 } 00674 }
| void HcalIsoTrkAnalyzer::beginJob | ( | const edm::EventSetup & | ) | [private, virtual] |
Reimplemented from edm::EDAnalyzer.
Definition at line 678 of file HcalIsoTrkAnalyzer.cc.
References cells, cells3x3, cellsPF, eventNumber, exampleP4, iEtaHit, input_to_L3, iPhiHit, rootFile, runNumber, targetE, and tree.
00679 { 00680 // MyL3Algo = new MinL3AlgoUniv<HcalDetId>(eventWeight); 00681 00682 input_to_L3.open("input_to_l3.txt"); 00683 00684 rootFile = new TFile("rootFile.root", "RECREATE"); 00685 tree = new TTree("hcalCalibTree", "Tree for IsoTrack Calibration"); 00686 // Do this in BeginJob() ----------------------------------------- 00687 cells = new TClonesArray("TCell", 10000); 00688 cells3x3 = new TRefArray; 00689 cellsPF = new TRefArray; 00690 exampleP4 = new TLorentzVector(); 00691 00692 tree->Branch("cells", &cells, 64000, 0); // works also without the last two arguments 00693 tree->Branch("targetE", &targetE, "targetE/F"); 00694 tree->Branch("exampleP4", "TLorentzVector", &exampleP4); 00695 tree->Branch("cells3x3", &cells3x3, 64000, 0); 00696 tree->Branch("cellsPF", &cellsPF, 64000, 0); 00697 tree->Branch("iEtaHit", &iEtaHit, "iEtaHit/I"); 00698 tree->Branch("iPhiHit", &iPhiHit, "iPhiHit/i"); 00699 tree->Branch("eventNumber", &eventNumber, "eventNumber/i"); 00700 tree->Branch("runNumber", &runNumber, "runNumber/i"); 00701 //_________________________________________________________________ 00702 00703 }
Reimplemented from edm::EDAnalyzer.
Definition at line 707 of file HcalIsoTrkAnalyzer.cc.
References cells, cells3x3, cellsPF, exampleP4, input_to_L3, and rootFile.
00707 { 00708 00709 /* 00710 // M+ 00711 // perform the calibration with given number of iterations 00712 cout<<" Begin of solution "<< EnergyVector.size() << " "<< EventMatrix.size()<< " "<< EventIds.size()<<endl; 00713 solution = MyL3Algo->iterate(EventMatrix, EventIds, EnergyVector, nIterations); 00714 cout<<" End of solution "<<endl; 00715 00716 // print the solution and make a few plots 00717 map<HcalDetId,float>::iterator ii; 00718 for (ii = solution.begin(); ii != solution.end(); ii++) 00719 { 00720 int curr_eta = ii->first.ieta(); 00721 int curr_phi = ii->first.iphi(); 00722 int curr_depth = ii->first.depth(); 00723 cout << "solution[eta=" << curr_eta << ",phi=" << curr_phi << ",subdet=" << ii->first.subdet() 00724 << ",depth=" << curr_depth << "] = " << ii->second 00725 << " Stat " << numbers[curr_eta+29][curr_phi] << " " << numbers2[curr_eta+29][curr_phi] 00726 << endl; 00727 } 00728 00729 for (ii = solution.begin(); ii != solution.end(); ii++) 00730 { 00731 int curr_eta = ii->first.ieta(); 00732 int curr_phi = ii->first.iphi(); 00733 int curr_depth = ii->first.depth(); 00734 cout << " "<< ii->first.subdet() << " " << curr_depth << " "<< curr_eta << " " << curr_phi << " "<< ii->second 00735 << endl; 00736 } 00737 */ 00738 00739 input_to_L3.close(); 00740 00741 rootFile->Write(); 00742 rootFile->Close(); 00743 00744 if (cells) delete cells; 00745 if (cells3x3) delete cells3x3; 00746 if (cellsPF) delete cellsPF; 00747 if (exampleP4) delete exampleP4; 00748 }
bool HcalIsoTrkAnalyzer::allowMissingInputs_ [private] |
Definition at line 133 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::associationConeSize_ [private] |
Definition at line 131 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
string HcalIsoTrkAnalyzer::AxB_ [private] |
Definition at line 132 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
UInt_t HcalIsoTrkAnalyzer::cell [private] |
Float_t HcalIsoTrkAnalyzer::cellEnergy [private] |
TClonesArray* HcalIsoTrkAnalyzer::cells [private] |
Definition at line 177 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), beginJob(), and endJob().
TRefArray* HcalIsoTrkAnalyzer::cells3x3 [private] |
Definition at line 178 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), beginJob(), and endJob().
TRefArray* HcalIsoTrkAnalyzer::cellsPF [private] |
Definition at line 179 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), beginJob(), and endJob().
char HcalIsoTrkAnalyzer::dirname[50] [private] |
Definition at line 165 of file HcalIsoTrkAnalyzer.cc.
double HcalIsoTrkAnalyzer::ecRHen[1500] [private] |
double HcalIsoTrkAnalyzer::ecRHeta[1500] [private] |
double HcalIsoTrkAnalyzer::ecRHphi[1500] [private] |
double HcalIsoTrkAnalyzer::eecal [private] |
double HcalIsoTrkAnalyzer::ehcal [private] |
InputTag HcalIsoTrkAnalyzer::eLabel_ [private] |
Definition at line 120 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
Float_t HcalIsoTrkAnalyzer::emEnergy [private] |
double HcalIsoTrkAnalyzer::energyECALmip [private] |
Definition at line 163 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::energyMaxIso [private] |
Definition at line 162 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::energyMinIso [private] |
Definition at line 162 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
vector<float> HcalIsoTrkAnalyzer::EnergyVector [private] |
vector<vector<HcalDetId> > HcalIsoTrkAnalyzer::EventIds [private] |
vector<vector<float> > HcalIsoTrkAnalyzer::EventMatrix [private] |
UInt_t HcalIsoTrkAnalyzer::eventNumber [private] |
float HcalIsoTrkAnalyzer::eventWeight [private] |
TLorentzVector* HcalIsoTrkAnalyzer::exampleP4 [private] |
Definition at line 185 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), beginJob(), and endJob().
const CaloGeometry* HcalIsoTrkAnalyzer::geo [private] |
InputTag HcalIsoTrkAnalyzer::hbheLabel_ [private] |
Definition at line 118 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::HCAL3x3[9] [private] |
double HcalIsoTrkAnalyzer::HCAL5x5[25] [private] |
double HcalIsoTrkAnalyzer::HcalAxBxDepthEnergy [private] |
int HcalIsoTrkAnalyzer::hcRHdepth[1500] [private] |
double HcalIsoTrkAnalyzer::hcRHen[1500] [private] |
double HcalIsoTrkAnalyzer::hcRHeta[1500] [private] |
int HcalIsoTrkAnalyzer::hcRHieta[1500] [private] |
int HcalIsoTrkAnalyzer::hcRHiphi[1500] [private] |
double HcalIsoTrkAnalyzer::hcRHphi[1500] [private] |
char HcalIsoTrkAnalyzer::hname[20] [private] |
Definition at line 166 of file HcalIsoTrkAnalyzer.cc.
InputTag HcalIsoTrkAnalyzer::hoLabel_ [private] |
Definition at line 119 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
int HcalIsoTrkAnalyzer::hoRHdepth[1500] [private] |
double HcalIsoTrkAnalyzer::hoRHen[1500] [private] |
double HcalIsoTrkAnalyzer::hoRHeta[1500] [private] |
int HcalIsoTrkAnalyzer::hoRHieta[1500] [private] |
int HcalIsoTrkAnalyzer::hoRHiphi[1500] [private] |
double HcalIsoTrkAnalyzer::hoRHphi[1500] [private] |
char HcalIsoTrkAnalyzer::htitle[80] [private] |
Definition at line 167 of file HcalIsoTrkAnalyzer.cc.
Int_t HcalIsoTrkAnalyzer::iEtaHit [private] |
ofstream HcalIsoTrkAnalyzer::input_to_L3 [private] |
Definition at line 189 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), beginJob(), and endJob().
UInt_t HcalIsoTrkAnalyzer::iPhiHit [private] |
std::string HcalIsoTrkAnalyzer::m_ebInstance [private] |
std::string HcalIsoTrkAnalyzer::m_ecalLabel [private] |
std::string HcalIsoTrkAnalyzer::m_eeInstance [private] |
std::string HcalIsoTrkAnalyzer::m_hcalLabel [private] |
int HcalIsoTrkAnalyzer::m_histoFlag [private] |
Definition at line 129 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
std::string HcalIsoTrkAnalyzer::m_inputTrackLabel [private] |
Definition at line 124 of file HcalIsoTrkAnalyzer.cc.
double HcalIsoTrkAnalyzer::MaxHhitEnergy [private] |
double HcalIsoTrkAnalyzer::maxPNear [private] |
Definition at line 163 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
int HcalIsoTrkAnalyzer::MinNTECHitsEndcap [private] |
Definition at line 160 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
int HcalIsoTrkAnalyzer::MinNTrackHitsBarrel [private] |
Definition at line 160 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
int HcalIsoTrkAnalyzer::nECRecHits [private] |
int HcalIsoTrkAnalyzer::nHCRecHits [private] |
int HcalIsoTrkAnalyzer::nHORecHits [private] |
int HcalIsoTrkAnalyzer::nIterations [private] |
UInt_t HcalIsoTrkAnalyzer::numberOfCells [private] |
int HcalIsoTrkAnalyzer::numbers[60][72] [private] |
int HcalIsoTrkAnalyzer::numbers2[60][72] [private] |
string HcalIsoTrkAnalyzer::outputFileName_ [private] |
Definition at line 115 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::ptrack [private] |
TFile* HcalIsoTrkAnalyzer::rootFile [private] |
UInt_t HcalIsoTrkAnalyzer::runNumber [private] |
double HcalIsoTrkAnalyzer::rvert [private] |
map<HcalDetId,float> HcalIsoTrkAnalyzer::solution [private] |
Definition at line 159 of file HcalIsoTrkAnalyzer.cc.
Float_t HcalIsoTrkAnalyzer::targetE [private] |
Definition at line 114 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::trackEta [private] |
InputTag HcalIsoTrkAnalyzer::trackLabel1_ [private] |
Definition at line 122 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
InputTag HcalIsoTrkAnalyzer::trackLabel_ [private] |
Definition at line 121 of file HcalIsoTrkAnalyzer.cc.
Referenced by analyze(), and HcalIsoTrkAnalyzer().
double HcalIsoTrkAnalyzer::trackPhi [private] |
double HcalIsoTrkAnalyzer::trackPt [private] |
TTree* HcalIsoTrkAnalyzer::tree [private] |
1.5.4