Jet-Study/src/JetUtils.cc

0001 #include "JetUtils.h"
0002 #include <math.h>
0003 #include <sstream>
0004 #include <fstream>
0005
0006 using std::stringstream;
0007 using std::cout;
0008 using std::endl;
0009 using std::make_pair;
0010
0011 vector<string> JetUtils::m_triggers = {"Clock"
0012                                       ,"ZDC South"
0013                                       ,"ZDC North"
0014                                       ,"ZDC Coincidence"
0015                                       ,"HCAL Singles/Coincidence"
0016                                       ,"Clock2"
0017                                       ,"unknown6"
0018                                       ,"unknown7"
0019                                       ,"MBD S >= 1"
0020                                       ,"MBD N >= 1"
0021                                       ,"MBD N&S >= 1"
0022                                       ,"MBD N&S >= 2"
0023                                       ,"MBD N&S >= 1, vtx < 10 cm"
0024                                       ,"MBD N&S >= 1, vtx < 30 cm"
0025                                       ,"MBD N&S >= 1, vtx < 60 cm"
0026                                       ,"HCAL Singles + MBD NS >= 1"
0027                                       ,"Jet 6 GeV + MBD NS >= 1"
0028                                       ,"Jet 8 GeV + MBD NS >= 1"
0029                                       ,"Jet 10 GeV + MBD NS >= 1"
0030                                       ,"Jet 12 GeV + MBD NS >= 1"
0031                                       ,"Jet 6 GeV"
0032                                       ,"Jet 8 GeV"
0033                                       ,"Jet 10 GeV"
0034                                       ,"Jet 12 GeV"
0035                                       ,"Photon 2 GeV+ MBD NS >= 1"
0036                                       ,"Photon 3 GeV + MBD NS >= 1"
0037                                       ,"Photon 4 GeV + MBD NS >= 1"
0038                                       ,"Photon 5 GeV + MBD NS >= 1"
0039                                       ,"Photon 2 GeV"
0040                                       ,"Photon 3 GeV"
0041                                       ,"Photon 4 GeV"
0042                                       ,"Photon 5 GeV"
0043                                       ,"Jet 6 GeV, MBD N&S >= 1, vtx < 10 cm"
0044                                       ,"Jet 8 GeV, MBD N&S >= 1, vtx < 10 cm"
0045                                       ,"Jet 10 GeV, MBD N&S >= 1, vtx < 10 cm"
0046                                       ,"Jet 12 GeV, MBD N&S >= 1, vtx < 10 cm"
0047                                       ,"Photon 3 GeV, MBD N&S >= 1, vtx < 10 cm"
0048                                       ,"Photon 4 GeV, MBD N&S >= 1, vtx < 10 cm"
0049                                       ,"Photon 5 GeV, MBD N&S >= 1, vtx < 10 cm"
0050                                       ,"unknown39"
0051                                       ,"MBD Laser"
0052                                       ,"None"
0053                                      };
0054
0055 int JetUtils::readEventList(const string &input, vector<pair<int, int>> &vec, int nEvents, bool verbose)
0056 {
0057   std::ifstream file(input);
0058
0059   // Check if the file was successfully opened
0060   if (!file.is_open())
0061   {
0062     cout << "Failed to open file: " << input << endl;
0063     return 1;
0064   }
0065
0066   string line;
0067
0068   // keep header
0069   std::getline(file, line);
0070
0071   // loop over each run
0072   while (std::getline(file, line))
0073   {
0074     std::istringstream lineStream(line);
0075
0076     int run;
0077     int event;
0078     char comma;
0079
0080     if (lineStream >> run >> comma >> event)
0081     {
0082       for (int i = event - nEvents; i <= event + nEvents; ++i)
0083       {
0084         vec.push_back(make_pair(run, i));
0085       }
0086     }
0087     else
0088     {
0089       cout << "Failed to parse line: " << line << endl;
0090       return 1;
0091     }
0092   }
0093
0094   // print the events of interest
0095   if (verbose)
0096   {
0097     for (auto p : vec)
0098     {
0099       cout << "Run: " << p.first << ", Event: " << p.second << endl;
0100     }
0101   }
0102
0103   // Close the file
0104   file.close();
0105
0106   return 0;
0107 }
0108
0109 bool JetUtils::failsLoEmFracETCut(float emFrac, float ET)
0110 {
0111   return emFrac < 0.1 && ET > (50*emFrac+20);
0112 }
0113
0114 bool JetUtils::failsHiEmFracETCut(float emFrac, float ET)
0115 {
0116   return emFrac > 0.9 && ET > (-50*emFrac+70);
0117 }
0118
0119 bool JetUtils::check_bad_jet_eta(float jet_eta, float zvtx, float jet_radius)
0120 {
0121   float emcal_mineta = get_emcal_mineta_zcorrected(zvtx);
0122   float emcal_maxeta = get_emcal_maxeta_zcorrected(zvtx);
0123   float ihcal_mineta = get_ihcal_mineta_zcorrected(zvtx);
0124   float ihcal_maxeta = get_ihcal_maxeta_zcorrected(zvtx);
0125   float ohcal_mineta = get_ohcal_mineta_zcorrected(zvtx);
0126   float ohcal_maxeta = get_ohcal_maxeta_zcorrected(zvtx);
0127   float minlimit = emcal_mineta;
0128   if (ihcal_mineta > minlimit) minlimit = ihcal_mineta;
0129   if (ohcal_mineta > minlimit) minlimit = ohcal_mineta;
0130   float maxlimit = emcal_maxeta;
0131   if (ihcal_maxeta < maxlimit) maxlimit = ihcal_maxeta;
0132   if (ohcal_maxeta < maxlimit) maxlimit = ohcal_maxeta;
0133   minlimit += jet_radius;
0134   maxlimit -= jet_radius;
0135   return jet_eta < minlimit || jet_eta > maxlimit;
0136 }
0137
0138 float JetUtils::get_ohcal_mineta_zcorrected(float zvtx)
0139 {
0140   float z = minz_OH - zvtx;
0141   float eta_zcorrected = asinh(z / (float) radius_OH);
0142   return eta_zcorrected;
0143 }
0144
0145 float JetUtils::get_ohcal_maxeta_zcorrected(float zvtx)
0146 {
0147   float z = maxz_OH - zvtx;
0148   float eta_zcorrected = asinh(z / (float) radius_OH);
0149   return eta_zcorrected;
0150 }
0151
0152 float JetUtils::get_ihcal_mineta_zcorrected(float zvtx)
0153 {
0154   float z = minz_IH - zvtx;
0155   float eta_zcorrected = asinh(z / (float) radius_IH);
0156   return eta_zcorrected;
0157 }
0158
0159 float JetUtils::get_ihcal_maxeta_zcorrected(float zvtx)
0160 {
0161   float z = maxz_IH - zvtx;
0162   float eta_zcorrected = asinh(z / (float) radius_IH);
0163   return eta_zcorrected;
0164 }
0165
0166 float JetUtils::get_emcal_mineta_zcorrected(float zvtx)
0167 {
0168   float z = minz_EM - zvtx;
0169   float eta_zcorrected = asinh(z / (float) radius_EM);
0170   return eta_zcorrected;
0171 }
0172
0173 float JetUtils::get_emcal_maxeta_zcorrected(float zvtx)
0174 {
0175   float z = maxz_EM - zvtx;
0176   float eta_zcorrected = asinh(z / (float) radius_EM);
0177   return eta_zcorrected;
0178 }
0179
0180 vector<string> JetUtils::split(const string &s, const char delimiter) {
0181     vector<string> result;
0182
0183     stringstream ss(s);
0184     string temp;
0185
0186     while(getline(ss,temp,delimiter)) {
0187         if(!temp.empty()) {
0188             result.push_back(temp);
0189         }
0190     }
0191
0192     return result;
0193 }