#include <RecoParticleFlow/PFBlockProducer/interface/PFSimParticleProducer.h>
Public Types | |
typedef edm::Handle < reco::PFRecTrackCollection > | TrackHandle |
Public Member Functions | |
virtual void | beginJob (const edm::EventSetup &c) |
void | getSimIDs (const TrackHandle &trackh, std::vector< unsigned > &recTrackSimID) |
PFSimParticleProducer (const edm::ParameterSet &) | |
virtual void | produce (edm::Event &, const edm::EventSetup &) |
~PFSimParticleProducer () | |
Private Attributes | |
edm::InputTag | inputTagEcalRecHitsEB_ |
edm::InputTag | inputTagEcalRecHitsEE_ |
edm::InputTag | inputTagFamosSimHits_ |
edm::InputTag | inputTagRecTracks_ |
edm::InputTag | inputTagSim_ |
module label for retrieving input simtrack and simvertex | |
bool | mctruthMatchingInfo_ |
FSimEvent * | mySimEvent |
edm::ParameterSet | particleFilter_ |
bool | processParticles_ |
process particles on/off | |
bool | verbose_ |
verbose ? |
Definition at line 38 of file PFSimParticleProducer.h.
Definition at line 49 of file PFSimParticleProducer.h.
PFSimParticleProducer::PFSimParticleProducer | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 67 of file PFSimParticleProducer.cc.
References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), inputTagEcalRecHitsEB_, inputTagEcalRecHitsEE_, inputTagFamosSimHits_, inputTagRecTracks_, inputTagSim_, mctruthMatchingInfo_, mySimEvent, particleFilter_, processParticles_, and verbose_.
00068 { 00069 00070 00071 processParticles_ = 00072 iConfig.getUntrackedParameter<bool>("process_Particles",true); 00073 00074 00075 inputTagSim_ 00076 = iConfig.getParameter<InputTag>("sim"); 00077 00078 //retrieving collections for MC Truth Matching 00079 00080 //modif-beg 00081 inputTagFamosSimHits_ 00082 = iConfig.getUntrackedParameter<InputTag>("famosSimHits"); 00083 mctruthMatchingInfo_ = 00084 iConfig.getUntrackedParameter<bool>("MCTruthMatchingInfo",false); 00085 //modif-end 00086 00087 inputTagRecTracks_ 00088 = iConfig.getParameter<InputTag>("RecTracks"); 00089 inputTagEcalRecHitsEB_ 00090 = iConfig.getParameter<InputTag>("ecalRecHitsEB"); 00091 inputTagEcalRecHitsEE_ 00092 = iConfig.getParameter<InputTag>("ecalRecHitsEE"); 00093 00094 verbose_ = 00095 iConfig.getUntrackedParameter<bool>("verbose",false); 00096 00097 00098 // register products 00099 produces<reco::PFSimParticleCollection>(); 00100 00101 particleFilter_ = iConfig.getParameter<ParameterSet> 00102 ( "ParticleFilter" ); 00103 00104 mySimEvent = new FSimEvent( particleFilter_ ); 00105 }
PFSimParticleProducer::~PFSimParticleProducer | ( | ) |
Definition at line 109 of file PFSimParticleProducer.cc.
References mySimEvent.
00110 { 00111 delete mySimEvent; 00112 }
void PFSimParticleProducer::beginJob | ( | const edm::EventSetup & | c | ) | [virtual] |
Reimplemented from edm::EDProducer.
Definition at line 116 of file PFSimParticleProducer.cc.
References edm::EventSetup::getData(), FBaseSimEvent::initializePdt(), ParticleTable::instance(), and mySimEvent.
00117 { 00118 00119 // init Particle data table (from Pythia) 00120 edm::ESHandle < HepPDT::ParticleDataTable > pdt; 00121 // edm::ESHandle < DefaultConfig::ParticleDataTable > pdt; 00122 es.getData(pdt); 00123 if ( !ParticleTable::instance() ) ParticleTable::instance(&(*pdt)); 00124 mySimEvent->initializePdt(&(*pdt)); 00125 }
void PFSimParticleProducer::getSimIDs | ( | const TrackHandle & | trackh, | |
std::vector< unsigned > & | recTrackSimID | |||
) |
Definition at line 520 of file PFSimParticleProducer.cc.
References i, edm::Handle< T >::isValid(), it, SiTrackerGSMatchedRecHit2D::simtrackId(), and reco::PFRecTrack::trackRef().
Referenced by produce().
00522 { 00523 00524 if( trackh.isValid() ) { 00525 // cout << "Size=" << trackh->size() << endl; 00526 00527 for(unsigned i=0;i<trackh->size(); i++) { 00528 00529 reco::PFRecTrackRef ref( trackh,i ); 00530 const reco::PFRecTrack& PFT = *ref; 00531 const reco::TrackRef trackref = PFT.trackRef(); 00532 00533 // double Pt = trackref->pt(); 00534 // double DPt = trackref->ptError(); 00535 // cout << " PFBlockProducer: PFrecTrack->Track Pt= " 00536 // << Pt << " DPt = " << DPt << endl; 00537 00538 trackingRecHit_iterator rhitbeg 00539 = trackref->recHitsBegin(); 00540 trackingRecHit_iterator rhitend 00541 = trackref->recHitsEnd(); 00542 for (trackingRecHit_iterator it = rhitbeg; 00543 it != rhitend; it++){ 00544 00545 if( it->get()->isValid() ){ 00546 00547 const SiTrackerGSMatchedRecHit2D * rechit 00548 = (const SiTrackerGSMatchedRecHit2D*) (it->get()); 00549 00550 // cout << "rechit" 00551 // << " corresponding simId " 00552 // << rechit->simtrackId() 00553 // << endl; 00554 00555 recTrackSimID.push_back( rechit->simtrackId() ); 00556 break; 00557 } 00558 }//loop track rechit 00559 }//loop recTracks 00560 }//track handle valid 00561 }
void PFSimParticleProducer::produce | ( | edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) | [virtual] |
Implements edm::EDProducer.
Definition at line 128 of file PFSimParticleProducer.cc.
References reco::PFTrajectoryPoint::BeamPipeOrEndVertex, FSimTrack::charge(), reco::PFTrajectoryPoint::ClosestApproach, FSimTrack::daughters(), dummy, reco::PFTrajectoryPoint::ECALEntrance, FSimTrack::ecalEntrance(), lat::endl(), FSimTrack::endVertex(), err, edm::EventID::event(), Exception, FSimEvent::fill(), edm::Event::getByLabel(), getSimIDs(), reco::PFTrajectoryPoint::HCALEntrance, FSimTrack::hcalEntrance(), i, edm::Event::id(), FSimTrack::id(), inputTagEcalRecHitsEB_, inputTagFamosSimHits_, inputTagRecTracks_, inputTagSim_, edm::Handle< T >::isValid(), it, FSimTrack::layer1Entrance(), FSimTrack::layer2Entrance(), LogDebug, mctruthMatchingInfo_, FSimTrack::momentum(), FSimTrack::mother(), mySimEvent, FSimTrack::noEndVertex(), FSimTrack::noMother(), FSimEvent::nTracks(), FSimTrack::onEcal(), FSimTrack::onHcal(), FSimTrack::onLayer1(), FSimTrack::onLayer2(), FSimVertex::position(), FBaseSimEvent::print(), processParticles_, edm::Handle< T >::product(), reco::PFTrajectoryPoint::PS1, reco::PFTrajectoryPoint::PS2, edm::Event::put(), edm::EventID::run(), FBaseSimEvent::track(), CoreSimTrack::type(), verbose_, FSimTrack::vertex(), RawParticle::x(), RawParticle::y(), and RawParticle::z().
00130 { 00131 00132 LogDebug("PFSimParticleProducer")<<"START event: "<<iEvent.id().event() 00133 <<" in run "<<iEvent.id().run()<<endl; 00134 00135 //MC Truth Matching only with Famos and UnFoldedMode option to true!! 00136 00137 //vector to store the trackIDs of rectracks corresponding 00138 //to the simulated particle. 00139 std::vector<unsigned> recTrackSimID; 00140 00141 //In order to know which simparticule contribute to 00142 //a given Ecal RecHit energy, we need to access 00143 //the PCAloHit from FastSim. 00144 00145 typedef std::pair<double, unsigned> hitSimID; 00146 typedef std::list< std::pair<double, unsigned> >::iterator ITM; 00147 std::vector< std::list <hitSimID> > caloHitsEBID(62000); 00148 std::vector<double> caloHitsEBTotE(62000,0.0); 00149 00150 if(mctruthMatchingInfo_){ 00151 00152 //getting the PCAloHit 00153 edm::Handle<edm::PCaloHitContainer> pcalohits; 00154 // bool found_phit 00155 // = iEvent.getByLabel("famosSimHits","EcalHitsEB", 00156 // pcalohits); 00157 //modif-beg 00158 bool found_phit 00159 = iEvent.getByLabel(inputTagFamosSimHits_,pcalohits); 00160 //modif-end 00161 00162 if(!found_phit) { 00163 ostringstream err; 00164 err<<"could not find pcaloHit "<<"famosSimHits:EcalHitsEB"; 00165 LogError("PFSimParticleProducer")<<err.str()<<endl; 00166 00167 throw cms::Exception( "MissingProduct", err.str()); 00168 } 00169 else { 00170 assert( pcalohits.isValid() ); 00171 00172 // cout << "PFSimParticleProducer: number of pcalohits=" 00173 // << pcalohits->size() << endl; 00174 00175 edm::PCaloHitContainer::const_iterator it 00176 = pcalohits.product()->begin(); 00177 edm::PCaloHitContainer::const_iterator itend 00178 = pcalohits.product()->end(); 00179 00180 //loop on the PCaloHit from FastSim Calorimetry 00181 for(;it!=itend;++it) 00182 { 00183 EBDetId detid(it->id()); 00184 00185 // cout << detid << " " << detid.rawId() 00186 // << " " << detid.hashedIndex() 00187 // << " " << it->energy() 00188 // << " " << it->id() 00189 // << " trackId=" 00190 // << it->geantTrackId() 00191 // << endl; 00192 00193 if(it->energy() > 0.0) { 00194 std::pair<double, unsigned> phitsimid 00195 = make_pair(it->energy(),it->geantTrackId()); 00196 caloHitsEBID[detid.hashedIndex()].push_back(phitsimid); 00197 caloHitsEBTotE[detid.hashedIndex()] 00198 += it->energy(); //summing pcalhit energy 00199 }//energy > 0 00200 00201 }//loop PcaloHits 00202 }//pcalohit handle access 00203 00204 //Retrieving the PFRecTrack collection for 00205 //Monte Carlo Truth Matching tool 00206 Handle< reco::PFRecTrackCollection > recTracks; 00207 try{ 00208 LogDebug("PFSimParticleProducer")<<"getting PFRecTracks"<<endl; 00209 iEvent.getByLabel(inputTagRecTracks_, recTracks); 00210 00211 } catch (cms::Exception& err) { 00212 LogError("PFSimParticleProducer")<<err 00213 <<" cannot get collection " 00214 <<"particleFlowBlock"<<":" 00215 <<"" 00216 <<endl; 00217 }//pfrectrack handle access 00218 00219 //getting the simID corresponding to 00220 //each of the PFRecTracks 00221 getSimIDs( recTracks, recTrackSimID ); 00222 00223 }//mctruthMatchingInfo_ //modif 00224 00225 // deal with true particles 00226 if( processParticles_) { 00227 00228 auto_ptr< reco::PFSimParticleCollection > 00229 pOutputPFSimParticleCollection(new reco::PFSimParticleCollection ); 00230 00231 Handle<vector<SimTrack> > simTracks; 00232 bool found = iEvent.getByLabel(inputTagSim_,simTracks); 00233 if(!found) { 00234 00235 ostringstream err; 00236 err<<"cannot find sim tracks: "<<inputTagSim_; 00237 LogError("PFSimParticleProducer")<<err.str()<<endl; 00238 00239 throw cms::Exception( "MissingProduct", err.str()); 00240 } 00241 00242 00243 00244 Handle<vector<SimVertex> > simVertices; 00245 found = iEvent.getByLabel(inputTagSim_,simVertices); 00246 if(!found) { 00247 LogError("PFSimParticleProducer") 00248 <<"cannot find sim vertices: "<<inputTagSim_<<endl; 00249 return; 00250 } 00251 00252 00253 00254 // for(unsigned it = 0; it<simTracks->size(); it++ ) { 00255 // cout<<"\t track "<< (*simTracks)[it]<<" " 00256 // <<(*simTracks)[it].momentum().vect().perp()<<" " 00257 // <<(*simTracks)[it].momentum().e()<<endl; 00258 // } 00259 00260 mySimEvent->fill( *simTracks, *simVertices ); 00261 00262 if(verbose_) 00263 mySimEvent->print(); 00264 00265 // const std::vector<FSimTrack>& fsimTracks = *(mySimEvent->tracks() ); 00266 for(unsigned i=0; i<mySimEvent->nTracks(); i++) { 00267 00268 const FSimTrack& fst = mySimEvent->track(i); 00269 00270 int motherId = -1; 00271 if( ! fst.noMother() ) // a mother exist 00272 motherId = fst.mother().id(); 00273 00274 //This is finding out the simID corresponding 00275 //to the recTrack 00276 // cout << "Particle " << i 00277 // << " " << fst.genpartIndex() 00278 // << " -------------------------------------" << endl; 00279 00280 //GETTING THE TRACK ID 00281 unsigned recTrackID = 99999; 00282 vector<unsigned> recHitContrib; //modif 00283 vector<double> recHitContribFrac; //modif 00284 00285 if(mctruthMatchingInfo_){ //modif 00286 00287 for(unsigned lo=0; lo<recTrackSimID.size(); 00288 lo++) { 00289 if( i == recTrackSimID[lo] ) { 00290 // cout << "Corresponding Rec Track " 00291 // << lo << endl; 00292 recTrackID = lo; 00293 }//match track 00294 }//loop rectrack 00295 // if( recTrackID == 99999 ) 00296 // cout << "Sim Track not reconstructed pT=" << 00297 // fst.momentum().pt() << endl; 00298 00299 // get the ecalBarrel rechits for MC truth matching tool 00300 edm::Handle<EcalRecHitCollection> rhcHandle; 00301 bool found = iEvent.getByLabel(inputTagEcalRecHitsEB_, 00302 rhcHandle); 00303 if(!found) { 00304 ostringstream err; 00305 err<<"could not find rechits "<< inputTagEcalRecHitsEB_; 00306 LogError("PFSimParticleProducer")<<err.str()<<endl; 00307 00308 throw cms::Exception( "MissingProduct", err.str()); 00309 } 00310 else { 00311 assert( rhcHandle.isValid() ); 00312 // cout << "PFSimParticleProducer: number of rechits=" 00313 // << rhcHandle->size() << endl; 00314 00315 EBRecHitCollection::const_iterator it_rh 00316 = rhcHandle.product()->begin(); 00317 EBRecHitCollection::const_iterator itend_rh 00318 = rhcHandle.product()->end(); 00319 00320 for(;it_rh!=itend_rh;++it_rh) 00321 { 00322 unsigned rhit_hi 00323 = EBDetId(it_rh->id()).hashedIndex(); 00324 EBDetId detid(it_rh->id()); 00325 // cout << detid << " " << detid.rawId() 00326 // << " " << detid.hashedIndex() 00327 // << " " << it_rh->energy() << endl; 00328 00329 ITM it_phit = caloHitsEBID[rhit_hi].begin(); 00330 ITM itend_phit = caloHitsEBID[rhit_hi].end(); 00331 for(;it_phit!=itend_phit;++it_phit) 00332 { 00333 if(i == it_phit->second) 00334 { 00335 //Alex (08/10/08) TO BE REMOVED, eliminating 00336 //duplicated rechits 00337 bool alreadyin = false; 00338 for( unsigned ihit = 0; ihit < recHitContrib.size(); 00339 ++ihit ) 00340 if(detid.rawId() == recHitContrib[ihit]) 00341 alreadyin = true; 00342 00343 if(!alreadyin){ 00344 double pcalofraction = 0.0; 00345 if(caloHitsEBTotE[rhit_hi] != 0.0) 00346 pcalofraction 00347 = (it_phit->first/caloHitsEBTotE[rhit_hi])*100.0; 00348 00349 //store info 00350 recHitContrib.push_back(it_rh->id()); 00351 recHitContribFrac.push_back(pcalofraction); 00352 }//selected rechits 00353 }//matching 00354 }//loop pcalohit 00355 00356 }//loop rechits 00357 00358 }//getting the rechits 00359 00360 }//mctruthMatchingInfo_ //modif 00361 00362 // cout << "This particule has " << recHitContrib.size() 00363 // << " rechit contribution" << endl; 00364 // for( unsigned ih = 0; ih < recHitContrib.size(); ++ih ) 00365 // cout << recHitContrib[ih] 00366 // << " f=" << recHitContribFrac[ih] << " "; 00367 // cout << endl; 00368 00369 reco::PFSimParticle particle( fst.charge(), 00370 fst.type(), 00371 fst.id(), 00372 motherId, 00373 fst.daughters(), 00374 recTrackID, 00375 recHitContrib, 00376 recHitContribFrac); 00377 00378 00379 const FSimVertex& originVtx = fst.vertex(); 00380 00381 math::XYZPoint posOrig( originVtx.position().x(), 00382 originVtx.position().y(), 00383 originVtx.position().z() ); 00384 00385 math::XYZTLorentzVector momOrig( fst.momentum().px(), 00386 fst.momentum().py(), 00387 fst.momentum().pz(), 00388 fst.momentum().e() ); 00389 reco::PFTrajectoryPoint 00390 pointOrig(-1, 00391 reco::PFTrajectoryPoint::ClosestApproach, 00392 posOrig, momOrig); 00393 00394 // point 0 is origin vertex 00395 particle.addPoint(pointOrig); 00396 00397 00398 if( ! fst.noEndVertex() ) { 00399 const FSimVertex& endVtx = fst.endVertex(); 00400 00401 math::XYZPoint posEnd( endVtx.position().x(), 00402 endVtx.position().y(), 00403 endVtx.position().z() ); 00404 // cout<<"end vertex : " 00405 // <<endVtx.position().x()<<" " 00406 // <<endVtx.position().y()<<endl; 00407 00408 math::XYZTLorentzVector momEnd; 00409 00410 reco::PFTrajectoryPoint 00411 pointEnd( -1, 00412 reco::PFTrajectoryPoint::BeamPipeOrEndVertex, 00413 posEnd, momEnd); 00414 00415 particle.addPoint(pointEnd); 00416 } 00417 else { // add a dummy point 00418 reco::PFTrajectoryPoint dummy; 00419 particle.addPoint(dummy); 00420 } 00421 00422 00423 if( fst.onLayer1() ) { // PS layer1 00424 const RawParticle& rp = fst.layer1Entrance(); 00425 00426 math::XYZPoint posLayer1( rp.x(), rp.y(), rp.z() ); 00427 math::XYZTLorentzVector momLayer1( rp.px(), rp.py(), rp.pz(), rp.e() ); 00428 reco::PFTrajectoryPoint layer1Pt(-1, reco::PFTrajectoryPoint::PS1, 00429 posLayer1, momLayer1); 00430 00431 particle.addPoint( layer1Pt ); 00432 00433 // extrapolate to cluster depth 00434 } 00435 else { // add a dummy point 00436 reco::PFTrajectoryPoint dummy; 00437 particle.addPoint(dummy); 00438 } 00439 00440 if( fst.onLayer2() ) { // PS layer2 00441 const RawParticle& rp = fst.layer2Entrance(); 00442 00443 math::XYZPoint posLayer2( rp.x(), rp.y(), rp.z() ); 00444 math::XYZTLorentzVector momLayer2( rp.px(), rp.py(), rp.pz(), rp.e() ); 00445 reco::PFTrajectoryPoint layer2Pt(-1, reco::PFTrajectoryPoint::PS2, 00446 posLayer2, momLayer2); 00447 00448 particle.addPoint( layer2Pt ); 00449 00450 // extrapolate to cluster depth 00451 } 00452 else { // add a dummy point 00453 reco::PFTrajectoryPoint dummy; 00454 particle.addPoint(dummy); 00455 } 00456 00457 if( fst.onEcal() ) { 00458 const RawParticle& rp = fst.ecalEntrance(); 00459 00460 math::XYZPoint posECAL( rp.x(), rp.y(), rp.z() ); 00461 math::XYZTLorentzVector momECAL( rp.px(), rp.py(), rp.pz(), rp.e() ); 00462 reco::PFTrajectoryPoint ecalPt(-1, 00463 reco::PFTrajectoryPoint::ECALEntrance, 00464 posECAL, momECAL); 00465 00466 particle.addPoint( ecalPt ); 00467 00468 // extrapolate to cluster depth 00469 } 00470 else { // add a dummy point 00471 reco::PFTrajectoryPoint dummy; 00472 particle.addPoint(dummy); 00473 } 00474 00475 // add a dummy point for ECAL Shower max 00476 reco::PFTrajectoryPoint dummy; 00477 particle.addPoint(dummy); 00478 00479 if( fst.onHcal() ) { 00480 00481 const RawParticle& rpin = fst.hcalEntrance(); 00482 00483 math::XYZPoint posHCALin( rpin.x(), rpin.y(), rpin.z() ); 00484 math::XYZTLorentzVector momHCALin( rpin.px(), rpin.py(), rpin.pz(), 00485 rpin.e() ); 00486 reco::PFTrajectoryPoint hcalPtin(-1, 00487 reco::PFTrajectoryPoint::HCALEntrance, 00488 posHCALin, momHCALin); 00489 00490 particle.addPoint( hcalPtin ); 00491 00492 00493 // const RawParticle& rpout = fst.hcalExit(); 00494 00495 // math::XYZPoint posHCALout( rpout.x(), rpout.y(), rpout.z() ); 00496 // math::XYZTLorentzVector momHCALout( rpout.px(), rpout.py(), rpout.pz(), 00497 // rpout.e() ); 00498 // reco::PFTrajectoryPoint 00499 // hcalPtout(0, reco::PFTrajectoryPoint::HCALEntrance, 00500 // posHCAL, momHCAL); 00501 00502 // particle.addPoint( hcalPtout ); 00503 } 00504 else { // add a dummy point 00505 reco::PFTrajectoryPoint dummy; 00506 particle.addPoint(dummy); 00507 } 00508 00509 pOutputPFSimParticleCollection->push_back( particle ); 00510 } 00511 00512 iEvent.put(pOutputPFSimParticleCollection); 00513 } 00514 00515 00516 LogDebug("PFSimParticleProducer")<<"STOP event: "<<iEvent.id().event() 00517 <<" in run "<<iEvent.id().run()<<endl; 00518 }
Definition at line 66 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
Definition at line 62 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
Definition at line 65 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
module label for retrieving input simtrack and simvertex
Definition at line 57 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
Definition at line 61 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
FSimEvent* PFSimParticleProducer::mySimEvent [private] |
Definition at line 72 of file PFSimParticleProducer.h.
Referenced by beginJob(), PFSimParticleProducer(), produce(), and ~PFSimParticleProducer().
bool PFSimParticleProducer::processParticles_ [private] |
process particles on/off
Definition at line 77 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().
bool PFSimParticleProducer::verbose_ [private] |
verbose ?
Definition at line 80 of file PFSimParticleProducer.h.
Referenced by PFSimParticleProducer(), and produce().