#ifndef __JRUNHISTOGRAMS__
#define __JRUNHISTOGRAMS__
/**
* \author rgruiz, adomi
*/
#include "JSupport/JSupport.hh"
#include "JSupport/JSupportToolkit.hh"
#include "JTrigger/JTriggerInterface.hh"
#include "JTrigger/JTriggerBits.hh"
#include "JDetector/JModuleRouter.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "km3net-dataformat/online/JDAQClock.hh"
#include "km3net-dataformat/online/JDAQ.hh"
#include "km3net-dataformat/online/JDAQTimeslice.hh"
#include "JGizmo/JGizmoToolkit.hh"
#include "JROOT/JManager.hh"
#include "TH1D.h"
#include "TH2D.h"
#include "TProfile2D.h"
#include "TAxis.h"
#include "TCanvas.h"
#include "TPaveText.h"
#include "TStyle.h"
#include "TString.h"
#include "TMath.h"
#include "TColor.h"
#include "TDirectory.h"
#include "TPRegexp.h"
#include "TObjArray.h"
#include "TObjString.h"
using namespace std ;
using namespace KM3NETDAQ ;
using namespace JLANG ;
using namespace JPP ;
using namespace JSUPPORT ;
/*
* Structure to store histograms obtained from the JDAQSummary TTree.
*/
struct SummaryHistograms{
TH1D* h_mean_fifo;
TH2D* h_fifo_per_dom;
TH2D* h_daq_status_per_dom;
TH2D* h_hrv_per_dom;
TH2D* h_rate_summary;
TH1D* h_pmt_rate_distribution;
TH1D* h_dom_rate_distribution;
/* One histogram for each DU */
JManager < string , TH2D >* m_mean_summary_rate;
JManager < string , TH1D >* m_mean_summary_rate_distribution;
/* One histogram for each module */
JManager < string , TH2D >* m_summary_rate_distribution;
SummaryHistograms(){
h_fifo_per_dom = NULL;
h_daq_status_per_dom = NULL;
h_hrv_per_dom = NULL;
h_rate_summary = NULL;
h_pmt_rate_distribution = NULL;
h_dom_rate_distribution = NULL;
m_mean_summary_rate = NULL;
m_mean_summary_rate_distribution = NULL;
m_summary_rate_distribution = NULL;
}
/*
* Initializes the histograms for summary slices
*/
void initialize(std::set<int> & du_ids , int modules_per_string){
h_fifo_per_dom = new TH2D ("h_fifo_per_dom", " FIFO ; String ; Floor ; Number of slices with FIFO almost full x number of PMTs ",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 , modules_per_string , 0.5 , modules_per_string + 0.5 );
h_daq_status_per_dom = new TH2D ("h_daq_status_per_dom",
" Fraction of wrong DAQ Status [%]; String ; Floor ; Fraction of slices with wrong DAQ status of packets [%]",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 , modules_per_string , 0.5 , modules_per_string + 0.5 );
h_hrv_per_dom = new TH2D ("h_hrv_per_dom", "HRV ; String ; Floor ; Number of slices x number of PMTs in HRV",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 , modules_per_string , 0.5 , modules_per_string + 0.5 );
h_rate_summary = new TH2D ("h_rate_summary", "Summary slices ; String ; Floor ; Mean rate over all summary slices [kHz]",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 , modules_per_string , 0.5 , modules_per_string + 0.5 );
h_pmt_rate_distribution = new TH1D ("h_pmt_rate_distribution", "PMT rate distribution from summary slices ; rate [kHz] ; Counts",
JDAQRate::getN(), JDAQRate::getData(1.0e-4));
h_dom_rate_distribution = new TH1D ("h_dom_rate_distribution", "DOM rate distribution from summary slices ; rate [kHz] ; Counts", 50 , log10(50) , log10(2000));
setLogarithmicX (h_dom_rate_distribution);
TH1D* h = new TH1D("%/h_mean_summary_rate_distribution", " ; rate [kHz] ; # PMTs", 40 , 0 , log10(1000));
setLogarithmicX (h);
m_mean_summary_rate_distribution = new JManager < string , TH1D > (h);
m_mean_summary_rate = new JManager < string , TH2D > (new TH2D("%/h_mean_summary_rate", " ; TDC Channel ; Floor ; rate [kHz]",
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5,
modules_per_string , 0.5 , 0.5 + modules_per_string));
TH2D* h_summary_rate_distribution = new TH2D ("%/h_pmt_rate_distributions_Summaryslice", "Summaryslice ; TDC channel ; rate [kHz] ; counts",
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5, 100 , -1 , log10(10000));
setLogarithmicY (h_summary_rate_distribution);
m_summary_rate_distribution = new JManager < string , TH2D > (h_summary_rate_distribution);
}
};
struct TimesliceHistograms{
int min_logdt;
int max_logdt;
int nbins_logdt;
int nbins_time;
/* One histogram per timeslice type. I decided not to use a JManager here because the range of each histogram could be different for each timeslice type. */
vector < TH2D* > h_dom_mean_rates;
/* One JManager per timeslice type. Each JManager hosts a histogram for each DU */
vector < JManager < string , TH2D >* > m_mean_ToT;
vector < JManager < string , TH1D >* > m_mean_ToT_distribution;
/* One JManager per timeslice type. Each manager hosts a histogram for an optical module. The key is expected to follow the pattern SXXFXX */
vector < JManager < string , TH2D >* > m_pmt_tot_distributions;
vector < JManager < string , TH2D >* > m_pmt_rate_distributions;
TimesliceHistograms():
min_logdt (0),
max_logdt (9),
nbins_logdt (150),
nbins_time (200)
{
int number_of_timeslice_types = JLength<JDAQTimesliceTypes_t>::value ;
h_dom_mean_rates .resize (number_of_timeslice_types , NULL);
m_mean_ToT .resize (number_of_timeslice_types , NULL);
m_mean_ToT_distribution .resize (number_of_timeslice_types , NULL);
m_pmt_tot_distributions .resize (number_of_timeslice_types , NULL);
m_pmt_rate_distributions .resize (number_of_timeslice_types , NULL);
}
/*
* Initializes histograms for a given timeslice type.
*
* \param du_ids The list of ids for the DUs in the detector
* \param modules_per_string the number of modules in a string
* \param ts_type Index of the timeslice types on the JDAQTimesliceTypes_t typelist.
* \param ts_name The name of the timeslice type
*/
void initialize(std::set<int> du_ids , int modules_per_string , int ts_type , std::string ts_name){
h_dom_mean_rates[ts_type] = new TH2D (MAKE_STRING ("h_mean_dom_rates_" + ts_name).c_str(),
MAKE_STRING (ts_name + " ; DU number ; Floor number ; time slice averaged rate [Hz]").c_str() ,
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 ,
modules_per_string , 0.5 , 0.5 + modules_per_string);
m_mean_ToT[ts_type] = new JManager < string , TH2D > (new TH2D (MAKE_STRING ("%/h_mean_ToT_" + ts_name).c_str(),
MAKE_STRING (ts_name + " ; TDC channel ; Floor number ; mean ToT [ns] ").c_str(),
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5 ,
modules_per_string , 0.5 , 0.5 + modules_per_string));
m_mean_ToT_distribution[ts_type] = new JManager < string , TH1D > (new TH1D (MAKE_STRING ("%/h_mean_ToT_distribution" + ts_name).c_str(),
MAKE_STRING (ts_name + " ; ToT [ns] ; # PMTS ").c_str(),
256, -0.5, 255.5));
TH2D* h = new TH2D (MAKE_STRING ("%_" + ts_name + "_2SToT").c_str(),
MAKE_STRING (ts_name + " ; TDC channel ; ToT [ns] ; counts").c_str(),
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5 , 256, -0.5, 255.5);
h->Sumw2();
m_pmt_tot_distributions[ts_type] = new JManager < string , TH2D > (h);
int NBinsY=0;
double Ymax=0;
if( ts_type == JIndexOf<JDAQTimesliceTypes_t , JDAQTimeslice>::value ){
NBinsY=100;
Ymax=25;
}else if( ts_type == JIndexOf<JDAQTimesliceTypes_t , JDAQTimesliceL1>::value ){
NBinsY=100;
Ymax=2;
}else if( ts_type == JIndexOf<JDAQTimesliceTypes_t , JDAQTimesliceSN>::value ){
NBinsY=10;
Ymax=0.1;
}
TH2D* h_pmt_rate_distributions = new TH2D (MAKE_STRING ("%/h_pmt_rate_distributions_" + ts_name).c_str(),
MAKE_STRING (ts_name + " ; TDC channel ; rate [kHz] ; counts ").c_str(),
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5 ,
NBinsY, 0, Ymax);
m_pmt_rate_distributions[ts_type] = new JManager < string , TH2D > (h_pmt_rate_distributions);
}
/*
* Fills the mean ToT as a function of the PMT and floor number for a given DU
*
* \param table table with the ToT distributions for each PMT in a module, for every timeslice type
* \param string The string number
* \param floor The floor number
*/
void Fill_mean_ToT_histograms(){
int i = 0 ;
for (typename vector < JManager < string , TH2D >* >::const_iterator it = m_pmt_tot_distributions.begin() ; it != m_pmt_tot_distributions.end() ; ++it , ++i){
if ((*it)){
for (typename JManager < string , TH2D >::const_iterator j = (*it) -> begin() ; j != (*it) -> end() ; ++j){
TString s (MAKE_STRING(j -> first).c_str());
TPRegexp r ("(\\w+)/(\\DU)(\\d+)/(F)(\\d+)");
TObjArray* o = r.MatchS(s);
int String = ((TObjString *)o->At(3))->GetString().Atoi();
int Floor = ((TObjString *)o->At(5))->GetString().Atoi();
for (int pmt = 1 ; pmt <= (j -> second) -> GetXaxis() -> GetNbins() ; pmt++){
(*m_mean_ToT[i])[MAKE_STRING("Detector/DU" + to_string(String))] -> Fill((j->second) -> GetXaxis() -> GetBinCenter(pmt) , Floor , (j -> second) -> ProjectionY ("" , pmt , pmt) -> GetMean () );
(*m_mean_ToT_distribution[i])[MAKE_STRING("Detector/DU" + to_string(String))] -> Fill((j -> second) -> ProjectionY ("" , pmt , pmt) -> GetMean () );
}
}
}
}
}
};
/*
* Structure to store histograms obtained from the JDAQEvent TTree
*/
struct TriggerHistograms{
TH1D* h_Trigger_bit_event;
TH1D* h_Trigger_bit_hit;
TH1D* h_Snapshot_hits;
TH1D* h_Triggered_hits;
TH1D* h_Triggered_hits_3dmuon;
TH2D* h_Triggered_hits_3dmuon_per_module;
TH1D* h_Triggered_over_Snapshot_hits;
TH1D* h_event_duration;
TH1D* h_Number_of_overlays;
TH2D* h_Snapshot_hits_per_module;
TH2D* h_Triggered_hits_per_module;
TH1D* h_pmt_distribution_triggered_hits;
TH1D* h_tot_distribution_triggered_hits;
TH1D* h_pmt_distribution_snapshot_hits;
TH1D* h_tot_distribution_snapshot_hits;
/* One histogram per DU*/
JManager < string , TH2D >* m_Snapshot_hits_per_pmt;
/*
* Constructor
*/
TriggerHistograms(){
h_Trigger_bit_event = NULL;
h_Trigger_bit_hit = NULL;
h_Snapshot_hits = NULL;
h_Triggered_hits = NULL;
h_Triggered_hits_3dmuon = NULL;
h_Triggered_hits_3dmuon_per_module = NULL;
h_Triggered_over_Snapshot_hits = NULL;
h_event_duration = NULL;
h_Number_of_overlays = NULL;
h_Snapshot_hits_per_module = NULL;
h_Triggered_hits_per_module = NULL;
m_Snapshot_hits_per_pmt = NULL;
h_pmt_distribution_triggered_hits = NULL;
h_tot_distribution_triggered_hits = NULL;
h_pmt_distribution_snapshot_hits = NULL;
h_tot_distribution_snapshot_hits = NULL;
}
/*
* Initializes the histograms.
* \param du_ids The list of du ids in the detector
* \param frame_index_range The range of frame indices
* \param modules_per_string The number of modules in a string.
*/
void initialize(std::set<int> & du_ids , int modules_per_string){
h_Trigger_bit_event = new TH1D ("h_Trigger_bit_event",
"Number of events as a function of trigger bit in event ; Trigger Bit ; Counts",
NUMBER_OF_TRIGGER_BITS , -0.5, NUMBER_OF_TRIGGER_BITS - 0.5 );
for (int i = 0 ; i < (int) NUMBER_OF_TRIGGER_BITS ; i++){
if (getTriggerName(i)){
h_Trigger_bit_event -> GetXaxis() -> SetBinLabel(i+1 , getTriggerName(i) );
}else{
h_Trigger_bit_event -> GetXaxis() -> SetBinLabel(i+1 , "" );
}
}
h_Trigger_bit_event -> GetXaxis() -> LabelsOption("v");
h_Trigger_bit_hit = new TH1D ("h_Trigger_bit_hit",
"Number of hits per event as a function of trigger bit in hit ; Trigger Bit ; #Events",
NUMBER_OF_TRIGGER_BITS, -0.5, NUMBER_OF_TRIGGER_BITS - 0.5 );
for (int i = 0 ; i < (int) NUMBER_OF_TRIGGER_BITS ; i++){
if (getTriggerName(i)){
h_Trigger_bit_hit -> GetXaxis() -> SetBinLabel(i+1 , getTriggerName(i) );
}else{
h_Trigger_bit_hit -> GetXaxis() -> SetBinLabel(i+1 , "" );
}
}
h_Trigger_bit_hit -> GetXaxis() -> LabelsOption("v");
h_Snapshot_hits = new TH1D ("h_Snapshot_hits", " ; Number of snapshot hits; Events ", 50, 0, 4 );
setLogarithmicX (h_Snapshot_hits);
h_Triggered_hits = new TH1D ("h_Triggered_hits", " ; Number of triggered hits; Events ", 50 , 0, 4 );
setLogarithmicX (h_Triggered_hits);
h_Triggered_hits_3dmuon = new TH1D ("h_Triggered_hits_3dmuon", " ; Number of triggered hits for JTRIGGER3DMUON; Events ", 50 , 0, 3 );
setLogarithmicX(h_Triggered_hits_3dmuon);
h_Triggered_hits_3dmuon_per_module = new TH2D ("h_Triggered_hits_3dmuon_per_module",
"Number of triggered hits for JTRIGGER3DMUON; String ; Floor ; Number of JTRIGGER3DMUON hits",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5,
modules_per_string , 0.5 , modules_per_string + 0.5 );
h_Triggered_over_Snapshot_hits = new TH1D ("h_Triggered_over_Snapshot_hits", " ; Triggered/Snapshot hits; Events", 100 , 0, 0.5 );
h_event_duration = new TH1D ("h_event_duration", " ; Event Duration [ns]; Events", 60 , 1, 6 );
setLogarithmicX(h_event_duration);
h_Number_of_overlays = new TH1D ("h_Number_of_overlays", " ; Number of overlays; Events ", 1000, -0.5, 1000 - 0.5 );
h_Snapshot_hits_per_module = new TH2D ("h_Snapshot_hits_per_module",
" ; String ; Floor ; Number of snapshot hits ",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 ,
modules_per_string , 0.5 , modules_per_string + 0.5 );
h_Triggered_hits_per_module = new TH2D ("h_Triggered_hits_per_module",
" ; String ; Floor ; Number of triggered hits",
*du_ids.rbegin() , 0.5 , *du_ids.rbegin() + 0.5 ,
modules_per_string , 0.5 , modules_per_string + 0.5 );
m_Snapshot_hits_per_pmt = new JManager < string , TH2D > ( new TH2D ("%/h_Snapshot_hits_per_pmt",
" ; TDC Channel ; Floor ; Number of snapshot hits",
NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5 ,
modules_per_string , 0.5 , modules_per_string + 0.5 ) );
h_pmt_distribution_triggered_hits = new TH1D ("h_pmt_distribution_triggered_hits",
" ; TDC Channel ; Counts [a.u.]", NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5);
h_tot_distribution_triggered_hits = new TH1D ("h_tot_distribution_triggered_hits", " ; ToT [ns] ; Counts [a.u.]", 256 , -0.5, 255.5);
h_pmt_distribution_snapshot_hits = new TH1D ("h_pmt_distribution_snapshot_hits", " ; TDC Channel ; Counts [a.u.]", NUMBER_OF_PMTS , -0.5 , NUMBER_OF_PMTS - 0.5);
h_tot_distribution_snapshot_hits = new TH1D ("h_tot_distribution_snapshot_hits", " ; ToT [ns] ; Counts [a.u.]", 256, -0.5, 255.5);
}
};
/*
* Class to manage the histograms produced by JRunAnalyzer.
*/
class JRA_Histograms{
public:
SummaryHistograms h_summary;
TimesliceHistograms h_timeslice;
TriggerHistograms h_trigger;
std::set <int> du_ids;
int modules_per_string;
JRA_Histograms(){}
JRA_Histograms (const JDetector& detector){
du_ids = getStringIDs(detector);
modules_per_string = JDETECTOR::getNumberOfModules (detector) / getNumberOfStrings (detector);
h_summary = SummaryHistograms();
h_timeslice = TimesliceHistograms();
h_trigger = TriggerHistograms();
}
/*
* Initializes summary slice histograms.
*
* \param range The range of values for the summary slice indices in the run. The range of some histograms is based on the number of time slices, and their indices.
*/
void initialize_summary_histograms (){
h_summary.initialize(du_ids , modules_per_string);
}
/*
* Initializes summary slice histograms.
*
* \param range The range of values for the summary slice indices in the run. The range of some histograms is based on the number of time slices, and their indices.
*/
template <class T>
void initialize_timeslice_histograms (){
const int index = JIndexOf<JDAQTimesliceTypes_t , T>::value ;
const string prefix = "KM3NETDAQ::JDAQ" ;
string ts_name = T::Class_Name();
string::size_type pos = ts_name.find(prefix);
if (pos != string::npos) ts_name.replace(ts_name.find(prefix) , prefix.length() , "");
h_timeslice.initialize(du_ids , modules_per_string , index , ts_name);
}
/*
* Initializes JDAQEvent histograms.
*/
void initialize_trigger_histograms (){
h_trigger.initialize (du_ids , modules_per_string);
}
/*
* Checks if the histograms in a table have been initialized. If yes, it writes them on the indicated directory of a root file
*
* \param f The root file
* \param dirname The directory where the histograms will be written
* \param table A vector of histograms.
*/
template <class T>
void Write_histogram_table_to_file(TFile & f , string dirname , vector < vector < T* > > table){
if(f.GetDirectory(dirname.c_str()) == 0) f.mkdir (dirname.c_str());
f.cd (dirname.c_str());
for (int i=0 ; i < (int)table.size(); i++){
for (int j=0 ; j< (int)table[i].size(); j++){
if (table[i][j]) table [i][j] -> Write();
}
}
}
/*
* Checks if the histograms in a table have been initialized. If yes, it writes them on the indicated directory of a root file
*
* \param f The root file
* \param dirname The directory where the histograms will be written
* \param table A vector of histograms.
*/
template <class T>
void Write_histogram_table_to_file(TFile & f , string dirname , vector < T* > table){
if(f.GetDirectory(dirname.c_str()) == 0) f.mkdir (dirname.c_str());
f.cd (dirname.c_str());
for (int i=0 ; i < (int)table.size(); i++){
if (table[i]) table[i] -> Write();
}
}
/*
* Checks if the histograms in a table have been initialized. If yes, it writes them on the indicated directory of a root file
*
* \param f The root file
* \param dirname The directory where the histograms will be written
* \param table A vector of histograms.
*/
template <class T , class V>
void Write_manager_to_file(TFile & f , string dirname , JManager < T , V >* table){
if(f.GetDirectory(dirname.c_str()) == 0) f.mkdir (dirname.c_str());
f.cd (dirname.c_str());
for (typename JManager < T , V >::const_iterator i = table -> begin() ; i != table -> end() ; ++i){
i -> second -> Write();
}
}
/*
* Replaces wildcard in manager objects titles by their keys.
* \param A manager.
* \param wc The wildcard
*/
template <class T , class V>
void Replace_wildcard_in_name(JManager < T , V >* manager , char wc = '%'){
for (typename JManager < T , V >::const_iterator i = manager -> begin() ; i != manager -> end() ; ++i){
if (i -> second -> GetTitle()){
std::string buffer = i -> second -> GetTitle();
string::size_type ipos = buffer.find(wc);
if (ipos!=std::string::npos){
ostringstream os;
os << i -> first ;
buffer.replace(ipos, 1, os.str());
i -> second -> SetTitle(buffer.c_str());
}
}
}
}
/*
* Writes the contents of a JManager into a file. Each object in each will be written in a directory specified by its key.
* \param f The root file
* \param A manager.
*/
template < class T , class V >
void Write_manager_in_key_dir(TFile & f ,JManager <T , V>* manager){
for (typename JManager < T , V >::const_iterator i = manager -> begin() ; i != manager -> end() ; ++i){
std::string fullpath = MAKE_STRING(i->second->GetName());
int pos = fullpath.rfind ('/');
std::string name = fullpath.substr (pos + 1);
std::string path = fullpath.substr (0 , pos);
if (f.GetDirectory(path.c_str()) == 0) f.mkdir (path.c_str());
f.cd(path.c_str());
i -> second -> SetName(name.c_str());
i -> second -> Write();
}
}
/*
* Writes the contents of a vector of JManager objects into a file. Each object in each will be written in a directory specified by its key.
* \param f The root file
* \param table the list of JManager objects.
*/
template < class T , class V >
void Write_manager_table_in_key_dir(TFile & f , vector < JManager < T , V >* > table){
for (typename vector < JManager < T , V >* >::const_iterator i = table.begin() ; i != table.end() ; ++i){
if ((*i)){
for (typename JManager < T , V >::const_iterator j = (*i) -> begin() ; j != (*i) -> end() ; ++j){
std::string fullpath = MAKE_STRING(j-> second -> GetName());
int pos = fullpath.rfind ('/');
std::string name = fullpath.substr (pos + 1);
std::string path = fullpath.substr (0 , pos);
if (f.GetDirectory(path.c_str()) == 0) f.mkdir (path.c_str());
f.cd(path.c_str());
j -> second -> SetName(name.c_str());
j -> second -> Write();
}
}
}
}
};
#endif