offline/script/get_calibration_factor.C

0001 #include <iostream>
0002 #include <fstream>
0003 #include <TFile.h>
0004 #include <TMath.h>
0005
0006 float stringtodouble(std::string str, bool& havenan) {
0007   if (str == "nan" || str == "NaN" || str == "NAN") {
0008     havenan = true;
0009     return 999999;
0010   } else {
0011     return std::stod(str);
0012   }
0013 }
0014
0015 void reading_simulation_data(std::ifstream& simulation_data, float hcal_sim_mpv[24][64], float hcal_sim_mpverror[24][64], float hcal_sim_chi2[24][64]) {
0016   int ieta, iphi;
0017   int ieta_check = 0;
0018   int iphi_check = 0;
0019   std::string mpv_str, mpv_error_str, chi2_str, drop1_str, drop2_str;
0020   std::string line;
0021   while (std::getline(simulation_data, line)) {
0022     std::istringstream iss(line);
0023
0024     if (!(iss >> ieta >> iphi)) {
0025       std::cerr << "Error reading ieta_check and iphi_check." << std::endl;
0026       return;
0027     }
0028     if (ieta != ieta_check || iphi != iphi_check) {
0029       std::cerr << "Error reading sim ieta_check and iphi_check. Mismatch: " << ieta << " - " << iphi << "    " << ieta_check << " - " << iphi_check << std::endl;
0030       return;
0031     }
0032     if (!(iss >> mpv_str >> mpv_error_str >> drop1_str >> drop2_str >> chi2_str)) {
0033       std::cerr << "Error reading data fields." << std::endl;
0034       return;
0035     }
0036
0037     bool havenan = false;
0038     hcal_sim_mpv[ieta_check][iphi_check] = stringtodouble(mpv_str, havenan);
0039     hcal_sim_mpv[ieta_check][iphi_check] /= 32.;// Convert to low gain.

0040     hcal_sim_mpverror[ieta_check][iphi_check] = stringtodouble(mpv_error_str, havenan);
0041     hcal_sim_mpverror[ieta_check][iphi_check] /= 32.;// Convert to low gain.

0042     hcal_sim_chi2[ieta_check][iphi_check] = stringtodouble(chi2_str, havenan);
0043     if (havenan) {
0044       hcal_sim_chi2[ieta_check][iphi_check] = 999999;
0045       std::cout << "Reading Nan from Tower " << ieta_check << " - " << iphi_check << std::endl;
0046     }
0047
0048     iphi_check = iphi+1;
0049     if (iphi_check == 64) {
0050       ieta_check = ieta+1;
0051       iphi_check = 0;
0052     }
0053   }
0054 }
0055
0056 void get_calibration_factor(int min_run, int max_run) {
0057   static const int n_etabin = 24;
0058   static const int n_phibin = 64;
0059
0060   // Read simulation.

0061   std::ifstream simulation_ohcaldata("script/fittingresult_ohcal.txt");
0062   std::ifstream simulation_ihcaldata("script/fittingresult_ihcal.txt");
0063   float ohcal_sim_mpv[n_etabin][n_phibin], ohcal_sim_mpverror[n_etabin][n_phibin], ohcal_sim_chi2[n_etabin][n_phibin];
0064   float ihcal_sim_mpv[n_etabin][n_phibin], ihcal_sim_mpverror[n_etabin][n_phibin], ihcal_sim_chi2[n_etabin][n_phibin];
0065   //std::cout << "Reading ohcal simulation data" << std::endl;

0066   reading_simulation_data(simulation_ohcaldata, ohcal_sim_mpv, ohcal_sim_mpverror, ohcal_sim_chi2);
0067   //std::cout << "Reading ihcal simulation data" << std::endl;

0068   reading_simulation_data(simulation_ihcaldata, ihcal_sim_mpv, ihcal_sim_mpverror, ihcal_sim_chi2);
0069
0070   // Read precalib.

0071   std::ifstream precalib_ohcaldata("script/ohcal_precalib.txt");
0072   std::ifstream precalib_ihcaldata("script/ihcal_precalib.txt");
0073   float ohcal_precalib[n_etabin][n_phibin], ihcal_precalib[n_etabin][n_phibin];
0074   int ieta_check, iphi_check;
0075   int ieta, iphi;
0076   std::string mpv_str, mpv_error_str, chi2_str, drop1_str, drop2_str;
0077   for (int ieta = 0; ieta < n_etabin; ++ieta) {
0078     for (int iphi = 0; iphi < n_phibin; ++iphi) {
0079       precalib_ohcaldata >> ieta_check >> iphi_check >> ohcal_precalib[ieta][iphi];
0080       if (ieta_check != ieta || iphi_check != iphi) {
0081         std::cerr << "Error reading ohcal precalib data" << std::endl;
0082         return;
0083       }
0084       precalib_ihcaldata >> ieta_check >> iphi_check >> ihcal_precalib[ieta][iphi];
0085       if (ieta_check != ieta || iphi_check != iphi) {
0086         std::cerr << "Error reading ihcal precalib data" << std::endl;
0087         return;
0088       }
0089     }
0090   }
0091
0092   // Read single fit result.

0093   float ohcal_data_mpv[n_etabin][n_phibin];
0094   std::ifstream ohcalfile("output/ohcal_mpv_" + std::to_string(min_run) + "_" + std::to_string(max_run) + ".txt");
0095   if (ohcalfile.is_open()) {
0096     std::string singleline;
0097     int single_ieta, single_iphi;
0098     std::string single_mpv_str;
0099     while (std::getline(ohcalfile, singleline)) {
0100       std::istringstream single_iss(singleline);
0101       single_iss >> single_ieta >> single_iphi >> single_mpv_str;
0102       bool havenan = false;
0103       ohcal_data_mpv[single_ieta][single_iphi] = stringtodouble(single_mpv_str, havenan);
0104     }
0105   }
0106   ohcalfile.close();
0107
0108   float ihcal_data_mpv[n_etabin][n_phibin];
0109   std::ifstream ihcalfile("output/ihcal_mpv_" + std::to_string(min_run) + "_" + std::to_string(max_run) + ".txt");
0110   if (ihcalfile.is_open()) {
0111     std::string singleline;
0112     int single_ieta, single_iphi;
0113     std::string single_mpv_str;
0114     while (std::getline(ihcalfile, singleline)) {
0115       std::istringstream single_iss(singleline);
0116       single_iss >> single_ieta >> single_iphi >> single_mpv_str;
0117       bool havenan = false;
0118       ihcal_data_mpv[single_ieta][single_iphi] = stringtodouble(single_mpv_str, havenan);
0119     }
0120   }
0121   ihcalfile.close();
0122
0123   // Calculate calibration factor.

0124   float ohcal_newcalib[n_etabin][n_phibin];
0125   ofstream ohcal_outputfile("output/ohcal_calibfactor_" + std::to_string(min_run) + "_" + std::to_string(max_run) + ".txt");
0126   for (int ieta = 0; ieta < n_etabin; ++ieta) {
0127     for (int iphi = 0; iphi < n_phibin; ++iphi) {
0128       // Convert to raw data.

0129       ohcal_newcalib[ieta][iphi] = ohcal_data_mpv[ieta][iphi]/(float)ohcal_precalib[ieta][iphi];
0130       // Convert to low gain.

0131       ohcal_newcalib[ieta][iphi] /= 32.;
0132       // Get ADC/GeV.

0133       ohcal_newcalib[ieta][iphi] /= ohcal_sim_mpv[ieta][iphi];
0134       // Get GeV/ADC.

0135       if (ohcal_newcalib[ieta][iphi] != 0) {
0136         ohcal_newcalib[ieta][iphi] = 1/(float)ohcal_newcalib[ieta][iphi];
0137       }
0138
0139       ohcal_outputfile << ieta << " " << iphi << " " << ohcal_newcalib[ieta][iphi] << std::endl;
0140     }
0141   }
0142   ohcal_outputfile.close();
0143
0144   float ihcal_newcalib[n_etabin][n_phibin];
0145   ofstream ihcal_outputfile("output/ihcal_calibfactor_" + std::to_string(min_run) + "_" + std::to_string(max_run) + ".txt");
0146   for (int ieta = 0; ieta < n_etabin; ++ieta) {
0147     for (int iphi = 0; iphi < n_phibin; ++iphi) {
0148       // Convert to raw data.

0149       ihcal_newcalib[ieta][iphi] = ihcal_data_mpv[ieta][iphi]/(float)ihcal_precalib[ieta][iphi];
0150       // Convert to low gain.

0151       ihcal_newcalib[ieta][iphi] /= 32.;
0152       // Get ADC/GeV.

0153       ihcal_newcalib[ieta][iphi] /= ihcal_sim_mpv[ieta][iphi];
0154       // Get GeV/ADC.

0155       if (ihcal_newcalib[ieta][iphi] != 0) {
0156         ihcal_newcalib[ieta][iphi] = 1/(float)ihcal_newcalib[ieta][iphi];
0157       }
0158
0159       ihcal_outputfile << ieta << " " << iphi << " " << ihcal_newcalib[ieta][iphi] << std::endl;
0160     }
0161   }
0162   ihcal_outputfile.close();
0163 }