hcal_tsc_cos/figureMaker/plot.C

0001 #include "sPhenixStyle.C"
0002
0003 void myText(Double_t x,Double_t y,Color_t color, const char *text, Double_t tsize=0.04) {
0004
0005   TLatex l; l.SetTextAlign(12); l.SetTextSize(tsize);
0006   l.SetNDC();
0007   l.SetTextColor(color);
0008   l.DrawLatex(x,y,text);
0009 }
0010
0011
0012 TH2F* DivideColumnsByAverage(TH2F* hist);
0013
0014
0015 void plot(){
0016
0017
0018 SetsPhenixStyle();
0019
0020 std::string file_path = "runList.txt";
0021
0022 std::ifstream file(file_path);
0023
0024 std::vector<int> run_numbers;
0025 std::string line;
0026 while (std::getline(file, line)) {
0027     std::istringstream iss(line);
0028     int run_number;
0029     if (iss >> run_number) {
0030         run_numbers.push_back(run_number);
0031     } else {
0032         std::cerr << "Warning: Invalid integer on line: " << line << std::endl;
0033     }
0034 }
0035
0036 file.close();
0037 if (run_numbers.empty()) {
0038     std::cerr << "Error: No valid run numbers found in the file." << std::endl;
0039     return 1;
0040 }
0041
0042 int first_run = run_numbers.front();
0043 int last_run = run_numbers.back();
0044
0045
0046 TFile* f_in = new TFile("../fitout_hcalin.root");
0047 TFile* f_out = new TFile("../fitout_hcalout.root");
0048
0049 TH2F* tsc_in = (TH2F*) f_in->Get("corrPat");
0050 TH2F* tsc_out = (TH2F*) f_out->Get("corrPat");
0051
0052
0053 TCanvas* c1 = new TCanvas("c1","c1",600,600);
0054 tsc_in->Draw("COLZ");
0055 tsc_in->SetXTitle("ieta");
0056 tsc_in->SetYTitle("iphi");
0057
0058 gPad->SetRightMargin(0.15);
0059
0060 myText(0.1,0.98,1,Form("%d-%d",first_run,last_run));
0061 myText(0.7,0.98,1,"#it{#bf{sPHENIX}} Internal");
0062 myText(0.5,0.98,1,"HCal IN");
0063
0064 c1->SaveAs("../plots/tsc_in.pdf");
0065
0066
0067 TCanvas* c2 = new TCanvas("c2","c2",600,600);
0068 tsc_out->Draw("COLZ");
0069 tsc_out->SetXTitle("ieta");
0070 tsc_out->SetYTitle("iphi");
0071
0072 gPad->SetRightMargin(0.15);
0073
0074 myText(0.1,0.98,1,Form("%d-%d",first_run,last_run));
0075 myText(0.7,0.98,1,"#it{#bf{sPHENIX}} Internal");
0076 myText(0.5,0.98,1,"HCal OUT");
0077
0078 c2->SaveAs("../plots/tsc_out.pdf");
0079
0080
0081 ////////////////////////////////////////
0082 // plot beam data output from calo valid
0083 int runCV = 23727;
0084 TFile* fin = new TFile(Form("../condor/combine_out/out_%d.root",runCV));
0085 TH2F* h_ohcal_hit = (TH2F*) fin->Get("h_ohcal_etaphi_wQA");
0086 TH2F* h_ihcal_hit = (TH2F*) fin->Get("h_ihcal_etaphi_wQA");
0087 //TH2F* h_ihcal_hit = (TH2F*) fin->Get("h_ihcal_etaphi_wQA");
0088
0089 TH2F*  h_ohcal_hit_avg = DivideColumnsByAverage(h_ohcal_hit);
0090 TH2F*  h_ihcal_hit_avg = DivideColumnsByAverage(h_ihcal_hit);
0091
0092 TCanvas* c3 = new TCanvas("c3","c3",600,600);
0093 h_ohcal_hit_avg->Draw("COLZ");
0094 h_ohcal_hit_avg->SetXTitle("ohcal ieta");
0095 h_ohcal_hit_avg->SetYTitle("ohcal iphi");
0096 h_ohcal_hit_avg->GetZaxis()->SetRangeUser(0.5,1.5);
0097
0098 gPad->SetRightMargin(0.15);
0099
0100 myText(0.1,0.98,1,Form("%d",runCV));
0101 myText(0.7,0.98,1,"#it{#bf{sPHENIX}} Internal");
0102 myText(0.5,0.98,1,"HCal OUT");
0103
0104 c3->SaveAs("../plots/hitsNorm_ohcal.pdf");
0105
0106
0107 TCanvas* c4 = new TCanvas("c4","c4",600,600);
0108 h_ihcal_hit_avg->Draw("COLZ");
0109 h_ihcal_hit_avg->SetXTitle("ihcal ieta");
0110 h_ihcal_hit_avg->SetYTitle("ihcal iphi");
0111 h_ihcal_hit_avg->GetZaxis()->SetRangeUser(0.8,1.2);
0112
0113 gPad->SetRightMargin(0.15);
0114
0115 myText(0.1,0.98,1,Form("%d",runCV));
0116 myText(0.7,0.98,1,"#it{#bf{sPHENIX}} Internal");
0117 myText(0.5,0.98,1,"HCal IN");
0118
0119 c4->SaveAs("../plots/hitsNorm_ihcal.pdf");
0120
0121
0122 }
0123
0124
0125
0126 TH2F* DivideColumnsByAverage(TH2F* hist) {
0127     int numBinsX = hist->GetNbinsX();
0128     int numBinsY = hist->GetNbinsY();
0129
0130     TH2F* dividedHistogram = new TH2F("dividedHistogram", "Histogram with Columns Divided by Average",
0131                                       numBinsX, hist->GetXaxis()->GetXmin(), hist->GetXaxis()->GetXmax(),
0132                                       numBinsY, hist->GetYaxis()->GetXmin(), hist->GetYaxis()->GetXmax());
0133
0134     for (int binX = 1; binX <= numBinsX; binX++) {
0135         double sum = 0.0;
0136         for (int binY = 1; binY <= numBinsY; binY++) {
0137             sum += hist->GetBinContent(binX, binY);
0138         }
0139         double average = (numBinsY > 0) ? sum / numBinsY : 0.0;
0140         for (int binY = 1; binY <= numBinsY; binY++) {
0141             double content = hist->GetBinContent(binX, binY);
0142             double dividedContent = (average != 0.0) ? content / average : 0.0;
0143             dividedHistogram->SetBinContent(binX, binY, dividedContent);
0144         }
0145     }
0146     return dividedHistogram;
0147 }