PhysicsNtuple/Hist/HistMan.cxx

Go to the documentation of this file.

// $Id: HistMan.cxx,v 1.12 2009/08/19 21:17:22 jyuko Exp $

// C++
#include <cmath>
#include <cstdlib>
#include <iostream>
#include <iomanip>

// MINOS
#include "Util/UtilString.h"

// ROOT
#include "TH1.h"
#include "TH2.h"
#include "TProfile.h"

// Local
#include "PhysicsNtuple/Default.h"

// Local
#include "HistMan.h"
#include "tinyxml.h"

//
// Initialize static member variables
//
Anp::HistMan *Anp::HistMan::fgInstance = 0;
Anp::Mutex Anp::HistMan::fMutex;

//
// Disable automatic addition of histograms to current TDirectory
//
//TH1::AddDirectory(kFALSE);

//------------------------------------------------------------------------------------------
// return (create for a first call) factory instance
//------------------------------------------------------------------------------------------
Anp::HistMan& Anp::HistMan::Instance()
{
   if(!fgInstance)
   {
      Anp::Lock<Mutex> lock(fMutex);
      
      if(!fgInstance)
      {
         fgInstance = new Anp::HistMan();
      }
   }
   
   return *fgInstance;
}

using namespace std;

//-----------------------------------------------------------------------------------------
Anp::HistMan::HistMan()
   :fNMiss(0),
    fKeyWidth(4)
{
}

//-----------------------------------------------------------------------------------------
Anp::HistMan::~HistMan()
{
}

//-----------------------------------------------------------------------------------------
bool Anp::HistMan::ReadFile(const string file)
{
   //
   // Parse xml file with histogram definitions. These are top level keys
   // which must be located within "histograms" object:
   //    keybase - creates multiple keys for one histogram
   //    dirbase - stores keys in "sub-directory"
   //    h1d - histogram definition
   //
   Anp::Lock<Mutex> lock(fMutex);

   TiXmlDocument doc(file.c_str());
   const bool loadok = doc.LoadFile();
   if(!loadok)
   {
      cerr << "HistMan::ReadFile - failed to load: " << file << endl;
      return false;
   }
   
   for(TiXmlNode *pChild = doc.FirstChild(); pChild != 0; pChild = pChild -> NextSibling()) 
   {
      if(!pChild)
      {
         cerr << "HistMan::ReadFile - null child..." << endl;
         continue;
      }

      if(pChild -> Type() != TiXmlNode::ELEMENT)
      {
         continue;
      }

      if(pChild -> Value() && strcmp(pChild -> Value(), "histograms"))
      {
         continue;
      }
      
      TiXmlElement* pHist = pChild -> ToElement();
      if(pHist)
      {
         ReadHistogramBlock(*pHist);
      }
      else
      {
         cerr << "HistMan::ReadFile - failed to cast root element for histograms" << endl;
         continue;
      }
   }

   return true;
}

//-----------------------------------------------------------------------------------------
bool Anp::HistMan::ReadHistogramBlock(TiXmlElement &block)
{
   InfoMap hmap;
   vector<vector<int> > basevec;
   vector<string> dirvec;

   for(TiXmlNode *pChild = block.FirstChild(); pChild != 0; pChild = pChild -> NextSibling()) 
   {
      if(!pChild)
      {
         cerr << "HistMan::ReadFile - null child..." << endl;
         continue;
      }

      if(pChild -> Type() != TiXmlNode::ELEMENT)
      {
         continue;
      }
      
      TiXmlElement* hist = pChild -> ToElement();
      if(!hist)
      {
         cerr << "HistMan::ReadFile - failed to cast xml node to element..." << endl;
         continue;
      }

      if(!hist -> Value())
      {
         cerr << "HistMan::ReadFile - node has no value" << endl;
         continue;
      }

      if(strcmp(hist -> Value(), "keybase") == 0 && hist -> Attribute("keys"))
      {
         const vector<int> ivec = GetIntVec(hist -> Attribute("keys"));
         if(!ivec.empty())
         {
            basevec.push_back(ivec);
         }
         
         continue;
      }

      if(strcmp(hist -> Value(), "dirbase") == 0 && hist -> Attribute("dirs"))
      {
         const vector<string> svec = GetStringVec(hist -> Attribute("dirs"));
         if(!svec.empty())
         {
            dirvec.insert(dirvec.end(), svec.begin(), svec.end());
         }

         continue;
      }

      HistInfo info;      
      if(strcmp(hist -> Value(), "h1d") != 0 || !ReadHistogram(*hist, info))
      {
         continue;
      }

      if(!hmap.insert(InfoMap::value_type(info.GetKey(), info)).second)
      {
         cerr << "HistMan::ReadFile - histogram key " << info.GetKey() << " already exists" << endl;
      }

      const int intkey = std::atoi(info.GetKey().c_str());
      if(intkey != 0)
      {
         const unsigned int width = ComputeWidth(intkey);
         fKeyWidth = max<int>(width, fKeyWidth);

         if(width != info.GetKey().size())
         {
           cerr << "HistMan::ReadFile - integer width and string size are not equal for key " <<info.GetKey()<<" "<<info.GetTitle()<< endl;
         }
      }
   }
   
   vector<int> allvec;
   for(vector<vector<int> >::const_iterator bit = basevec.begin(); bit != basevec.end(); ++bit)
   {
      const vector<int> prevvec = allvec;
      const vector<int> &currvec = *bit;
      
      for(vector<int>::const_iterator cit = currvec.begin(); cit != currvec.end(); ++cit)
      {
         if(prevvec.empty())
         {
            allvec.push_back(*cit);
         }
         else
         {
            for(vector<int>::const_iterator pit = prevvec.begin(); pit != prevvec.end(); ++pit)
            {
               allvec.push_back(*cit + *pit);
            }
         } 
      }
   }

   InfoMap nmap = hmap;

   for(vector<int>::const_iterator it = allvec.begin(); it != allvec.end(); ++it)
   {
      const int base = *it;

      for(InfoMap::const_iterator hit = hmap.begin(); hit != hmap.end(); ++hit)
      {
         const int oldkey = std::atoi((hit -> first).c_str());
         if(oldkey == 0)
         {
            continue;
         }

         const int newkey = base + oldkey;
         HistInfo info = hit -> second;

         stringstream strkey;
         strkey << newkey;
         info.fKey = strkey.str();

         if(!nmap.insert(InfoMap::value_type(info.GetKey(), info)).second)
         {
            cerr << "HistMan::ReadFile - histogram key " << info.GetKey() << " already exists" << endl;
         }

         const unsigned int width = ComputeWidth(newkey);
         fKeyWidth = max<int>(width, fKeyWidth);

         if(width != info.GetKey().size())
         {
            cerr << "HistMan::ReadFile - integer width and string size are not equal" << endl;
         }
      }
   }

   if(dirvec.empty())
   {
      dirvec.push_back(".");
   }

   for(vector<string>::const_iterator dit = dirvec.begin(); dit != dirvec.end(); ++dit)
   {
      InfoMap &imap = fDirMap[*dit];

      for(InfoMap::const_iterator nit = nmap.begin(); nit != nmap.end(); ++nit)
      {
         if(!imap.insert(InfoMap::value_type(*nit)).second)
         {
            cerr << "HistMan::ReadFile - failed to add " << nit -> first << " to " << *dit << endl;
         }
      }
   }

   return true;
}

//-----------------------------------------------------------------------------------------
bool Anp::HistMan::ReadHistogram(const TiXmlElement &hist, HistInfo &info)
{
   if(hist.Value() && strcmp(hist.Value(), "h1d") != 0)
   {
      return false;
   }

   AxisInfo xaxis, yaxis;
 
   bool axis1dok = false, axis2dok = false;
   for(const TiXmlNode *pChild = hist.FirstChild(); pChild != 0; pChild = pChild -> NextSibling()) 
   {
      if(!pChild)
      {
         cerr << "Null child..." << endl;
         continue;
      }

      if(pChild -> Type() != TiXmlNode::ELEMENT)
      {
         continue;
      }
      
      const TiXmlElement* axis = pChild -> ToElement();
      if(!axis)
      {
         cerr << "Failed to cast xml node to element..." << endl;
         continue;
      }
      
      if(!axis -> Value())
      {
         cerr << "Node has no value" << endl;
         continue;
      }

      if(strcmp(axis -> Value(), "axis1d") == 0 && ReadAxis(*axis, xaxis))
      {
         axis1dok = true;
      }
      if(strcmp(axis -> Value(), "axis2d") == 0 && ReadAxis(*axis, yaxis))
      {
         axis2dok = true;
      }
   }

   if(!axis1dok)
   {
      cerr << "Failed to read 1d axis info" << endl;
      return false;
   }
   else
   {
      info.fXaxis = xaxis;
   }

   if(axis2dok)
   {
      info.fYaxis = yaxis;
   }

   const char *key = hist.Attribute("key");
   if(key)
   {
      info.fKey = string(key);
   }
   else
   {
      cerr << "HistMan::ReadHistogram() - key string is empty" << endl;
      return false;
   }   

   const char *name = hist.Attribute("name");
   if(name)
   {
      info.fName = name;
   }

   const char *title = hist.Attribute("title");
   if(title)
   {
      info.fTitle = title;
   }

   const char *axis = hist.Attribute("axis");
   if(axis)
   {
      info.fXaxis.fTitle = axis;
   }

   return true;
}

//-----------------------------------------------------------------------------------------
bool Anp::HistMan::ReadAxis(const TiXmlElement &axis, AxisInfo &info)
{
   bool read_bin = false, read_bins = false;
 
   int bin = 0;
   float min = 0.0, max = 0.0;

   if(axis.QueryIntAttribute("bin", &bin)   == TIXML_SUCCESS &&
      axis.QueryFloatAttribute("min", &min) == TIXML_SUCCESS &&
      axis.QueryFloatAttribute("max", &max) == TIXML_SUCCESS)
   {
      read_bin = true;
   }
   
   const char *xtitle = axis.Attribute("title");
   if(xtitle)
   {
      info.fTitle = xtitle;
   }
   
   vector<double> bvec;
   
   const char *bins = axis.Attribute("bins");
   if(bins)
   {
      vector<string> datavec;
      UtilString::StringTok(datavec, bins, ", ");
      
      for(vector<string>::const_iterator sit = datavec.begin(); sit != datavec.end(); ++sit)
      {
         if(sit -> size() > 0)
         {
           bvec.push_back(std::atof(sit -> c_str()));
         }
      }

      std::sort(bvec.begin(), bvec.end());

      //for(unsigned int i = 0; i < bvec.size(); ++i)
      //{
      //cout << "bin = " << i << " edge = " << bvec[i] << endl;
      //}

      if(!bvec.empty())
      {
         read_bins = true;
      }
   }
   
   if(read_bin)
   {
      info.fNbins = bin;
      info.fMin = min;
      info.fMax = max;
      info.fValid = true;
   }
   else if(read_bins)
   {
      info.fBins = bvec;
   }
  
   return true;
}

//-----------------------------------------------------------------------------------------      
bool Anp::HistMan::KeyExists(const int key, const string &dir) const
{
   return HistMan::KeyExists(HistMan::Convert(key), dir);
}

//-----------------------------------------------------------------------------------------      
bool Anp::HistMan::KeyExists(const string &key, const string &dir) const
{
   Anp::Lock<Mutex> lock(fMutex);
   
   DirIter dit = fDirMap.find(dir);
   if(dit == fDirMap.end())
   {
      return false;
   }

   const InfoMap &imap = dit -> second;

   return (imap.find(key) != imap.end());
}

//-----------------------------------------------------------------------------------------      
const Anp::Hist1d<double> Anp::HistMan::CreateHist1d(const int key, const string &dir)
{
   return HistMan::CreateHist1d(HistMan::Convert(key), dir);
}

//-----------------------------------------------------------------------------------------      
const Anp::Hist1d<double> Anp::HistMan::CreateHist1d(const string &key, const string &dir)
{
   Anp::Lock<Mutex> lock(fMutex);

   DirIter dit = fDirMap.find(dir);
   if(dit == fDirMap.end())
   {
      ++fNMiss;
      return Hist1d<double>(0, 0.0, 0.0);
   }

   const InfoMap &imap = dit -> second;

   InfoIter hit = imap.find(key);
   if(hit == imap.end())
   {
      ++fNMiss;
      return Hist1d<double>(0, 0.0, 0.0);
   }

   return HistMan::CreateHist1d(hit -> second);
}

//-----------------------------------------------------------------------------------------      
TH1* Anp::HistMan::CreateTH1(const int key, const string &dir)
{
   return HistMan::CreateTH1(HistMan::Convert(key), dir);
}

//-----------------------------------------------------------------------------------------      
TH1* Anp::HistMan::CreateTH1(const string &key, const string &dir_)
{
   Anp::Lock<Mutex> lock(fMutex);
   
   //
   // Extract option
   //
   string opt;
   string dir = dir_;
   if     (dir.find("profile") != string::npos) opt += "profile";
   else if(dir.find("float")   != string::npos) opt += "float";
   
   if(!opt.empty() && dir.find(opt.c_str()) != string::npos)
   {
      dir.erase(dir.find(opt.c_str()), opt.size());
   }

   if(dir.empty())
   {
      dir = ".";
   }

   DirIter dit = fDirMap.find(dir);
   if(dit == fDirMap.end())
   {
      ++fNMiss;
      return 0;
   }

   const InfoMap &imap = dit -> second;

   InfoIter hit = imap.find(key);
   if(hit == imap.end())
   {
      ++fNMiss;
      return 0;
   }

   return HistMan::CreateTH1(hit -> second, opt);
}

//-----------------------------------------------------------------------------------------      
TH2* Anp::HistMan::CreateTH2(const int key, const string &dir)
{
   return HistMan::CreateTH2(Convert(key), dir);
}

//-----------------------------------------------------------------------------------------      
TH2* Anp::HistMan::CreateTH2(const string &key, const string &dir_)
{
   Anp::Lock<Mutex> lock(fMutex);

   //
   // Extract option
   //
   string opt;
   string dir = dir_;
   if     (dir.find("double") != string::npos) opt += "double";
   
   if(!opt.empty() && dir.find(opt.c_str()) != string::npos)
   {
      dir.erase(dir.find(opt.c_str()), opt.size());
   }

   if(dir.empty())
   {
      dir = ".";
   }


   DirIter dit = fDirMap.find(dir);
   if(dit == fDirMap.end())
   {
      ++fNMiss;
      return 0;
   }

   const InfoMap &imap = dit -> second;

   InfoIter hit = imap.find(key);
   if(hit == imap.end())
   {
      ++fNMiss;
      return 0;
   }

   return HistMan::CreateTH2(hit -> second, opt);
}

//-----------------------------------------------------------------------------------------      
TH2* Anp::HistMan::GetTH2(const int xkey, const int ykey, const string &dir)
{
   return HistMan::GetTH2(Convert(xkey), HistMan::Convert(ykey), dir);
}

//-----------------------------------------------------------------------------------------      
TH2* Anp::HistMan::GetTH2(const string &xkey, const string &ykey, const string &dir)
{
   Anp::Lock<Mutex> lock(fMutex);

   DirIter dit = fDirMap.find(dir);
   if(dit == fDirMap.end())
   {
      ++fNMiss;
      return 0;
   }

   const InfoMap &imap = dit -> second;

   InfoIter xhit = imap.find(xkey);
   if(xhit == imap.end())
   {
      ++fNMiss;
      return 0;
   }

   InfoIter yhit = imap.find(ykey);
   if(yhit == imap.end())
   {
      ++fNMiss;
      return 0;
   }

   const HistInfo &xinfo = xhit -> second;
   const HistInfo &yinfo = yhit -> second;

   HistInfo info;

   info.fKey = "keys_" + xinfo.GetKey() + "_and_" + yinfo.GetKey();
   info.fXaxis = xinfo.GetXaxis();
   info.fYaxis = yinfo.GetXaxis();
   info.fName  = info.GetKey();
   info.fTitle = info.GetKey();

   if(xinfo.GetYaxis().Valid())
   {
      info.fXaxis = xinfo.GetYaxis();
   }
   if(yinfo.GetYaxis().Valid())
   {
      info.fYaxis = yinfo.GetYaxis();
   }

   return HistMan::CreateTH2(info);
}

//-----------------------------------------------------------------------------------------      
const Anp::Hist1d<double> Anp::HistMan::CreateHist1d(const HistInfo &info)
{
   const AxisInfo &xaxis = info.GetXaxis();

   if(!xaxis.GetBins().empty())
   {
      return Hist1d<double>(xaxis.GetBins());
   }

   return Hist1d<double>(xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax());      
}

//-----------------------------------------------------------------------------------------      
TH1* Anp::HistMan::CreateTH1(const HistInfo &info, const string &opt)
{
   string strkey = info.GetKey();
   string name   = info.GetName();
   string title  = info.GetTitle();

   const int intkey = std::atoi(strkey.c_str());   
   if(intkey != 0)
   {
      strkey = "key_" + Convert(intkey, fKeyWidth);
   }
   if(name.empty())
   {
      name = strkey;
   }
   if(title.empty())
   {
      title = info.GetName();
   }

   const AxisInfo &xaxis = info.GetXaxis();
   const AxisInfo &yaxis = info.GetYaxis();

   TH1 *h = 0;
   if(xaxis.GetBins().empty())
   {
      Anp::Lock<Anp::Mutex> lock(Anp::GetROOTMutex());

      if(opt.find("profile") != string::npos)
      {
         h = new TProfile(name.c_str(), title.c_str(), 
                          xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax());  
      }
      else if(opt.find("float") != string::npos)
      {
         h = new TH1F(name.c_str(), title.c_str(), 
                      xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax());      
      }
      else
      {
         h = new TH1D(name.c_str(), title.c_str(),
                      xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax());      
      }
   }
   else
   {
      Anp::Lock<Anp::Mutex> lock(Anp::GetROOTMutex());

      const vector<double> &bvec = xaxis.GetBins();

      double *darray = new double[bvec.size()];
      
      for(unsigned int ibin = 0; ibin < bvec.size(); ++ibin)
      {
         darray[ibin] = bvec[ibin];
      }
      
      if(opt.find("profile") != string::npos)
      {
         h = new TProfile(name.c_str(), title.c_str(), bvec.size() - 1, darray);
      }
      else if(opt.find("float") != string::npos)
      {
         h = new TH1F(name.c_str(), title.c_str(),  bvec.size() - 1, darray);
      }
      else
      {
         h = new TH1D(name.c_str(), title.c_str(), bvec.size() - 1, darray);
      }
      
      delete [] darray;      
   }
   
   //
   // Set NULL TDirectory
   //
   Anp::SetDir(h, 0);

   if(xaxis.GetTitle().size() > 0)
   {
      h -> GetXaxis() -> SetTitle(xaxis.GetTitle().c_str());
      h -> GetXaxis() -> CenterTitle();
   }
   if(yaxis.GetTitle().size() > 0)
   {
      h -> GetYaxis() -> SetTitle(yaxis.GetTitle().c_str());
      h -> GetYaxis() -> CenterTitle();
   }

   return h;
}

//-----------------------------------------------------------------------------------------      
TH2* Anp::HistMan::CreateTH2(const HistInfo &info, const string &opt)
{
   string strkey = info.GetKey();
   string name   = info.GetName();
   string title  = info.GetTitle();

   const int intkey = std::atoi(strkey.c_str());   
   if(intkey != 0)
   {
      strkey = Convert(intkey, fKeyWidth);
   }
   if(name.empty())
   {
      name = strkey;
   }
   if(title.empty())
   {
      title = info.GetName();
   }

   const AxisInfo &xaxis = info.GetXaxis();
   const AxisInfo &yaxis = info.GetYaxis();

   TH2 *h = 0;
   if(!name.empty())
   {
      Anp::Lock<Anp::Mutex> lock(Anp::GetROOTMutex());

      if(opt.find("double") != string::npos)
        {
          h = new TH2D(name.c_str(), title.c_str(),
                   xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax(),
                   yaxis.GetNbins(), yaxis.GetMin(), yaxis.GetMax());   
          
        }
      else
        {
          h = new TH2F(name.c_str(), title.c_str(),
                   xaxis.GetNbins(), xaxis.GetMin(), xaxis.GetMax(),
                   yaxis.GetNbins(), yaxis.GetMin(), yaxis.GetMax());      
        }

      if(xaxis.GetTitle().size() > 0)
      {
         h -> GetXaxis() -> SetTitle(xaxis.GetTitle().c_str());
         h -> GetXaxis() -> CenterTitle();
      }
      if(yaxis.GetTitle().size() > 0)
      {
         h -> GetYaxis() -> SetTitle(yaxis.GetTitle().c_str());
         h -> GetYaxis() -> CenterTitle();
      }
   }

   //
   // Set NULL TDirectory
   //
   Anp::SetDir(h, 0);

   return h;
}

//-----------------------------------------------------------------------------------------      
const vector<int> Anp::HistMan::GetIntVec(const string &list) const
{
   vector<string> namelist;
   UtilString::StringTok(namelist, list, ", ");
   
   vector<int> ivec;
   for(vector<string>::const_iterator sit = namelist.begin(); sit != namelist.end(); ++sit)
   {
      const int newbase = std::atoi(sit -> c_str());
      if(newbase != 0)
      {
         ivec.push_back(newbase);
      }
   }
   
   return ivec;
}

//-----------------------------------------------------------------------------------------      
const vector<string> Anp::HistMan::GetStringVec(const string &list) const
{
   vector<string> namelist;
   UtilString::StringTok(namelist, list, ", ");

   vector<string> nvec;
   for(vector<string>::const_iterator sit = namelist.begin(); sit != namelist.end(); ++sit)
   {
      if(!sit -> empty())
      {
         nvec.push_back(*sit);
      }
   }
   
   return nvec;
}

//-----------------------------------------------------------------------------------------      
const string Anp::HistMan::Convert(const int key, const int width) const
{
   stringstream strkey;
   if(width == 0)
   {
      strkey << setw(ComputeWidth(key)) << key;
   }
   else
   {
      strkey << setfill('0') << setw(width) << key;
   }

   return strkey.str();
}

//-----------------------------------------------------------------------------------------      
unsigned short Anp::HistMan::ComputeWidth(int key) const
{
   //
   // Compute key width
   //
   return static_cast<unsigned int>(std::floor(std::log10(double(std::abs(key) + 1)))) + 1;
}

//-----------------------------------------------------------------------------------------      
unsigned int Anp::HistMan::NMiss() const
{
   Anp::Lock<Mutex> lock(fMutex);
   return fNMiss;
}

//-----------------------------------------------------------------------------------------      
void Anp::HistMan::Print(ostream& o) const
{
   o << "HistMan::Print" << endl;
   for(DirIter dit = fDirMap.begin(); dit != fDirMap.end(); ++dit)
   {
      o << "Directory = " << dit -> first << endl; 
      const InfoMap &imap = dit -> second;

      for(InfoIter it = imap.begin(); it != imap.end(); ++it)
      {
         it -> second.Print(o);
      }
   }
}

//-----------------------------------------------------------------------------------------      
Anp::HistInfo::HistInfo()
   :fKey(),
    fName(),
    fTitle(),
    fXaxis(),
    fYaxis()
{
}

//-----------------------------------------------------------------------------------------      
Anp::HistInfo::~HistInfo()
{
}

//-----------------------------------------------------------------------------------------      
const std::string& Anp::HistInfo::GetName() const 
{
   return fName;
}

//-----------------------------------------------------------------------------------------      
const std::string& Anp::HistInfo::GetTitle() const
{
   return fTitle;
}

//-----------------------------------------------------------------------------------------      
const Anp::AxisInfo& Anp::HistInfo::GetXaxis() const
{
   return fXaxis;
}

//-----------------------------------------------------------------------------------------      
const Anp::AxisInfo& Anp::HistInfo::GetYaxis() const
{
   return fYaxis;
}

//-----------------------------------------------------------------------------------------      
const std::string& Anp::HistInfo::GetKey() const
{
   return fKey;
}

//-----------------------------------------------------------------------------------------      
bool Anp::HistInfo::operator<(const HistInfo &rhs) const
{
   if(fKey < rhs.fKey)
   {
      return true;
   }

   return false;
}

//-----------------------------------------------------------------------------------------
bool Anp::HistInfo::operator==(const HistInfo &rhs) const
{
   if(fKey == rhs.fKey)
   {
      return true;
   }

   return false;
}

//----------------------------------------------------------------------
void Anp::HistInfo::Print(ostream& o) const
{
   o << "Printing HistInfo: key = " << GetKey() << std::endl
     << "   name = " << GetName()
     << ", title = " << GetTitle()
     << "  " << fXaxis
     << "  " << fYaxis;
}

//----------------------------------------------------------------------
std::ostream& Anp::operator<<(ostream& o, const Anp::HistInfo &self)
{
   self.Print(o);
   return o;
}

//-----------------------------------------------------------------------------------------      
Anp::AxisInfo::AxisInfo()
   :fNbins(0),
    fMin(0),
    fMax(0),
    fValid(false),
    fTitle()
{
}

//-----------------------------------------------------------------------------------------      
Anp::AxisInfo::~AxisInfo()
{
}

//-----------------------------------------------------------------------------------------      
const std::string& Anp::AxisInfo::GetTitle() const
{
   return fTitle;
}

//-----------------------------------------------------------------------------------------      
double Anp::AxisInfo::GetMin() const
{
   return fMin;
}

//-----------------------------------------------------------------------------------------      
double Anp::AxisInfo::GetMax() const
{
   return fMax;
}

//-----------------------------------------------------------------------------------------      
int Anp::AxisInfo::GetNbins() const
{
   return fNbins;
}

//-----------------------------------------------------------------------------------------      
const std::vector<double>& Anp::AxisInfo::GetBins() const
{
   return fBins;
}

//----------------------------------------------------------------------
bool Anp::AxisInfo::Valid() const
{
   return fValid;
}

//----------------------------------------------------------------------
void Anp::AxisInfo::Print(std::ostream& o) const
{
   o << "AxisInfo (" << fTitle << "): (nbin, min, max) = "
     << "(" << GetNbins() << ", " << GetMin() << ", " << GetMax() << ")" << std::endl;
}

//----------------------------------------------------------------------
ostream& Anp::operator<<(std::ostream& o, const Anp::AxisInfo &self)
{
   self.Print(o);
   return o;
}

//----------------------------------------------------------------------
TH1* Anp::MakeTH1(const string &dir, TDirectory *dir_, const string &key)
{
   TH1 *h = Anp::HistMan::Instance().CreateTH1(key, dir);
   if(h)
   {
      Anp::SetDir(h, dir_);
   }
   
   return h;
}

//----------------------------------------------------------------------
TH2* Anp::MakeTH2(const string &dir, TDirectory *dir_, const string &key)
{
   TH2 *h = Anp::HistMan::Instance().CreateTH2(key, dir);
   if(h)
   {
      Anp::SetDir(h, dir_);
   }
   
   return h;
}

---