Go to the documentation of this file.00001
00002
00003
00004
00005
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #include <memory>
00023
00024
00025 #include "FWCore/Framework/interface/Frameworkfwd.h"
00026 #include "FWCore/Framework/interface/EDProducer.h"
00027
00028 #include "FWCore/Framework/interface/Event.h"
00029 #include "FWCore/Framework/interface/MakerMacros.h"
00030
00031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00032
00033 #include "FWCore/Framework/interface/ESHandle.h"
00034 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
00035 #include "DataFormats/DetId/interface/DetIdCollection.h"
00036
00037 #include "CondFormats/HcalObjects/interface/HcalChannelQuality.h"
00038 #include "CondFormats/DataRecord/interface/HcalChannelQualityRcd.h"
00039
00040 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputer.h"
00041 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputerRcd.h"
00042
00043
00044
00045
00046
00047 class HcalHitSelection : public edm::EDProducer {
00048 public:
00049 explicit HcalHitSelection(const edm::ParameterSet&);
00050 ~HcalHitSelection();
00051
00052 private:
00053 virtual void beginJob() ;
00054 virtual void produce(edm::Event&, const edm::EventSetup&);
00055 virtual void endJob() ;
00056
00057 edm::InputTag hbheTag,hoTag,hfTag;
00058 std::vector<edm::InputTag> interestingDetIdCollections;
00059
00060
00061 edm::ESHandle<HcalChannelQuality> theHcalChStatus;
00062 edm::ESHandle<HcalSeverityLevelComputer> theHcalSevLvlComputer;
00063 std::set<DetId> toBeKept;
00064 template <typename CollectionType> void skim( const edm::Handle<CollectionType> & input, std::auto_ptr<CollectionType> & output);
00065
00066
00067 };
00068
00069 template <class CollectionType> void HcalHitSelection::skim( const edm::Handle<CollectionType> & input, std::auto_ptr<CollectionType> & output){
00070 output->reserve(input->size());
00071 typename CollectionType::const_iterator begin=input->begin();
00072 typename CollectionType::const_iterator end=input->end();
00073 typename CollectionType::const_iterator hit=begin;
00074
00075 for (;hit!=end;++hit){
00076
00077 const DetId & id = hit->detid();
00078 const uint32_t & recHitFlag = hit->flags();
00079
00080
00081 const uint32_t & dbStatusFlag = theHcalChStatus->getValues(id)->getValue();
00082 int severityLevel = theHcalSevLvlComputer->getSeverityLevel(id, recHitFlag, dbStatusFlag);
00083
00084 if (severityLevel!=0){
00085 output->push_back(*hit);
00086 }else{
00087
00088 if (toBeKept.find(id)!=toBeKept.end())
00089 output->push_back(*hit);
00090 }
00091 }
00092 }
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106 HcalHitSelection::HcalHitSelection(const edm::ParameterSet& iConfig)
00107 {
00108 hbheTag=iConfig.getParameter<edm::InputTag>("hbheTag");
00109 hfTag=iConfig.getParameter<edm::InputTag>("hfTag");
00110 hoTag=iConfig.getParameter<edm::InputTag>("hoTag");
00111
00112 interestingDetIdCollections = iConfig.getParameter< std::vector<edm::InputTag> >("interestingDetIds");
00113
00114 produces<HBHERecHitCollection>(hbheTag.label());
00115 produces<HFRecHitCollection>(hfTag.label());
00116 produces<HORecHitCollection>(hoTag.label());
00117
00118 }
00119
00120
00121 HcalHitSelection::~HcalHitSelection()
00122 {
00123
00124
00125
00126
00127 }
00128
00129
00130
00131
00132
00133
00134
00135 void
00136 HcalHitSelection::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
00137 {
00138 iSetup.get<HcalChannelQualityRcd>().get(theHcalChStatus);
00139 iSetup.get<HcalSeverityLevelComputerRcd>().get(theHcalSevLvlComputer);
00140
00141 edm::Handle<HBHERecHitCollection> hbhe;
00142 edm::Handle<HFRecHitCollection> hf;
00143 edm::Handle<HORecHitCollection> ho;
00144
00145 iEvent.getByLabel(hbheTag,hbhe);
00146 iEvent.getByLabel(hfTag,hf);
00147 iEvent.getByLabel(hoTag,ho);
00148
00149 toBeKept.clear();
00150 edm::Handle<DetIdCollection > detId;
00151 for( unsigned int t = 0; t < interestingDetIdCollections.size(); ++t )
00152 {
00153 iEvent.getByLabel(interestingDetIdCollections[t],detId);
00154 if (!detId.isValid())
00155 continue;
00156 toBeKept.insert(detId->begin(),detId->end());
00157 }
00158
00159 std::auto_ptr<HBHERecHitCollection> hbhe_out(new HBHERecHitCollection());
00160 skim(hbhe,hbhe_out);
00161 iEvent.put(hbhe_out,hbheTag.label());
00162
00163 std::auto_ptr<HFRecHitCollection> hf_out(new HFRecHitCollection());
00164 skim(hf,hf_out);
00165 iEvent.put(hf_out,hfTag.label());
00166
00167 std::auto_ptr<HORecHitCollection> ho_out(new HORecHitCollection());
00168 skim(ho,ho_out);
00169 iEvent.put(ho_out,hoTag.label());
00170
00171
00172
00173 }
00174
00175
00176 void
00177 HcalHitSelection::beginJob()
00178 {
00179 }
00180
00181
00182 void
00183 HcalHitSelection::endJob() {
00184 }
00185
00186
00187 DEFINE_FWK_MODULE(HcalHitSelection);