AlCaECALRecHitReducer.cc

Go to the documentation of this file.

#include "Calibration/EcalAlCaRecoProducers/plugins/AlCaECALRecHitReducer.h"
//#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
#include "DataFormats/EgammaCandidates/interface/Photon.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "DataFormats/EgammaReco/interface/BasicCluster.h"
#include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
#include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
#include "FWCore/Framework/interface/ESHandle.h"

//#define ALLrecHits
//#define DEBUG

//#define QUICK -> if commented loop over the recHits of the SC and add them to the list of recHits to be saved
//                 comment it if you want a faster module but be sure the window is large enough


AlCaECALRecHitReducer::AlCaECALRecHitReducer(const edm::ParameterSet& iConfig) {
  ebRecHitsToken_ = consumes<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>("ebRecHitsLabel"));
  eeRecHitsToken_ = consumes<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>("eeRecHitsLabel"));
  //  esRecHitsLabel_ = iConfig.getParameter< edm::InputTag > ("esRecHitsLabel");

  std::vector<edm::InputTag> srcLabels = iConfig.getParameter<std::vector<edm::InputTag> >("srcLabels");
  for (auto inputTag = srcLabels.begin(); inputTag != srcLabels.end(); ++inputTag) {
    eleViewTokens_.push_back(consumes<edm::View<reco::RecoCandidate> >(*inputTag));
  }

  //eleViewToken_ = consumes<edm::View <reco::RecoCandidate> > (iConfig.getParameter< edm::InputTag > ("electronLabel"));
  photonToken_ = consumes<reco::PhotonCollection>(iConfig.getParameter<edm::InputTag>("photonLabel"));
  EESuperClusterToken_ =
      consumes<reco::SuperClusterCollection>(iConfig.getParameter<edm::InputTag>("EESuperClusterCollection"));

  minEta_highEtaSC_ = iConfig.getParameter<double>("minEta_highEtaSC");

  alcaBarrelHitsCollection_ = iConfig.getParameter<std::string>("alcaBarrelHitCollection");
  alcaEndcapHitsCollection_ = iConfig.getParameter<std::string>("alcaEndcapHitCollection");
  alcaCaloClusterCollection_ = iConfig.getParameter<std::string>("alcaCaloClusterCollection");

  //  alcaPreshowerHitsCollection_ = iConfig.getParameter<std::string>("alcaPreshowerHitCollection");

  etaSize_ = iConfig.getParameter<int>("etaSize");
  phiSize_ = iConfig.getParameter<int>("phiSize");
  // FIXME: minimum size of etaSize_ and phiSize_
  if (phiSize_ % 2 == 0 || etaSize_ % 2 == 0)
    edm::LogError("AlCaECALRecHitReducerError") << "Size of eta/phi should be odd numbers";

  //  esNstrips_  = iConfig.getParameter<int> ("esNstrips");
  //  esNcolumns_ = iConfig.getParameter<int> ("esNcolumns");

  //register your products
  produces<EBRecHitCollection>(alcaBarrelHitsCollection_);
  produces<EERecHitCollection>(alcaEndcapHitsCollection_);
  produces<reco::CaloClusterCollection>(alcaCaloClusterCollection_);
  //  produces< ESRecHitCollection > (alcaPreshowerHitsCollection_) ;
}

AlCaECALRecHitReducer::~AlCaECALRecHitReducer() {}

// ------------ method called to produce the data  ------------
void AlCaECALRecHitReducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
  using namespace edm;
  //using namespace std;

  edm::ESHandle<CaloTopology> theCaloTopology;
  iSetup.get<CaloTopologyRecord>().get(theCaloTopology);
  const CaloTopology* caloTopology = theCaloTopology.product();

  // Get Photons
  Handle<reco::PhotonCollection> phoHandle;
  iEvent.getByToken(photonToken_, phoHandle);

  // get RecHits
  Handle<EBRecHitCollection> barrelRecHitsHandle;
  iEvent.getByToken(ebRecHitsToken_, barrelRecHitsHandle);
  const EBRecHitCollection* barrelHitsCollection = barrelRecHitsHandle.product();

  // get RecHits
  Handle<EERecHitCollection> endcapRecHitsHandle;
  iEvent.getByToken(eeRecHitsToken_, endcapRecHitsHandle);
  const EERecHitCollection* endcapHitsCollection = endcapRecHitsHandle.product();

  //   // get ES RecHits
  //   Handle<ESRecHitCollection> preshowerRecHitsHandle;
  //   iEvent.getByToken(esRecHitsToken_,preshowerRecHitsHandle);

  //   const ESRecHitCollection * preshowerHitsCollection = 0 ;
  //   if (preshowerIsFull)
  //     preshowerHitsCollection = preshowerRecHitsHandle.product () ;

  //   // make a vector to store the used ES rechits:
  //   set<ESDetId> used_strips;
  //   used_strips.clear();

  // for Z->ele+SC
  Handle<reco::SuperClusterCollection> EESCHandle;
  iEvent.getByToken(EESuperClusterToken_, EESCHandle);

  //Create empty output collections
  auto miniEBRecHitCollection = std::make_unique<EBRecHitCollection>();
  auto miniEERecHitCollection = std::make_unique<EERecHitCollection>();
  //  auto miniESRecHitCollection = std::make_unique<ESRecHitCollection>();

  std::set<DetId> reducedRecHit_EBmap;
  std::set<DetId> reducedRecHit_EEmap;

  //  std::set< edm::Ref<reco::CaloCluster> > reducedCaloClusters_map;

  auto reducedCaloClusterCollection = std::make_unique<reco::CaloClusterCollection>();

  //Photons:
#ifdef shervin
  for (reco::PhotonCollection::const_iterator phoIt = phoHandle->begin(); phoIt != phoHandle->end(); phoIt++) {
    const reco::SuperCluster& sc = *(phoIt->superCluster());

    if (phoIt->isEB()) {
      AddMiniRecHitCollection(sc, reducedRecHit_EBmap, caloTopology);
    } else {  // endcap
      AddMiniRecHitCollection(sc, reducedRecHit_EEmap, caloTopology);
    }  // end of endcap

    reco::CaloCluster_iterator it = sc.clustersBegin();
    reco::CaloCluster_iterator itend = sc.clustersEnd();
    for (; it != itend; ++it) {
      reco::CaloCluster caloClus(**it);
      reducedCaloClusterCollection->push_back(caloClus);
    }
  }
#endif

  Handle<edm::View<reco::RecoCandidate> > eleViewHandle;
  for (auto iToken = eleViewTokens_.begin(); iToken != eleViewTokens_.end(); iToken++) {
    iEvent.getByToken(*iToken, eleViewHandle);

    //Electrons:
    for (auto eleIt = eleViewHandle->begin(); eleIt != eleViewHandle->end(); eleIt++) {
      const reco::SuperCluster& sc = *(eleIt->superCluster());

      if (fabs(sc.eta()) < 1.479) {
        AddMiniRecHitCollection(sc, reducedRecHit_EBmap, caloTopology);
      } else {  // endcap
        AddMiniRecHitCollection(sc, reducedRecHit_EEmap, caloTopology);
      }  // end of endcap

      reco::CaloCluster_iterator it = sc.clustersBegin();
      reco::CaloCluster_iterator itend = sc.clustersEnd();
      for (; it != itend; ++it) {
        reco::CaloCluster caloClus(**it);
        reducedCaloClusterCollection->push_back(caloClus);
      }
    }
  }

  //saving recHits for highEta SCs for highEta studies
  for (reco::SuperClusterCollection::const_iterator SC_iter = EESCHandle->begin(); SC_iter != EESCHandle->end();
       SC_iter++) {
    if (fabs(SC_iter->eta()) < minEta_highEtaSC_)
      continue;
    AddMiniRecHitCollection(*SC_iter, reducedRecHit_EEmap, caloTopology);

    const reco::SuperCluster& sc = *(SC_iter);
    reco::CaloCluster_iterator it = sc.clustersBegin();
    reco::CaloCluster_iterator itend = sc.clustersEnd();
    for (; it != itend; ++it) {
      reco::CaloCluster caloClus(**it);
      reducedCaloClusterCollection->push_back(caloClus);
    }
  }

  //------------------------------ fill the alcareco reduced recHit collection
  for (std::set<DetId>::const_iterator itr = reducedRecHit_EBmap.begin(); itr != reducedRecHit_EBmap.end(); itr++) {
    if (barrelHitsCollection->find(*itr) != barrelHitsCollection->end())
      miniEBRecHitCollection->push_back(*(barrelHitsCollection->find(*itr)));
  }

  for (std::set<DetId>::const_iterator itr = reducedRecHit_EEmap.begin(); itr != reducedRecHit_EEmap.end(); itr++) {
    if (endcapHitsCollection->find(*itr) != endcapHitsCollection->end())
      miniEERecHitCollection->push_back(*(endcapHitsCollection->find(*itr)));
  }

  //--------------------------------------- Put selected information in the event
  iEvent.put(std::move(miniEBRecHitCollection), alcaBarrelHitsCollection_);
  iEvent.put(std::move(miniEERecHitCollection), alcaEndcapHitsCollection_);
  //  iEvent.put(std::move(miniESRecHitCollection),alcaPreshowerHitsCollection_ );
  iEvent.put(std::move(reducedCaloClusterCollection), alcaCaloClusterCollection_);
}

void AlCaECALRecHitReducer::AddMiniRecHitCollection(const reco::SuperCluster& sc,
                                                    std::set<DetId>& reducedRecHitMap,
                                                    const CaloTopology* caloTopology) {
  DetId seed = (sc.seed()->seed());
  int phiSize = phiSize_, etaSize = etaSize_;
  if (seed.subdetId() != EcalBarrel) {  // if not EB, take a square window
    etaSize = std::max(phiSize_, etaSize_);
    phiSize = etaSize;
  }

  std::vector<DetId> recHit_window = caloTopology->getWindow(seed, phiSize, etaSize);
  for (unsigned int i = 0; i < recHit_window.size(); i++) {
    reducedRecHitMap.insert(recHit_window[i]);
  }

  const std::vector<std::pair<DetId, float> >& scHits = sc.hitsAndFractions();
  for (std::vector<std::pair<DetId, float> >::const_iterator scHit_itr = scHits.begin(); scHit_itr != scHits.end();
       scHit_itr++) {
    // the map fills just one time (avoiding double insert of recHits)
    reducedRecHitMap.insert(scHit_itr->first);
  }

  return;
}

DEFINE_FWK_MODULE(AlCaECALRecHitReducer);