LHCOpticalFunctionsSet.cc

Go to the documentation of this file.

// Original Author:  Jan Kašpar

#include "CondFormats/CTPPSReadoutObjects/interface/LHCOpticalFunctionsSet.h"

#include "FWCore/Utilities/interface/Exception.h"

#include "TFile.h"
#include "TGraph.h"

//----------------------------------------------------------------------------------------------------

LHCOpticalFunctionsSet::LHCOpticalFunctionsSet(const std::string &fileName, const std::string &directoryName, double z) :
  m_z(z)
{
  TFile *f_in = TFile::Open(fileName.c_str());
  if (f_in == nullptr)
    throw cms::Exception("LHCOpticalFunctionsSet") << "Cannot open file " << fileName << ".";

  std::vector<TGraph*> graphs(nFunctions);
  for (unsigned int fi = 0; fi < nFunctions; ++fi)
  {
    std::string tag;
    if (fi == evx) tag = "v_x";
    else if (fi == eLx) tag = "L_x";
    else if (fi == e14) tag = "E_14";
    else if (fi == exd) tag = "x_D";
    else if (fi == evpx) tag = "vp_x";
    else if (fi == eLpx) tag = "Lp_x";
    else if (fi == e24) tag = "E_24";
    else if (fi == expd) tag = "xp_D";
    else if (fi == e32) tag = "E_32";
    else if (fi == evy) tag = "v_y";
    else if (fi == eLy) tag = "L_y";
    else if (fi == eyd) tag = "y_D";
    else if (fi == e42) tag = "E_42";
    else if (fi == evpy) tag = "vp_y";
    else if (fi == eLpy) tag = "Lp_y";
    else if (fi == eypd) tag = "yp_D";
    else
      throw cms::Exception("LHCOpticalFunctionsSet") << "Invalid tag for optical functions: \"" << fi << "\"";

    std::string objPath = directoryName + "/g_" + tag + "_vs_xi";
    auto gr_obj = dynamic_cast<TGraph*>( f_in->Get(objPath.c_str()) );
    if (!gr_obj)
      throw cms::Exception("LHCOpticalFunctionsSet") << "Cannot load object " << objPath << " from file " << fileName << ".";

    graphs[fi] = gr_obj;
  }

  const unsigned int num_xi_vals = graphs[0]->GetN();
  m_xi_values.resize(num_xi_vals);

  m_fcn_values.resize(nFunctions);

  for (unsigned int fi = 0; fi < nFunctions; ++fi)
    m_fcn_values[fi].resize(num_xi_vals);

  for (unsigned int pi = 0; pi < num_xi_vals; ++pi)
  {
    const double xi = graphs[0]->GetX()[pi];
    m_xi_values[pi] = xi;

    for (unsigned int fi = 0; fi < m_fcn_values.size(); ++fi)
      m_fcn_values[fi][pi] = graphs[fi]->Eval(xi);
  }

  delete f_in;
}