JEnergyCorrection.hh

Go to the documentation of this file.

#ifndef __JFIT__JENERGYCORRECTION__
#define __JFIT__JENERGYCORRECTION__

#include <istream>
#include <ostream>
#include <sstream>
#include <fstream>
#include <iostream>
#include <cmath>
#include <cctype>

#include "TFile.h"
#include "TString.h"
#include "TFormula.h"

#include "JLang/JException.hh"
#include "JLang/JPointer.hh"
#include "Jeep/JeepToolkit.hh"


/**
 * \author mdejong
 */

namespace JFIT {}
namespace JPP { using namespace JFIT; }

namespace JFIT {

  using JLANG::JPointer;
  using JLANG::JProtocolException;
  using JLANG::JFileOpenException;


  /**
   * Auxiliary class for correction of energy determined by JEnergy.cc.
   * Note that the correction is applied to <tt>10log(E)</tt> with <tt>E</tt> in GeV.
   */
  class JEnergyCorrection :
    public JPointer<const TFormula>
  {
  public:
    /**
     * Default constructor.
     * No correction is applied to the energy.
     */
    JEnergyCorrection()
    {}
    

    /**
     * Constructor.
     *
     * \param  formula       formula
     */
    JEnergyCorrection(const std::string& formula)
    {
      setFormula(formula);
    }
    

    /**
     * Constructor.
     *
     * \param  fcn           pointer to ROOT formula
     */
    JEnergyCorrection(const TFormula* fcn) :
      JPointer<const TFormula>(fcn)
    {}


    /**
     * Get formula.
     *
     * \return               formula
     */
    TString getFormula() const
    {
      if (is_valid())
        return get()->GetExpFormula();
      else
        return TString();
    }


    /**
     * Set formula.
     *
     * The input can be:
     *   - TFormula compatible expression;
     *   - name of ROOT  file (extension <tt>.root</tt>);
     *   - name of ASCII file (extension <tt>.txt</tt>);
     *
     * \param  formula       formula
     */
    void setFormula(const std::string& formula)
    {
      if (formula != "") {

        TString buffer(formula.c_str());
      
        if (buffer.EndsWith(".root") ||
            buffer.EndsWith(".txt")) {
          
          load(buffer);
          
        } else {

          reset(new TFormula(JEnergyCorrection::getName(), buffer));
        }
      }
    }


    /**
     * Get corrected energy.
     *
     * \param  E             energy [GeV]
     * \return               corrected energy [GeV]
     */
    double operator()(const double E) const
    {
      if (is_valid()) {

        const double y = get()->Eval(log10(E));

        return pow(10.0, y);
      }

      return E;
    }


    /**
     * Get name of energy correction formula.
     *
     * \return               name of formula
     */
    static const char* getName() 
    {
      return "energy_correction";
    }


    /**
     * Load formula from file.
     *
     * Supported file formats:
     *   - ROOT  file (extension <tt>.root</tt>) containing TFormula object with name JEnergyCorrection::getName();
     *   - ASCII file (extension <tt>.txt</tt>)  containing TFormula compatible expression;
     *
     * Note that the method JEEP::getFullFilename() is used to search for a possible location of
     * a file with the given name in the LD_LIBRARY_PATH environment variable.
     *
     * \param  file_name     file name
     */
    void load(const char* file_name)
    {
      using namespace std;
      using namespace JLANG;
      using namespace JEEP;

      TString buffer(getFullFilename(LD_LIBRARY_PATH, file_name).c_str());

      if        (buffer.EndsWith(".root")) {

        TFile in(buffer, "READ");
        
        if (in.IsOpen()) {

          in.GetObject(JEnergyCorrection::getName(), getReference());
      
          in.Close();

        } else {

          THROW(JFileOpenException, "Error opening file " << buffer);
        }

      } else if (buffer.EndsWith(".txt")) {

        ifstream in(buffer);

        if (in) {

          buffer.ReadFile(in);

          // remove control characters

          for (Ssiz_t i = 0; i != buffer.Length(); ++i) {
            if (iscntrl(buffer[i])) {
              buffer[i] = ' ';
            }
          }
          
          in.close();

          reset(new TFormula(JEnergyCorrection::getName(), buffer));

        } else {
          
          THROW(JFileOpenException, "Error opening file " << buffer);
        }

      } else {

        throw JProtocolException("Protocol not defined.");
      }
    }


    /**
     * Store formula to file.
     *
     * Supported file formats:
     *   - ROOT  file (extension <tt>.root</tt>) containing TFormula object with name JEnergyCorrection::getName;
     *   - ASCII file (extension <tt>.txt</tt>)  containing TFormula compatible expression;
     *
     * \param  file_name     file name
     */
    void store(const char* file_name)
    {
      using namespace std;
      using namespace JLANG;
      using namespace JEEP;

      if (is_valid()) {

        const TString buffer(file_name);

        if        (buffer.EndsWith(".root")) {

          TFile out(file_name, "RECREATE");

          out.WriteObject(get(), JEnergyCorrection::getName());
          
          out.Write();
          out.Close();

        } else if (buffer.EndsWith(".txt")) {

          ofstream out(file_name);

          out << *this;

          out.close();
        
        } else {

          throw JProtocolException("Protocol not defined.");
        }
      }
    }
    

    /**
     * Read energy correction from input.
     *
     * In case a file name is specified, the method load() is used
     * to read the energy correction from the corresponding file.
     *
     * \param  in            input stream
     * \param  object        energy correction
     * \return               input stream
     */
    friend inline std::istream& operator>>(std::istream& in, JEnergyCorrection& object)
    {
      std::string buffer;

      getline(in, buffer);

      object.setFormula(buffer);

      return in;
    }


    /**
     * Write energy correction to output.
     *
     * \param  out           output stream
     * \param  object        energy correction
     * \return               output stream
     */
    friend inline std::ostream& operator<<(std::ostream& out, const JEnergyCorrection& object)
    {
      return out << object.getFormula();
    }
  };
}

#endif