CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HLTCaloObjInRegionsProducer.h
Go to the documentation of this file.
1 #ifndef RecoEgamma_EgammaHLTProducers_HLTCaloObjInRegionsProducer_h_
2 #define RecoEgamma_EgammaHLTProducers_HLTCaloObjInRegionsProducer_h_
3 
4 #include <memory>
5 
16 
25 
30 
32 
36 
39 
40 /**************************************************************
41 / purpose: enable filtering of calo objects in eta/phi or deltaR
42 / regions around generic objects
43 /
44 / operation : accepts all objects with
45 / (dEta <dEtaMax && dPhi < dPhiMax) || dR < dRMax
46 / so the OR of a rectangular region and cone region
47 ****************************************************************/
48 
49 //this is a struct which contains all the eta/phi regions
50 //from which to filter the calo objs
52 private:
53  float centreEta_;
54  float centrePhi_;
55  float maxDeltaR2_;
56  float maxDEta_;
57  float maxDPhi_;
58  public:
59  EtaPhiRegion(float iEta,float iPhi,float iDR,float iDEta,float iDPhi):
60  centreEta_(iEta),centrePhi_(iPhi),maxDeltaR2_(iDR*iDR),
61  maxDEta_(iDEta),maxDPhi_(iDPhi){}
63  bool operator()(float eta,float phi)const{
64  return reco::deltaR2(eta,phi,centreEta_,centrePhi_)<maxDeltaR2_ ||
65  (std::abs(eta-centreEta_)<maxDEta_ && std::abs(reco::deltaPhi(phi,centrePhi_))<maxDPhi_);}
66 
67 };
68 
70 public:
72  virtual ~EtaPhiRegionDataBase() = default;
73  virtual void getEtaPhiRegions(const edm::Event&,std::vector<EtaPhiRegion>&)const=0;
74 };
75 
76 
77 //this class stores the tokens to access the objects around which we wish to filter
78 //it makes a vector of EtaPhiRegions which are then used to filter the CaloObjs
79 template<typename T1> class EtaPhiRegionData : public EtaPhiRegionDataBase{
80 private:
81  float minEt_;
82  float maxEt_;
83  float maxDeltaR_;
84  float maxDEta_;
85  float maxDPhi_;
87 public:
89  minEt_(para.getParameter<double>("minEt")),
90  maxEt_(para.getParameter<double>("maxEt")),
91  maxDeltaR_(para.getParameter<double>("maxDeltaR")),
92  maxDEta_(para.getParameter<double>("maxDEta")),
93  maxDPhi_(para.getParameter<double>("maxDPhi")),
94  token_(consumesColl.consumes<T1>(para.getParameter<edm::InputTag>("inputColl"))){}
95 
96  void getEtaPhiRegions(const edm::Event&,std::vector<EtaPhiRegion>&)const override;
97  template<typename T2>static typename T2::const_iterator beginIt(const T2& coll){return coll.begin();}
98  template<typename T2>static typename T2::const_iterator endIt(const T2& coll){return coll.end();}
99  template<typename T2>static typename BXVector<T2>::const_iterator beginIt(const BXVector<T2>& coll){return coll.begin(0);}
100  template<typename T2>static typename BXVector<T2>::const_iterator endIt(const BXVector<T2>& coll){return coll.end(0);}
101 
102 };
103 
104 
105 template<typename CaloObjType,
106  typename CaloObjCollType =edm::SortedCollection<CaloObjType> >
108 
109 
110  public:
111 
112 
113 
116 
117  void produce(edm::Event&, const edm::EventSetup&) override;
118  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
119 
120  private:
121  EtaPhiRegionDataBase* createEtaPhiRegionData(const std::string&,const edm::ParameterSet&,edm::ConsumesCollector &&); //calling function owns this
122  static std::unique_ptr<CaloObjCollType>
123  makeFilteredColl(const edm::Handle<CaloObjCollType>& inputColl,
124  const edm::ESHandle<CaloGeometry>& caloGeomHandle,
125  const std::vector<EtaPhiRegion>& regions);
126  static bool validIDForGeom(const DetId& id);
127  std::vector<std::string> outputProductNames_;
128  std::vector<edm::InputTag> inputCollTags_;
129  std::vector<edm::EDGetTokenT<CaloObjCollType>> inputTokens_;
130  std::vector<std::unique_ptr<EtaPhiRegionDataBase>> etaPhiRegionData_;
131 };
132 
133 
134 template<typename CaloObjType,typename CaloObjCollType>
136 {
137  const std::vector<edm::ParameterSet> etaPhiRegions = para.getParameter<std::vector<edm::ParameterSet>>("etaPhiRegions");
138  for(auto& pset : etaPhiRegions){
139  const std::string type=pset.getParameter<std::string>("type");
140  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
141  }
142 
143  outputProductNames_=para.getParameter<std::vector<std::string>>("outputProductNames");
144  inputCollTags_=para.getParameter<std::vector<edm::InputTag>>("inputCollTags");
145  if(outputProductNames_.size()!=inputCollTags_.size()){
146  throw cms::Exception("InvalidConfiguration") <<" error outputProductNames and inputCollTags must be the same size, they are "<<outputProductNames_.size()<<" vs "<<inputCollTags_.size();
147  }
148  for (unsigned int collNr=0; collNr<inputCollTags_.size(); collNr++) {
149  inputTokens_.push_back(consumes<CaloObjCollType>(inputCollTags_[collNr]));
150  produces<CaloObjCollType> (outputProductNames_[collNr]);
151  }
152 }
153 
154 template<typename CaloObjType,typename CaloObjCollType>
156 {
158  std::vector<std::string> outputProductNames;
159  outputProductNames.push_back("EcalRegionalRecHitsEB");
160  desc.add<std::vector<std::string>>("outputProductNames", outputProductNames);
161  std::vector<edm::InputTag> inputColls;
162  inputColls.push_back(edm::InputTag("hltHcalDigis"));
163  desc.add<std::vector<edm::InputTag>>("inputCollTags", inputColls);
164  std::vector<edm::ParameterSet> etaPhiRegions;
165 
166  edm::ParameterSet ecalCandPSet;
167  ecalCandPSet.addParameter<std::string>("type","RecoEcalCandidate");
168  ecalCandPSet.addParameter<double>("minEt",-1);
169  ecalCandPSet.addParameter<double>("maxEt",-1);
170  ecalCandPSet.addParameter<double>("maxDeltaR",0.5);
171  ecalCandPSet.addParameter<double>("maxDEta",0.);
172  ecalCandPSet.addParameter<double>("maxDPhi",0.);
173  ecalCandPSet.addParameter<edm::InputTag>("inputColl",edm::InputTag("hltEgammaCandidates"));
174  etaPhiRegions.push_back(ecalCandPSet);
175 
176 
177  edm::ParameterSetDescription etaPhiRegionDesc;
178  etaPhiRegionDesc.add<std::string>("type");
179  etaPhiRegionDesc.add<double>("minEt");
180  etaPhiRegionDesc.add<double>("maxEt");
181  etaPhiRegionDesc.add<double>("maxDeltaR");
182  etaPhiRegionDesc.add<double>("maxDEta");
183  etaPhiRegionDesc.add<double>("maxDPhi");
184  etaPhiRegionDesc.add<edm::InputTag>("inputColl");
185  desc.addVPSet("etaPhiRegions",etaPhiRegionDesc,etaPhiRegions);
186 
188 }
189 
190 
191 template<typename CaloObjType,typename CaloObjCollType>
193 
194  // get the collection geometry:
195  edm::ESHandle<CaloGeometry> caloGeomHandle;
196  setup.get<CaloGeometryRecord>().get(caloGeomHandle);
197 
198  std::vector<EtaPhiRegion> regions;
199  std::for_each(etaPhiRegionData_.begin(),etaPhiRegionData_.end(),
200  [&event,&regions](const std::unique_ptr<EtaPhiRegionDataBase>& input)
201  {input->getEtaPhiRegions(event,regions);}
202  );
203 
204  for(size_t inputCollNr=0;inputCollNr<inputTokens_.size();inputCollNr++){
206  event.getByToken(inputTokens_[inputCollNr],inputColl);
207 
208  if (!(inputColl.isValid())) {
209  edm::LogError("ProductNotFound")<< "could not get a handle on the "<<typeid(CaloObjCollType).name() <<" named "<< inputCollTags_[inputCollNr].encode() << std::endl;
210  continue;
211  }
212  auto outputColl = makeFilteredColl(inputColl,caloGeomHandle,regions);
213  event.put(std::move(outputColl),outputProductNames_[inputCollNr]);
214  }
215 }
216 
217 template<typename CaloObjType,typename CaloObjCollType>
218 std::unique_ptr<CaloObjCollType>
221  const edm::ESHandle<CaloGeometry>& caloGeomHandle,
222  const std::vector<EtaPhiRegion>& regions)
223 {
224 
225  auto outputColl = std::make_unique<CaloObjCollType>();
226  if(!inputColl->empty()){
227  const CaloSubdetectorGeometry* subDetGeom=caloGeomHandle->getSubdetectorGeometry(inputColl->begin()->id());
228  if(!regions.empty()){
229  for(const CaloObjType& obj : *inputColl){
230  auto objGeom = subDetGeom->getGeometry(obj.id());
231  if(objGeom==nullptr){
232  //wondering what to do here
233  //something is very very wrong
234  //given HLT should never crash or throw, decided to log an error
235  //update: so turns out HCAL can pass through calibration channels in QIE11 so for that module, its an expected behaviour
236  //so we check if the ID is valid
237  if(validIDForGeom(obj.id())) {
238  edm::LogError("HLTCaloObjInRegionsProducer") << "for an object of type "<<typeid(CaloObjType).name()<<" the geometry returned null for id "<<DetId(obj.id()).rawId()<<" with initial ID "<<DetId(inputColl->begin()->id()).rawId()<<" in HLTCaloObjsInRegion, this shouldnt be possible and something has gone wrong, auto accepting hit";
239  }
240  outputColl->push_back(obj);
241  continue;
242  }
243  float eta = objGeom->getPosition().eta();
244  float phi = objGeom->getPosition().phi();
245 
246  for(const auto& region : regions){
247  if(region(eta,phi)) {
248  outputColl->push_back(obj);
249  break;
250  }
251  }
252  }
253  }//end check of empty regions
254  }//end check of empty rec-hits
255  return outputColl;
256 
257 }
258 
259 //tells us if an ID should have a valid geometry
260 //it assumes that all IDs do except those specifically mentioned
261 //HCAL for example have laser calibs in the digi collection so
262 //so we have to ensure that HCAL is HB,HE or HO
263 template<typename CaloObjType,typename CaloObjCollType>
266 {
267  if(id.det()==DetId::Hcal){
268  if(id.subdetId() == HcalSubdetector::HcalEmpty ||
269  id.subdetId() == HcalSubdetector::HcalOther){
270  return false;
271  }
272  }
273  return true;
274 }
275 
276 
277 
278 
279 template<typename CaloObjType,typename CaloObjCollType>
283  edm::ConsumesCollector && consumesColl)
284 {
285  if(type=="L1EGamma"){
286  return new EtaPhiRegionData<l1t::EGammaBxCollection>(para,consumesColl);
287  }else if(type=="L1Jet"){
288  return new EtaPhiRegionData<l1t::JetBxCollection>(para,consumesColl);
289  }else if(type=="L1Muon"){
290  return new EtaPhiRegionData<l1t::MuonBxCollection>(para,consumesColl);
291  }else if(type=="L1Tau"){
292  return new EtaPhiRegionData<l1t::TauBxCollection>(para,consumesColl);
293  }else if(type=="RecoEcalCandidate"){
294  return new EtaPhiRegionData<reco::RecoEcalCandidateCollection>(para,consumesColl);
295  }else if(type=="RecoChargedCandidate"){
296  return new EtaPhiRegionData<reco::RecoChargedCandidateCollection>(para,consumesColl);
297  }else if(type=="Electron"){
298  return new EtaPhiRegionData<reco::Electron>(para,consumesColl);
299  }else{
300  //this is a major issue and could lead to rather subtle efficiency losses, so if its incorrectly configured, we're aborting the job!
301  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;
302  }
303 
304 }
305 
306 template<typename CandCollType>
307 void EtaPhiRegionData<CandCollType>::getEtaPhiRegions(const edm::Event& event,std::vector<EtaPhiRegion>&regions)const
308 {
310  event.getByToken(token_,cands);
311 
312  for(auto candIt = beginIt(*cands);candIt!=endIt(*cands);++candIt){
313  if(candIt->et() >= minEt_ && (maxEt_<0 || candIt->et() < maxEt_)){
314  regions.push_back(EtaPhiRegion(candIt->eta(),candIt->phi(),
315  maxDeltaR_,maxDEta_,maxDPhi_));
316  }
317 
318  }
319 }
320 
321 
322 
323 #endif
324 
325 
const_iterator end(int bx) const
type
Definition: HCALResponse.h:21
T getParameter(std::string const &) const
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:44
ParameterDescriptionBase * addVPSet(U const &iLabel, ParameterSetDescription const &validator, std::vector< ParameterSet > const &defaults)
std::vector< edm::EDGetTokenT< CaloObjCollType > > inputTokens_
std::vector< edm::InputTag > inputCollTags_
EtaPhiRegionDataBase * createEtaPhiRegionData(const std::string &, const edm::ParameterSet &, edm::ConsumesCollector &&)
ProductID id() const
Definition: HandleBase.cc:15
def setup(process, global_tag, zero_tesla=False)
Definition: GeneralSetup.py:1
static BXVector< T2 >::const_iterator beginIt(const BXVector< T2 > &coll)
std::string defaultModuleLabel()
static T2::const_iterator beginIt(const T2 &coll)
static bool validIDForGeom(const DetId &id)
static std::string const input
Definition: EdmProvDump.cc:44
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
std::vector< std::string > outputProductNames_
static std::unique_ptr< CaloObjCollType > makeFilteredColl(const edm::Handle< CaloObjCollType > &inputColl, const edm::ESHandle< CaloGeometry > &caloGeomHandle, const std::vector< EtaPhiRegion > &regions)
static T2::const_iterator endIt(const T2 &coll)
void addParameter(std::string const &name, T const &value)
Definition: ParameterSet.h:144
bool operator()(float eta, float phi) const
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
EtaPhiRegionData(const edm::ParameterSet &para, edm::ConsumesCollector &consumesColl)
SubDetector subDetGeom[21]
ParameterDescriptionBase * add(U const &iLabel, T const &value)
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
bool isValid() const
Definition: HandleBase.h:74
HLTCaloObjInRegionsProducer(const edm::ParameterSet &ps)
double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:22
Definition: DetId.h:18
JetCorrectorParametersCollection coll
Definition: classes.h:10
edm::EDGetTokenT< T1 > token_
virtual std::shared_ptr< const CaloCellGeometry > getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
const T & get() const
Definition: EventSetup.h:58
void add(std::string const &label, ParameterSetDescription const &psetDescription)
T1 deltaR2(T1 eta1, T2 phi1, T3 eta2, T4 phi2)
Definition: deltaR.h:36
HLT enums.
static BXVector< T2 >::const_iterator endIt(const BXVector< T2 > &coll)
void getEtaPhiRegions(const edm::Event &, std::vector< EtaPhiRegion > &) const override
const_iterator begin(int bx) const
EtaPhiRegion(float iEta, float iPhi, float iDR, float iDEta, float iDPhi)
std::vector< std::unique_ptr< EtaPhiRegionDataBase > > etaPhiRegionData_
def move(src, dest)
Definition: eostools.py:510
Definition: event.py:1
void produce(edm::Event &, const edm::EventSetup &) override