JPlot2D.cc

Go to the documentation of this file.

#include <string>
#include <iostream>
#include <vector>
#include <set>
#include <cmath>

#include "TROOT.h"
#include "TFile.h"
#include "TClass.h"
#include "TApplication.h"
#include "TCanvas.h"
#include "TKey.h"
#include "TStyle.h"
#include "TAttMarker.h"
#include "TAttLine.h"
#include "TH2.h"
#include "TH3.h"
#include "TH2D.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TGraph2D.h"
#include "TGraph2DErrors.h"
#include "TEllipse.h"
#include "TLine.h"
#include "TF2.h"
#include "TString.h"
#include "TRegexp.h"
#include "TText.h"

#include "JTools/JRange.hh"
#include "JLang/JEquals.hh"
#include "JROOT/JStyle.hh"
#include "JROOT/JCanvas.hh"
#include "JROOT/JMarkerAttributes.hh"
#include "JGizmo/JRootObjectID.hh"
#include "JGizmo/JRootObject.hh"
#include "JGizmo/JGizmoToolkit.hh"

#include "Jeep/JProperties.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

namespace {

  /**
   * Set x of object to logarithmic.
   *
   * \param  object         pointer to object
   * \return                true if set; else false
   */  
  template<class T>
  inline bool setLogX(TObject* object)
  {
    using namespace JPP;

    T* p = dynamic_cast<T*>(object);

    if (p != NULL) {

      setLogarithmicX(p);

      return true;
    }

    return false;
  }
  
  /**
   * Set y of object to logarithmic.
   *
   * \param  object         pointer to object
   * \return                true if set; else false
   */  
  template<class T>
  inline bool setLogY(TObject* object)
  {
    using namespace JPP;

    T* p = dynamic_cast<T*>(object);

    if (p != NULL) {

      setLogarithmicY(p);

      return true;
    }

    return false;
  }

  const std::string MASTER   = "__H__"; //!< Name of prototype

  const char* const JName_t  = "?";     //!< Draw histogram name  as title of plot
  const char* const JTitle_t = "%";     //!< Draw histogram title as title of plot


  /**
   * Auxiliary data structure for master from TH2 or TGraph.
   */
  struct JMaster :
    public JLANG::JEquals<JMaster, TObject*>
  {
    JMaster() :
      h2(NULL),
      g2(NULL)
    {}

    bool equals(const JMaster& master) const
    {
      return (this->h2 == master.h2 &&
              this->g2 == master.g2);
    }

    bool equals(const TObject* p) const
    {
      return (this->h2 == dynamic_cast<const TH2*>(p)      &&
              this->g2 == dynamic_cast<const TGraph2D*>(p));
    }

    JMaster* operator->()
    {
      return this;
    }

    void operator=(TObject* p)
    {
      h2 = dynamic_cast<TH2*>     (p);
      g2 = dynamic_cast<TGraph2D*>(p);
    }

    const char* GetName() const
    {
      return (h2 != NULL ? h2->GetName() :
              g2 != NULL ? g2->GetName() :
              NULL);
    }

    void SetTitle(const char* title)
    {
      if (h2 != NULL) { h2->SetTitle(title); }
      if (g2 != NULL) { g2->SetTitle(title); }
    }

    void SetStats(const Bool_t stats)
    {
      if (h2 != NULL) { h2->SetStats(stats); }
    }

    void SetMaximum(Double_t maximum)
    {
      if (h2 != NULL) { h2->SetMaximum(maximum); }
      if (g2 != NULL) { g2->SetMaximum(maximum); }
    }

    void SetMinimum(Double_t minimum)
    {
      if (h2 != NULL) { h2->SetMinimum(minimum); }
      if (g2 != NULL) { g2->SetMinimum(minimum); }
    }

    void SetNdivisions(Int_t n, Option_t *axis)
    {
      if (h2 != NULL) { h2->SetNdivisions(n, axis); }
    }

    TAxis* GetXaxis()
    {
      return (h2 != NULL ? h2->GetXaxis() :
              g2 != NULL ? g2->GetXaxis() :
              NULL);
    }

    TAxis* GetYaxis()
    {
      return (h2 != NULL ? h2->GetYaxis() :
              g2 != NULL ? g2->GetYaxis() :
              NULL);
    }

    TAxis* GetZaxis()
    {
      return (h2 != NULL ? h2->GetZaxis() :
              g2 != NULL ? g2->GetZaxis() :
              NULL);
    }

    void Draw(Option_t *option="")
    {
      if (h2 != NULL) { h2->Draw(option); }
      if (g2 != NULL) { g2->Draw(option); }
    }

    void setLogarithmicX()
    {
      if (h2 != NULL) { JGIZMO::setLogarithmicX(h2); }
      if (g2 != NULL) { JGIZMO::setLogarithmicX(g2); }
    }

    void setLogarithmicY()
    {
      if (h2 != NULL) { JGIZMO::setLogarithmicY(h2); }
      if (g2 != NULL) { JGIZMO::setLogarithmicY(g2); }
    }

  private:
    TH2*      h2;
    TGraph2D* g2;
  };

  inline void setLogarithmicX(JMaster& master) { master.setLogarithmicX(); }
  inline void setLogarithmicY(JMaster& master) { master.setLogarithmicY(); }
}


/**
 * \file
 * General purpose plot program for 2D ROOT objects.
 * The option <tt>-f</tt> corresponds to <tt><file name>:<object name></tt>.
 * \author mdejong
 */
int main(int argc, char **argv)
{
  using namespace std;
  using namespace JPP;

  typedef JRange<double>  JRange_t;

  vector<JRootObjectID>   inputFile;
  string                  outputFile;
  JCanvas                 canvas;
  JStyle::JParameters     parameters;
  int                     stats;
  JRange_t                X;
  JRange_t                Y;
  JRange_t                Z;
  JCounter                logx;
  JCounter                logy;
  bool                    logz;
  string                  project;
  string                  xLabel;
  string                  yLabel;
  string                  zLabel;
  double                  markerSize;
  string                  option;
  set<char>               grid;
  bool                    batch;
  string                  title;
  map<string, int>        Ndivisions;
  int                     palette;
  int                     debug;
 
  try {

    JProperties properties = parameters.getProperties();

    JParser<> zap("General purpose plot program for 2D ROOT objects.");
    
    zap['f'] = make_field(inputFile,  "<input file>:<object name>");
    zap['o'] = make_field(outputFile, "graphics output")                    = "";
    zap['w'] = make_field(canvas,     "size of canvas <nx>x<ny> [pixels]")  = JCanvas(500, 500);
    zap['@'] = make_field(properties)                                       = JPARSER::initialised();
    zap['s'] = make_field(stats)                                            = -1;
    zap['x'] = make_field(X,          "x-abscissa range")                   = JRange_t::DEFAULT_RANGE();
    zap['y'] = make_field(Y,          "y-abscissa range")                   = JRange_t::DEFAULT_RANGE();
    zap['z'] = make_field(Z,          "ordinate range")                     = JRange_t::DEFAULT_RANGE();
    zap['X'] = make_field(logx,       "logarithmic x-axis (-XX log10 axis)");
    zap['Y'] = make_field(logy,       "logarithmic y-axis (-YY log10 axis)");
    zap['Z'] = make_field(logz,       "logarithmic z-axis");
    zap['P'] = make_field(project,    "projection")                         = "", "xy", "yx", "xz", "zx", "yz", "zy";
    zap['>'] = make_field(xLabel,     "x-axis label")                       = "";
    zap['<'] = make_field(yLabel,     "y-axis label")                       = "";
    zap['^'] = make_field(zLabel,     "z-axis label")                       = "";
    zap['S'] = make_field(markerSize, "marker size")                        = 1.0;
    zap['O'] = make_field(option,     "plotting option")                    = "";
    zap['G'] = make_field(grid,       "grid lines [X][Y]")                  = JPARSER::initialised();
    zap['B'] = make_field(batch,      "batch processing");
    zap['T'] = make_field(title,      "graphics title ("
                          << "\"" << JName_t  << "\" -> ROOT name; "
                          << "\"" << JTitle_t << "\" -> ROOT title)")       = "KM3NeT preliminary";
    zap['N'] = make_field(Ndivisions, "axis divisioning (e.g. \"X 505\")")  = JPARSER::initialised();
    zap['p'] = make_field(palette,    "palette")                            = -1;
    zap['d'] = make_field(debug)      = 0;
     
    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }


  gROOT->SetBatch(batch);

  TApplication* tp = new TApplication("user", NULL, NULL);
  TCanvas*      cv = new TCanvas("c1", "c1", canvas.x, canvas.y);

  unique_ptr<TStyle> gStyle(new JStyle("gplot", cv->GetWw(), cv->GetWh(), parameters));

  if (palette != -1) {
    gStyle->SetPalette(palette);
  }

  gROOT->SetStyle("gplot");
  gROOT->ForceStyle();


  cv->SetFillStyle(4000);
  cv->SetFillColor(kWhite);
  cv->Divide(1,1);
  cv->cd(1);

  JMarkerAttributes::getInstance().setMarkerSize(markerSize);

  Double_t xmin = numeric_limits<double>::max();
  Double_t xmax = numeric_limits<double>::lowest();
  Double_t ymin = numeric_limits<double>::max();
  Double_t ymax = numeric_limits<double>::lowest();
  Double_t zmin = numeric_limits<double>::max();
  Double_t zmax = numeric_limits<double>::lowest();

  vector<JRootObject> listOfObjects;

  //TH2* master = NULL;
  JMaster master;
  
  for (vector<JRootObjectID>::const_iterator input = inputFile.begin(); input != inputFile.end(); ++input) {

    DEBUG("Input: " << *input << endl);

    TDirectory* dir = getDirectory(*input);

    if (dir == NULL) {
      ERROR("File: " << input->getFullFilename() << " not opened." << endl);
      continue;
    }

    const TRegexp regexp(input->getObjectName());

    TIter iter(dir->GetListOfKeys());

    for (TKey* key; (key = (TKey*) iter.Next()) != NULL; ) {

      const TString tag(key->GetName());

      DEBUG("Key: " << tag << " match = " << tag.Contains(regexp) << endl);
        
      // option match

      if (tag.Contains(regexp) && isTObject(key)) {

        if        (title == JName_t) {
          title = key->GetName();
        } else if (title == JTitle_t) {
          title = key->GetTitle();
        }

        JRootObject object(key->ReadObj());

        try {
          
          TH3& h3 = dynamic_cast<TH3&>(*object);

          object = h3.Project3D(project.c_str());
        }
        catch(exception&) {}

        try {
          
          TH2& h2 = dynamic_cast<TH2&>(*object);

          h2.SetStats(stats != -1);
          
          xmin = min(xmin, h2.GetXaxis()->GetXmin());
          xmax = max(xmax, h2.GetXaxis()->GetXmax());
          ymin = min(ymin, h2.GetYaxis()->GetXmin());
          ymax = max(ymax, h2.GetYaxis()->GetXmax());
          zmin = min(zmin, logz ? h2.GetMinimum(0.0) : h2.GetMinimum());
          zmax = max(zmax, h2.GetMaximum());
        }
        catch(exception&) {}

        try {
          
          TGraph& g1 = dynamic_cast<TGraph&>(*object);
          
          for (Int_t i = 0; i != g1.GetN(); ++i) {
            xmin = min(xmin, g1.GetX()[i] - numeric_limits<float>::epsilon());
            xmax = max(xmax, g1.GetX()[i] + numeric_limits<float>::epsilon());
            ymin = min(ymin, g1.GetY()[i] - numeric_limits<float>::epsilon());
            ymax = max(ymax, g1.GetY()[i] + numeric_limits<float>::epsilon());
          }

          int ng = 0;
          
          for (vector<JRootObject>::iterator i = listOfObjects.begin(); i != listOfObjects.end(); ++i) {
            if (dynamic_cast<TGraph*>(i->get()) != NULL) {
              ++ng;
            }
          }

          static_cast<TAttMarker&>(g1) = JMarkerAttributes::getInstance().get(ng);
        }
        catch(exception&) {}
        
        try {
          
          TGraph2D& g2 = dynamic_cast<TGraph2D&>(*object);
          
          for (Int_t i = 0; i != g2.GetN(); ++i) {
            if (!logz || g2.GetZ()[i] > 0.0) {
              xmin = min(xmin, g2.GetX()[i]);
              xmax = max(xmax, g2.GetX()[i]);
              ymin = min(ymin, g2.GetY()[i]);
              ymax = max(ymax, g2.GetY()[i]);
              zmin = min(zmin, g2.GetZ()[i]);
              zmax = max(zmax, g2.GetZ()[i]);
            }
          }

          if (Z.is_valid()) {
            g2.SetMinimum(Z.getLowerLimit());
            g2.SetMaximum(Z.getUpperLimit());
          }
        }
        catch(exception&) {}
        
        try {
          
          TF2& f2 = dynamic_cast<TF2&>(*object);

          f2.SetLineColor(kRed);
          f2.SetTitle((title +  ";" + xLabel + ";" + yLabel + ";" + zLabel).c_str());

          double __xmin;
          double __xmax;
          double __ymin;
          double __ymax;

          f2.GetRange(__xmin, __ymin, __xmax, __ymax);

          xmin = min(xmin, __xmin);
          xmax = max(xmax, __xmax);
          ymin = min(ymin, __ymin);
          ymax = max(ymax, __ymax);
          zmin = min(zmin, f2.GetMinimum());
          zmax = max(zmax, f2.GetMaximum());
        }
        catch(exception&) {}

        // label

        for (TString buffer[] = { object.getLabel(), input->getFilename().c_str(), "" }, *i = buffer; *i != ""; ++i) {

          *i = (*i)(TRegexp("\\[.*\\]"));
            
          if ((*i).Length() > 2) {
            object.setLabel((*i)(1, (*i).Length() - 2));
          }
        }

        if (dynamic_cast<TH2*>     (object.get())  !=  NULL ||
            dynamic_cast<TGraph*>  (object.get())  !=  NULL ||
            dynamic_cast<TGraph2D*>(object.get())  !=  NULL ||
            dynamic_cast<TEllipse*>(object.get())  !=  NULL ||
            dynamic_cast<TLine*>   (object.get())  !=  NULL ||
            dynamic_cast<TText*>   (object.get())  !=  NULL ||
            dynamic_cast<TF2*>     (object.get())  !=  NULL) {
                  
          DEBUG("Add object: " << tag << " with label <" << object.getLabel() << ">" << endl);

          if (master == NULL) {
            //master = dynamic_cast<TH2*>(object.get());
            master = object.get();
          }

          listOfObjects.push_back(object);
          
        } else {
          ERROR("For other objects than 2D histograms, use JPlot1D" << endl);
        }
      }
    }
  }

  if (listOfObjects.empty()) {
    ERROR("Nothing to draw." << endl);
  }

  // plot frame

  cv->cd(1);

  if (option.find("COLZ") != string::npos ||
      option.find("colz") != string::npos) {
    gPad->SetRightMargin(0.20);
  }

  if (X.is_valid()) {
    xmin = X.getLowerLimit();
    xmax = X.getUpperLimit();
  }

  if (Y.is_valid()) {
    ymin = Y.getLowerLimit();
    ymax = Y.getUpperLimit();
  }
  
  if (Z.is_valid()) {
    zmin = Z.getLowerLimit();
    zmax = Z.getUpperLimit();
  } else if (zmax > zmin) {
    setRange(zmin, zmax, logz);
  }

  if (!listOfObjects.empty()) {

    if (master == NULL) {

      master = new TH2D(MASTER.c_str(), NULL, 
                        100, xmin, xmax,
                        100, ymin, ymax);

      master->SetStats(kFALSE);
    }
  }
      
  if (master == NULL) {

    TText* p = new TText(0.5, 0.5, "No data");

    p->SetTextAlign(21);
    p->SetTextAngle(45);
    p->Draw();

  } else {

    if (logx) { gPad->SetLogx(); }
    if (logy) { gPad->SetLogy(); }
    if (logz) { gPad->SetLogz(); }

    master->GetXaxis()->SetRangeUser(xmin, xmax);
    master->GetYaxis()->SetRangeUser(ymin, ymax);

    if (logx > 1) { setLogarithmicX(master); }
    if (logy > 1) { setLogarithmicY(master); }

    master->SetTitle(title.c_str());

    master->SetMinimum(zmin);
    master->SetMaximum(zmax);

    if (xLabel != "") { master->GetXaxis()->SetTitle(xLabel.c_str()); master->GetXaxis()->CenterTitle(true); }
    if (yLabel != "") { master->GetYaxis()->SetTitle(yLabel.c_str()); master->GetYaxis()->CenterTitle(true); }
    if (zLabel != "") { master->GetZaxis()->SetTitle(zLabel.c_str()); master->GetZaxis()->CenterTitle(true); }
  
    master->GetXaxis()->SetMoreLogLabels((logx == 1 && log10(xmax/xmin) < 2) ||
                                         (logx >  1 &&       xmax-xmin  < 2));
    master->GetYaxis()->SetMoreLogLabels((logy == 1 && log10(ymax/ymin) < 2) ||
                                         (logy >  1 &&       ymax-ymin  < 2));

    for (map<string, int>::const_iterator i = Ndivisions.begin();  i != Ndivisions.end(); ++i) {
      master->SetNdivisions(i->second, i->first.c_str());
    }
    
    DEBUG("Draw " << master->GetName() << ' ' << option << endl);

    master->Draw(option.c_str());
  }

  if (logx > 1 || logy > 1) {
    for (vector<JRootObject>::iterator i = listOfObjects.begin(); i != listOfObjects.end(); ++i) {
      if (dynamic_cast<TH2*>(i->get()) != master) {
        if (logx > 1) { 
          if      (setLogX<TH2>           (i->get())) {}
          else if (setLogX<TF2>           (i->get())) {}
          else if (setLogX<TGraph2DErrors>(i->get())) {}
          else if (setLogX<TGraph2D>      (i->get())) {}
          else if (setLogX<TGraphErrors>  (i->get())) {}
          else if (setLogX<TGraph>        (i->get())) {}
          else if (setLogX<TLine>         (i->get())) {}
          else if (setLogX<TEllipse>      (i->get())) {}
        }
        if (logy > 1) { 
          if      (setLogY<TH2>           (i->get())) {}
          else if (setLogY<TF2>           (i->get())) {}
          else if (setLogY<TGraph2DErrors>(i->get())) {}
          else if (setLogY<TGraph2D>      (i->get())) {}
          else if (setLogY<TGraphErrors>  (i->get())) {}
          else if (setLogY<TGraph>        (i->get())) {}
          else if (setLogY<TLine>         (i->get())) {}
          else if (setLogY<TEllipse>      (i->get())) {}
        }
      }
    }
  }

  if (grid.count('x') || grid.count('X')) { gPad->SetGridx(); }
  if (grid.count('y') || grid.count('Y')) { gPad->SetGridy(); }

  if (stats != -1)
    gStyle->SetOptStat(stats);
  else
    gStyle->SetOptFit(kFALSE);
  

  for (vector<JRootObject>::iterator i = listOfObjects.begin(); i != listOfObjects.end(); ++i) {

    if (i->get() != master) {

      DEBUG("Draw " << (*i)->GetName() << ' ' << (*i)->GetTitle() << endl);

      string buffer(option);

      buffer += "SAME";
    
      try {
      
        TH2& h2 = dynamic_cast<TH2&>(*(*i));

        h2.SetMinimum(zmin);
        h2.SetMaximum(zmax);
      }
      catch(exception&) {}

      try {
      
        TGraph& g1 = dynamic_cast<TGraph&>(*(*i));

        buffer = "P";
      }
      catch(exception&) {}

      try {
      
        TGraph2D& g2 = dynamic_cast<TGraph2D&>(*(*i));

        g2.SetMinimum(zmin);
        g2.SetMaximum(zmax);
      }
      catch(exception&) {}

      try {
      
        TF2& f2 = dynamic_cast<TF2&>(*(*i));

        f2.SetNpx(1000);
        f2.SetNpy(1000);
        /*
        f2.SetRange(xmin, ymin, zmin,
                    xmax, ymax, zmax);
        */
      }
      catch(exception&) {}

      (*i)->Draw(buffer.c_str());
    }
  }

  //gPad->RedrawAxis();
  
  cv->Update();

  if (outputFile != "") {
    cv->SaveAs(outputFile.c_str());
  }

  if (!batch) {
    tp->Run();
  }

  return (master != NULL ? 0 : 1);
}