Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-08 05:47:14

0001 /*
0002 UncleanSCRecoveryProducer:
0003 ^^^^^^^^^^^^^^^^^^^^^^^^^^
0004 
0005 takes the collections of flagged clean and unclean only SC 
0006 (this is the output of UnifiedSCCollectionProducer) and
0007 recovers the original collection of unclean SC.
0008 
0009 18 Aug 2010
0010 Nikolaos Rompotis and Chris Seez  - Imperial College London
0011 many thanks to David Wardrope, Shahram Rahatlou and Federico Ferri
0012 */
0013 
0014 #include "DataFormats/CaloRecHit/interface/CaloCluster.h"
0015 #include "DataFormats/Common/interface/Handle.h"
0016 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0017 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
0018 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0019 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0020 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
0021 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0022 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
0023 #include "FWCore/Framework/interface/ESHandle.h"
0024 #include "FWCore/Framework/interface/Event.h"
0025 #include "FWCore/Framework/interface/EventSetup.h"
0026 #include "FWCore/Framework/interface/global/EDProducer.h"
0027 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0028 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0029 
0030 #include <iostream>
0031 #include <memory>
0032 #include <vector>
0033 
0034 class UncleanSCRecoveryProducer : public edm::global::EDProducer<> {
0035 public:
0036   UncleanSCRecoveryProducer(const edm::ParameterSet& ps);
0037 
0038   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0039 
0040 private:
0041   // the clean collection
0042   const edm::EDGetTokenT<reco::BasicClusterCollection> cleanBcCollection_;
0043   const edm::EDGetTokenT<reco::SuperClusterCollection> cleanScCollection_;
0044   // the uncleaned collection
0045   const edm::EDGetTokenT<reco::BasicClusterCollection> uncleanBcCollection_;
0046   const edm::EDGetTokenT<reco::SuperClusterCollection> uncleanScCollection_;
0047   // the names of the products to be produced:
0048   const std::string bcCollection_;
0049   const std::string scCollection_;
0050 };
0051 
0052 #include "FWCore/Framework/interface/MakerMacros.h"
0053 DEFINE_FWK_MODULE(UncleanSCRecoveryProducer);
0054 
0055 UncleanSCRecoveryProducer::UncleanSCRecoveryProducer(const edm::ParameterSet& ps)
0056     : cleanBcCollection_(consumes<reco::BasicClusterCollection>(ps.getParameter<edm::InputTag>("cleanBcCollection"))),
0057       cleanScCollection_(consumes<reco::SuperClusterCollection>(ps.getParameter<edm::InputTag>("cleanScCollection"))),
0058       uncleanBcCollection_(
0059           consumes<reco::BasicClusterCollection>(ps.getParameter<edm::InputTag>("uncleanBcCollection"))),
0060       uncleanScCollection_(
0061           consumes<reco::SuperClusterCollection>(ps.getParameter<edm::InputTag>("uncleanScCollection"))),
0062       bcCollection_(ps.getParameter<std::string>("bcCollection")),
0063       scCollection_(ps.getParameter<std::string>("scCollection")) {
0064   // the products:
0065   produces<reco::BasicClusterCollection>(bcCollection_);
0066   produces<reco::SuperClusterCollection>(scCollection_);
0067 }
0068 
0069 void UncleanSCRecoveryProducer::produce(edm::StreamID, edm::Event& evt, const edm::EventSetup& es) const {
0070   // __________________________________________________________________________
0071   //
0072   // cluster collections:
0073   edm::Handle<reco::BasicClusterCollection> pCleanBC;
0074   edm::Handle<reco::SuperClusterCollection> pCleanSC;
0075   //
0076   edm::Handle<reco::BasicClusterCollection> pUncleanBC;
0077   edm::Handle<reco::SuperClusterCollection> pUncleanSC;
0078   // clean collections ________________________________________________________
0079   evt.getByToken(cleanScCollection_, pCleanSC);
0080   const reco::SuperClusterCollection cleanSC = *(pCleanSC.product());
0081 
0082   // unclean collections ______________________________________________________
0083   evt.getByToken(uncleanBcCollection_, pUncleanBC);
0084   const reco::BasicClusterCollection uncleanBC = *(pUncleanBC.product());
0085   //
0086   evt.getByToken(uncleanScCollection_, pUncleanSC);
0087   const reco::SuperClusterCollection uncleanSC = *(pUncleanSC.product());
0088   int uncleanSize = pUncleanSC->size();
0089   int cleanSize = pCleanSC->size();
0090 
0091   LogTrace("EcalCleaning") << "Size of Clean Collection: " << cleanSize << ", uncleanSize: " << uncleanSize;
0092 
0093   // collections are all taken now ____________________________________________
0094   //
0095   // the collections to be produced ___________________________________________
0096   reco::BasicClusterCollection basicClusters;
0097   reco::SuperClusterCollection superClusters;
0098   //
0099   //
0100   // collect all the basic clusters of the SC that belong to the unclean
0101   // collection and put them into the basicClusters vector
0102   // keep the information of which SC they belong to
0103   //
0104   // loop over the unclean sc: all SC will join the new collection
0105   std::vector<std::pair<int, int> > basicClusterOwner;  // counting all
0106 
0107   std::vector<DetId> scUncleanSeedDetId;  // counting the unclean
0108   for (int isc = 0; isc < uncleanSize; ++isc) {
0109     const reco::SuperCluster unsc = uncleanSC[isc];
0110     scUncleanSeedDetId.push_back(unsc.seed()->seed());
0111     reco::CaloCluster_iterator bciter = unsc.clustersBegin();
0112     for (; bciter != unsc.clustersEnd(); ++bciter) {
0113       // the basic clusters
0114       basicClusters.push_back(**bciter);
0115       // index of the unclean SC
0116       basicClusterOwner.push_back(std::make_pair(isc, 0));
0117     }
0118   }
0119   // loop over the clean: only the ones which are in common with the unclean
0120   // are taken into account
0121 
0122   std::vector<DetId> scCleanSeedDetId;  // counting the clean
0123   std::vector<int> isToBeKept;
0124   for (int isc = 0; isc < cleanSize; ++isc) {
0125     reco::SuperClusterRef cscRef(pCleanSC, isc);
0126     scCleanSeedDetId.push_back(cscRef->seed()->seed());
0127     for (reco::CaloCluster_iterator bciter = cscRef->clustersBegin(); bciter != cscRef->clustersEnd(); ++bciter) {
0128       // the basic clusters
0129       basicClusters.push_back(**bciter);
0130       // index of the clean SC
0131       basicClusterOwner.push_back(std::make_pair(isc, 1));
0132     }
0133     if (cscRef->isInUnclean())
0134       isToBeKept.push_back(1);
0135     else
0136       isToBeKept.push_back(0);
0137   }
0138   //
0139   // now export the basic clusters into the event and get them back
0140   auto basicClusters_p = std::make_unique<reco::BasicClusterCollection>();
0141   basicClusters_p->assign(basicClusters.begin(), basicClusters.end());
0142   edm::OrphanHandle<reco::BasicClusterCollection> bccHandle = evt.put(std::move(basicClusters_p), bcCollection_);
0143   if (!(bccHandle.isValid())) {
0144     edm::LogWarning("MissingInput") << "could not handle the new BasicClusters!";
0145     return;
0146   }
0147   reco::BasicClusterCollection basicClustersProd = *bccHandle;
0148 
0149   LogTrace("EcalCleaning") << "Got the BasicClusters from the event again";
0150   int bcSize = bccHandle->size();
0151   //
0152   // now we can create the SC collection
0153   //
0154   // run over the unclean SC: all to be kept here
0155   for (int isc = 0; isc < uncleanSize; ++isc) {
0156     reco::CaloClusterPtrVector clusterPtrVector;
0157     // the seed is the basic cluster with the highest energy
0158     reco::CaloClusterPtr seed;
0159     for (int jbc = 0; jbc < bcSize; ++jbc) {
0160       std::pair<int, int> theBcOwner = basicClusterOwner[jbc];
0161       if (theBcOwner.first == isc && theBcOwner.second == 0) {
0162         reco::CaloClusterPtr currentClu = reco::CaloClusterPtr(bccHandle, jbc);
0163         clusterPtrVector.push_back(currentClu);
0164         if (scUncleanSeedDetId[isc] == currentClu->seed()) {
0165           seed = currentClu;
0166         }
0167       }
0168     }
0169     const reco::SuperCluster unsc = uncleanSC[isc];
0170     reco::SuperCluster newSC(unsc.energy(), unsc.position(), seed, clusterPtrVector);
0171     newSC.setFlags(reco::CaloCluster::uncleanOnly);
0172     superClusters.push_back(newSC);
0173   }
0174   // run over the clean SC: only those who are in common between the
0175   // clean and unclean collection are kept
0176   for (int isc = 0; isc < cleanSize; ++isc) {
0177     reco::SuperClusterRef cscRef(pCleanSC, isc);
0178     if (not cscRef->isInUnclean())
0179       continue;
0180     reco::CaloClusterPtrVector clusterPtrVector;
0181     // the seed is the basic cluster with the highest energy
0182     reco::CaloClusterPtr seed;
0183     for (int jbc = 0; jbc < bcSize; ++jbc) {
0184       std::pair<int, int> theBcOwner = basicClusterOwner[jbc];
0185       if (theBcOwner.first == isc && theBcOwner.second == 1) {
0186         reco::CaloClusterPtr currentClu = reco::CaloClusterPtr(bccHandle, jbc);
0187         clusterPtrVector.push_back(currentClu);
0188         if (scCleanSeedDetId[isc] == currentClu->seed()) {
0189           seed = currentClu;
0190         }
0191       }
0192     }
0193     reco::SuperCluster newSC(cscRef->energy(), cscRef->position(), seed, clusterPtrVector);
0194     newSC.setFlags(reco::CaloCluster::common);
0195     superClusters.push_back(newSC);
0196   }
0197 
0198   auto superClusters_p = std::make_unique<reco::SuperClusterCollection>();
0199   superClusters_p->assign(superClusters.begin(), superClusters.end());
0200 
0201   evt.put(std::move(superClusters_p), scCollection_);
0202 
0203   LogTrace("EcalCleaning") << "Clusters (Basic/Super) added to the Event! :-)";
0204 
0205   // ----- debugging ----------
0206   // print the new collection SC quantities
0207   // print out the clean collection SC
0208   LogTrace("EcalCleaning") << "Clean Collection SC ";
0209   for (int i = 0; i < cleanSize; ++i) {
0210     const reco::SuperCluster csc = cleanSC[i];
0211     LogTrace("EcalCleaning") << " >>> clean    #" << i << "; Energy: " << csc.energy() << " eta: " << csc.eta()
0212                              << " sc seed detid: " << csc.seed()->seed().rawId();
0213   }
0214   // the unclean SC
0215   LogTrace("EcalCleaning") << "Unclean Collection SC ";
0216   for (int i = 0; i < uncleanSize; ++i) {
0217     const reco::SuperCluster usc = uncleanSC[i];
0218     LogTrace("EcalCleaning") << " >>> unclean  #" << i << "; Energy: " << usc.energy() << " eta: " << usc.eta()
0219                              << " sc seed detid: " << usc.seed()->seed().rawId();
0220   }
0221   // the new collection
0222   LogTrace("EcalCleaning") << "The new SC clean collection with size " << superClusters.size();
0223   for (unsigned int i = 0; i < superClusters.size(); ++i) {
0224     const reco::SuperCluster nsc = superClusters[i];
0225     LogTrace("EcalCleaning") << " >>> newSC    #" << i << "; Energy: " << nsc.energy() << " eta: " << nsc.eta()
0226                              << " isClean=" << nsc.isInClean() << " isUnclean=" << nsc.isInUnclean()
0227                              << " sc seed detid: " << nsc.seed()->seed().rawId();
0228   }
0229 }