JProject2D.cc

Go to the documentation of this file.

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

#include "TROOT.h"
#include "TFile.h"
#include "TKey.h"
#include "TH2.h"
#include "TProfile2D.h"
#include "TString.h"
#include "TRegexp.h"

#include "JTools/JRange.hh"

#include "JGizmo/JRootObjectID.hh"
#include "JGizmo/JGizmoToolkit.hh"

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


/**
 * \file
 * Auxiliary program to project 2D histograms.
 * \author mdejong
 */
int main(int argc, char **argv)
{
  using namespace std;
  using namespace JPP;

  typedef JRange<double>  JRange_t;
  
  vector<JRootObjectID> inputFile;
  string                outputFile;
  vector<JRange_t>      X;
  vector<JRange_t>      Y;
  char                  project;
  bool                  overflow;
  string                format;
  int                   debug;

  try {

    JParser<> zap("Auxiliary program to project 2D histograms.");

    zap['f'] = make_field(inputFile,  "<input file>:<object name>");
    zap['P'] = make_field(project,    "projection")                  = ' ', 'x', 'X', 'y', 'Y';
    zap['o'] = make_field(outputFile, "ROOT file with histogram(s)") = "project.root";
    zap['x'] = make_field(X,          "x-abscissa ranges")           = JPARSER::initialised();
    zap['y'] = make_field(Y,          "y-abscissa ranges")           = JPARSER::initialised();
    zap['O'] = make_field(overflow);
    zap['F'] = make_field(format,     "format, e.g. \"%s %i\" or \"%s %f %f\"") = "";
    zap['d'] = make_field(debug)      = 1;

    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }

  const bool px = (project == 'x' || project == 'X' || !Y.empty());  // projection on x-axis
  const bool py = (project == 'y' || project == 'Y' || !X.empty());  // projection on y-axis

  if (px == py) {
    FATAL("Invalid operation: "
          << (px ? "" : "no") << " X projection " << " and " 
          << (py ? "" : "no") << " X projection " << endl);
  }

  if (format == "") {

    format = "%s_";
    
    if (px) { format += "px"; }
    if (py) { format += "py"; }

    if (X.empty() && Y.empty()) 
      format += "[%i]";
    else
      format += "[%f,%f]";
  }
  
  vector<TObject*> listOfObjects;

  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)) {
        
        TH2* h2 = dynamic_cast<TH2*>(key->ReadObj());

        if (h2 != NULL) {
          
          if        (px) {

            if (Y.empty()) {
              
              for (Int_t i = (overflow ? 0 : 1); i <= h2->GetYaxis()->GetNbins() + (overflow ? 1 : 0); ++i) {
                listOfObjects.push_back(h2->ProjectionX(TString::Format(format.c_str(), h2->GetName(), i),
                                                        i,
                                                        i));
              }

            } else {

              for (Int_t i = 0; i != (Int_t) Y.size(); ++i) {
                listOfObjects.push_back(h2->ProjectionX(TString::Format(format.c_str(), h2->GetName(), Y[i].getLowerLimit(), Y[i].getUpperLimit()),
                                                        h2->GetYaxis()->FindBin(Y[i].getLowerLimit()),
                                                        h2->GetYaxis()->FindBin(Y[i].getUpperLimit()) - 1));
              }
            }
            
          } else if (py) {

            if (X.empty()) {
              
              for (Int_t i = (overflow ? 0 : 1); i <= h2->GetXaxis()->GetNbins() + (overflow ? 1 : 0); ++i) {
                listOfObjects.push_back(h2->ProjectionY(TString::Format(format.c_str(), h2->GetName(), i), i, i));
              }

            } else {

              for (Int_t i = 0; i != (Int_t) X.size(); ++i) {
                listOfObjects.push_back(h2->ProjectionY(TString::Format(format.c_str(), h2->GetName(), X[i].getLowerLimit(), X[i].getUpperLimit()),
                                                        h2->GetXaxis()->FindBin(X[i].getLowerLimit()),
                                                        h2->GetXaxis()->FindBin(X[i].getUpperLimit()) - 1));
              }
            }
          }
        }
      }
    }
  }

  if (!listOfObjects.empty()) {

    TFile out(outputFile.c_str(), "recreate");
    
    for (vector<TObject*>::const_iterator i = listOfObjects.begin(); i != listOfObjects.end(); ++i) {
      (*i)->Write();
    }
      
    out.Write();
    out.Close();
  }
}