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/JSinglePointer.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/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
}


/**
 * \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;
  typedef pair<string, int>  division_type;

  vector<JRootObjectID>   inputFile;
  string                  outputFile;
  JCanvas                 canvas;
  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;
  vector<division_type>   Ndivisions;
  int                     palette;
  int                     debug;
 
  try {

    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['s'] = make_field(stats)                                            = -1;
    zap['x'] = make_field(X,          "x-abscissa range")                   = JRange_t();
    zap['y'] = make_field(Y,          "y-abscissa range")                   = JRange_t();
    zap['z'] = make_field(Z,          "ordinate range")                     = JRange_t();
    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,      "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);

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

  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;
  
  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 != JRange_t()) {
            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&) {}

        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<TF2*>     (object.get())  !=  NULL) {
                  
          DEBUG("Add object: " << tag << " with label " << object.getLabel() << endl);

          if (master == NULL) {
            master = dynamic_cast<TH2*>(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 != JRange_t()) {    
    xmin = X.getLowerLimit();
    xmax = X.getUpperLimit();
  }

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

  if (master == NULL) {

    if        (X != JRange_t() && 
               Y != JRange_t()) {

      master = new TH2D(MASTER.c_str(), NULL, 
                        100, X.getLowerLimit(), X.getUpperLimit(),
                        100, Y.getLowerLimit(), Y.getUpperLimit());

      master->SetStats(kFALSE);

    } else if (xmin < xmax && ymin < ymax) {

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

      master->SetStats(kFALSE);

    } else {

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

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

  if (master != NULL) {

    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 (vector<division_type>::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(MAKE_CSTRING(option << (MASTER == master->GetName() ? "AXIS" : "")));
  }

  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) {

    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);
}