File indexing completed on 2023-03-17 10:50:05
0001 #include <string>
0002 #include <iostream>
0003 #include <vector>
0004 #include <array>
0005 #include <cstdint>
0006
0007 #include "FWCore/Framework/interface/stream/EDAnalyzer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012
0013 #include "DataFormats/Common/interface/Handle.h"
0014 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0015 #include "DataFormats/HcalRecHit/interface/CaloRecHitAuxSetter.h"
0016
0017 namespace {
0018 std::ostream& operator<<(std::ostream& s, const HFQIE10Info& i) {
0019 s << i.id() << ": " << i.energy() << " GeV"
0020 << ", t= " << i.timeRising() << " to " << i.timeFalling() << " ns";
0021 return s;
0022 }
0023
0024 std::ostream& operator<<(std::ostream& s, const HFPreRecHit& hit) {
0025 s << "{ ";
0026 const HFQIE10Info* i = hit.getHFQIE10Info(0);
0027 if (i) {
0028 s << *i;
0029 }
0030 s << " }, ";
0031 s << "{ ";
0032 i = hit.getHFQIE10Info(1);
0033 if (i) {
0034 s << *i;
0035 }
0036 s << " }";
0037 return s;
0038 }
0039
0040 template <std::size_t N>
0041 void printBits(std::ostream& s, const std::array<uint32_t, N>& allbits, const std::vector<int>& bits) {
0042 const int maxbit = N * 32;
0043 const unsigned len = bits.size();
0044 for (unsigned i = 0; i < len; ++i) {
0045 const int bitnum = bits[i];
0046 if (bitnum >= 0 && bitnum < maxbit) {
0047 const unsigned ibit = bitnum % 32;
0048 const bool bit = (allbits[bitnum / 32] & (1U << ibit)) >> ibit;
0049 s << bit;
0050 } else
0051 s << '-';
0052 }
0053 }
0054
0055 void printRecHitAuxInfo(std::ostream& s, const HFPreRecHit& i, const std::vector<int>& bits, bool) {}
0056
0057 void printRecHitAuxInfo(std::ostream& s, const HBHERecHit& i, const std::vector<int>& bits, const bool plan1) {
0058 if (plan1 && i.isMerged()) {
0059
0060 std::vector<HcalDetId> ids;
0061 i.getMergedIds(&ids);
0062 const unsigned n = ids.size();
0063 s << "; merged " << n << ": ";
0064 for (unsigned j = 0; j < n; ++j) {
0065 if (j)
0066 s << ", ";
0067 s << ids[j];
0068 }
0069 }
0070 if (!bits.empty()) {
0071 std::array<uint32_t, 4> allbits;
0072 allbits[0] = i.flags();
0073 allbits[1] = i.aux();
0074 allbits[2] = i.auxHBHE();
0075 allbits[3] = i.auxPhase1();
0076 s << "; bits: ";
0077 printBits(s, allbits, bits);
0078 }
0079
0080
0081 s << "; tdc:";
0082 const uint32_t auxTDC = i.auxTDC();
0083 if (auxTDC) {
0084 const unsigned six_bits_mask = 0x3f;
0085 for (unsigned ts = 0; ts < 5; ++ts)
0086 s << ' ' << CaloRecHitAuxSetter::getField(auxTDC, six_bits_mask, ts * 6);
0087 } else {
0088 s << " none";
0089 }
0090 }
0091
0092 void printRecHitAuxInfo(std::ostream& s, const HFRecHit& i, const std::vector<int>& bits, bool) {
0093 if (!bits.empty()) {
0094 std::array<uint32_t, 3> allbits;
0095 allbits[0] = i.flags();
0096 allbits[1] = i.aux();
0097 allbits[2] = i.getAuxHF();
0098 s << "; bits: ";
0099 printBits(s, allbits, bits);
0100 }
0101 }
0102 }
0103
0104 using namespace std;
0105
0106 class HcalRecHitDump : public edm::stream::EDAnalyzer<> {
0107 public:
0108 explicit HcalRecHitDump(edm::ParameterSet const& conf);
0109 virtual void analyze(edm::Event const& e, edm::EventSetup const& c) override;
0110
0111 private:
0112 string hbhePrefix_;
0113 string hfPrefix_;
0114 string hfprePrefix_;
0115 std::vector<int> bits_;
0116 bool printPlan1Info_;
0117
0118 edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0119 edm::EDGetTokenT<HFRecHitCollection> tok_hf_;
0120 edm::EDGetTokenT<HFPreRecHitCollection> tok_prehf_;
0121
0122 unsigned long long counter_;
0123
0124 template <class Collection, class Token>
0125 void analyzeT(edm::Event const& e,
0126 const Token& tok,
0127 const char* name,
0128 const string& prefix,
0129 const bool printPlan1Info = false) const {
0130 cout << prefix << " rechit dump " << counter_ << endl;
0131
0132 edm::Handle<Collection> coll;
0133 bool found = false;
0134 try {
0135 e.getByToken(tok, coll);
0136 found = true;
0137 } catch (...) {
0138 cout << prefix << " Error: no " << name << " rechit data" << endl;
0139 }
0140 if (found) {
0141 for (typename Collection::const_iterator j = coll->begin(); j != coll->end(); ++j) {
0142 cout << prefix << *j;
0143 printRecHitAuxInfo(cout, *j, bits_, printPlan1Info);
0144 cout << endl;
0145 }
0146 }
0147 }
0148 };
0149
0150 HcalRecHitDump::HcalRecHitDump(edm::ParameterSet const& conf)
0151 : hbhePrefix_(conf.getUntrackedParameter<string>("hbhePrefix", "")),
0152 hfPrefix_(conf.getUntrackedParameter<string>("hfPrefix", "")),
0153 hfprePrefix_(conf.getUntrackedParameter<string>("hfprePrefix", "")),
0154 bits_(conf.getUntrackedParameter<std::vector<int> >("bits")),
0155 printPlan1Info_(conf.getUntrackedParameter<bool>("printPlan1Info", false)),
0156 counter_(0) {
0157 if (!hbhePrefix_.empty())
0158 tok_hbhe_ = consumes<HBHERecHitCollection>(conf.getParameter<edm::InputTag>("tagHBHE"));
0159 if (!hfPrefix_.empty())
0160 tok_hf_ = consumes<HFRecHitCollection>(conf.getParameter<edm::InputTag>("tagHF"));
0161 if (!hfprePrefix_.empty())
0162 tok_prehf_ = consumes<HFPreRecHitCollection>(conf.getParameter<edm::InputTag>("tagPreHF"));
0163 }
0164
0165 void HcalRecHitDump::analyze(edm::Event const& e, edm::EventSetup const& c) {
0166 if (!hbhePrefix_.empty())
0167 analyzeT<HBHERecHitCollection>(e, tok_hbhe_, "HBHE", hbhePrefix_, printPlan1Info_);
0168 if (!hfPrefix_.empty())
0169 analyzeT<HFRecHitCollection>(e, tok_hf_, "HF", hfPrefix_);
0170 if (!hfprePrefix_.empty())
0171 analyzeT<HFPreRecHitCollection>(e, tok_prehf_, "PreHF", hfprePrefix_);
0172 ++counter_;
0173 }
0174
0175 DEFINE_FWK_MODULE(HcalRecHitDump);