Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:58:53

0001 
0002 {
0003   #include <map>
0004 
0005   gStyle -> SetOptStat(0);
0006 //  gStyle->SetPalette(1);
0007   gStyle->SetPadGridX(1);
0008   gStyle->SetPadGridY(1);
0009 
0010 TString imgpath("~/afs/public_html/pfcorrs/v02/");
0011 
0012  TFile* fcorrs = new TFile("~/nobackup/arch/hcalCorrsFile_resp01.root","OPEN");
0013  TProfile* respcorrs = (TProfile*)fcorrs->Get("corrs1");
0014 
0015  TFile* tf = new TFile("~/nobackup/arch/hcalCorrPFv13_35cm.root","OPEN");
0016  //TFile* tf = new TFile("~/afs/arch/hcalCorrPFv2_26.2cm.root","OPEN");
0017 //TFile* tf = new TFile("./HcalCorrPF.root","OPEN");
0018  ofstream new_pfcorrs("newPFcorrs.txt");
0019  ofstream new_pfcorrs_subnoise("newPFcorrsNoiseSubtracted.txt");
0020 
0021 // ofstream new_pfcorrs("newPFcorrs26.2cm.txt");
0022 // ofstream new_pfcorrs_subnoise("newPFcorrsNoiseSubtracted26.2cm.txt");
0023   
0024 
0025 /*
0026  TString imgpath("~/afs/public_html/validation/pfcorrs/30cm/");
0027  TFile* tf = new TFile("./HcalCorrPFcone30cm.root","OPEN");    
0028  ofstream new_pfcorrs("newPFcorrs30cm.txt");
0029  ofstream new_pfcorrs_subnoise("newPFcorrsNoiseSubtracted30cm.txt");
0030 */
0031 /*  
0032  TString imgpath("~/afs/public_html/validation/pfcorrs/40cm/");
0033  TFile* tf = new TFile("./HcalCorrPFcone40cm.root","OPEN");
0034  ofstream new_pfcorrs("newPFcorrs40cm.txt");
0035  ofstream new_pfcorrs_subnoise("newPFcorrsNoiseSubtracted.txt");
0036   */
0037 
0038 /*
0039  TString imgpath("~/afs/public_html/validation/pfcorrs/50cm/");
0040  TFile* tf = new TFile("./HcalCorrPFcone50cm.root","OPEN");
0041  ofstream new_pfcorrs("newPFcorrs50cm.txt");
0042  ofstream new_pfcorrs_subnoise("newPFcorrsNoiseSubtracted50cm.txt");
0043 */
0044 
0045 /*
0046  TProfile* p1 = (TProfile*)tf->Get("hcalRecoAnalyzer/enHcal");  
0047  TProfile* p2 = (TProfile*)tf->Get("hcalRecoAnalyzer/enHcalNoise");  
0048  TProfile* p3 = (TProfile*)tf->Get("hcalRecoAnalyzer/nCells");
0049  TProfile* p4 = (TProfile*)tf->Get("hcalRecoAnalyzer/nCellsNoise");
0050  TProfile* p5 = (TProfile*)tf->Get("hcalRecoAnalyzer/enEcal");
0051 */
0052 
0053  TCut tr_cut = "eParticle>45 && eParticle<55";  label = new TText(0.03,0.2, "50 GeV");
0054   
0055  label -> SetNDC(); label->SetTextAlign(11); label->SetTextSize(0.05); label->SetTextAngle(90); label->SetTextColor(kRed);
0056 
0057   TCut trkQual = "";
0058   //  TCut trkQual = "(trkQual[0]==1 && abs(etaParticle)<2.4) || abs(etaParticle)>2.4";
0059   TCut mip_cut = "eECAL09cm<1.";
0060     TCut hit_dist = "delR<16";
0061   TCut mc_dist = "";
0062   // TCut mc_dist = "(abs(etaParticle)<2.4 && delRmc[0]<2) || abs(etaParticle)>2.4";
0063   //TCut low_resp_cut = "";
0064   TCut low_resp_cut = "eHcalCone/eParticle>0.2";
0065   TCut maxPNearBy = "";
0066   //  TCut neutral_iso_cut ="";
0067   TCut neutral_iso_cut = "(abs(iEta)<=13&&(eECAL40cm-eECAL09cm)<6.8)||(abs(iEta)==14&&(eECAL40cm-eECAL09cm)<6.6)||(abs(iEta)==15&&(eECAL40cm-eECAL09cm)<6.6)||(abs(iEta)==16&&(eECAL40cm-eECAL09cm)<9.8)||(abs(iEta)==17&&(eECAL40cm-eECAL09cm)<10.4)||(abs(iEta)==18&&(eECAL40cm-eECAL09cm)<9.8)||(abs(iEta)==19&&(eECAL40cm-eECAL09cm)<11.0)||(abs(iEta)==20&&(eECAL40cm-eECAL09cm)<12.3)||(abs(iEta)==21&&(eECAL40cm-eECAL09cm)<13.6)||(abs(iEta)==22&&(eECAL40cm-eECAL09cm)<15.2)||(abs(iEta)==23&&(eECAL40cm-eECAL09cm)<15.4)||(abs(iEta)==24&&(eECAL40cm-eECAL09cm)<16.3)||(abs(iEta)==25&&(eECAL40cm-eECAL09cm)<16.1)||(abs(iEta)==26&&(eECAL40cm-eECAL09cm)<15.4)||(abs(iEta)==27&&(eECAL40cm-eECAL09cm)<15.4)||abs(iEta)>27";
0068   
0069   TCut selection = trkQual && neutral_iso_cut && tr_cut && mip_cut && hit_dist && mc_dist && maxPNearBy && low_resp_cut;
0070   
0071 
0072   TTree* pftree = (TTree*)tf->Get("hcalPFcorrs/pfTree");
0073 //  pftree -> Draw("eHcalCone/eParticle:iEta>>p1(82, -40.5, 41.5)", "", "prof goff");
0074   pftree -> Draw("eHcalCone/eParticle:iEta>>p1(82, -40.5, 41.5)", selection, "prof goff");
0075   pftree -> Draw("(eHcalCone-eHcalConeNoise)/eParticle:iEta>>p2(82, -40.5, 41.5)", selection, "prof goff");
0076   pftree -> Draw("UsedCells:iEta>>p3(82, -40.5, 41.5)", selection, "prof goff");
0077   pftree -> Draw("UsedCellsNoise:iEta>>p4(82, -40.5, 41.5)", selection, "prof goff");
0078   pftree -> Draw("eECAL:iEta>>p5(82, -40.5, 41.5)", selection, "prof goff");
0079 
0080 TProfile* diff1 = new TProfile("diff1", "(oldvalue - new)/oldvalue", 84,-42,42);
0081 TProfile* diff2 = new TProfile("diff2", "(pfcorr - pfSubNoise)/pfcorr", 84,-42,42);
0082 
0083 
0084 TProfile* corrs2 = new TProfile("corrs2", "PF corrs", 84,-42,42);
0085 TProfile* corrs3 = new TProfile("corrs3", "PF corrs", 84,-42,42);
0086 TProfile* corrs4 = new TProfile("corrs4", "PF corrs", 84,-42,42);
0087 
0088 //ifstream old_pfcorrs("../data/HcalPFCorrs_v1.03_mc.txt");
0089 //ifstream old_pfcorrs("../data/HcalPFCorrs_v2.00_mc.txt");
0090 ifstream old_pfcorrs("../data/response_corrections.txt");
0091 
0092 
0093   //FILE *new_pfcorrs = fopen("new_pfcorrs.txt", "w+"); 
0094   Int_t   iEta;
0095   UInt_t  iPhi;
0096   Int_t  depth;
0097   //TString sdName;
0098   string sdName;
0099   UInt_t  detId;
0100   Float_t value;
0101 //  HcalSubdetector sd;
0102   map<Int_t, Float_t> CorrValues;
0103   map<Int_t, Float_t> CorrValuesSubNoise;
0104   
0105 for (Int_t i=0; i<=82; i++)
0106 {
0107 
0108   if (p1->GetBinContent(i)!=0) CorrValues[p1->GetBinCenter(i)]= 1./p1->GetBinContent(i);
0109   if (p2->GetBinContent(i)!=0) CorrValuesSubNoise[p2->GetBinCenter(i)]= 1./p2->GetBinContent(i);
0110 
0111   cout<<"bin: "<< p1->GetBinCenter(i)<<"   response: "<< p1->GetBinContent(i)<<"   corrs: "<< CorrValues[p1->GetBinCenter(i)]<<"   subnoise corrs: "<< CorrValuesSubNoise[p1->GetBinCenter(i)]<<endl;  
0112 
0113 
0114 //if(depth>0) {
0115 
0116 //diff1-> Fill(iEta, (value-CorrValues[iEta])/value);
0117 //   diff2-> Fill(i, (CorrValues[i] - CorrValuesSubNoise[i])/CorrValues[i]);
0118 
0119    //corrs1-> Fill(iEta, value);
0120  }  
0121  // }
0122 
0123 for (Int_t i=-41; i<=41; i++)
0124 {
0125   corrs2-> Fill(i, CorrValues[i]);
0126   corrs3-> Fill(i, CorrValuesSubNoise[i]);
0127 
0128  }
0129 
0130 
0131 
0132 TCanvas* c1 = new TCanvas("c1","all",0,0,350,350);
0133   c1-> cd();
0134   p1 -> Draw("");
0135   // p1 -> SetMaximum(60.);
0136   p1 -> SetXTitle("iEta");
0137   c1->SaveAs(imgpath+"a01.png");  
0138 
0139 
0140   p2 -> Draw("");
0141   p2 -> SetXTitle("iEta");
0142    c1->SaveAs(imgpath+"a02.png");  
0143 
0144 
0145   p3 -> Draw("");
0146   p3 -> SetXTitle("iEta");
0147    c1->SaveAs(imgpath+"a03.png");  
0148 
0149   p4 -> Draw("");
0150   p4 -> SetXTitle("iEta");
0151    c1->SaveAs(imgpath+"a04.png");  
0152 
0153 
0154   p5 -> Draw("");
0155   p5 -> SetXTitle("iEta");
0156   c1->SaveAs(imgpath+"a05.png");  
0157   
0158   
0159   respcorrs -> Draw("");  
0160   respcorrs -> SetMaximum(1.6);
0161   respcorrs -> SetMinimum(0.6);
0162   corrs2 -> Draw("same");
0163   corrs2 -> SetLineColor(kBlue);  
0164   c1->SaveAs(imgpath+"c01.png");  
0165 
0166   respcorrs -> Draw("");  
0167   corrs3 -> Draw("same");
0168   respcorrs -> SetMaximum(1.6);
0169   respcorrs -> SetMinimum(0.6);
0170   //corrs2 -> SetLineColor(kBlue);  
0171   corrs3 -> SetLineColor(kRed+1);  
0172   c1->SaveAs(imgpath+"c02.png");  
0173 
0174   //diff1 -> Draw("");  
0175   //c1->SaveAs(imgpath+"c03.png");
0176   
0177   //diff2 -> Draw("");  
0178   //c1->SaveAs(imgpath+"c04.png");  
0179   
0180 
0181   tf->Close();
0182 
0183 
0184 
0185   std::string line;
0186   
0187   while (getline(old_pfcorrs, line)) 
0188     {
0189       if(!line.size() || line[0]=='#') 
0190     {
0191       //fprintf(new_pfcorrs, "%1s%16s%16s%16s%16s%9s%11s\n", "#", "eta", "phi", "depth", "det", "value", "DetId");   
0192       new_pfcorrs<<"#             eta             phi           depth        det           value       DetId"<<endl;
0193       continue;
0194     }
0195       
0196       std::istringstream linestream(line);
0197       linestream >> iEta >> iPhi >> depth >> sdName >> value >> hex >> detId;
0198       // corrs1-> Fill(abs(iEta), value);
0199       
0200       //cout<<" I'm HERE"<<iEta<<iPhi<<endl;      
0201 
0202 
0203       //format:
0204       //               -1               1               1              HB    0.97635     42004081
0205       
0206       //  fprintf(new_pfcorrs, "%17i%16i%16i%16s%9.5f%11X\n", iEta, iPhi, depth, sdName, value, detId);     
0207       
0208 
0209       new_pfcorrs.width(17);
0210       new_pfcorrs<<dec<<iEta;
0211       new_pfcorrs.width(16);
0212       new_pfcorrs<<iPhi;
0213       new_pfcorrs.width(16);
0214       new_pfcorrs<<depth;
0215       new_pfcorrs.width(11);
0216       new_pfcorrs<<sdName;
0217       new_pfcorrs.width(16);
0218       new_pfcorrs.setf(ios::fixed, ios::floatfield);
0219       new_pfcorrs.precision(5); 
0220       if (sdName=="HO" || depth==-99) {new_pfcorrs<<value;} else new_pfcorrs<<CorrValues[iEta];
0221       new_pfcorrs.width(13);
0222       new_pfcorrs.setf(ios::uppercase);
0223       new_pfcorrs<<hex<<detId<<endl;
0224 
0225       
0226       new_pfcorrs_subnoise.width(17);
0227       new_pfcorrs_subnoise<<dec<<iEta;
0228       new_pfcorrs_subnoise.width(16);
0229       new_pfcorrs_subnoise<<iPhi;
0230       new_pfcorrs_subnoise.width(16);
0231       new_pfcorrs_subnoise<<depth;
0232       new_pfcorrs_subnoise.width(11);
0233       new_pfcorrs_subnoise<<sdName;
0234       new_pfcorrs_subnoise.width(16);
0235       new_pfcorrs_subnoise.setf(ios::fixed, ios::floatfield);
0236       new_pfcorrs_subnoise.precision(5); 
0237       if (sdName=="HO" || depth==-99) {new_pfcorrs_subnoise<<value;} else new_pfcorrs_subnoise<<CorrValues[iEta];
0238       new_pfcorrs_subnoise.width(13);
0239       new_pfcorrs_subnoise.setf(ios::uppercase);
0240       new_pfcorrs_subnoise<<hex<<detId<<endl;
0241 
0242     }
0243 }
0244