Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:34:14

0001 // -*- C++ -*-
0002 //
0003 // Package:    CommonTools/RecoAlgos
0004 // Class:      PFClusterMatchedToPhotonsSelector
0005 //
0006 /**\class PFClusterMatchedToPhotonsSelector PFClusterMatchedToPhotonsSelector.cc CommonTools/RecoAlgos/plugins/PFClusterMatchedToPhotonsSelector.cc
0007 
0008  Description: Matches ECAL PF clusters to photons that do not convert
0009 
0010  Implementation:
0011 
0012 */
0013 //
0014 // Original Author:  RCLSA
0015 //         Created:  Wed, 22 Mar 2017 18:01:40 GMT
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 #include <iostream>
0022 
0023 // user include files
0024 #include "FWCore/Framework/interface/Frameworkfwd.h"
0025 #include "FWCore/Framework/interface/stream/EDProducer.h"
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0029 
0030 #include "DataFormats/Common/interface/ValueMap.h"
0031 #include "DataFormats/Common/interface/Ref.h"
0032 #include "DataFormats/Math/interface/deltaR.h"
0033 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0034 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
0035 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0036 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
0037 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0038 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
0039 
0040 typedef reco::PFCluster::EEtoPSAssociation::value_type EEPSPair;
0041 bool sortByKey(const EEPSPair& a, const EEPSPair& b) { return a.first < b.first; }
0042 
0043 class PFClusterMatchedToPhotonsSelector : public edm::stream::EDProducer<> {
0044 public:
0045   PFClusterMatchedToPhotonsSelector(const edm::ParameterSet&);
0046   static void fillDescriptions(edm::ConfigurationDescriptions&);
0047 
0048 private:
0049   void produce(edm::Event&, const edm::EventSetup&) override;
0050 
0051   edm::EDGetTokenT<reco::GenParticleCollection> genParticleToken_;  // genParticles
0052   edm::EDGetTokenT<reco::PFClusterCollection> particleFlowClusterECALToken_;
0053   edm::EDGetTokenT<reco::PFCluster::EEtoPSAssociation> associationToken_;
0054   edm::EDGetTokenT<TrackingParticleCollection> trackingParticleToken_;
0055   double matchMaxDR2_;
0056   double matchMaxDEDR2_;
0057 
0058   double volumeZ_EB_;
0059   double volumeRadius_EB_;
0060   double volumeZ_EE_;
0061 };
0062 
0063 PFClusterMatchedToPhotonsSelector::PFClusterMatchedToPhotonsSelector(const edm::ParameterSet& iConfig) {
0064   //now do what ever initialization is needed
0065   particleFlowClusterECALToken_ =
0066       consumes<reco::PFClusterCollection>(iConfig.getParameter<edm::InputTag>("pfClustersTag"));
0067   associationToken_ =
0068       consumes<reco::PFCluster::EEtoPSAssociation>(iConfig.getParameter<edm::InputTag>("pfClustersTag"));
0069   trackingParticleToken_ =
0070       consumes<TrackingParticleCollection>(iConfig.getParameter<edm::InputTag>("trackingParticleTag"));
0071   genParticleToken_ = consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("genParticleTag"));
0072 
0073   matchMaxDR2_ = iConfig.getParameter<double>("maxDR2");
0074   matchMaxDEDR2_ = iConfig.getParameter<double>("maxDEDR2");
0075   volumeZ_EB_ = iConfig.getParameter<double>("volumeZ_EB");
0076   volumeRadius_EB_ = iConfig.getParameter<double>("volumeRadius_EB");
0077   volumeZ_EE_ = iConfig.getParameter<double>("volumeZ_EE");
0078 
0079   produces<reco::PFClusterCollection>();
0080   produces<reco::PFCluster::EEtoPSAssociation>();
0081   produces<edm::ValueMap<reco::GenParticleRef> >();
0082 }
0083 
0084 void PFClusterMatchedToPhotonsSelector::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0085   edm::ParameterSetDescription desc;
0086   desc.add<edm::InputTag>("pfClustersTag", edm::InputTag("particleFlowClusterECAL"));
0087   desc.add<edm::InputTag>("trackingParticleTag", edm::InputTag("mix", "MergedTrackTruth"));
0088   desc.add<edm::InputTag>("genParticleTag", edm::InputTag("genParticles"));
0089   desc.add<double>("maxDR2", 0.1 * 0.1);
0090   desc.add<double>("maxDEDR2", 0.5 * 0.5);
0091   desc.add<double>("volumeZ_EB", 304.5);
0092   desc.add<double>("volumeRadius_EB", 123.8);
0093   desc.add<double>("volumeZ_EE", 317.0);
0094   descriptions.add("pfClusterMatchedToPhotonsSelector", desc);
0095 }
0096 
0097 void PFClusterMatchedToPhotonsSelector::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0098   edm::Handle<reco::PFClusterCollection> particleFlowClusterECALHandle_;
0099   edm::Handle<reco::PFCluster::EEtoPSAssociation> associationHandle_;
0100   edm::Handle<TrackingParticleCollection> trackingParticleHandle_;
0101   edm::Handle<reco::GenParticleCollection> genParticleHandle_;
0102   iEvent.getByToken(particleFlowClusterECALToken_, particleFlowClusterECALHandle_);
0103   iEvent.getByToken(trackingParticleToken_, trackingParticleHandle_);
0104   iEvent.getByToken(genParticleToken_, genParticleHandle_);
0105   iEvent.getByToken(associationToken_, associationHandle_);
0106 
0107   std::unique_ptr<reco::PFClusterCollection> out = std::make_unique<reco::PFClusterCollection>();
0108   std::unique_ptr<reco::PFCluster::EEtoPSAssociation> association_out =
0109       std::make_unique<reco::PFCluster::EEtoPSAssociation>();
0110   std::unique_ptr<edm::ValueMap<reco::GenParticleRef> > genmatching_out =
0111       std::make_unique<edm::ValueMap<reco::GenParticleRef> >();
0112 
0113   std::vector<reco::GenParticleRef> genmatching;
0114 
0115   size_t iN(0);
0116   for (size_t iP = 0; iP < particleFlowClusterECALHandle_->size(); iP++) {
0117     auto&& pfCluster = particleFlowClusterECALHandle_->at(iP);
0118     bool isMatched = false;
0119     reco::GenParticleRef::key_type matchedKey;
0120 
0121     // Preselect PFclusters
0122     if (pfCluster.energy() <= 0) {
0123       continue;
0124     }
0125 
0126     for (auto&& trackingParticle : *trackingParticleHandle_) {
0127       if (trackingParticle.pdgId() != 22)
0128         continue;
0129       if (trackingParticle.status() != 1)
0130         continue;
0131       matchedKey = trackingParticle.genParticles().at(0).key();
0132       float dR2 = reco::deltaR2(trackingParticle, pfCluster.position());
0133       if (dR2 > matchMaxDR2_)
0134         continue;
0135       float dE = 1. - trackingParticle.genParticles().at(0)->energy() / pfCluster.energy();
0136       if ((dR2 + dE * dE) > matchMaxDEDR2_)
0137         continue;
0138 
0139       bool isConversion = false;
0140       for (auto&& vertRef : trackingParticle.decayVertices()) {
0141         if (vertRef->position().rho() > volumeRadius_EB_ && std::abs(vertRef->position().z()) < volumeZ_EB_)
0142           continue;
0143         if (std::abs(vertRef->position().z()) > volumeZ_EE_)
0144           continue;
0145 
0146         for (auto&& tpRef : vertRef->daughterTracks()) {
0147           if (std::abs(tpRef->pdgId()) == 11)
0148             isConversion = true;
0149           break;
0150         }
0151         if (isConversion)
0152           break;
0153       }
0154       if (isConversion)
0155         continue;
0156 
0157       isMatched = true;
0158       break;
0159     }
0160 
0161     if (isMatched) {
0162       out->push_back(pfCluster);
0163       for (size_t i = 0; i < associationHandle_.product()->size(); i++) {
0164         if (associationHandle_.product()->at(i).first == iP) {
0165           association_out->push_back(std::make_pair(iN, associationHandle_.product()->at(i).second));
0166         }
0167       }
0168       genmatching.push_back(edm::Ref<reco::GenParticleCollection>(genParticleHandle_, matchedKey));
0169     }
0170   }
0171 
0172   std::sort(association_out->begin(), association_out->end(), sortByKey);
0173   edm::OrphanHandle<reco::PFClusterCollection> pfClusterHandle = iEvent.put(std::move(out));
0174   iEvent.put(std::move(association_out));
0175 
0176   edm::ValueMap<reco::GenParticleRef>::Filler mapFiller(*genmatching_out);
0177   mapFiller.insert(pfClusterHandle, genmatching.begin(), genmatching.end());
0178   mapFiller.fill();
0179   iEvent.put(std::move(genmatching_out));
0180 }
0181 
0182 //define this as a plug-in
0183 DEFINE_FWK_MODULE(PFClusterMatchedToPhotonsSelector);