Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:26:37

0001 #include "FWCore/Framework/interface/Event.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "DataFormats/Common/interface/ValueMap.h"
0007 #include "RecoParticleFlow/PFProducer/test/PFSuperClusterReader.h"
0008 
0009 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0010 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0011 #include "DataFormats/EgammaReco/interface/PreshowerCluster.h"
0012 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0013 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0014 #include "DataFormats/ParticleFlowReco/interface/PFBlockElement.h"
0015 #include "DataFormats/ParticleFlowReco/interface/PFBlockFwd.h"
0016 #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
0017 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0018 
0019 PFSuperClusterReader::PFSuperClusterReader(const edm::ParameterSet& iConfig) {
0020   inputTagGSFTracks_ = iConfig.getParameter<edm::InputTag>("GSFTracks");
0021   inputTagValueMapSC_ = iConfig.getParameter<edm::InputTag>("SuperClusterRefMap");
0022   inputTagValueMapMVA_ = iConfig.getParameter<edm::InputTag>("MVAMap");
0023   inputTagPFCandidates_ = iConfig.getParameter<edm::InputTag>("PFCandidate");
0024 }
0025 
0026 PFSuperClusterReader::~PFSuperClusterReader() { ; }
0027 
0028 void PFSuperClusterReader::beginRun(edm::Run const&, edm::EventSetup const&) { ; }
0029 
0030 void PFSuperClusterReader::analyze(const edm::Event& iEvent, const edm::EventSetup& c) {
0031   edm::Handle<reco::GsfTrackCollection> gsfTracksH;
0032   bool found = iEvent.getByLabel(inputTagGSFTracks_, gsfTracksH);
0033   if (!found) {
0034     std::ostringstream err;
0035     err << " cannot get GsfTracks: " << inputTagGSFTracks_ << std::endl;
0036     edm::LogError("PFSuperClusterReader") << err.str();
0037     throw cms::Exception("MissingProduct", err.str());
0038   }
0039 
0040   edm::Handle<reco::PFCandidateCollection> pfCandidatesH;
0041   found = iEvent.getByLabel(inputTagPFCandidates_, pfCandidatesH);
0042   if (!found) {
0043     std::ostringstream err;
0044     err << " cannot get PFCandidates: " << inputTagPFCandidates_ << std::endl;
0045     edm::LogError("PFSuperClusterReader") << err.str();
0046     throw cms::Exception("MissingProduct", err.str());
0047   }
0048 
0049   edm::Handle<edm::ValueMap<reco::SuperClusterRef> > pfClusterTracksH;
0050   found = iEvent.getByLabel(inputTagValueMapSC_, pfClusterTracksH);
0051   if (!found) {
0052     std::ostringstream err;
0053     err << " cannot get SuperClusterRef Map: " << inputTagValueMapSC_ << std::endl;
0054     edm::LogError("PFSuperClusterReader") << err.str();
0055     throw cms::Exception("MissingProduct", err.str());
0056   }
0057 
0058   edm::Handle<edm::ValueMap<float> > pfMVAH;
0059   found = iEvent.getByLabel(inputTagValueMapMVA_, pfMVAH);
0060   if (!found) {
0061     std::ostringstream err;
0062     err << " cannot get MVA Map: " << inputTagValueMapSC_ << std::endl;
0063     edm::LogError("PFSuperClusterReader") << err.str();
0064     throw cms::Exception("MissingProduct", err.str());
0065   }
0066 
0067   const edm::ValueMap<reco::SuperClusterRef>& mySCValueMap = *pfClusterTracksH;
0068   const edm::ValueMap<float>& myMVAValueMap = *pfMVAH;
0069 
0070   unsigned ngsf = gsfTracksH->size();
0071   for (unsigned igsf = 0; igsf < ngsf; ++igsf) {
0072     reco::GsfTrackRef theTrackRef(gsfTracksH, igsf);
0073     if (mySCValueMap[theTrackRef].isNull()) {
0074       continue;
0075     }
0076     const reco::SuperCluster& mySuperCluster(*(mySCValueMap[theTrackRef]));
0077     float mva = myMVAValueMap[theTrackRef];
0078     float eta = mySuperCluster.position().eta();
0079     float et = mySuperCluster.energy() * sin(mySuperCluster.position().theta());
0080     std::cout << " Super Cluster energy, eta, Et , EtaWidth, PhiWidth " << mySuperCluster.energy() << " ";
0081     std::cout << eta << " " << et << " " << mySuperCluster.etaWidth() << " " << mySuperCluster.phiWidth() << std::endl;
0082 
0083     const reco::PFCandidate* myPFCandidate = findPFCandidate(pfCandidatesH.product(), theTrackRef);
0084 
0085     if (mySuperCluster.seed().isNull()) {
0086       continue;
0087     }
0088     std::cout << " Super Cluster seed energy " << mySuperCluster.seed()->energy() << std::endl;
0089     std::cout << " Preshower contribution " << mySuperCluster.preshowerEnergy() << std::endl;
0090     std::cout << " MVA value " << mva << std::endl;
0091     std::cout << " List of basic clusters " << std::endl;
0092     reco::CaloCluster_iterator it = mySuperCluster.clustersBegin();
0093     reco::CaloCluster_iterator it_end = mySuperCluster.clustersEnd();
0094     float etotbasic = 0;
0095     for (; it != it_end; ++it) {
0096       std::cout << " Basic cluster " << (*it)->energy() << std::endl;
0097       etotbasic += (*it)->energy();
0098     }
0099     it = mySuperCluster.preshowerClustersBegin();
0100     it_end = mySuperCluster.preshowerClustersEnd();
0101     for (; it != it_end; ++it) {
0102       std::cout << " Preshower cluster " << (*it)->energy() << std::endl;
0103     }
0104 
0105     std::cout << " Comparison with PFCandidate : Energy " << myPFCandidate->ecalEnergy()
0106               << " SC : " << mySuperCluster.energy() << std::endl;
0107     std::cout << " Sum of Basic clusters :" << etotbasic;
0108     std::cout << " Calibrated preshower energy : " << mySuperCluster.preshowerEnergy();
0109     etotbasic += mySuperCluster.preshowerEnergy();
0110     std::cout << " Basic Clusters + PS :" << etotbasic << std::endl;
0111     std::cout << " Summary " << mySuperCluster.energy() << " " << eta << " " << et;
0112     std::cout << " " << mySuperCluster.preshowerEnergy() << " " << mva << std::endl;
0113     std::cout << std::endl;
0114   }
0115 }
0116 
0117 const reco::PFCandidate* PFSuperClusterReader::findPFCandidate(const reco::PFCandidateCollection* coll,
0118                                                                const reco::GsfTrackRef& ref) {
0119   const reco::PFCandidate* result = 0;
0120   unsigned ncand = coll->size();
0121   for (unsigned icand = 0; icand < ncand; ++icand) {
0122     if (!(*coll)[icand].gsfTrackRef().isNull() && (*coll)[icand].gsfTrackRef() == ref) {
0123       result = &((*coll)[icand]);
0124       return result;
0125     }
0126   }
0127   return result;
0128 }
0129 
0130 DEFINE_FWK_MODULE(PFSuperClusterReader);