HCAL Digi Filtering Module : 80X

RecoEgamma/EgammaHLTProducers/interface/HLTCaloObjInRegionsProducer.h

@@ -0,0 +1,306 @@
+#ifndef RecoEgamma_EgammaHLTProducers_HLTCaloObjInRegionsProducer_h_
+#define RecoEgamma_EgammaHLTProducers_HLTCaloObjInRegionsProducer_h_
+
+#include <memory>
+
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/EventSetup.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+#include "FWCore/Framework/interface/ConsumesCollector.h"
+#include "FWCore/Framework/interface/ESHandle.h"
+
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+
+// Reco candidates (might not need)
+#include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
+#include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
+#include "DataFormats/RecoCandidate/interface/RecoEcalCandidateFwd.h"
+#include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
+#include "DataFormats/EgammaCandidates/interface/Electron.h"
+#include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
+
+
+// Geometry and topology
+#include "Geometry/Records/interface/CaloGeometryRecord.h"
+#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
+#include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
+#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
+
+#include "RecoEcal/EgammaCoreTools/interface/EcalEtaPhiRegion.h"
+
+#include "DataFormats/L1Trigger/interface/EGamma.h"
+#include "DataFormats/L1Trigger/interface/Jet.h"
+#include "DataFormats/L1Trigger/interface/Muon.h"
+
+#include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
+#include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
+#include "DataFormats/HcalDigi/interface/HBHEDataFrame.h"
+
+#include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
+
+#include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
+
+/**************************************************************
+/ purpose: enable filtering of calo objects in eta/phi or deltaR
+/          regions around generic objects 
+/
+/ operation : accepts all objects with
+/             (dEta <dEtaMax  && dPhi < dPhiMax) || dR < dRMax
+/             so the OR of a rectangular region and cone region
+****************************************************************/
+
+//this is a struct which contains all the eta/phi regions
+//from which to filter the calo objs
+class EtaPhiRegion{
+private:
+    float centreEta_;
+    float centrePhi_;
+    float maxDeltaR2_;
+    float maxDEta_;
+    float maxDPhi_;
+  public:
+    EtaPhiRegion(float iEta,float iPhi,float iDR,float iDEta,float iDPhi):
+      centreEta_(iEta),centrePhi_(iPhi),maxDeltaR2_(iDR*iDR),
+      maxDEta_(iDEta),maxDPhi_(iDPhi){}
+    ~EtaPhiRegion(){}
+  bool operator()(float eta,float phi)const{
+    return reco::deltaR2(eta,phi,centreEta_,centrePhi_)<maxDeltaR2_ ||
+      (std::abs(eta-centreEta_)<maxDEta_ && std::abs(reco::deltaPhi(phi,centrePhi_)<maxDPhi_));}
+
+};
+
+class EtaPhiRegionDataBase {
+public:
+  EtaPhiRegionDataBase(){}
+  virtual void getEtaPhiRegions(const edm::Event&,std::vector<EtaPhiRegion>&)const=0;
+};
+
+
+//this class stores the tokens to access the objects around which we wish to filter
+//it makes a vector of EtaPhiRegions which are then used to filter the CaloObjs
+template<typename T1> class EtaPhiRegionData : public EtaPhiRegionDataBase{
+private:
+  float minEt_;
+  float maxEt_;
+  float maxDeltaR_;
+  float maxDEta_;
+  float maxDPhi_;
+  edm::EDGetTokenT<T1> token_;
+public:
+  EtaPhiRegionData(const edm::ParameterSet& para,edm::ConsumesCollector & consumesColl):
+    minEt_(para.getParameter<double>("minEt")),
+    maxEt_(para.getParameter<double>("maxEt")),
+    maxDeltaR_(para.getParameter<double>("maxDeltaR")),
+    maxDEta_(para.getParameter<double>("maxDEta")),
+    maxDPhi_(para.getParameter<double>("maxDPhi")),
+    token_(consumesColl.consumes<T1>(para.getParameter<edm::InputTag>("inputColl"))){}
+
+  void getEtaPhiRegions(const edm::Event&,std::vector<EtaPhiRegion>&)const override;
+  template<typename T2>static typename T2::const_iterator beginIt(const T2& coll){return coll.begin();}
+  template<typename T2>static typename T2::const_iterator endIt(const T2& coll){return coll.end();}
+  template<typename T2>static typename BXVector<T2>::const_iterator beginIt(const BXVector<T2>& coll){return coll.begin(0);}
+  template<typename T2>static typename BXVector<T2>::const_iterator endIt(const BXVector<T2>& coll){return coll.end(0);}
+
+};
+
+
+template<typename CaloObjType,
+	 typename CaloObjCollType =edm::SortedCollection<CaloObjType> >
+class HLTCaloObjInRegionsProducer : public edm::stream::EDProducer<> {
+
+
+ public:
+
+
+
+  HLTCaloObjInRegionsProducer(const edm::ParameterSet& ps);
+  ~HLTCaloObjInRegionsProducer(){}
+
+  void produce(edm::Event&, const edm::EventSetup&) override;
+  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
+
+ private:
+  EtaPhiRegionDataBase* createEtaPhiRegionData(const std::string&,const edm::ParameterSet&,edm::ConsumesCollector &&); //calling function owns this
+  static std::auto_ptr<CaloObjCollType> 
+  makeFilteredColl(const edm::Handle<CaloObjCollType>& inputColl,
+		   const edm::ESHandle<CaloGeometry>& caloGeomHandle,
+		   const std::vector<EtaPhiRegion>& regions);
+  std::vector<std::string> outputProductNames_;
+  std::vector<edm::InputTag> inputCollTags_;
+  std::vector<edm::EDGetTokenT<CaloObjCollType>> inputTokens_;
+  std::vector<std::unique_ptr<EtaPhiRegionDataBase>> etaPhiRegionData_;
+};
+
+
+template<typename CaloObjType,typename CaloObjCollType> 
+HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>::HLTCaloObjInRegionsProducer(const edm::ParameterSet& para)
+{
+  const std::vector<edm::ParameterSet> etaPhiRegions = para.getParameter<std::vector<edm::ParameterSet>>("etaPhiRegions");
+  for(auto& pset : etaPhiRegions){
+    const std::string type=pset.getParameter<std::string>("type");
+    etaPhiRegionData_.emplace_back(createEtaPhiRegionData(type,pset,consumesCollector())); //meh I was going to use a factory but it was going to be overly complex for my needs
+  }
+
+  outputProductNames_=para.getParameter<std::vector<std::string>>("outputProductNames");
+  inputCollTags_=para.getParameter<std::vector<edm::InputTag>>("inputCollTags");
+  for (unsigned int collNr=0; collNr<inputCollTags_.size(); collNr++) { 
+    inputTokens_.push_back(consumes<CaloObjCollType>(inputCollTags_[collNr]));
+    produces<CaloObjCollType> (outputProductNames_[collNr]);
+  }
+}
+
+template<typename CaloObjType,typename CaloObjCollType> 
+void HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>::fillDescriptions(edm::ConfigurationDescriptions& descriptions)
+{
+  edm::ParameterSetDescription desc;
+  std::vector<std::string> outputProductNames;
+  outputProductNames.push_back("EcalRegionalRecHitsEB");
+  desc.add<std::vector<std::string>>("outputProductNames", outputProductNames);
+  std::vector<edm::InputTag> inputColls;
+  inputColls.push_back(edm::InputTag("hltHcalDigis"));
+  desc.add<std::vector<edm::InputTag>>("inputCollTags", inputColls);
+  std::vector<edm::ParameterSet> etaPhiRegions;
+
+  edm::ParameterSet ecalCandPSet;
+  ecalCandPSet.addParameter<std::string>("type","RecoEcalCandidate");
+  ecalCandPSet.addParameter<double>("minEt",-1);
+  ecalCandPSet.addParameter<double>("maxEt",-1);
+  ecalCandPSet.addParameter<double>("maxDeltaR",0.5);
+  ecalCandPSet.addParameter<double>("maxDEta",0.);
+  ecalCandPSet.addParameter<double>("maxDPhi",0.);
+  ecalCandPSet.addParameter<edm::InputTag>("inputColl",edm::InputTag("hltEgammaCandidates"));
+  etaPhiRegions.push_back(ecalCandPSet);
+
+
+  edm::ParameterSetDescription etaPhiRegionDesc;
+  etaPhiRegionDesc.add<std::string>("type");
+  etaPhiRegionDesc.add<double>("minEt");
+  etaPhiRegionDesc.add<double>("maxEt");
+  etaPhiRegionDesc.add<double>("maxDeltaR");
+  etaPhiRegionDesc.add<double>("maxDEta");
+  etaPhiRegionDesc.add<double>("maxDPhi");
+  etaPhiRegionDesc.add<edm::InputTag>("inputColl");
+  desc.addVPSet("etaPhiRegions",etaPhiRegionDesc,etaPhiRegions);
+
+  descriptions.add(defaultModuleLabel<HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>>(), desc); 
+}
+
+
+template<typename CaloObjType,typename CaloObjCollType>
+void HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>::produce(edm::Event& event, const edm::EventSetup& setup) {
+
+  // get the collection geometry:
+  edm::ESHandle<CaloGeometry> caloGeomHandle;
+  setup.get<CaloGeometryRecord>().get(caloGeomHandle);
+
+  std::vector<EtaPhiRegion> regions;
+  std::for_each(etaPhiRegionData_.begin(),etaPhiRegionData_.end(),
+		[&event,&regions](const std::unique_ptr<EtaPhiRegionDataBase>& input)
+		{input->getEtaPhiRegions(event,regions);}
+		);
+
+  for(size_t inputCollNr=0;inputCollNr<inputTokens_.size();inputCollNr++){
+    edm::Handle<CaloObjCollType> inputColl;
+    event.getByToken(inputTokens_[inputCollNr],inputColl);
+
+    if (!(inputColl.isValid())) {
+      edm::LogError("ProductNotFound")<< "could not get a handle on the "<<typeid(CaloObjCollType).name() <<" named "<< inputCollTags_[inputCollNr].encode() << std::endl;
+      continue;
+    }
+    auto outputColl = makeFilteredColl(inputColl,caloGeomHandle,regions);
+    event.put(outputColl,outputProductNames_[inputCollNr]);
+  }
+}
+
+template<typename CaloObjType,typename CaloObjCollType>
+std::auto_ptr<CaloObjCollType> 
+HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>::
+makeFilteredColl(const edm::Handle<CaloObjCollType>& inputColl,
+		 const edm::ESHandle<CaloGeometry>& caloGeomHandle,
+		 const std::vector<EtaPhiRegion>& regions)
+{  
+
+  std::auto_ptr<CaloObjCollType> outputColl(new CaloObjCollType);
+  if(!inputColl->empty()){
+    const CaloSubdetectorGeometry* subDetGeom=caloGeomHandle->getSubdetectorGeometry(inputColl->front().id());
+    if(!regions.empty()){
+      for(const CaloObjType& obj : *inputColl){
+	const CaloCellGeometry* objGeom = subDetGeom->getGeometry(obj.id());
+	if(objGeom==nullptr){
+	  //wondering what to do here
+	  //something is very very wrong
+	  //given HLT should never crash or throw, decided to log an error and 
+	  edm::LogError("HLTCaloObjInRegionsProducer") << "for an object of type "<<typeid(CaloObjType).name()<<" the geometry returned null for id "<<obj.id()<<" in HLTCaloObjsInRegion, this shouldnt be possible and something has gone wrong, auto accepting hit";
+	  outputColl->push_back(obj);
+	}
+	float eta = objGeom->getPosition().eta();
+	float phi = objGeom->getPosition().phi();
+
+	for(const auto& region : regions){
+	  if(region(eta,phi)) {
+	    outputColl->push_back(obj);
+	    break;
+	  }
+	}
+      }
+    }//end check of empty regions
+  }//end check of empty rec-hits
+  return outputColl;
+
+}
+
+
+
+
+
+template<typename CaloObjType,typename CaloObjCollType> 
+EtaPhiRegionDataBase* 
+HLTCaloObjInRegionsProducer<CaloObjType,CaloObjCollType>::
+createEtaPhiRegionData(const std::string& type,const edm::ParameterSet& para,
+		       edm::ConsumesCollector && consumesColl)
+{
+  if(type=="L1EGamma"){
+    return new EtaPhiRegionData<l1t::EGammaBxCollection>(para,consumesColl);
+  }else if(type=="L1Jet"){
+    return new EtaPhiRegionData<l1t::JetBxCollection>(para,consumesColl);
+  }else if(type=="L1Muon"){
+    return new EtaPhiRegionData<l1t::MuonBxCollection>(para,consumesColl);
+  }else if(type=="L1Tau"){
+    return new EtaPhiRegionData<l1t::TauBxCollection>(para,consumesColl);
+  }else if(type=="RecoEcalCandidate"){
+    return new EtaPhiRegionData<reco::RecoEcalCandidateCollection>(para,consumesColl);
+  }else if(type=="RecoChargedCandiate"){
+    return new EtaPhiRegionData<reco::RecoChargedCandidateCollection>(para,consumesColl);
+  }else if(type=="Electron"){
+    return new EtaPhiRegionData<reco::Electron>(para,consumesColl);
+  }else{
+    //this is a major issue and could lead to rather subtle efficiency losses, so if its incorrectly configured, we're aborting the job!
+    throw cms::Exception("InvalidConfig") << " type "<<type<<" is not recognised, this means the rec-hit you think you are keeping may not be and you should fix this error as it can lead to hard to find efficiency loses"<<std::endl;
+  }
+
+}
+
+template<typename CandCollType>
+void EtaPhiRegionData<CandCollType>::getEtaPhiRegions(const edm::Event& event,std::vector<EtaPhiRegion>&regions)const
+{
+  edm::Handle<CandCollType> cands;
+  event.getByToken(token_,cands);
+
+  for(auto candIt = beginIt(*cands);candIt!=endIt(*cands);++candIt){
+    if(candIt->et() >= minEt_ && (maxEt_<0 || candIt->et() < maxEt_)){
+       regions.push_back(EtaPhiRegion(candIt->eta(),candIt->phi(),
+				      maxDeltaR_,maxDEta_,maxDPhi_));
+    }
+
+  }
+}
+
+
+typedef HLTCaloObjInRegionsProducer<HBHEDataFrame> HLTHcalDigisInRegionsProducer;
+DEFINE_FWK_MODULE(HLTHcalDigisInRegionsProducer);
+
+#endif
+
+

RecoEgamma/EgammaHLTProducers/src/HLTCaloObjInRegionsProducer.cc

@@ -0,0 +1 @@
+#include "RecoEgamma/EgammaHLTProducers/interface/HLTCaloObjInRegionsProducer.h"