#include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "TH2D.h" #include "TProfile.h" #include "JDAQ/JDAQSummarysliceIO.hh" #include "JDetector/JDetectorSupportkit.hh" #include "JTrigger/JTriggerToolkit.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JSupport.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Example program to histogram KM3NETDAQ::JDAQSummaryslice. * \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; bool correct; int debug; try { JParser<> zap("Example program to histogram summary data."); zap['f'] = make_field(inputFile); zap['o'] = make_field(outputFile) = "summaryslice.root"; zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['c'] = make_field(correct); zap['d'] = make_field(debug) = 2; zap(argc, argv); } catch(const exception& error) { FATAL(error.what() << endl); } const double factor = 1.0e-3; // [kHz] TFile out(outputFile.c_str(), "recreate"); TH1D h0("h0", NULL, JDAQRate::getN(), JDAQRate::getData(factor)); TProfile h1("h1", NULL, 32,-0.5, 31.5); TProfile h2("h2", NULL, 32,-0.5, 31.5); TH2D hu("hu", NULL, 201, -0.5, +200.5, 201, -0.5, +200.5); TH2D hv("hv", NULL, 51, -0.5, +50.5, 200, 0.0, 50.0); TH2D hw("hw", NULL, 51, -0.5, +50.5, 200, 0.0, 50.0); TH2D hx("hx", NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5, 200, 0.0, 50.0); TH2D hy("hy", NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5, 200, 0.0, 50.0); while (inputFile.hasNext()) { STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl); JDAQSummaryslice* summaryslice = inputFile.next(); DEBUG("Summary: " << setw(8) << inputFile.getCounter() << ' ' << setw(8) << summaryslice->getRunNumber() << ' ' << setw(8) << summaryslice->getFrameIndex() << ' ' << setw(6) << summaryslice->size() << endl); const JDetectorAddressMap& demo = getDetectorAddressMap(summaryslice->getDetectorID()); for (JDAQSummaryslice::const_iterator frame = summaryslice->begin(); frame != summaryslice->end(); ++frame) { const JModuleAddressMap& memo = demo.get(frame->getModuleID()); int N[2] = { 0 }; double R[2] = { 0.0 }; int lower[2] = { 0 }; int upper[2] = { 0 }; for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) { const int index = (!frame->testHighRateVeto(pmt) && !frame->testFIFOStatus(pmt) ? 0 : 1); N[index] += 1; R[index] += correct ? getRate(*frame, pmt, factor) : frame->getRate(pmt, factor); if (memo.getAddressTranslator(pmt).ring <= 'D') lower[index] += 1; else upper[index] += 1; h0.Fill(correct ? getRate(*frame, pmt, factor) : frame->getRate(pmt, factor), frame->getWeight(pmt, factor)); h1.Fill((Double_t) pmt, (frame->testHighRateVeto(pmt) ? 1.0 : 0.0)); h2.Fill((Double_t) pmt, (frame->testFIFOStatus (pmt) ? 1.0 : 0.0)); } hu.Fill((double) frame->getUDPMaximalSequenceNumber(), (double) frame->getUDPNumberOfReceivedPackets()); if (N[0] != 0) { hv.Fill((double) frame->getUDPMaximalSequenceNumber(), R[0] / N[0]); } if (N[1] != 0) { hw.Fill((double) frame->getUDPMaximalSequenceNumber(), R[1] / N[1]); } if (N[0] != 0) { hx.Fill((double) N[1], R[0] / N[0]); if (lower[0] != 0 && upper[0] != 0) { hy.Fill((double) N[1], R[0] / N[0]); } } const bool status = (frame->getUDPNumberOfReceivedPackets() == frame->getUDPMaximalSequenceNumber() + 1 && frame->hasUDPTrailer()); h1.Fill((Double_t) 31, (frame->testWhiteRabbitStatus() ? 1.0 : 0.0)); h2.Fill((Double_t) 31, (status ? 1.0 : 0.0)); } } STATUS(endl); out.Write(); out.Close(); }