hcal_calib_year2/genFinalCalib/genFinalCalib.C

0001 #include"tsc_cos_merge.C"
0002
0003 #include <TSystem.h>
0004
0005 #include <format>
0006 #include <fstream>
0007 #include <set>
0008 #include <sstream>
0009 #include <string>
0010
0011 struct TowerData {
0012     int tower_ieta;
0013     int tower_iphi;
0014     float calib_factor;
0015     int run_range_low;
0016     int run_range_high;
0017 };
0018 void loadCSV(const std::string& filename, std::vector<TowerData>& data);
0019 std::vector<std::vector<float>> getTowersForRun(const std::vector<TowerData>& data, int run_number);
0020
0021
0022 //make sure TS_list has relative or full paths to actual root files
0023 void genFinalCalib(const std::string &runs = "cos_runs.txt")
0024 {
0025
0026   ///////////////////////////////////
0027   // Load the half hieght data from the CSV file
0028   ///////////////////////////////////////
0029   std::vector<TowerData> hh_data;
0030   loadCSV("halfHieghtList.csv", hh_data);
0031
0032
0033   /////////////////////////////////////////////////
0034   // Get cosmic calibrations
0035   ///////////////////////////////////////
0036   std::ifstream runlist(runs);
0037   std::string line;
0038   std::vector<int> cosmic_runnumbers;
0039   while(getline(runlist,line))
0040   {
0041     cosmic_runnumbers.push_back(stoi(line));
0042
0043     //std::string oh_out = "fitResults/oh_cosmicsMerged_fitted_" + line + ".root";
0044     std::string oh_out = "../tsc_cos_comb/fitResults/avgCorr.root";
0045 //    std::string ih_out = "fitResults/ih_cosmicsMerged_fitted_" + line + ".root";
0046
0047     std::string cosmic_CDB_file_oh = std::format("/sphenix/u/bseidlitz/work/macros/calibrations/calo/hcal_towerSlope_y2/cosmicCalibFiles/ohcal_cosmic_calibration_{}.root",line);
0048
0049     tsc_cos_merge(oh_out,cosmic_CDB_file_oh,std::format("cdbFiles/ohcal_cdb_tsc_cos_calib_{}.root",line),0);
0050
0051   }
0052
0053
0054   ///////////////////////////////////////
0055   // temperature correction run-by-run
0056   ///////////////////////////////////////
0057   std::ifstream runlistBeam("runList.txt");
0058   std::set<int> tempset;
0059   while(getline(runlistBeam,line)){
0060     int runnumber = stoi(line);
0061     tempset.insert(runnumber);
0062   }
0063   std::vector<int> runnumbersBeam(tempset.begin(),tempset.end());
0064
0065   TFile* ftemp = new TFile("../fileQA/output/anaOut.root");
0066
0067   TH1* h_runnumbers_temp = (TH1F*) ftemp->Get("runNumberHist");
0068   TH1* h_temp_run = (TH1F*) ftemp->Get("h_temp_run");
0069
0070
0071   ///////////////////////////
0072   // get cosmic run temps
0073   ///////////////////////////
0074   // get the temperature of the cosmics runs
0075   // could be improved, right now it uses the cosest beam run
0076   std::vector<float> cosmic_temps;
0077   for (int cos_run : cosmic_runnumbers){
0078     int minDelta = 1000000;
0079 //    int minRun = 0;
0080     float cos_temp = 0;
0081     for( int it=0; it<h_runnumbers_temp->GetNbinsX(); it++){
0082       int delta = std::abs(cos_run - h_runnumbers_temp->GetBinContent(it));
0083       if (delta < minDelta  && h_temp_run->GetBinContent(it) > 10){
0084 //        minRun = h_runnumbers_temp->GetBinContent(it);
0085         minDelta = delta;
0086         cos_temp =h_temp_run->GetBinContent(it);
0087       }
0088     }
0089     cosmic_temps.push_back(cos_temp);
0090     std::cout << "cosmics run " << cos_run << "  temp="  << cos_temp << " size=" << cosmic_temps.size() << std::endl;
0091   }
0092
0093
0094   ///////////////////////////////////////
0095   // Generate the run-by-run calibrations
0096   ///////////////////////////////////////
0097   // This follows the logic laid out in the note
0098   // cosmics calibration is the base
0099   // the tsc is the used to adjust the cosmic calib
0100   // then the temperture difference between the cosmic
0101   // run and the run of interest to correct the
0102   // adjusted cosmics for the rinal calibration
0103   for(int irb=0; irb<runnumbersBeam.size(); irb++)
0104   {
0105     int runnumber = runnumbersBeam[irb];
0106
0107     int runnumber_last;
0108     if (irb+1  < runnumbersBeam.size() ) runnumber_last = runnumbersBeam[irb+1]-1;
0109     else runnumber_last = 53880;
0110
0111     std::cout << "----------------------------------" << std::endl;
0112     std::cout << "starting run " << runnumber << std::endl;
0113
0114     int asoCosmic = 0;
0115     int minCosDelta = 1000000;
0116     float asoCosmicTemp = 0;
0117     for(int ic=0; ic< cosmic_runnumbers.size(); ic++){
0118        int cos_run = cosmic_runnumbers[ic];
0119        int delta = runnumber - cos_run;
0120        if (minCosDelta  > delta && delta > 0){
0121          asoCosmic     = cos_run;
0122          asoCosmicTemp =  cosmic_temps[ic];
0123          minCosDelta = delta;
0124        }
0125     }
0126     std::cout << "closeset cosmic run " << asoCosmic << std::endl;
0127
0128
0129     for( int it=0; it<h_runnumbers_temp->GetNbinsX(); it++){
0130       if(runnumber == h_runnumbers_temp->GetBinContent(it)){
0131          float temp = h_temp_run->GetBinContent(it);
0132          float alpha = 0.045;
0133          float deltaT =  temp - asoCosmicTemp;
0134
0135          std::cout << "of this run = " << temp << "   temp of cosmic run =" << asoCosmicTemp << std::endl;
0136
0137          float correction = 1+ deltaT*alpha;
0138          if (temp < 10) continue;
0139          std::cout << "correction=" << correction << std::endl;
0140
0141          // std::string cdbFileName = std::format("cdbFiles/cdb_ohcal_tsc_temp_cosmic_{}.root",runnumber);
0142          //std::string cosmic_CDB_file_oh = std::format("/sphenix/u/bseidlitz/work/macros/calibrations/calo/hcal_towerSlope_y2/cosmicCalibFiles/ohcal_cosmic_calibration_{}.root",asoCosmic);
0143          std::string cosmic_CDB_file_oh = std::format("cdbFiles/ohcal_cdb_tsc_cos_calib_{}.root",asoCosmic);
0144
0145           // gen calibration file with cosmics, tsc and temp correction
0146           // but no half hieght corrections
0147          //genCdbCorr(correction,cosmic_CDB_file_oh,cdbFileName,0);
0148
0149          // with half hieght correction
0150          std::vector<std::vector<float>> hh_towersForRun = getTowersForRun(hh_data, runnumber);
0151          std::string cdbFileName_hh = std::format("final_ohcal_calib_y2_namceChange/cdb_ohcal_tsc_temp_cosmic_hh_{}_{}.root",runnumber,runnumber_last);
0152          genCdbCorr_HH(correction,cosmic_CDB_file_oh,cdbFileName_hh,0,hh_towersForRun);
0153       }
0154     }
0155   }
0156
0157
0158
0159 }//end macro
0160
0161
0162
0163 void loadCSV(const std::string& filename, std::vector<TowerData>& data) {
0164     std::ifstream file(filename);
0165     std::string line;
0166     if(!file.is_open())
0167     {
0168       std::cout << std::endl << "ccan't find file" << std::endl;
0169       gSystem->Exit(1);
0170     }
0171
0172     while (getline(file, line)) {
0173         std::stringstream ss(line);
0174         TowerData towerData{};
0175
0176         // Parse the line and fill the struct
0177         char comma;  // To consume the commas between values
0178         ss >> towerData.tower_ieta >> comma
0179            >> towerData.tower_iphi >> comma
0180            >> towerData.calib_factor >> comma
0181            >> towerData.run_range_low >> comma
0182            >> towerData.run_range_high;
0183
0184         // Store the parsed data
0185         data.push_back(towerData);
0186     }
0187 }
0188
0189
0190 std::vector<std::vector<float>> getTowersForRun(const std::vector<TowerData>& data, int run_number)
0191 {
0192     std::vector<std::vector<float>> towersForRun;
0193
0194     // Loop through all the entries and check if the run_number is within the range
0195     for (const auto& towerData : data) {
0196         if (run_number > towerData.run_range_low && run_number < towerData.run_range_high) {
0197             // If it matches the range, store the tower data
0198             towersForRun.push_back({static_cast<float>(towerData.tower_ieta),
0199                                     static_cast<float>(towerData.tower_iphi),
0200                                     towerData.calib_factor});
0201         }
0202     }
0203
0204     return towersForRun;
0205 }
0206