JTriggerMonitor.cc

Go to the documentation of this file.

#include <string>
#include <iostream>
#include <iomanip>

#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"

#include "km3net-dataformat/online/JDAQ.hh"
#include "km3net-dataformat/online/JDAQTriggerMask.hh"
#include "JDAQ/JDAQEventIO.hh"
#include "JDAQ/JDAQTimesliceIO.hh"
#include "JDAQ/JDAQSummarysliceIO.hh"

#include "JTrigger/JTriggerParameters.hh"
#include "JTrigger/JTriggerBits.hh"

#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JTriggerParametersSupportkit.hh"
#include "JROOT/JRootToolkit.hh"

#include "Jeep/JPrint.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"


static const int WIDTH = 20;


/**
 * Print histogram parameters.
 *
 * \param  h1       histogram
 * \param  out      output stream
 */
inline void print(const TH1& h1, std::ostream& out)
{
  using namespace std;

  out << setw(WIDTH) << left << h1.GetName() << ' '
      << FIXED(8,0)  << h1.GetEntries()      << ' '
      << FIXED(8,2)  << h1.GetMean()         << ' '
      << FIXED(8,2)  << h1.GetRMS()          << ' '
      << endl;
}


/**
 * \file
 *
 * Auxiliary program to histogram and print trigger statistics.
 * \author mdejong
 */
int main(int argc, char **argv)
{
  using namespace std;
  using namespace JPP;
  using namespace KM3NETDAQ;

  JMultipleFileScanner<>     inputFile;
  JLimit_t& numberOfEvents = inputFile.getLimit();
  string       outputFile;
  int          qaqc;
  int          debug;

  try {

    JParser<> zap("Auxiliary program to histogram and print trigger statistics.");

    zap['f'] = make_field(inputFile);
    zap['o'] = make_field(outputFile)          = "";
    zap['n'] = make_field(numberOfEvents)      = JLimit::max();
    zap['Q'] = make_field(qaqc)                = 0;
    zap['d'] = make_field(debug)               = 1;

    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }


  cout.tie(&cerr);


  TH1D h1("Event    ", NULL, 100, 0.0, 1.0e1);
  TH1D h2("Summary  ", NULL, 100, 0.0, 1.0e1);

  TH1D hm("Trigger  mask ", NULL, NUMBER_OF_TRIGGER_BITS, -0.5, NUMBER_OF_TRIGGER_BITS - 0.5);
  TH1D hn("Trigger  hits ", NULL, 1000, -0.5,  1000 - 0.5);
  TH1D hs("Snapshot hits ", NULL, 1000, -0.5, 10000 - 0.5);
  TH1D hr("PMT rate [kHz]", NULL,  100,  0.0,    50.0);


  const JTriggerParameters trigger_parameters = getTriggerParameters(inputFile);

  NOTICE(trigger_parameters << endl);


  for (JMultipleFileScanner<JDAQEvent> input = inputFile; input.hasNext(); ) {

    STATUS("entry: " << setw(10) << input.getCounter() << '\r'); DEBUG(endl);

    const JDAQEvent* object = input.next();

    h1.Fill((Double_t) getSizeof(*object) * 1e-6);

    hn.Fill((Double_t) object->size<JDAQTriggeredHit>());
    hs.Fill((Double_t) object->size<JDAQSnapshotHit> ());
    
    for (unsigned int i = 0; i != NUMBER_OF_TRIGGER_BITS; ++i) {
      if (object->hasTriggerBit(i)) {
        hm.Fill((Double_t) i);
      }
    }

    {
      unsigned int n = 0;

      typedef JDAQSnapshotHit JHit_t;
      
      for (JDAQEvent::const_iterator<JHit_t> i = object->begin<JHit_t>(); i != object->end<JHit_t>(); ++i) {
        
        if (object->getTriggerMask(*i) != 0) {
          ++n;
        }
      }
      
      if (object->size<JHit_t>() != 0 && n != object->size<JDAQTriggeredHit>()) {
        ERROR("Number of snapshot hits with trigger mask " << n << " != " << object->size<JDAQTriggeredHit>() << endl);
      }
    }
  }


  long long int numberOfSummaryslices = 0;

  for (JMultipleFileScanner<JDAQSummaryslice> input = inputFile; input.hasNext(); ++numberOfSummaryslices) {

    STATUS("entry: " << setw(10) << input.getCounter() << '\r'); DEBUG(endl);

    const JDAQSummaryslice* object = input.next();

    h2.Fill((Double_t) getSizeof(*object) * 1e-6);

    for (JDAQSummaryslice::const_iterator i = object->begin(); i != object->end(); ++i) {
      for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
        hr.Fill(i->getRate(pmt) * 1.0e-3);
      }
    }
  }


  if (debug >= JEEP::status_t) {
    print(h1, cout);
    print(h2, cout);
    print(hn, cout);
    print(hs, cout);
    print(hr, cout);
  }

  NOTICE(endl);
  NOTICE("Background rate estimate from snapshot hits " << setprecision(2) 
         << (hs.GetMean() - hn.GetMean()) * numberOfSummaryslices / (hr.GetEntries() * 2e-6 * trigger_parameters.TMaxEvent_ns)
         << " [kHz]" << endl);
  NOTICE(endl);


  const double T = numberOfSummaryslices * 1.0e-9 * getFrameTime();    // [s]

  NOTICE(setw(WIDTH) << left << "Total run duration (based on time slices) [s] " << FIXED(7,1) << T << endl);


  for (Int_t i = 1; i <= hm.GetNbinsX(); ++i) {

    const JTriggerbit_t x = (JTriggerbit_t) hm.GetBinCenter (i);
    const Double_t      y = (Double_t)      hm.GetBinContent(i);

    const char* const   name = getTriggerName(x);

    if (name != NULL || y != 0.0) {

      NOTICE("Trigger"
             << "[" << setw(2) << x << "]"                      << ' '
             << setw(16) << left << (name != NULL ? name : "?") << ' '
             << FIXED(7,0) << y);

      if (numberOfSummaryslices != 0) {
        NOTICE(FIXED(7,2) << y / T << " [Hz]");
      }

      NOTICE(endl);
    }
  }
  
  NOTICE("Total  "
         << setw(4)  << " " << ' '
         << setw(16) << " " << ' '
         << FIXED(7,0) << h1.GetEntries());

  if (numberOfSummaryslices != 0) {
    NOTICE(FIXED(7,2) << h1.GetEntries() / T << " [Hz]");
  }

  NOTICE(endl);

  if (outputFile != "") {

    TFile out(outputFile.c_str(), "recreate");

    out << h1 << h2;
    out << hm << hn << hs << hr;

    out.Write();
    out.Close();
  }

  
  for (Int_t i = 1; i <= hm.GetNbinsX(); ++i) {

    const JTriggerbit_t x = (JTriggerbit_t) hm.GetBinCenter (i);
    const Double_t      y = (Double_t)      hm.GetBinContent(i);

    const char* const name = getTriggerName(x);

    if (name != NULL) {
      QAQC(' ' << FIXED(7,0) << y);
    }
  }

  QAQC(endl);

  return 0;
}