NCUtility.cxx

Go to the documentation of this file.

#include "NCUtility.h"

#include "TH1.h"
#include "TH2.h"
#include "TH3.h"

#include "TStopwatch.h"

#include "MessageService/MsgService.h"

#include <cassert>

#include "sys/stat.h"

CVSID("$Id: NCUtility.cxx,v 1.7 2010/01/20 15:41:32 rodriges Exp $");

TH1* NC::Utility::CloneFast(const TH1* h)
{
  const TString className = h->ClassName();

  if(className == "TH1F") return new TH1F(*((TH1F*)h));
  if(className == "TH1D") return new TH1D(*((TH1D*)h));
  if(className == "TH2F") return new TH2F(*((TH2F*)h));
  if(className == "TH2D") return new TH2D(*((TH2D*)h));
  if(className == "TH3F") return new TH3F(*((TH3F*)h));
  if(className == "TH3D") return new TH3D(*((TH3D*)h));

  MSG("NCUtility", Msg::kWarning) << "Histogram type "
                                  << className
                                  << " is not handled."
                                  << " Falling back to Clone()"
                                  << " which is slow." << endl;
  return (TH1*)h->Clone();
}


inline void GetArray(const TH1* h, float*& fArray, double*& dArray)
{
  fArray = 0; dArray = 0;

  const TString className = h->ClassName();

  if(className == "TH1F") fArray = ((TH1F*)h)->fArray;
  if(className == "TH2F") fArray = ((TH2F*)h)->fArray;
  if(className == "TH3F") fArray = ((TH3F*)h)->fArray;
  if(className == "TH1D") dArray = ((TH1D*)h)->fArray;
  if(className == "TH2D") dArray = ((TH2D*)h)->fArray;
  if(className == "TH3D") dArray = ((TH3D*)h)->fArray;
}


int NC::Utility::GetArraySize(const TH1* h)
{
  // Horrible hack to get at the number of entries in the underlying
  // representation of a TH1 subclass.
  struct CellsGetter: public TH1
  {
    int NCells()
    {
      return fNcells;
    }
  };

  return ((CellsGetter*)h)->NCells();
}


void NC::Utility::AddFast(TH1* a, const TH1* b, double c)
{
  if(c == 0) return;

  assert(a->ClassName() == b->ClassName());

  const int N = GetArraySize(a);
  assert(GetArraySize(b) == N);

  float* fArrayA = 0;
  float* fArrayB = 0;
  double* dArrayA = 0;
  double* dArrayB = 0;

  GetArray(a, fArrayA, dArrayA);
  GetArray(b, fArrayB, dArrayB);

  if(!(fArrayA && fArrayB) && !(dArrayA && dArrayB)){
    MSG("NCUtility", Msg::kWarning) << "Histogram type "
                                    << a->ClassName()
                                    << " is not handled."
                                    << " Falling back to Add()"
                                    << " which is slow." << endl;
    a->Add(b, c);
    return;
  }

  if(fArrayA && fArrayB) for(int n = 0; n < N; ++n) fArrayA[n] += fArrayB[n]*c;
  if(dArrayA && dArrayB) for(int n = 0; n < N; ++n) dArrayA[n] += dArrayB[n]*c;
}


void NC::Utility::MultiplyFast(TH1* a, const TH1* b, double c)
{
  assert(a->ClassName() == b->ClassName());

  const int N = GetArraySize(a);
  assert(GetArraySize(b) == N);

  float* fArrayA = 0;
  float* fArrayB = 0;
  double* dArrayA = 0;
  double* dArrayB = 0;

  GetArray(a, fArrayA, dArrayA);
  GetArray(b, fArrayB, dArrayB);

  if(!(fArrayA && fArrayB) && !(dArrayA && dArrayB)){
    MSG("NCUtility", Msg::kWarning) << "Histogram type "
                                    << a->ClassName()
                                    << " is not handled."
                                    << " Falling back to Multiply()"
                                    << " which is slow." << endl;
    a->Multiply(b);
    a->Scale(c);
    return;
  }

  if(fArrayA && fArrayB) for(int n = 0; n < N; ++n) fArrayA[n] *= fArrayB[n]*c;
  if(dArrayA && dArrayB) for(int n = 0; n < N; ++n) dArrayA[n] *= dArrayB[n]*c;
}


void NC::Utility::DivideFast(TH1* a, const TH1* b, double c)
{
  assert(a->ClassName() == b->ClassName());

  const int N = GetArraySize(a);
  assert(GetArraySize(b) == N);

  float* fArrayA = 0;
  float* fArrayB = 0;
  double* dArrayA = 0;
  double* dArrayB = 0;

  GetArray(a, fArrayA, dArrayA);
  GetArray(b, fArrayB, dArrayB);

  if(!(fArrayA && fArrayB) && !(dArrayA && dArrayB)){
    MSG("NCUtility", Msg::kWarning) << "Histogram type "
                                    << a->ClassName()
                                    << " is not handled."
                                    << " Falling back to Divide()"
                                    << " which is slow." << endl;
    a->Divide(b);
    a->Scale(1/c);
    return;
  }

  if(fArrayA && fArrayB){
    for(int n = 0; n < N; ++n){
      const double denom = fArrayB[n]*c;
      if(denom) fArrayA[n] /= denom; else fArrayA[n] = 0;
    }
  }

  if(dArrayA && dArrayB){
    for(int n = 0; n < N; ++n){
      const double denom = dArrayB[n]*c;
      if(denom) dArrayA[n] /= denom; else dArrayA[n] = 0;
    }
  }
}


vector<TString> NC::Utility::ParseStringList(const char* sl)
{
  TString strList = sl;

  vector<TString> ret;

  strList.Remove(TString::kBoth, ' '); // strip whitespace

  if(strList == "") return ret;

  int space = strList.First(' ');

  // If there aren't any spaces then we want 'head' to be the whole string
  if(space < 0) space = strList.Length();

  // Slicing operator takes start and length
  const TString head = strList(0, space);
  const TString tail = strList(space+1, strList.Length()-space);

  ret.push_back(head);

  // Recurse to parse the remainder of the string
  const vector<TString> rest = ParseStringList(tail);
  ret.insert(ret.end(), rest.begin(), rest.end());

  return ret;
}


vector<int> NC::Utility::ParseNumberList(const char* nl)
{
  const vector<TString> sl = ParseStringList(nl);

  vector<int> ret;

  for(unsigned int n = 0; n < sl.size(); ++n){
    assert(sl[n].IsDigit());
    ret.push_back(sl[n].Atoi());
  }

  return ret;
}


void NC::Utility::ReportProgress(double est_frac, TStopwatch& sw)
{
  // Check if we're outputing to a file. If so don't bother showing off
  // with the progress bar, it won't work.
  struct stat buf;
  fstat(fileno(stdout), &buf);
  const bool isFile = (buf.st_mode & S_IFREG) || (buf.st_mode & S_IFIFO);

  const double elapsed = sw.RealTime();
  sw.Start(kFALSE);

  if(est_frac < 0) est_frac = 0;
  if(est_frac > 1) est_frac = 1;

  const int remain = est_frac ? int(elapsed*(1/est_frac-1)) : 0;

  const TString txt = TString::
    Format("%d%% complete, approx %dh%dm%ds remain",
           int(100*est_frac), remain/3600, (remain/60)%60, remain%60);

  if(isFile){
    MSG("NCUtility", Msg::kInfo) << txt << std::endl;
  }
  else{
    TString bar;
    const int extra = 50-txt.Length();
    const int before = extra/2;
    const int after = extra-before;

    string plain_bar;
    for(int n = 0; n < before; ++n) plain_bar += " ";
    plain_bar += txt;
    for(int n = 0; n < after; ++n) plain_bar += " ";

    const char kColor_Reset[]   = { 0x1B, '[', '0', 'm', 0 };
    const char kColor_BgRed[]   = { 0x1B, '[', '4', '1', 'm', 0 };
    const char kColor_BgGreen[] = { 0x1B, '[', '4', '2', 'm', 0 };

    int fracInt = int(est_frac*50);
    if(fracInt < 0) fracInt = 0;
    if(fracInt > 50) fracInt = 50;
    bar += kColor_BgGreen;
    for(int n = 0; n < 50; ++n){
      if(n >= fracInt) bar += kColor_BgRed;
      bar += plain_bar[n];
    }
    bar += kColor_Reset;

    MSG("NCUtility", Msg::kInfo) << bar <<  "\r" << std::flush;
  }
}

---