RateMonitor_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: RateMonitor
3 // Module Type: analyzer
4 // File: RateMonitor_module.cc
5 //
6 // Generated at Fri Jan 24 13:19:10 2014 by Matthew Tamsett using artmod
7 // from cetpkgsupport v1_04_02.
8 ////////////////////////////////////////////////////////////////////////
9 
17 #include "fhiclcpp/ParameterSet.h"
18 
24 
26 #include <iomanip>
27 #include <stdlib.h>
28 #include "TH1.h"
29 
30 //---------------------------------------------------------------------
31 namespace novaddt {
32  class RateMonitor;
33 }
34 //---------------------------------------------------------------------
36  public:
37  explicit RateMonitor(fhicl::ParameterSet const & p);
38  virtual ~RateMonitor();
39  void analyze(art::Event const & e) override;
40  void endJob() override;
41  double convertTDCtoSeconds(uint64_t tdc);
42  uint64_t quantiseToMicroSlice(uint64_t TDC);
43  private:
44  double average(std::vector<double> const& vec) const;
45 
46  // Declare member data here.
47  // FHiCL configurations
48  std::string _sliceModuleLabel; ///< label of module making the TDC sorted hits
49  std::string _sliceInstanceLabel; ///< instance label making the TDC sorted hits
50  std::vector< std::string > _triggerLabels; ///< vector of labels of trigger modules
51  std::string _quantisation; ///< Method of microslice quantisation, one of: minimum, maximum, stochastic
54  bool _hit_rates;
55  bool _verbose;
56  // Counters
57  uint64_t _timeSum = 0;
58  unsigned int _nEvents = 0;
59  std::map <std::string, int> d_nEvents;
60  std::map <std::string, uint64_t> d_sumTime;
61  std::map <std::string, uint64_t> d_quantisedSumTime;
62 
63  std::map <std::string, std::vector<double>> hr_ratios;
64 
65  // Histograms
66  uint64_t running_sumTime = 0;
67  std::map <std::string, int> d_running_nEvents;
68  std::map <std::string, uint64_t> d_running_sumTime;
69  std::map <std::string, uint64_t> d_running_quantisedSumTime;
70  std::map <std::string, std::vector<double> > d_rates;
71  std::map <std::string, std::vector<double> > d_acceptances;
72  std::map <std::string, std::vector<double> > d_quantisedAcceptances;
73 };
74 //---------------------------------------------------------------------
76  : EDAnalyzer(p)
77  , _sliceModuleLabel (p.get<std::string >("SliceModuleLabel"))
78  , _sliceInstanceLabel(p.get<std::string >("SliceInstanceLabel"))
79  , _triggerLabels(p.get<std::vector<std::string>>("TriggerLabels"))
80  , _quantisation(p.get<std::string>("QuantisationMethod"))
81  , _writeHistograms(p.get<bool>("WriteHistograms"))
82  , _average_over_n_events(p.get<int>("AverageOverNEvents"))
83  , _hit_rates(p.get<bool>("HitRates"))
84  , _verbose(p.get<bool>("Verbose"))
85 
86 {
87  std::cout << "--- novaddt::RateMonitor instantiate" << std::endl;
88  std::cout << "\t SliceModuleLabel: " << _sliceModuleLabel << std::endl;
89  std::cout << "\t SliceInstanceLabel: " << _sliceInstanceLabel << std::endl;
90  std::cout << "\t Quantisation method: " << _quantisation << std::endl;
91  std::cout << "\t Monitoring rates for the following triggers:\n";
92  for (auto t : _triggerLabels){
93  std::cout << "\t - " << t << std::endl;
94  } // end of loop on trigger labels
95  std::cout << "\t Write histograms: " << _writeHistograms << std::endl;
96  std::cout << "\t Averaged over n events " << _average_over_n_events << std::endl;
97  std::cout << "\t Verbose: " << _verbose << std::endl;
98 
99  // Initalise maps
100  for (auto t : _triggerLabels){
101  d_nEvents[t] = 0;
102  d_sumTime[t] = 0;
103  // histogram stuff
104  d_running_nEvents[t] = 0;
105  d_running_sumTime[t] = 0;
107  } // end of loop on trigger labels
108 
109  // check that the quantisation is configured
110  assert( (_quantisation == "minimum") ||
111  (_quantisation == "maximum") ||
112  (_quantisation == "stochastic")
113  );
114 }
115 //---------------------------------------------------------------------
117 {
118  // Clean up dynamic memory and other resources here.
119 }
120 //---------------------------------------------------------------------
122 {
123  // Implementation of required member function here.
124  if(_verbose) std::cout << "--- novaddt::RateMonitor analyze. Event: "
125  << event.id().event()
126  << std::endl;
127  _nEvents++;
128 
129  // Get the hits
131  event.getByLabel(_sliceModuleLabel, _sliceInstanceLabel, hits);
132  if(hits.failedToGet()){
133  if(_verbose) std::cout << "No slices found, skipping event" << std::endl;
134  _nEvents--;
135  return;
136  }
137 
138  if(_verbose) std::cout << "\tgot: " << hits->size() << " hits\n";
139  uint64_t min_TDC = 0;
140  uint64_t max_TDC = 0;
141  uint64_t d_TDC = 0;
142  if (hits->size() >0){
143  min_TDC = hits->at(0).TDC().val;
144  max_TDC = hits->at(hits->size()-1).TDC().val;
145  assert( max_TDC > min_TDC ); // sanity check
146  d_TDC = max_TDC - min_TDC;
147  }
148  _timeSum = _timeSum + d_TDC;
149 
150  if(_verbose) std::cout << "\ttime range: " << min_TDC
151  << " to " << max_TDC
152  << ", dTDC: " << d_TDC
153  << std::endl;
154 
155  // Loop over the trigger decisions
156  for (auto t : _triggerLabels){
157  art::Handle<std::vector<TriggerDecision>> trigger_decisions;
158  event.getByLabel(t, trigger_decisions);
159  if (trigger_decisions.failedToGet()){
160  if (_verbose) std::cout << "\tno trigger decision objects found for trigger module " << t << std::endl;
161  continue;
162  }
163  if(_verbose) std::cout << "\t" << t
164  << ":\t" << trigger_decisions->size()
165  << " triggers\n";
166 
167  if(trigger_decisions->size() == 0) continue;
168  assert(hits->size() >0); // sanity check - if there are no hits, how can there be a trigger pass?
169 
170  for(auto td : *trigger_decisions)
171  {
172  if(_verbose)
173  std::cout << "\t\t" << td.start()
174  << " to " << td.start()+td.duration()
175  << ", dTDC: "<< td.duration()
176  << std::endl;
177 
178  d_nEvents[t]++;
179  d_sumTime[t] = d_sumTime[t]+td.duration();
181 
182  if(_writeHistograms){
183  d_running_nEvents[t]++;
184  d_running_sumTime[t] = d_running_sumTime[t]+td.duration();
186  } // end of if on histogram filling
187 
188  if (_hit_rates)
189  {
190  // Let's also output the number of hits in this time window
191  auto start = std::lower_bound(
192  hits->begin(), hits->end(),
193  TDC(td.start()), CompareDAQHit<TDC>());
194 
195  auto finish = std::lower_bound
196  (start, hits->end(),
197  TDC(td.start() + td.duration()), CompareDAQHit<TDC>());
198 
199  if (start != hits->end() && finish != hits->end())
200  {
201  // 1. calculate average cosmic rate
202  double TDC_to_second = 1.0 / 64.0e6;
203  double d_TDC_seconds = d_TDC * TDC_to_second;
204  double cosmic_rate = hits->size() / d_TDC_seconds;
205  if (_verbose)
206  std::cout << "\t\t Average Cosmic Rate = "
207  << cosmic_rate << " hits/s" << std::endl;
208 
209  // 2. calculate rate in trigger window
210  double dt = (finish->TDC() - start->TDC()) * TDC_to_second;
211  double rate = (finish - start) / dt;
212  if (_verbose)
213  std::cout << "\t\t Rate in Trigger Window = "
214  << rate << " hits/s" << std::endl;
215 
216  // 3. combine into a 50 us averaged rate
217  double avg_rate =
218  (dt * rate + (50.0e-6 - dt) * cosmic_rate) / 50.0e-6;
219  if (_verbose)
220  std::cout << "\t\t 50 us Averaged Rate = "
221  << avg_rate << " hits/s" << std::endl;
222 
223  // 4. report the ratio of 3./1.
224  double hr_ratio = avg_rate / cosmic_rate;
225  if (_verbose)
226  std::cout << "\t\t 50 us Averaged Rate / "
227  << "Average Cosmic Rate = "
228  << hr_ratio << std::endl;
229 
230  hr_ratios[t].push_back(hr_ratio);
231  }
232  } // if _hit_rates
233  } // end of loop on decisions
234  } // end of loop on trigger labels
235 
236  // Convert running sums to rates/acceptances
237  if(_writeHistograms){
239 
240  // average over some number of events
241  if ( (event.id().event() !=0) && ( (event.id().event()%_average_over_n_events) == 0 ) ){
242 
243  if(_verbose) std::cout << "\tFilling histograms\n";
244 
245  for (auto t : _triggerLabels){
247  d_acceptances[t].push_back( (double(d_running_sumTime[t]) / double(running_sumTime) ) * 100 );
248  d_quantisedAcceptances[t].push_back( (double(d_running_quantisedSumTime[t]) / double(running_sumTime) ) * 100 );
249  } // end of loop on trigger labels
250 
251  // reset
252  running_sumTime = 0;
253  for (auto t : _triggerLabels){
254  d_running_nEvents[t] = 0;
255  d_running_sumTime[t] = 0;
257  } // end of loop on trigger labels
258  } // end on if on average over some number of events
259  } // end of if on histogram filling
260 }
261 //---------------------------------------------------------------------
263 {
264  // Implementation of optional member function here.
265  std::cout << "--- novaddt::RateMonitor endJob" << std::endl;
266  std::cout << "\tNumber of events: " << _nEvents << std::endl;
267  std::cout << "\tTotal TDC sum: " << _timeSum << std::endl;
268  std::cout << "\tTotal time [s]: " << convertTDCtoSeconds(_timeSum) << std::endl
269  << std::endl;
270 
271  // print formatted table
273  << std::setw(20) << "Trigger"
274  << std::right
275  << " | " << std::setw(12) << "Passes"
276  << " | " << std::setw(12) << "Sum TDC"
277  << " | " << std::setw(12) << "Rate [Hz]"
278  << " | " << std::setw(14) << "Acceptance [%]"
279  << " | " << std::setw(24) << "Quantised Acceptance [%]"
280  << " | " << std::setw(20) << "50 us Hit Rate Ratio"
281  << std::endl;
282 
283  int width = 132;
284  for (int i = 0; i < width; i++) std::cout << "-";
285  std::cout << std::endl;
286 
287  for (auto t : _triggerLabels){
288  double rate = double(d_nEvents[t]) / convertTDCtoSeconds(_timeSum);
289  double acceptance = (double(d_sumTime[t]) / double(_timeSum) ) * 100;
290  double q_acceptance = (double(d_quantisedSumTime[t]) / double(_timeSum) ) * 100;
291 
292  double hr_avg = 0.0;
293  if (_hit_rates)
294  hr_avg = average(hr_ratios[t]);
295 
297  << std::setw(20) << t
298  << std::right
299  << " | " << std::setw(12) << d_nEvents[t]
300  << " | " << std::setw(12) << d_sumTime[t]
301  << " | " << std::setw(12) << rate
302  << " | " << std::setw(14) << acceptance
303  << " | " << std::setw(24) << q_acceptance
304  << " | " << std::setw(20) << hr_avg
305  << std::endl;
306  } // end of loop on trigger labels
307 
308  for (int i = 0; i < width; i++) std::cout << "-";
309  std::cout << std::endl;
310 
311  // Histogram filling
312  if(_writeHistograms){
314  std::cout << "\tWriting histograms\n";
315  TH1F* h_sumTime;
316  h_sumTime = tfs->make<TH1F>( (std::string("_sum_time")).c_str(),
317  (std::string(";x;")+std::string(" sum time [s]")).c_str(),
318  1,0,1);
319  h_sumTime->SetBinContent(1,convertTDCtoSeconds(_timeSum));
320 
321  std::map <std::string, TH1F* > d_h_total_events;
322  std::map <std::string, std::vector<TH1F*> > d_h_rates;
323  std::map <std::string, std::vector<TH1F*> > d_h_acceptances;
324  std::map <std::string, std::vector<TH1F*> > d_h_quantisedAcceptances;
325 
326  for (auto t : _triggerLabels){
327  assert( d_rates[t].size() == d_acceptances[t].size() );
329  int n_bins = d_rates[t].size();
330  float low_edge = -0.5;
331  float high_edge = float(n_bins)-0.5;
332 
333  // declare histogram
334  d_h_total_events[t] = tfs->make<TH1F>( (t+std::string("_total_events")).c_str(),
335  (std::string(";x;")+t+std::string(" events")).c_str(),
336  1,0,1);
337  d_h_total_events[t]->SetBinContent(1,d_nEvents[t]);
338 
339  d_h_rates[t].push_back(tfs->make<TH1F>( (t+std::string("_rate")).c_str(),
340  (std::string(";Event index;")+t+std::string(" rate [Hz]")).c_str(),
341  n_bins, low_edge, high_edge) );
342  d_h_acceptances[t].push_back(tfs->make<TH1F>( (t+std::string("_acceptance")).c_str(),
343  (std::string(";Event index;")+t+std::string(" acceptance [%]")).c_str(),
344  n_bins, low_edge, high_edge) );
345 
346  d_h_quantisedAcceptances[t].push_back(tfs->make<TH1F>( (t+std::string("_quantised_acceptance")).c_str(),
347  (std::string(";Event index;")+t+std::string(" quantised acceptance [%]")).c_str(),
348  n_bins, low_edge, high_edge) );
349  // fill histograms
350  for (int i=0;i<n_bins;i++){
351  d_h_rates[t][d_h_rates[t].size()-1]->SetBinContent(i+1,d_rates[t][i]);
352  d_h_acceptances[t][d_h_rates[t].size()-1]->SetBinContent(i+1,d_acceptances[t][i]);
353  d_h_quantisedAcceptances[t][d_h_rates[t].size()-1]->SetBinContent(i+1,d_quantisedAcceptances[t][i]);
354  } // end of loop on n bins
355  } // end of loop on trigger labels
356  } // end of if on write histograms
357 }
358 //---------------------------------------------------------------------
360 {
361  // TDC clock ticks at 64 MHz
362  double one_tdc = 1. / 64e6;
363  return double(TDC) * one_tdc;
364 }
365 //---------------------------------------------------------------------
367 {
368  // Microslices are 50 us long = 3,200 TDC ticks.
369  // We can only store trigger decisions in quanta of MicroSlices,
370  // this function will provide this quantisation.
371 
372  // Minimum - assumes that the start of trigger decision window and microslice boundary allign
373  if(_quantisation == "minimum") return ((TDC / 3200)+1)*3200;
374 
375  // Maximum - assumes that the centre of trigger decision window and microslice boundary allign
376  if(_quantisation == "maximum") return ((TDC / 3200)+2)*3200;
377 
378  // Stochastic - add random offset then quantise
379  if(_quantisation == "stochastic"){
380  uint64_t offset = uint64_t(rand() % 3200);
381  return (((TDC+offset) / 3200)+1)*3200;
382  }
383 
384  return 0;
385 }
386 //---------------------------------------------------------------------
387 double novaddt::RateMonitor::average(std::vector<double> const& vec) const
388 {
389  double sum = 0.0;
390 
391  for (auto const& v : vec)
392  sum += v;
393 
394  return sum / vec.size();
395 }
std::vector< std::string > _triggerLabels
vector of labels of trigger modules
constexpr auto const & right(const_AssnsIter< L, R, D, Dir > const &a, const_AssnsIter< L, R, D, Dir > const &b)
Definition: AssnsIter.h:112
std::string _sliceInstanceLabel
instance label making the TDC sorted hits
value_type val
Definition: BaseProducts.h:34
std::map< std::string, std::vector< double > > d_acceptances
void analyze(art::Event const &e) override
const char * p
Definition: xmltok.h:285
DEFINE_ART_MODULE(TestTMapFile)
std::string _quantisation
Method of microslice quantisation, one of: minimum, maximum, stochastic.
double average(std::vector< double > const &vec) const
Definition: Cand.cxx:23
std::map< std::string, int > d_running_nEvents
void hits()
Definition: readHits.C:15
RateMonitor(fhicl::ParameterSet const &p)
std::string _sliceModuleLabel
label of module making the TDC sorted hits
std::map< std::string, std::vector< double > > d_quantisedAcceptances
double convertTDCtoSeconds(uint64_t tdc)
std::map< std::string, int > d_nEvents
std::map< std::string, uint64_t > d_quantisedSumTime
std::map< std::string, std::vector< double > > d_rates
Eigen::VectorXd vec
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
OStream cout
Definition: OStream.cxx:6
std::map< std::string, uint64_t > d_running_sumTime
uint64_t quantiseToMicroSlice(uint64_t TDC)
EventNumber_t event() const
Definition: EventID.h:116
constexpr auto const & left(const_AssnsIter< L, R, D, Dir > const &a, const_AssnsIter< L, R, D, Dir > const &b)
Definition: AssnsIter.h:104
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
std::map< std::string, uint64_t > d_sumTime
T * make(ARGS...args) const
std::map< std::string, std::vector< double > > hr_ratios
assert(nhit_max >=nhit_nbins)
std::map< std::string, uint64_t > d_running_quantisedSumTime
Double_t sum
Definition: plot.C:31
Float_t e
Definition: plot.C:35
EventID id() const
Definition: Event.h:56
bool failedToGet() const
Definition: Handle.h:196