JPythia.cc

Go to the documentation of this file.

 
#include <string>
#include <iostream>
#include <iomanip>
 
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
 
#include "JSirene/JPythia.hh"
#include "JGizmo/JManager.hh"
 
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"
 
 
/**
 * \file
 * Auxiliary program histogram plot EM-equivalent energy for hadrons.
 * \author mdejong
 */
int main(int argc, char **argv)
{
  using namespace std;
  using namespace JPP;
 
  string   outputFile;
  double   epsilon;
  bool     debug;
 
  try {
 
    JParser<> zap("Auxiliary program histogram plot EM-equivalent energy for hadrons."\
                  "\nNote that the results have a small offset to effectively be visualised.");
 
    zap['o'] = make_field(outputFile)  = "pythia.root";
    zap['e'] = make_field(epsilon)     = 0.004;
    zap['d'] = make_field(debug);
    
    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }
 
 
 
  JGIZMO::JManager<int, TH1D> zmap(new TH1D("E[%]",  NULL, 500, -0.8, 9.0));
 
  const int buffer[] = {
    TRACK_TYPE_PHOTON,
    TRACK_TYPE_NEUTRAL_PION,
    TRACK_TYPE_NEUTRAL_ANTIPION,
    TRACK_TYPE_ELECTRON,
    TRACK_TYPE_ANTIELECTRON,
    TRACK_TYPE_CHARGED_PION_PLUS,
    TRACK_TYPE_CHARGED_PION_MINUS,
    -1 };
 
  map<int, double> Y;
 
  for (const int *type = buffer; *type != -1; ++type) {
    Y[*type] = epsilon * (type - buffer);
  }
  
  for (int i = 1; i <= zmap->GetNbinsX(); ++i) {
 
    const double x = zmap->GetBinCenter(i);
    const double E = pow(10.0, x);
 
    for (const int *type = buffer; *type != -1; ++type) {
      zmap[*type]->SetBinContent(i, pythia(*type, E)/E + Y[*type]);
    }
  }
 
  zmap.Write(outputFile.c_str());
}