rtang/ANASEN_analysis
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

modified: Armory/HistPlotter.h

Browse Source 
modified:   TrackRecon.C made changes to account for frequent crashes due to the number of histograms. New changes allow for flsuhing of data to histograms every n events
	modified:   run_27Al.sh
This commit is contained in:
Vignesh Sitaraman 2026-07-31 14:05:18 -04:00
parent a292695a4d
commit 14999546b7
3 changed files with 279 additions and 163 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

194
Armory/HistPlotter.h
View File

@ -20,10 +20,15 @@
#include <set> #include <set>
#include <TGraphErrors.h> #include <TGraphErrors.h>
class HistPlotter { class HistPlotter
{
private: private:
long long barrier_count, barrier_limit; // meant to keep track of how often to call FillN() on histograms long long barrier_count, barrier_limit; // meant to keep track of how often to call FillN() on histograms
enum {TFILE, TMEMFILE} filetype; enum
{
TFILE,
TMEMFILE
} filetype;
std::map<std::string, TObject *> oMap; //!< Maps std::string to all TH1, TH2 objects in the class std::map<std::string, TObject *> oMap; //!< Maps std::string to all TH1, TH2 objects in the class
std::unordered_map<std::string, TObject *> cutsMap; //!< Maps std::string to TCutG objects held by the class std::unordered_map<std::string, TObject *> cutsMap; //!< Maps std::string to TCutG objects held by the class
std::set<std::string> folderList; //!< List of all folder names used to nest objects std::set<std::string> folderList; //!< List of all folder names used to nest objects
@ -35,22 +40,32 @@ private:
std::map<std::string, std::vector<double>> onedimcache; std::map<std::string, std::vector<double>> onedimcache;
std::map<std::string, std::pair<std::vector<double>, std::vector<double>>> twodimcache; std::map<std::string, std::pair<std::vector<double>, std::vector<double>>> twodimcache;
inline void FillN_All_Histograms(); inline void FillN_All_Histograms();
public: public:
HistPlotter(std::string outfile, std::string type); HistPlotter(std::string outfile, std::string type);
inline void FlushToDisk(int integral); //!< Writes all objects to file before closing, nesting objects in folders as is found necessary inline void FlushToDisk(int integral); //!< Writes all objects to file before closing, nesting objects in folders as is found necessary
inline void PrintObjects(); //!< Dump objects to std::cout for inspection inline void PrintObjects(); //!< Dump objects to std::cout for inspection
inline void ReadCuts(std::string); inline void ReadCuts(std::string);
inline TCutG* FindCut(std::string cut) { inline TCutG *FindCut(std::string cut)
{
return static_cast<TCutG *>(cutsMap.at(cut)); return static_cast<TCutG *>(cutsMap.at(cut));
} }
inline void set_barrier_limit(long long limit) { barrier_limit = limit; } inline void set_barrier_limit(long long limit) { barrier_limit = limit; }
inline void barrier_increment() { inline void barrier_increment()
{
barrier_count++; barrier_count++;
if(barrier_count == barrier_limit) { if (barrier_count == barrier_limit)
{
FillN_All_Histograms(); FillN_All_Histograms();
barrier_count = 0; barrier_count = 0;
} }
} }
inline void force_flush_caches()
{
FillN_All_Histograms();
barrier_count = 0;
}
/*! \fn void FindCut() /*! \fn void FindCut()
\brief \brief
- Searches for a cut by name 'cut' in the internal list of cuts 'cutsMap'. Ugly fails (via unresolved at()) if such a cut isn't found. - Searches for a cut by name 'cut' in the internal list of cuts 'cutsMap'. Ugly fails (via unresolved at()) if such a cut isn't found.
@ -58,9 +73,11 @@ public:
\return Pointer to the TCutG object that matches the name. Very useful to use this as plotter.FindCut("protonbarrelpid")->IsInside(deltaE, E) for instance. \return Pointer to the TCutG object that matches the name. Very useful to use this as plotter.FindCut("protonbarrelpid")->IsInside(deltaE, E) for instance.
*/ */
inline void SetNewTitle(std::string name, std::string title) { inline void SetNewTitle(std::string name, std::string title)
{
auto result = oMap.find(name); // result is an iterator auto result = oMap.find(name); // result is an iterator
if(result==oMap.end()) return; //no warnings, could be changed in future if (result == oMap.end())
return; // no warnings, could be changed in future
else else
static_cast<TNamed *>(oMap.at(name))->SetTitle(title.c_str()); // set new title static_cast<TNamed *>(oMap.at(name))->SetTitle(title.c_str()); // set new title
} }
@ -68,27 +85,31 @@ public:
// Smart functions that create a new histogram if it doesn't exist. // Smart functions that create a new histogram if it doesn't exist.
inline void FillGraph(const std::string &name, float valuex, float valuey, float errx = 0, float erry = 0); inline void FillGraph(const std::string &name, float valuex, float valuey, float errx = 0, float erry = 0);
inline void Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value); inline void Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value);
inline void Fill2D(const std::string& name,int nbinsx, float xlow, float xhigh inline void Fill2D(const std::string &name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey);
,int nbinsy, float ylow, float yhigh, float valuex, float valuey);
inline void Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value, const std::string &folder); inline void Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value, const std::string &folder);
inline void Fill2D(const std::string& name,int nbinsx, float xlow, float xhigh inline void Fill2D(const std::string &name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey, const std::string &folder);
,int nbinsy, float ylow, float yhigh, float valuex, float valuey, const std::string& folder);
// TObject* findObject(std::string key); // TObject* findObject(std::string key);
}; };
HistPlotter::HistPlotter(std::string outfile, std::string type="") { HistPlotter::HistPlotter(std::string outfile, std::string type = "")
{
/*! /*!
\brief Constructor. Opens a TFile instance with the specified filename \brief Constructor. Opens a TFile instance with the specified filename
\param outfile : std::string that holds the desired output ROOT filename \param outfile : std::string that holds the desired output ROOT filename
\return None \return None
*/ */
if(type=="" || type == "TFILE") { if (type == "" || type == "TFILE")
{
ofile = new TFile(outfile.c_str(), "recreate"); ofile = new TFile(outfile.c_str(), "recreate");
filetype = TFILE; filetype = TFILE;
} else if(type =="TMEMFILE") { }
else if (type == "TMEMFILE")
{
omfile = new TMemFile(outfile.c_str(), "recreate"); omfile = new TMemFile(outfile.c_str(), "recreate");
filetype = TMEMFILE; filetype = TMEMFILE;
} else { }
else
{
std::cout << "Unknown type " << type << " specified for HistPlotter (use \"TFILE\" or \"TMEMFILE\"), using default \"TFILE\" " << std::endl; std::cout << "Unknown type " << type << " specified for HistPlotter (use \"TFILE\" or \"TMEMFILE\"), using default \"TFILE\" " << std::endl;
ofile = new TFile(outfile.c_str(), "recreate"); ofile = new TFile(outfile.c_str(), "recreate");
filetype = TFILE; filetype = TFILE;
@ -97,16 +118,21 @@ HistPlotter::HistPlotter(std::string outfile, std::string type="") {
barrier_limit = 1000; barrier_limit = 1000;
} }
void HistPlotter::FillN_All_Histograms() { void HistPlotter::FillN_All_Histograms()
for(auto it=oMap.begin(); it!=oMap.end(); it++ ) { {
for (auto it = oMap.begin(); it != oMap.end(); it++)
{
// it->first is std::string 'name', it->second is the TObject // it->first is std::string 'name', it->second is the TObject
if(it->second->InheritsFrom("TH1F")) { if (it->second->InheritsFrom("TH1F"))
{
// FillN(size, array-of-doubles, array-of-weights); //we set array-of-weights to (1,1,1,.. (size) // FillN(size, array-of-doubles, array-of-weights); //we set array-of-weights to (1,1,1,.. (size)
static_cast<TH1F *>(it->second)->FillN(onedimcache[it->first].size(), // size static_cast<TH1F *>(it->second)->FillN(onedimcache[it->first].size(), // size
onedimcache[it->first].data(), // array onedimcache[it->first].data(), // array
std::vector<double>(onedimcache[it->first].size(), 1.0).data()); // weight of ones std::vector<double>(onedimcache[it->first].size(), 1.0).data()); // weight of ones
onedimcache[it->first].clear(); onedimcache[it->first].clear();
} else if(it->second->InheritsFrom("TH2F")) { }
else if (it->second->InheritsFrom("TH2F"))
{
// FillN(size, array-of-doubles, array-of-weights); //we set array-of-weights to (1,1,1,.. (size)) // FillN(size, array-of-doubles, array-of-weights); //we set array-of-weights to (1,1,1,.. (size))
static_cast<TH2F *>(it->second)->FillN(twodimcache[it->first].first.size(), // size static_cast<TH2F *>(it->second)->FillN(twodimcache[it->first].first.size(), // size
twodimcache[it->first].first.data(), // x array twodimcache[it->first].first.data(), // x array
@ -119,30 +145,37 @@ void HistPlotter::FillN_All_Histograms() {
std::cout << "." << std::endl; std::cout << "." << std::endl;
} }
void HistPlotter::FlushToDisk(int min_integral=0) { void HistPlotter::FlushToDisk(int min_integral = 0)
{
/*! \fn void FlushToDisk() /*! \fn void FlushToDisk()
\brief Function that can be used at any point to exit smoothly by saving all ROOT objects in memory \brief Function that can be used at any point to exit smoothly by saving all ROOT objects in memory
to the output file before closing it. Obeys the binding of histograms to separate folders, if so specified. to the output file before closing it. Obeys the binding of histograms to separate folders, if so specified.
\return No return -- void \return No return -- void
*/ */
if(filetype==TMEMFILE && omfile) { if (filetype == TMEMFILE && omfile)
{
std::cout << "Not flushing a TMemfile .. exiting .." << std::endl; std::cout << "Not flushing a TMemfile .. exiting .." << std::endl;
delete omfile; delete omfile;
return; return;
} }
if(ofile->IsZombie() || !ofile) { if (ofile->IsZombie() || !ofile)
{
std::cerr << "Output file is zombie, finishing up without writing to disk!" << std::endl; std::cerr << "Output file is zombie, finishing up without writing to disk!" << std::endl;
return; return;
} }
FillN_All_Histograms(); FillN_All_Histograms();
for(auto it=oMap.begin(); it!=oMap.end(); it++ ) { for (auto it = oMap.begin(); it != oMap.end(); it++)
{
// omap maps: name(first) to object address(second). // omap maps: name(first) to object address(second).
// foldersForObjects maps: object address(first) to foldername(second) // foldersForObjects maps: object address(first) to foldername(second)
auto result = foldersForObjects.find(it->second); // returns <TObject* histogram,std::string foldername> pair if found auto result = foldersForObjects.find(it->second); // returns <TObject* histogram,std::string foldername> pair if found
if(result!=foldersForObjects.end()) { //we try to create folder if needed and cd to it if (result != foldersForObjects.end())
{ // we try to create folder if needed and cd to it
ofile->mkdir(result->second.c_str(), "", kTRUE); // args: name, title, returnExistingDirectory ofile->mkdir(result->second.c_str(), "", kTRUE); // args: name, title, returnExistingDirectory
ofile->cd(result->second.c_str()); ofile->cd(result->second.c_str());
} else { }
else
{
ofile->cd(); // toplevel for all default histograms. Default setting ofile->cd(); // toplevel for all default histograms. Default setting
} }
if (((TH1F *)it->second)->Integral() > min_integral) if (((TH1F *)it->second)->Integral() > min_integral)
@ -152,7 +185,8 @@ void HistPlotter::FlushToDisk(int min_integral=0) {
// Create a directory for all cuts, and save all cuts in them // Create a directory for all cuts, and save all cuts in them
ofile->mkdir("gCUTS", "", kTRUE); ofile->mkdir("gCUTS", "", kTRUE);
ofile->cd("gCUTS"); ofile->cd("gCUTS");
for(auto it=cutsMap.begin(); it!=cutsMap.end(); it++) { for (auto it = cutsMap.begin(); it != cutsMap.end(); it++)
{
(static_cast<TNamed *>(it->second))->SetName(it->first.c_str()); (static_cast<TNamed *>(it->second))->SetName(it->first.c_str());
it->second->Write(); it->second->Write();
} }
@ -160,7 +194,8 @@ void HistPlotter::FlushToDisk(int min_integral=0) {
std::cout << "Wrote " << oMap.size() << " histograms to TFile " << std::string(ofile->GetName()) << std::endl; std::cout << "Wrote " << oMap.size() << " histograms to TFile " << std::string(ofile->GetName()) << std::endl;
} }
void HistPlotter::FillGraph(const std::string& name, float valuex, float valuey, float errx, float erry) { void HistPlotter::FillGraph(const std::string &name, float valuex, float valuey, float errx, float erry)
{
/*! \fn void FillGraph() /*! \fn void FillGraph()
\brief \brief
- Creates a TGraphError in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey - Creates a TGraphError in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey
@ -174,12 +209,14 @@ void HistPlotter::FillGraph(const std::string& name, float valuex, float valuey,
\return No return void \return No return void
*/ */
auto result = oMap.find(name); auto result = oMap.find(name);
if(result==oMap.end()) { if (result == oMap.end())
{
TGraphErrors *tempG = new TGraphErrors(); TGraphErrors *tempG = new TGraphErrors();
tempG->SetName(name.c_str()); tempG->SetName(name.c_str());
oMap.insert(std::make_pair(name, static_cast<TObject *>(tempG))); oMap.insert(std::make_pair(name, static_cast<TObject *>(tempG)));
} }
if(!oMap.at(name)->InheritsFrom("TGraphErrors")) { if (!oMap.at(name)->InheritsFrom("TGraphErrors"))
{
std::cerr << "Object " << name << " refers to something other than a TGraph*, not filling it hence!" << std::endl; std::cerr << "Object " << name << " refers to something other than a TGraph*, not filling it hence!" << std::endl;
std::cerr << "Abort.." << std::endl; std::cerr << "Abort.." << std::endl;
FlushToDisk(); FlushToDisk();
@ -188,7 +225,8 @@ void HistPlotter::FillGraph(const std::string& name, float valuex, float valuey,
// static_cast<TGraphErrors*>(oMap.at(name))->AddPointError(valuex,valuey,errx,erry); // static_cast<TGraphErrors*>(oMap.at(name))->AddPointError(valuex,valuey,errx,erry);
} }
void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float xhigh, float value) { void HistPlotter::Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value)
{
/*! \fn void Fill1D() /*! \fn void Fill1D()
\brief \brief
- Creates a TH1F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey - Creates a TH1F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey
@ -202,17 +240,21 @@ void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float
\return No return void \return No return void
*/ */
auto result = oMap.find(name); // result is an iterator auto result = oMap.find(name); // result is an iterator
if(result==oMap.end()) { if (result == oMap.end())
{
TH1F *temp1D = new TH1F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh); TH1F *temp1D = new TH1F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh);
oMap.insert(std::make_pair(name, static_cast<TObject *>(temp1D))); oMap.insert(std::make_pair(name, static_cast<TObject *>(temp1D)));
onedimcache.insert(std::make_pair(name, std::vector<double>())); onedimcache.insert(std::make_pair(name, std::vector<double>()));
onedimcache[name].reserve(16384); onedimcache[name].reserve(16384);
} else if(foldersForObjects.find(oMap.at(name))!=foldersForObjects.end()) { //shouldn't have a folder associated with it }
else if (foldersForObjects.find(oMap.at(name)) != foldersForObjects.end())
{ // shouldn't have a folder associated with it
std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in toplevel .." << std::endl; std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in toplevel .." << std::endl;
} }
// Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue // Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue
if(!oMap.at(name)->InheritsFrom("TH1F")) { if (!oMap.at(name)->InheritsFrom("TH1F"))
{
std::cerr << "Object " << name << " refers to something other than a TH1*, not filling it hence!" << std::endl; std::cerr << "Object " << name << " refers to something other than a TH1*, not filling it hence!" << std::endl;
std::cerr << "Abort.." << std::endl; std::cerr << "Abort.." << std::endl;
FlushToDisk(); FlushToDisk();
@ -222,7 +264,8 @@ void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float
// static_cast<TH1F*>(oMap.at(name))->Fill(value); // static_cast<TH1F*>(oMap.at(name))->Fill(value);
} }
void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float xhigh, float value, const std::string& foldername) { void HistPlotter::Fill1D(const std::string &name, int nbinsx, float xlow, float xhigh, float value, const std::string &foldername)
{
/*! \fn void Fill1D() /*! \fn void Fill1D()
\brief \brief
- Creates a TH1F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey - Creates a TH1F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey
@ -239,28 +282,37 @@ void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float
*/ */
auto result = oMap.find(name); // result is an iterator auto result = oMap.find(name); // result is an iterator
if(result==oMap.end()) { if (result == oMap.end())
{
TH1F *temp1D = new TH1F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh); TH1F *temp1D = new TH1F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh);
oMap.insert(std::make_pair(name, static_cast<TObject *>(temp1D))); oMap.insert(std::make_pair(name, static_cast<TObject *>(temp1D)));
onedimcache.insert(std::make_pair(name, std::vector<double>())); onedimcache.insert(std::make_pair(name, std::vector<double>()));
onedimcache[name].reserve(16384); onedimcache[name].reserve(16384);
if(foldername!="") { if (foldername != "")
if(folderList.find(foldername)==folderList.end()) { {
if (folderList.find(foldername) == folderList.end())
{
folderList.insert(foldername); folderList.insert(foldername);
} }
foldersForObjects.insert(std::make_pair(static_cast<TObject *>(temp1D), foldername)); foldersForObjects.insert(std::make_pair(static_cast<TObject *>(temp1D), foldername));
} }
} else { }
else
{
// object is present in map, but we enforce unique names // object is present in map, but we enforce unique names
// it must already have a folder attached to it // it must already have a folder attached to it
if(foldersForObjects.find(oMap.at(name))==foldersForObjects.end()) { if (foldersForObjects.find(oMap.at(name)) == foldersForObjects.end())
{
std::cerr << "Object " << name << " already registered at toplevel, choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl; std::cerr << "Object " << name << " already registered at toplevel, choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl;
} else if(foldersForObjects[oMap[name]]!=foldername) { }
else if (foldersForObjects[oMap[name]] != foldername)
{
std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl; std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl;
} }
} }
// Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue // Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue
if(!oMap.at(name)->InheritsFrom("TH1F")) { if (!oMap.at(name)->InheritsFrom("TH1F"))
{
std::cerr << "Object " << name << " refers to something other than a TH1*, not filling it hence!" << std::endl; std::cerr << "Object " << name << " refers to something other than a TH1*, not filling it hence!" << std::endl;
std::cerr << "Abort.." << std::endl; std::cerr << "Abort.." << std::endl;
FlushToDisk(); FlushToDisk();
@ -270,7 +322,8 @@ void HistPlotter::Fill1D(const std::string& name, int nbinsx, float xlow, float
// static_cast<TH1F*>(oMap.at(name))->Fill(value); // static_cast<TH1F*>(oMap.at(name))->Fill(value);
} }
void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey) { void HistPlotter::Fill2D(const std::string &name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey)
{
/*! \fn void Fill2D() /*! \fn void Fill2D()
\brief \brief
- Creates a TH2F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey - Creates a TH2F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey
@ -288,18 +341,22 @@ void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float
*/ */
auto result = oMap.find(name); // result is an iterator auto result = oMap.find(name); // result is an iterator
if(result==oMap.end()) { if (result == oMap.end())
{
TH2F *temp2D = new TH2F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh, nbinsy, ylow, yhigh); TH2F *temp2D = new TH2F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh, nbinsy, ylow, yhigh);
oMap.insert(std::make_pair(name, static_cast<TObject *>(temp2D))); oMap.insert(std::make_pair(name, static_cast<TObject *>(temp2D)));
twodimcache.insert(std::make_pair(name, std::make_pair(std::vector<double>(), std::vector<double>()))); twodimcache.insert(std::make_pair(name, std::make_pair(std::vector<double>(), std::vector<double>())));
twodimcache[name].first.reserve(16384); twodimcache[name].first.reserve(16384);
twodimcache[name].second.reserve(16384); twodimcache[name].second.reserve(16384);
} else if(foldersForObjects.find(oMap.at(name))!=foldersForObjects.end()) { //shouldn't have a folder associated with it }
else if (foldersForObjects.find(oMap.at(name)) != foldersForObjects.end())
{ // shouldn't have a folder associated with it
std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in toplevel .." << std::endl; std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in toplevel .." << std::endl;
} }
// Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue // Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue
if(!oMap.at(name)->InheritsFrom("TH2F")) { if (!oMap.at(name)->InheritsFrom("TH2F"))
{
std::cerr << "Object " << name << " refers to something other than a TH2*, not filling it hence!" << std::endl; std::cerr << "Object " << name << " refers to something other than a TH2*, not filling it hence!" << std::endl;
std::cerr << "Abort.." << std::endl; std::cerr << "Abort.." << std::endl;
FlushToDisk(); FlushToDisk();
@ -310,7 +367,8 @@ void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float
// static_cast<TH2F*>(oMap.at(name))->Fill(valuex,valuey); // static_cast<TH2F*>(oMap.at(name))->Fill(valuex,valuey);
} }
void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey, const std::string& foldername) { void HistPlotter::Fill2D(const std::string &name, int nbinsx, float xlow, float xhigh, int nbinsy, float ylow, float yhigh, float valuex, float valuey, const std::string &foldername)
{
/*! \fn void Fill2D() /*! \fn void Fill2D()
\brief \brief
- Creates a TH2F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey - Creates a TH2F in memory with name 'name' if it doesn't exist, and fills it with valuex, valuey
@ -331,30 +389,39 @@ void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float
*/ */
auto result = oMap.find(name); // result is an iterator auto result = oMap.find(name); // result is an iterator
if(result==oMap.end()) { if (result == oMap.end())
{
TH2F *temp2D = new TH2F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh, nbinsy, ylow, yhigh); TH2F *temp2D = new TH2F(name.c_str(), name.c_str(), nbinsx, xlow, xhigh, nbinsy, ylow, yhigh);
oMap.insert(std::make_pair(name, static_cast<TObject *>(temp2D))); oMap.insert(std::make_pair(name, static_cast<TObject *>(temp2D)));
twodimcache.insert(std::make_pair(name, std::make_pair(std::vector<double>(), std::vector<double>()))); twodimcache.insert(std::make_pair(name, std::make_pair(std::vector<double>(), std::vector<double>())));
twodimcache[name].first.reserve(16384); twodimcache[name].first.reserve(16384);
twodimcache[name].second.reserve(16384); twodimcache[name].second.reserve(16384);
if(foldername!="") { if (foldername != "")
if(folderList.find(foldername)==folderList.end()) { {
if (folderList.find(foldername) == folderList.end())
{
folderList.insert(foldername); folderList.insert(foldername);
} }
foldersForObjects.insert(std::make_pair(static_cast<TObject *>(temp2D), foldername)); foldersForObjects.insert(std::make_pair(static_cast<TObject *>(temp2D), foldername));
} }
} else { }
else
{
// object is present in map, but we enforce unique names // object is present in map, but we enforce unique names
// it must already have a folder attached to it // it must already have a folder attached to it
if(foldersForObjects.find(oMap.at(name))==foldersForObjects.end()) { if (foldersForObjects.find(oMap.at(name)) == foldersForObjects.end())
{
std::cerr << "Object " << name << " already registered at toplevel, choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl; std::cerr << "Object " << name << " already registered at toplevel, choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl;
} else if(foldersForObjects[oMap.at(name)]!=foldername) { }
else if (foldersForObjects[oMap.at(name)] != foldername)
{
std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl; std::cerr << "Object " << name << " already registered at " << foldersForObjects[oMap[name]] << ", choose a different name for the histogram to be stored in " << foldername << " folder .." << std::endl;
} }
} }
// Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue // Check if the string 'name' maps to a 1D hist. If there's any other object by this name raise issue
if(!oMap.at(name)->InheritsFrom("TH2F")) { if (!oMap.at(name)->InheritsFrom("TH2F"))
{
std::cerr << "Object " << name << " refers to something other than a TH2*, not filling it hence!" << std::endl; std::cerr << "Object " << name << " refers to something other than a TH2*, not filling it hence!" << std::endl;
std::cerr << "Abort.." << std::endl; std::cerr << "Abort.." << std::endl;
FlushToDisk(); FlushToDisk();
@ -365,7 +432,8 @@ void HistPlotter::Fill2D(const std::string& name, int nbinsx, float xlow, float
// static_cast<TH2F*>(oMap.at(name))->Fill(valuex,valuey); // static_cast<TH2F*>(oMap.at(name))->Fill(valuex,valuey);
} }
void HistPlotter::ReadCuts(std::string filename) { void HistPlotter::ReadCuts(std::string filename)
{
/*! \fn void ReadCuts() /*! \fn void ReadCuts()
\brief Reads a list of cuts from a file. The file must have the format below, two columns \brief Reads a list of cuts from a file. The file must have the format below, two columns
- Column#1 - path to a file that contains a single TCutG object named "CUTG", the default name in ROOT. - Column#1 - path to a file that contains a single TCutG object named "CUTG", the default name in ROOT.
@ -377,15 +445,18 @@ void HistPlotter::ReadCuts(std::string filename) {
std::ifstream infile; std::ifstream infile;
infile.open(filename); infile.open(filename);
std::string cutfilename, cutname; std::string cutfilename, cutname;
for(std::string line; std::getline(infile, line); ) { for (std::string line; std::getline(infile, line);)
{
if (line.size() != 0 && line[0] == '#') if (line.size() != 0 && line[0] == '#')
; // don't do anything with '#' lines ; // don't do anything with '#' lines
else { else
{
std::stringstream ss(line); std::stringstream ss(line);
ss >> cutfilename >> cutname; ss >> cutfilename >> cutname;
TFile f(cutfilename.c_str()); TFile f(cutfilename.c_str());
if(f.IsZombie()) { if (f.IsZombie())
{
std::cerr << "Cannot open cutfile " << cutfilename << " .. skipping.." << std::endl; std::cerr << "Cannot open cutfile " << cutfilename << " .. skipping.." << std::endl;
continue; continue;
} }
@ -397,7 +468,8 @@ void HistPlotter::ReadCuts(std::string filename) {
infile.close(); infile.close();
} }
void HistPlotter::PrintObjects() { void HistPlotter::PrintObjects()
{
/* /*
void PrintObjects() void PrintObjects()
Prints the contents of the unordered_maps oMap and cutsMap to facilitate debugging Prints the contents of the unordered_maps oMap and cutsMap to facilitate debugging
@ -405,10 +477,12 @@ void HistPlotter::PrintObjects() {
*/ */
std::cout << "Type | Name " << std::endl; std::cout << "Type | Name " << std::endl;
std::cout << "---- | --------------------- " << std::endl; std::cout << "---- | --------------------- " << std::endl;
for(auto it=oMap.begin(); it!=oMap.end(); it++ ) { for (auto it = oMap.begin(); it != oMap.end(); it++)
{
std::cout << it->second->ClassName() << " | " << it->first << std::endl; std::cout << it->second->ClassName() << " | " << it->first << std::endl;
} }
for(auto it=cutsMap.begin(); it!=cutsMap.end(); it++ ) { for (auto it = cutsMap.begin(); it != cutsMap.end(); it++)
{
std::cout << it->second->ClassName() << " | " << it->first << std::endl; std::cout << it->second->ClassName() << " | " << it->first << std::endl;
} }
std::cout << "---- | --------------------- " << std::endl; std::cout << "---- | --------------------- " << std::endl;

44
TrackRecon.C
View File

@ -56,7 +56,9 @@ bool process_alpha_proton_scattering = false,
// --- Geometry, Calibration, & Model Variables --- // --- Geometry, Calibration, & Model Variables ---
double source_vertex = 53.0, double source_vertex = 53.0,
z_entrance = -174.3 - 9.7 - 100.0, // z_entrance = -174.3 - 9.7 - 100.0,
z_entrance = -174.3 - 9.7 - 270.0, // new measurement of the chamber length puts the chamber at
// 1175mm instead of 1105, plus some part of the window actually lies outside the chamber
dither_sigma = 8.0, dither_sigma = 8.0,
dither_sigma_c0 = 16.0, dither_sigma_c0 = 16.0,
cathode_gain = 1.0, cathode_gain = 1.0,
@ -511,6 +513,8 @@ void TrackRecon::Begin(TTree * /*tree*/)
else else
plotter = new HistPlotter("Analyzer_SX3.root", "TFILE"); plotter = new HistPlotter("Analyzer_SX3.root", "TFILE");
plotter->set_barrier_limit(getenv("FLUSH_BARRIER") ? std::atoll(getenv("FLUSH_BARRIER")) : 50000);
if (getenv("reactiondata")) if (getenv("reactiondata"))
{ {
reactiondata = std::atoi(getenv("reactiondata")); reactiondata = std::atoi(getenv("reactiondata"));
@ -1008,8 +1012,46 @@ inline void pcEnergyCalibrationAccumulateProton(const std::vector<Event> &PC_Eve
} }
} }
// Reads VmRSS (resident memory, MB) for this process from /proc/self/status.
// Returns -1.0 if unavailable (e.g. non-Linux) so callers can skip the check.
inline double currentRSS_MB()
{
std::ifstream statusFile("/proc/self/status");
std::string line;
while (std::getline(statusFile, line))
{
if (line.compare(0, 6, "VmRSS:") == 0)
{
std::istringstream iss(line.substr(6));
double kb = -1.0;
iss >> kb;
return kb > 0.0 ? kb / 1024.0 : -1.0;
}
}
return -1.0;
}
Bool_t TrackRecon::Process(Long64_t entry) Bool_t TrackRecon::Process(Long64_t entry)
{ {
static const double maxRSS_MB = getenv("MAX_RSS_MB") ? std::atof(getenv("MAX_RSS_MB")) : 0.0;
static const Long64_t checkStride = getenv("MEMCHECK_STRIDE") ? std::atoll(getenv("MEMCHECK_STRIDE")) : 5000;
static Long64_t processedCount = 0;
++processedCount;
if (maxRSS_MB > 0.0 && (processedCount % checkStride == 0))
{
double rss = currentRSS_MB();
if (rss > 0.0 && rss > maxRSS_MB)
{
std::cout << "MAX_RSS_MB (" << maxRSS_MB << ") exceeded (RSS=" << rss
<< " MB) at entry " << entry << " -- forcing a cache flush and continuing." << std::endl;
plotter->force_flush_caches();
}
}
plotter->barrier_increment();
hitPos.Clear(); hitPos.Clear();
qqqenergy = -1; qqqenergy = -1;
qqqtimestamp = -1; qqqtimestamp = -1;

2
run_27Al.sh
View File

@ -38,7 +38,7 @@ process_run() {
export -f process_run export -f process_run
echo "Starting parallel processing..." echo "Starting parallel processing..."
time parallel --bar -j 6 process_run ::: {24..41} # time parallel --bar -j 6 process_run ::: {24..41}
time parallel --bar -j 8 process_run ::: {50..59} time parallel --bar -j 8 process_run ::: {50..59}
# time parallel --bar -j 4 process_run ::: 62 63 66 67 68 # time parallel --bar -j 4 process_run ::: 62 63 66 67 68
# time parallel --bar -j 1 process_run ::: 73 # time parallel --bar -j 1 process_run ::: 73