PPSObjects/src/LHCOpticalFunctionsSet.cc

0001 // Original Author:  Jan Kašpar
0002
0003 #include "CondFormats/PPSObjects/interface/LHCOpticalFunctionsSet.h"
0004
0005 #include "FWCore/Utilities/interface/Exception.h"
0006
0007 #include "TFile.h"
0008 #include "TGraph.h"
0009
0010 //----------------------------------------------------------------------------------------------------
0011
0012 LHCOpticalFunctionsSet::LHCOpticalFunctionsSet(const std::string &fileName, const std::string &directoryName, double z)
0013     : m_z(z) {
0014   TFile *f_in = TFile::Open(fileName.c_str());
0015   if (f_in == nullptr)
0016     throw cms::Exception("LHCOpticalFunctionsSet") << "Cannot open file " << fileName << ".";
0017
0018   std::vector<TGraph *> graphs(nFunctions);
0019   for (unsigned int fi = 0; fi < nFunctions; ++fi) {
0020     std::string tag;
0021     if (fi == evx)
0022       tag = "v_x";
0023     else if (fi == eLx)
0024       tag = "L_x";
0025     else if (fi == e14)
0026       tag = "E_14";
0027     else if (fi == exd)
0028       tag = "x_D";
0029     else if (fi == evpx)
0030       tag = "vp_x";
0031     else if (fi == eLpx)
0032       tag = "Lp_x";
0033     else if (fi == e24)
0034       tag = "E_24";
0035     else if (fi == expd)
0036       tag = "xp_D";
0037     else if (fi == e32)
0038       tag = "E_32";
0039     else if (fi == evy)
0040       tag = "v_y";
0041     else if (fi == eLy)
0042       tag = "L_y";
0043     else if (fi == eyd)
0044       tag = "y_D";
0045     else if (fi == e42)
0046       tag = "E_42";
0047     else if (fi == evpy)
0048       tag = "vp_y";
0049     else if (fi == eLpy)
0050       tag = "Lp_y";
0051     else if (fi == eypd)
0052       tag = "yp_D";
0053     else
0054       throw cms::Exception("LHCOpticalFunctionsSet") << "Invalid tag for optical functions: \"" << fi << "\"";
0055
0056     std::string objPath = directoryName + "/g_" + tag + "_vs_xi";
0057     auto gr_obj = dynamic_cast<TGraph *>(f_in->Get(objPath.c_str()));
0058     if (!gr_obj)
0059       throw cms::Exception("LHCOpticalFunctionsSet")
0060           << "Cannot load object " << objPath << " from file " << fileName << ".";
0061
0062     graphs[fi] = gr_obj;
0063   }
0064
0065   const unsigned int num_xi_vals = graphs[0]->GetN();
0066   m_xi_values.resize(num_xi_vals);
0067
0068   m_fcn_values.resize(nFunctions);
0069
0070   for (unsigned int fi = 0; fi < nFunctions; ++fi)
0071     m_fcn_values[fi].resize(num_xi_vals);
0072
0073   for (unsigned int pi = 0; pi < num_xi_vals; ++pi) {
0074     const double xi = graphs[0]->GetX()[pi];
0075     m_xi_values[pi] = xi;
0076
0077     for (unsigned int fi = 0; fi < m_fcn_values.size(); ++fi)
0078       m_fcn_values[fi][pi] = graphs[fi]->Eval(xi);
0079   }
0080
0081   delete f_in;
0082 }