CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
L1ABCDebugger.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: SiStripTools
4 // Class: L1ABCDebugger
5 //
13 //
14 // Original Author: Andrea Venturi
15 // Created: Tue Jul 19 11:56:00 CEST 2009
16 //
17 //
18 
19 // system include files
20 #include <memory>
21 
22 // user include files
23 #include "TH2F.h"
24 #include "TProfile.h"
26 
27 #include <vector>
28 #include <string>
29 
32 
36 
39 
41 
43 //
44 // class decleration
45 //
46 
48 public:
49  explicit L1ABCDebugger(const edm::ParameterSet&);
50  ~L1ABCDebugger() override;
51 
52 private:
53  void beginJob() override;
54  void analyze(const edm::Event&, const edm::EventSetup&) override;
55  void beginRun(const edm::Run&, const edm::EventSetup&) override;
56  void endJob() override;
57 
58  // ----------member data ---------------------------
59 
61  const unsigned int m_maxLS;
62  const unsigned int m_LSfrac;
63 
65 
66  TH2F** m_hoffsets;
67  TProfile** m_horboffvsorb;
68  TProfile** m_hbxoffvsorb;
69 };
70 
71 //
72 // constants, enums and typedefs
73 //
74 
75 //
76 // static data member definitions
77 //
78 
79 //
80 // constructors and destructor
81 //
83  : m_l1abccollectionToken(
84  consumes<L1AcceptBunchCrossingCollection>(iConfig.getParameter<edm::InputTag>("l1ABCCollection"))),
85  m_maxLS(iConfig.getUntrackedParameter<unsigned int>("maxLSBeforeRebin", 250)),
86  m_LSfrac(iConfig.getUntrackedParameter<unsigned int>("startingLSFraction", 16)),
87  m_rhm(consumesCollector()) {
88  //now do what ever initialization is needed
89 
91  "offsets", "Orbit vs BX offsets between SCAL and Event", 2 * 3564 + 1, -3564.5, 3564.5, 201, -100.5, 100.5);
93  m_rhm.makeTProfile("orboffvsorb", "SCAL Orbit offset vs orbit number", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
95  m_rhm.makeTProfile("bxoffvsorb", "SCAL BX offset vs orbit number", m_LSfrac * m_maxLS, 0, m_maxLS * 262144);
96 }
97 
99  // do anything here that needs to be done at desctruction time
100  // (e.g. close files, deallocate resources etc.)
101 }
102 
103 //
104 // member functions
105 //
106 
107 // ------------ method called to for each event ------------
109  using namespace edm;
110 
112  iEvent.getByToken(m_l1abccollectionToken, pIn);
113 
114  // offset computation
115  for (L1AcceptBunchCrossingCollection::const_iterator l1abc = pIn->begin(); l1abc != pIn->end(); ++l1abc) {
116  if (l1abc->l1AcceptOffset() == 0) {
117  if (m_hoffsets && *m_hoffsets)
118  (*m_hoffsets)
119  ->Fill((int)l1abc->bunchCrossing() - (int)iEvent.bunchCrossing(),
120  (long long)l1abc->orbitNumber() - (long long)iEvent.orbitNumber());
122  (*m_horboffvsorb)->Fill(iEvent.orbitNumber(), (long long)l1abc->orbitNumber() - (long long)iEvent.orbitNumber());
124  (*m_hbxoffvsorb)->Fill(iEvent.orbitNumber(), (int)l1abc->bunchCrossing() - (int)iEvent.bunchCrossing());
125  }
126  }
127 
128  // dump of L1ABC collection
129 
130  edm::LogInfo("L1ABCDebug") << "Dump of L1AcceptBunchCrossing Collection for event in orbit " << iEvent.orbitNumber()
131  << " and BX " << iEvent.bunchCrossing();
132 
133  for (L1AcceptBunchCrossingCollection::const_iterator l1abc = pIn->begin(); l1abc != pIn->end(); ++l1abc) {
134  edm::LogVerbatim("L1ABCDebug") << *l1abc;
135  }
136 }
137 
138 void L1ABCDebugger::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
139  m_rhm.beginRun(iRun);
140 
141  if (m_hoffsets && *m_hoffsets) {
142  (*m_hoffsets)->GetXaxis()->SetTitle("#Delta BX (SCAL-Event)");
143  (*m_hoffsets)->GetYaxis()->SetTitle("#Delta orbit (SCAL-Event)");
144  }
145  if (m_horboffvsorb && *m_horboffvsorb) {
146  (*m_horboffvsorb)->GetXaxis()->SetTitle("Orbit");
147  (*m_horboffvsorb)->GetYaxis()->SetTitle("#Delta orbit (SCAL-Event)");
148  (*m_horboffvsorb)->SetCanExtend(TH1::kXaxis);
149  }
150  if (m_hbxoffvsorb && *m_hbxoffvsorb) {
151  (*m_hbxoffvsorb)->GetXaxis()->SetTitle("Orbit");
152  (*m_hbxoffvsorb)->GetYaxis()->SetTitle("#Delta BX (SCAL-Event)");
153  (*m_hbxoffvsorb)->SetCanExtend(TH1::kXaxis);
154  }
155 }
156 // ------------ method called once each job just before starting event loop ------------
158 
159 // ------------ method called once each job just after ending the event loop ------------
161 
162 //define this as a plug-in
Log< level::Info, true > LogVerbatim
TH2F ** m_hoffsets
const unsigned int m_maxLS
~L1ABCDebugger() override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
RunHistogramManager m_rhm
int bunchCrossing() const
Definition: EventBase.h:64
edm::EDGetTokenT< L1AcceptBunchCrossingCollection > m_l1abccollectionToken
void analyze(const edm::Event &, const edm::EventSetup &) override
const unsigned int m_LSfrac
int iEvent
Definition: GenABIO.cc:224
m_rhm(consumesCollector())
TProfile ** m_hbxoffvsorb
TProfile ** makeTProfile(const char *name, const char *title, const unsigned int nbinx, const double xmin, const double xmax)
L1ABCDebugger(const edm::ParameterSet &)
std::vector< L1AcceptBunchCrossing > L1AcceptBunchCrossingCollection
int orbitNumber() const
Definition: EventBase.h:65
void endJob() override
Log< level::Info, false > LogInfo
void beginJob() override
void beginRun(const edm::Run &iRun)
TProfile ** m_horboffvsorb
TH2F ** makeTH2F(const char *name, const char *title, const unsigned int nbinx, const double xmin, const double xmax, const unsigned int nbiny, const double ymin, const double ymax)
Definition: Run.h:45
void beginRun(const edm::Run &, const edm::EventSetup &) override