#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 & 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::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 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::value ){ NBinsY=100; Ymax=25; }else if( ts_type == JIndexOf::value ){ NBinsY=100; Ymax=2; }else if( ts_type == JIndexOf::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 & 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 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 void initialize_timeslice_histograms (){ const int index = JIndexOf::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 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 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 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 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 * 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