Functions | Variables
attenuation_calibration_comparisons.C File Reference
#include "HistogramAttr.h"
#include "Header.h"

Go to the source code of this file.

Functions

 attenuation_calibration_comparisons ()
 
 uncorrected_PE (char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString hist1, TString hist2)
 
void uncorrected_PE_tricellvstrajectory (char *c1, TString det, TString mcdata, TString analysis1, TString analysis2, TString hist1, TString hist2)
 
void thresholdshadowing (char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString var)
 
void corrected_PE (char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString plane, TString cell)
 
void calibrationconstants (char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2)
 

Variables

bool uncorrectedPEplots = false
 
bool thresholdshadowingplots = false
 
bool correctedPEplots = true
 
bool calibrationconstantsplots = false
 
bool nd = true
 
bool fd = false
 
bool mc = false
 
bool data = true
 
bool tricellvstrajectory = false
 
bool onetoone = true
 
bool manytoone = false
 
bool nddatamuoncatcher = true
 

Function Documentation

attenuation_calibration_comparisons ( )

Definition at line 40 of file attenuation_calibration_comparisons.C.

References demo5::c1, demo5::c2, chisquared::c3, calibrationconstants(), calibrationconstantsplots, getBrightness::cell, corrected_PE(), correctedPEplots, fillBadChanDBTables::det, fd, hist1, hist2, mc, nd, nddatamuoncatcher, onetoone, NDAPDHVSetting::plane, thresholdshadowing(), thresholdshadowingplots, tricellvstrajectory, uncorrected_PE(), uncorrected_PE_tricellvstrajectory(), uncorrectedPEplots, and PandAna.Demos.tute_pid_validation::var.

40  {
41  gStyle->SetOptStat(" ");
42  //-----------------------------------------------------------------------------------------------------
43  //set of all .root files for uncorrected PE/cm vs. w plots
44  char* c1 = "/nova/ana/users/prabhjot/gridjob/Work_with_Adam_and_Brian/nd/NDPCHitData_S150224/mainbody_and_muoncatcher/mainboby/ConsolidateViewsMode_false/Fromgrid_aftercosmiccalib/histofiles/calibNDPCHitData_hist_mainbody_ConsolidateViewsMode_false_S150224.root"; // base file
45  char* c2 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/Period2_DDActivity/mainbody/calib_hist_2ndAna_fulldetector_nddata_mainbody_period2ddactivity.root"; // 2nd file
46 
47  //.root file for uncorrected PE/cm for tricell vs trajectory study
48  char* c3 = "";
49 
50  //set of .root files for threshold and shadowing corrections
51  char* cthreshold1 = ""; // base file
52  char* cthreshold2 = ""; // 2nd file
53 
54  //set of .root files to plot attenuation fitting plots
55  char* cattenuation1 = "/nova/app/users/prabhjot/nova_test_S150224/AttenuationCalibration/nd/pchitdata/fromgrid/mainbody_and_muoncatcher/mainboby/ConsolidateViewsMode_false/Atten_Fit_Results.root"; // base file
56  char* cattenuation2 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/Period2_DDActivity/mainbody/FitAttenuation/Atten_Fit_Results.root"; // 2nd file
57 
58  //set of .root files to plot attenuation fitting plots for the ND data muon catcher
59  char* cattenuation_muoncatcher1 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/Period2_DDActivity/muoncatcher/FitAttenuation/Atten_Fit_Results.root"; // base file
60  char* cattenuation_muoncatcher2 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/Epoch3b_DDActivity/muoncatcher/FitAttenuation/Atten_Fit_Results.root"; // 2nd file
61 
62  //set of attenuation constants .txt files for attenuation constants comparisons
63  char* cconstants1 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/comparison_plots/calib_atten_consts_uniq_sort_1st_ana.txt"; // base file
64  char* cconstants2 = "/nova/ana/users/prabhjot/AttenuationCalibration_2ndAnalysis/nddata/comparison_plots/calib_atten_consts_uniq_sort_period2.txt"; // 2nd file
65  //---------------------------------------------------------------------------------------------------
66  //---------------------------------------------------------------------------------------------------
67 
68  TString hist1; // base histogram from the base file. You want to compare other histograms wrt base histogram.
69  TString hist2; // histogram from the second file that you want to compare with the base histogram.
70  TString det; // name of the detector.
71  TString mcdata; // labeling for the Monte Carlo or the data .
72  TString analysis1; // labeling for the base file.
73  TString analysis2; // labelling for the second file.
74  TString var; // variable labelling on the x axis. It is W (cm).
75  TString plane; // plane number.
76  TString cell; // cell number.
77 
78  analysis1 = "2nd_analysis_period2";
79  analysis2 = "2nd_analysis_epoch3b";
80 
81  hist1 = "h_WPE_corr_xy";
82  hist2 = "h_WPE_corr_xy"; //_corr_xy_ for the tricell, _corr_traj_ for the trajectory
83 
84  var = "W";
85  plane = "003";
86  cell = "023";
87 
88  //ND MC two cases comparison
89  if(nd && onetoone){
90  det = "nd";
91  if(mc) mcdata = "mc";
92  if(data) mcdata = "data";
93 
94  if(tricellvstrajectory && uncorrectedPEplots) uncorrected_PE_tricellvstrajectory(c3, det, mcdata, analysis1, analysis2, hist1, hist2);
95  else if(uncorrectedPEplots) uncorrected_PE(c1, c2, det, mcdata, analysis1, analysis2, hist1, hist2);
96 
97  if(mc && thresholdshadowingplots) thresholdshadowing(cthreshold1, cthreshold2, det, mcdata, analysis1, analysis2, var);
98 
99  if(data && nddatamuoncatcher && correctedPEplots)corrected_PE(cattenuation_muoncatcher1, cattenuation_muoncatcher2, det, mcdata, analysis1, analysis2, plane, cell);
100  else if(correctedPEplots) corrected_PE(cattenuation1, cattenuation2, det, mcdata, analysis1, analysis2, plane, cell);
101 
102  if(calibrationconstantsplots) calibrationconstants(cconstants1, cconstants2, det, mcdata, analysis1, analysis2);
103  }//end of if(nd && mc && onetoone)
104 
105  //FD two case comparison
106  if(fd && onetoone){
107  det = "fd";
108  if(mc) mcdata = "mc";
109  if(data) mcdata = "data";
110 
111  if(tricellvstrajectory && uncorrectedPEplots) uncorrected_PE_tricellvstrajectory(c3, det, mcdata, analysis1, analysis2, hist1, hist2);
112  else if(uncorrectedPEplots) uncorrected_PE(c1, c2, det, mcdata, analysis1, analysis2, hist1, hist2);
113 
114  if(mc && thresholdshadowingplots) thresholdshadowing(cthreshold1, cthreshold2, det, mcdata, analysis1, analysis2, var);
115 
116  if(correctedPEplots) corrected_PE(cattenuation1, cattenuation2, det, mcdata, analysis1, analysis2, plane, cell);
117 
118  if(calibrationconstantsplots) calibrationconstants(cconstants1, cconstants2, det, mcdata, analysis1, analysis2);
119  }//end of if(fd && mc && onetoone)
120 }//attenuation_calibration_comparisons()
TH1D * hist2
Definition: plotHisto.C:12
void corrected_PE(char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString plane, TString cell)
const XML_Char const XML_Char * data
Definition: expat.h:268
c2
Definition: demo5.py:33
uncorrected_PE(char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString hist1, TString hist2)
TH1D * hist1
Definition: plotHisto.C:9
void thresholdshadowing(char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString var)
void uncorrected_PE_tricellvstrajectory(char *c1, TString det, TString mcdata, TString analysis1, TString analysis2, TString hist1, TString hist2)
c1
Definition: demo5.py:24
void calibrationconstants(char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2)
void calibrationconstants ( char *  c1,
char *  c2,
TString  det,
TString  mcdata,
TString  analysis1,
TString  analysis2 
)

Definition at line 419 of file attenuation_calibration_comparisons.C.

References Munits::c_2, channels, om::cout, allTimeWatchdog::endl, fclose(), check_time_usage::float, submit_syst::fout, HistogramAttr1D(), HistogramAttr2D(), make_mec_shifts_plots::legend, and confusionMatrixTree::out.

Referenced by attenuation_calibration_comparisons().

419  {
420  ifstream in_1;
421  ifstream in_2;
422 
423  ofstream out;
424 
425  FILE* fout = fopen("uncalibrated_channels_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+".csv","w");
426  fprintf(fout, "%s%50s%s\n", analysis1.Data(), " ", analysis2.Data());
427  fprintf(fout, "channel,coeff_exp,atten_length,background,chisq, %9s, channel,coeff_exp,atten_length,background,chisq\n", " ");
428 
429  in_1.open(c1);
430  in_2.open(c2);
431 
432 
433  //x=coeff_exp, y=atten_length, z=background, a=edgx_low, b=edgx_high, c=coeff_low
434  Float_t x_1, y_1, z_1, a_1, b_1, c_1, d_1, e_1, g_1;
435  Float_t x_2, y_2, z_2, a_2, b_2, c_2, d_2, e_2, g_2;
436  Float_t chisq = 0.2;
437 
438  Int_t channels;
439  Int_t totalchannels1 = 0;
440  Int_t totalchannels2 = 0;
441  Int_t calibratedchannels1 = 0;
442  Int_t calibratedchannels2 = 0;
443  Int_t f_1;
444  Int_t f_2;
445  Int_t bins1;
446 
447  Double_t limit1D = 3.;
448  Double_t limit2D = 1000.;
449 
450  if(det=="nd" && mcdata=="mc" ) { channels = 3063; }
451  if(det=="nd" && mcdata=="data" ) { channels = 1003063; }
452  if(det=="fd" && mcdata=="mc" ) { channels = 1383; }
453  if(det=="fd" && mcdata=="data" ) { channels = 895383; }
454 
455 
456  TH1F* hcoeffexp = new TH1F("hcoeffexp", "coeffexp difference distribution", 120, -limit1D, +limit1D );
457  TH1F* hattenlength = new TH1F("hattenlength", "attenlength difference distribution", 600, -limit1D, +limit1D );
458  TH1F* hbackground = new TH1F("hbackground", "background difference distribution", 120, -limit1D, +limit1D );
459  TH1F* hchisq = new TH1F("hchisq", "chisq difference distribution", 600, -limit1D, +limit1D );
460 
461  TH2F* h2coeffexp = new TH2F("h2coeffexp", " " , 2*(Int_t)limit2D, -limit2D, limit2D, 2*(Int_t)limit2D, -limit2D, limit2D );
462  TH2F* h2attenlength = new TH2F("h2attenlength"," " , 2*(Int_t)limit2D, -limit2D, limit2D, 2*(Int_t)limit2D, -limit2D, limit2D );
463  TH2F* h2background = new TH2F("h2background", " " , 2*(Int_t)limit2D, -limit2D, limit2D, 2*(Int_t)limit2D, -limit2D, limit2D );
464  TH2F* h2chisq = new TH2F("h2chisq"," ", 200,0.,2., 200,0.,2. );
465 
466  while(1){
467  in_1 >> x_1 >> y_1 >> z_1 >> a_1 >> b_1 >> c_1 >> d_1 >> e_1 >> f_1 >> g_1;
468  in_2 >> x_2 >> y_2 >> z_2 >> a_2 >> b_2 >> c_2 >> d_2 >> e_2 >> f_2 >> g_2;
469  if(!in_1.good()) break;
470  if(!in_2.good()) break;
471 
472  if(f_1<=channels && e_1 < chisq && e_2 < chisq){
473  hcoeffexp->Fill((x_1 - x_2)/x_1);
474  hattenlength->Fill((y_1 - y_2)/y_1);
475  hbackground->Fill((z_1 - z_2)/z_1);
476 
477  h2coeffexp->Fill(x_2, x_1);
478  h2attenlength->Fill(y_2, y_1);
479  h2background->Fill(z_2, z_1);
480  }//end of loop over channels and chisq conditions.
481  if(f_1<=channels){
482  totalchannels1++;
483  totalchannels2++;
484  if(e_1<chisq) calibratedchannels1++;
485  if(e_2<chisq) calibratedchannels2++;
486 
487  hchisq->Fill((e_1 - e_2)/e_1);
488  h2chisq->Fill(e_2, e_1);
489 
490  if(e_1>=chisq || e_2>=chisq){
491  fprintf(fout, "%i,%g,%g,%g,%g,%20s,%i,%g,%g,%g,%g\n",f_1,x_1,y_1,z_1,e_1, " ", f_2,x_2,y_2,z_2,e_2);
492  out << f_1 << setw(22) << x_1 << setw(22) << y_1 << setw(22) << z_1 << setw(22) << e_1 << f_2 << setw(22) << x_2 << setw(22) << y_2 << setw(22) << z_2 << setw(22) << e_2 << endl;
493  }
494  }//end of only channels condition for chisq plots. For chisq plots we have to consider all channels.
495  }//end of while condition.
496 
497  cout << " " << endl;
498  cout << "Total channels in "+analysis1+" = " << totalchannels1 << " and calibrated channels = " << calibratedchannels1 << " = " << ((float)calibratedchannels1/(float)totalchannels1)*100. << " %." << endl;
499  cout << "Total channels in "+analysis2+" = " << totalchannels2 << " and calibrated channels = " << calibratedchannels2 << " = " << ((float)calibratedchannels2/(float)totalchannels2)*100. << " %." << endl;
500  cout << " " << endl;
501 
502  //coeffexp plots
503  new TCanvas;
504  gPad->SetGrid();
505  HistogramAttr1D(hcoeffexp, "#Deltacoeffexp/coeffexp", "# of channels", 0.5*hcoeffexp->GetXaxis()->GetBinCenter(1), 0.5*hcoeffexp->GetXaxis()->GetBinCenter(hcoeffexp->GetXaxis()->GetNbins()), 0., 7e5, kBlack);
506  hcoeffexp->Draw();
507  TLegend *legend = new TLegend(0.00, 0.60, 0.47, 0.93);
508  legend->SetBorderSize(0);
509  legend->SetFillStyle(0);
510  legend->SetTextSize(0.05);
511  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
512  legend->AddEntry(hcoeffexp,"#splitline{"+analysis1+"}{vs "+analysis2+"}","P");
513  legend->Draw();
514  gPad->Print("Images/Delta_coeffexp_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
515 
516  new TCanvas;
517  gPad->SetGrid();
518  gPad->SetLogz();
519  HistogramAttr2D(h2coeffexp, "coeffexp_{"+analysis2+"}", "coeffexp_{"+analysis1+"}", " ", 0, 100, 0, 100);
520  h2coeffexp->Draw("colz");
521  gPad->Print("Images/coeffexp2D_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
522 
523  //attenlength plots
524  new TCanvas;
525  gPad->SetGrid();
526  HistogramAttr1D(hattenlength, "#Deltaattenlength/attenlength", "# of channels", 0.5*hattenlength->GetXaxis()->GetBinCenter(1), 0.5*hattenlength->GetXaxis()->GetBinCenter(hattenlength->GetXaxis()->GetNbins()), 0., 166500, kBlack);
527  hattenlength->Draw();
528  TLegend *legend = new TLegend(0.00, 0.60, 0.47, 0.93);
529  legend->SetBorderSize(0);
530  legend->SetFillStyle(0);
531  legend->SetTextSize(0.05);
532  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
533  legend->AddEntry(hattenlength,"#splitline{"+analysis1+"}{vs "+analysis2+"}","P");
534  legend->Draw();
535  gPad->Print("Images/Delta_attenlength_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
536 
537  new TCanvas;
538  gPad->SetGrid();
539  gPad->SetLogz();
540  HistogramAttr2D(h2attenlength, "attenlength_{"+analysis2+"} (cm)", "attenlength_{"+analysis1+"} (cm)", " ", 0, 1000, 0, 1000);
541  h2attenlength->Draw("colz");
542  gPad->Print("Images/attenlength2D_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
543 
544  //backgrounds plots
545  new TCanvas;
546  gPad->SetGrid();
547  HistogramAttr1D(hbackground, "#Deltabackground/background", "# of channels", 0.5*hbackground->GetXaxis()->GetBinCenter(1), 0.5*hbackground->GetXaxis()->GetBinCenter(hbackground->GetXaxis()->GetNbins()), 0., 3e4, kBlack);
548  hbackground->Draw();
549  TLegend *legend = new TLegend(0.00, 0.60, 0.47, 0.93);
550  legend->SetBorderSize(0);
551  legend->SetFillStyle(0);
552  legend->SetTextSize(0.05);
553  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
554  legend->AddEntry(hbackground,"#splitline{"+analysis1+"}{vs "+analysis2+"}","P");
555  legend->Draw();
556  gPad->Print("Images/Delta_background_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
557 
558  new TCanvas;
559  gPad->SetGrid();
560  gPad->SetLogz();
561  HistogramAttr2D(h2background, "background_{"+analysis2+"}", "background_{"+analysis1+"}", " ", -50, 50, -50, 50);
562  h2background->Draw("colz");
563  gPad->Print("Images/background2D_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
564 
565  //chisq plots
566  new TCanvas;
567  gPad->SetGrid();
568  HistogramAttr1D(hchisq, "#Deltachisq/chisq", "# of channels", 0.5*hchisq->GetXaxis()->GetBinCenter(1), 0.5*hchisq->GetXaxis()->GetBinCenter(hchisq->GetXaxis()->GetNbins()), 0., 4200., kBlack);
569  hchisq->Draw();
570  TLegend *legend = new TLegend(0.00, 0.60, 0.47, 0.93);
571  legend->SetBorderSize(0);
572  legend->SetFillStyle(0);
573  legend->SetTextSize(0.05);
574  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
575  legend->AddEntry(hchisq,"#splitline{"+analysis1+"}{vs "+analysis2+"}","P");
576  legend->Draw();
577  gPad->Print("Images/Delta_chisq_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
578 
579  new TCanvas;
580  gPad->SetGrid();
581  gPad->SetLogz();
582  HistogramAttr2D(h2chisq, "#chi^{2}_{"+analysis2+"}", "#chi^{2}_{"+analysis1+"}", " ", 0, 0.5, 0, 0.5);
583  h2chisq->Draw("colz");
584  gPad->Print("Images/chisq2D_"+analysis1+"_vs_"+analysis2+"_"+det+"_"+mcdata+".pdf");
585 
586  in_1.close();
587  in_2.close();
588  out.close();
589  fclose(fout);
590 }//end of calibrationconstants()
static constexpr Double_t c_2
Definition: Munits.h:245
c2
Definition: demo5.py:33
fclose(fg1)
std::map< ToFCounter, std::vector< unsigned int > > channels
OStream cout
Definition: OStream.cxx:6
ifstream in
Definition: comparison.C:7
void HistogramAttr2D(TH2 *h, char *titlx_x, char *titlx_y, char *titlx_z, Double_t binlowx, Double_t binhighx, Double_t binlowy, Double_t binhighy)
c1
Definition: demo5.py:24
HistogramAttr1D(TH1 *h, char *title_x, char *title_y, Double_t binlowx, Double_t binhighx, Double_t binlowy, Double_t binhighy, Color_t lcolor)
void corrected_PE ( char *  c1,
char *  c2,
TString  det,
TString  mcdata,
TString  analysis1,
TString  analysis2,
TString  plane,
TString  cell 
)

Definition at line 248 of file attenuation_calibration_comparisons.C.

References om::cout, allTimeWatchdog::endl, f1, f2, h1, h2, HistogramAttr1D(), make_mec_shifts_plots::legend, nddatamuoncatcher, Ratiofit(), and Ratiohist().

Referenced by attenuation_calibration_comparisons().

248  {
249 
250  if(nddatamuoncatcher){
251  if(plane.Atoi() > 3 || plane.Atoi() < 2){
252  std::cout << "It is a ND muon catcher data file and has only 002, 003 planes. Plane " << plane << " is out of range. " << std::endl;
253  abort();
254  }
255  if(cell.Atoi() > 95 || cell.Atoi() < 0){
256  std::cout << "It is a Data file. Cell " << cell << " is out of range. " << std::endl;
257  abort();
258  }
259  }
260  if(det.Data()=="nd"){
261  if(mcdata.Data()=="mc"){
262  if(plane.Atoi() > 3 || plane.Atoi() < 0){
263  std::cout << "It is a Monte Carlo file. Plane " << plane << " out of range. " << std::endl;
264  abort();
265  }
266  if(cell.Atoi() > 95 || cell.Atoi() < 0){
267  std::cout << "It is a Monte Carlo file. Cell " << cell << " out of range. " << std::endl;
268  abort();
269  }
270  }//end of MC condition
271  if(mcdata.Data()=="data"){
272  if(plane.Atoi() > 191 || plane.Atoi() < 0){
273  std::cout << "It is a ND data file and has only 191 main body planes. Plane " << plane << " out of range. " << std::endl;
274  abort();
275  }
276  if(cell.Atoi() > 95 || cell.Atoi() < 0){
277  std::cout << "It is a Data file. Cell " << cell << " out of range. " << std::endl;
278  abort();
279  }
280  }//end of data condition
281  }//end of the ND condition
282 
283  if(det.Data()=="fd"){
284  if(mcdata.Data()=="mc"){
285  if(plane.Atoi() > 1 || plane.Atoi() < 0){
286  std::cout << "It is a Monte Carlo file. Plane " << plane << " out of range. " << std::endl;
287  abort();
288  }
289  if(cell.Atoi() > 383 || cell.Atoi() < 0){
290  std::cout << "It is a Monte Carlo file. Cell " << cell << " out of range. " << std::endl;
291  abort();
292  }
293  }//end of MC condition
294  else if(mcdata.Data()=="data"){
295  if(plane.Atoi() > 895 || plane.Atoi() < 0){
296  std::cout << "It is a Data file. Plane " << plane << " out of range. " << std::endl;
297  abort();
298  }
299  if(cell.Atoi() > 383 || cell.Atoi() < 0){
300  std::cout << "It is a Data file. Cell " << cell << " out of range. " << std::endl;
301  abort();
302  }
303  }//end of data condition
304  }//end of FD condition
305 
306  TFile* f1 = new TFile(c1, "read");
307  TFile* f2 = new TFile(c2, "read");
308 
309  TF1* fit1;
310  TF1* fit2;
311 
312  f1->cd();
313  TH1F* h1 = (TH1F*)f1->Get("plane_"+plane+"/Atten_Fit_"+plane+"_"+cell);
314  Color_t color1 = kBlack;
315  if(mcdata=="mc"){
316  if(plane.Atoi()==0) fit1 = (TF1*)h1->FindObject("fitX");
317  else if(plane.Atoi()==1) fit1 = (TF1*)h1->FindObject("fitY");
318  else if(plane.Atoi()==2) fit1 = (TF1*)h1->FindObject("fitMuX");
319  else if(plane.Atoi()==3) fit1 = (TF1*)h1->FindObject("fitMuY");
320  if(fit1) fit1->SetLineColor(color1);
321  }//end of MC condition
322  else if(mcdata=="data"){
323  if(nddatamuoncatcher){
324  if(plane.Atoi()==2) fit1 = (TF1*)h1->FindObject("fitMuX");
325  else if(plane.Atoi()==3) fit1 = (TF1*)h1->FindObject("fitMuY");
326  if(fit1) fit1->SetLineColor(color1);
327  }//end of muon catcher condition
328  else if(plane.Atoi()%2==0) fit1 = (TF1*)h1->FindObject("fitY");
329  else if(plane.Atoi()%2!=0) fit1 = (TF1*)h1->FindObject("fitX");
330  if(fit1) fit1->SetLineColor(color1);
331  }//end of data condition
332 
333  f2->cd();
334  TH1F* h2 = (TH1F*)f2->Get("plane_"+plane+"/Atten_Fit_"+plane+"_"+cell);
335  Color_t color2 = kRed;
336  if(mcdata=="mc"){
337  if(plane.Atoi()==0) fit2 = (TF1*)h2->FindObject("fitX");
338  else if(plane.Atoi()==1) fit2 = (TF1*)h2->FindObject("fitY");
339  else if(plane.Atoi()==2) fit2 = (TF1*)h2->FindObject("fitMuX");
340  else if(plane.Atoi()==3) fit2 = (TF1*)h2->FindObject("fitMuY");
341  if(fit2) fit2->SetLineColor(color2);
342  }//end of MC condition
343  else if(mcdata=="data"){
344  if(nddatamuoncatcher){
345  if(plane.Atoi()==2) fit2 = (TF1*)h2->FindObject("fitMuX");
346  else if(plane.Atoi()==3) fit2 = (TF1*)h2->FindObject("fitMuY");
347  if(fit2) fit2->SetLineColor(color2);
348  }//end of muon catcher condition
349  else if(plane.Atoi()%2==0) fit2 = (TF1*)h2->FindObject("fitY");
350  else if(plane.Atoi()%2!=0) fit2 = (TF1*)h2->FindObject("fitX");
351  if(fit2) fit2->SetLineColor(color2);
352  }//end of data condition
353 
354  new TCanvas;
355  gPad->SetGrid();
356  h1->Draw();
357  HistogramAttr1D(h1, "W (cm)", "Mean PE/cm", h1->GetXaxis()->GetBinCenter(1), h1->GetXaxis()->GetBinCenter(h1->GetXaxis()->GetNbins()), 0., 150., color1);
358  h1->SetLineColor(color1);
359  h1->SetMarkerColor(color1);
360  h1->SetMarkerStyle(24);
361  h2->Draw("same");
362  h2->SetLineColor(color2);
363  h2->SetMarkerColor(color2);
364  h2->SetMarkerStyle(24);
365 
366  TLegend *legend = new TLegend(.12,.62,.64,.90);
367  legend->SetBorderSize(0);
368  legend->SetFillStyle(0);
369  legend->SetTextSize(0.04);
370  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
371  legend->AddEntry(h1, analysis1);
372  legend->AddEntry(h2, analysis2);
373  legend->Draw();
374  gPad->Print("Images/Atten_Fit_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_"+plane+"_"+cell+".pdf");
375 
376  //ratio of data points
377  new TCanvas;
378  gPad->SetGrid();
379  const TString hrationame = Form("hratio%s_%s_%s_%s", det.Data(), mcdata.Data(), plane.Data(), cell.Data());
380  TH1D* hratio = new TH1D(hrationame, "",h1->GetXaxis()->GetNbins(), h1->GetBinCenter(0), h1->GetBinCenter(h1->GetXaxis()->GetNbins()));
381  Ratiohist(h1, h2, hratio, kRed, 24, "W (cm)", "Ratio of Mean PE/cm data points", 0.3, 1.0);
382 
383  TLegend *legend = new TLegend(.12,.62,.64,.90);
384  legend->SetBorderSize(0);
385  legend->SetFillStyle(0);
386  legend->SetTextSize(0.04);
387  legend->AddEntry((TObject*)0, det+" "+mcdata+" : Plane "+plane+" Cell "+cell," ");
388  legend->AddEntry(hratio,analysis2+" to "+analysis1);
389  legend->Draw();
390  gPad->Print("Images/Atten_Fit_ratio_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_plane"+plane+"_cell"+cell+".pdf");
391 
392  //ratio of fitting curve
393  new TCanvas;
394  gPad->SetGrid();
395  const TString hratiofitname = Form("hratiofit%s_%s_%s_%s", det.Data(), mcdata.Data(), plane.Data(), cell.Data());
396  if(det=="nd"){
397  float x_limit = 250.;
398  int x_bins = 2*x_limit;
399  }
400  else if(det=="fd"){
401  float x_limit = 800.;
402  int x_bins = 2*x_limit;
403  }
404  TH1D* hfit = new TH1D("hfit", "", x_bins, -x_limit, +x_limit);
405  HistogramAttr1D(hfit, "W (cm)", "Ratio of Mean PE/cm fitting curves", -x_limit, +x_limit, 0.55, 0.68, kWhite);
406  hfit->Draw();
407  TH1D* hratiofit = new TH1D(hratiofitname, "", x_bins, -x_limit, +x_limit);
408  Ratiofit(fit1, fit2, hratiofit, plane.Data(), cell.Data(), color2, 24, "W (cm)", "Ratio of Mean PE/cm fitting curves", 0.55, 0.68);
409  TLegend *legend = new TLegend(.12,.62,.64,.90);
410  legend->SetBorderSize(0);
411  legend->SetFillStyle(0);
412  legend->SetTextSize(0.04);
413  legend->AddEntry((TObject*)0, det+" "+mcdata+" : Plane "+plane+" Cell "+cell," ");
414  legend->AddEntry(hratiofit,analysis2+" to "+analysis1);
415  legend->Draw();
416  gPad->Print("Images/Atten_Fit_fitratio_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_plane"+plane+"_cell"+cell+".pdf");
417 }//end of corrected_PE()
Float_t f2
c2
Definition: demo5.py:33
void Ratiohist(TH1 *h1, TH1 *h2, TH1 *hratio, Color_t color, Style_t mstyle, TString labelx, TString labely, Double_t lowy, Double_t highy)
Float_t f1
OStream cout
Definition: OStream.cxx:6
TH1F * h2
Definition: plot.C:45
TH1F * h1
void Ratiofit(TF1 *f1, TF1 *f2, TH1 *hratio, TString plane, TString cell, Color_t color, Style_t mstyle, TString labelx, TString labely, Double_t lowy, Double_t highy)
c1
Definition: demo5.py:24
HistogramAttr1D(TH1 *h, char *title_x, char *title_y, Double_t binlowx, Double_t binhighx, Double_t binlowy, Double_t binhighy, Color_t lcolor)
void thresholdshadowing ( char *  c1,
char *  c2,
TString  det,
TString  mcdata,
TString  analysis1,
TString  analysis2,
TString  var 
)

Definition at line 191 of file attenuation_calibration_comparisons.C.

References f1, f2, HistogramAttr1D(), make_mec_shifts_plots::legend, Ratiohist(), and POTSpillRate::view.

Referenced by attenuation_calibration_comparisons().

191  {
192  TFile* f1 = new TFile(c1, "read");
193  TFile* f2 = new TFile(c2, "read");
194  for(int view = 0; view < 2; ++view){
195  const TString viewstr = view ? "y" : "x";
196  const TString viewStr = view ? "Y" : "X";
197 
198  f1->cd();
199  TH1D* correction_1 = (TH1D*)f1->Get("a"+viewstr+"view");
200  correction_1->SetLineColor(kBlack);
201 
202  f2->cd();
203  TH1D* correction_2 = (TH1D*)f2->Get("a"+viewstr+"view");
204  correction_2->SetLineColor(kRed);
205 
206  new TCanvas;
207  gPad->SetGrid();
208  if(correction_1->FindObject("stats")){
209  TPaveStats* pave1 = (TPaveStats*)correction_1->FindObject("stats");
210  pave1->Delete();
211  }
212  if(correction_2->FindObject("stats")){
213  TPaveStats* pave2 = (TPaveStats*)correction_2->FindObject("stats");
214  pave2->Delete();
215  }
216  correction_1->Draw();
217  correction_2->Draw("same");
218  if(var=="W") TString labelx = "W (cm)";
219  HistogramAttr1D(correction_1, labelx.Data(), "correction factor", correction_1->GetXaxis()->GetBinCenter(1), correction_1->GetXaxis()->GetBinCenter(correction_1->GetXaxis()->GetNbins()), 0.5, 2.0, kBlack);
220 
221  TLegend *legend = new TLegend(.12,.62,.64,.90);
222  legend->SetBorderSize(0);
223  legend->SetFillStyle(0);
224  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
225  legend->AddEntry(correction_1, analysis1+" "+viewStr+" view");
226  legend->AddEntry(correction_2, analysis2+" "+viewStr+" view");
227  legend->Draw();
228  gPad->Print("Images/thresholdshadowing_vs_"+var+"_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_"+viewstr+"view.pdf");
229 
230  //ratios
231  new TCanvas;
232  gPad->SetGrid();
233  const TString hrationame = Form("hratio%s_%s_%s_%s", det.Data(), mcdata.Data(), viewstr.Data(), var.Data());
234  TH1D* hratio = new TH1D(hrationame, "",correction_1->GetXaxis()->GetNbins(), correction_1->GetBinCenter(0), correction_1->GetBinCenter(correction_1->GetXaxis()->GetNbins()));
235  Ratiohist(correction_1, correction_2, hratio, kRed, 24, labelx.Data(), "Correction factor ratio", 1.8, 1.2);
236 
237  TLegend *legend = new TLegend(.12,.62,.64,.90);
238  legend->SetBorderSize(0);
239  legend->SetFillStyle(0);
240  legend->SetTextSize(0.04);
241  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
242  legend->AddEntry(hratio,analysis2+" to "+analysis1+" "+viewStr+" view");
243  legend->Draw();
244  gPad->Print("Images/thresholdshadowinge_vs_"+var+"_ratio_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_"+viewstr+"view.pdf");
245  }//end of loop over views
246 }//end of thresholdshadowing()
Float_t f2
c2
Definition: demo5.py:33
void Ratiohist(TH1 *h1, TH1 *h2, TH1 *hratio, Color_t color, Style_t mstyle, TString labelx, TString labely, Double_t lowy, Double_t highy)
Float_t f1
c1
Definition: demo5.py:24
HistogramAttr1D(TH1 *h, char *title_x, char *title_y, Double_t binlowx, Double_t binhighx, Double_t binlowy, Double_t binhighy, Color_t lcolor)
uncorrected_PE ( char *  c1,
char *  c2,
TString  det,
TString  mcdata,
TString  analysis1,
TString  analysis2,
TString  hist1,
TString  hist2 
)

Definition at line 124 of file attenuation_calibration_comparisons.C.

References f1, f2, HistogramAttr1D(), make_mec_shifts_plots::legend, nddatamuoncatcher, Ratiohist(), and POTSpillRate::view.

Referenced by attenuation_calibration_comparisons(), and uncorrected_PE_tricellvstrajectory().

124  {
125  TFile* f1 = new TFile(c1, "read");
126  TFile* f2 = new TFile(c2, "read");
127  TString viewstr;
128  TString viewStr;
129  for(int view = 0; view < 2; ++view){
130  if(nddatamuoncatcher) viewstr = view ? "muy" : "mux";
131  else if(!nddatamuoncatcher) viewstr = view ? "y" : "x";
132  if(nddatamuoncatcher) viewStr = view ? "muY" : "muX";
133  else if(!nddatamuoncatcher) viewStr = view ? "Y" : "X";
134 
135  f1->cd();
136  if(f1->GetListOfKeys()->Contains("calib"))
137  TH2F* h2D_1 = (TH2F*)f1->Get("calib/"+hist1+"_in"+viewStr);
138  if(f1->GetListOfKeys()->Contains("make"))
139  TH2F* h2D_1 = (TH2F*)f1->Get("make/"+hist1+"_in"+viewStr);
140  TProfile* ProfileX_1 = h2D_1->ProfileX();
141  if(det=="nd" && mcdata=="data") ProfileX_1->Rebin();
142 
143  f2->cd();
144  if(f2->GetListOfKeys()->Contains("calib"))
145  TH2F* h2D_2 = (TH2F*)f2->Get("calib/"+hist2+"_in"+viewStr);
146  else if(f2->GetListOfKeys()->Contains("make"))
147  TH2F* h2D_2 = (TH2F*)f2->Get("make/"+hist2+"_in"+viewStr);
148  TProfile* ProfileX_2 = h2D_2->ProfileX();
149  if(det=="nd" && mcdata=="data") ProfileX_2->Rebin();
150 
151  new TCanvas;
152  gPad->SetGrid();
153  HistogramAttr1D(ProfileX_1, "W (cm)", "Uncorrected PE/cm", ProfileX_1->GetXaxis()->GetBinCenter(1), ProfileX_1->GetXaxis()->GetBinCenter(ProfileX_1->GetXaxis()->GetNbins()), 15, 80, kBlue);
154  ProfileX_1->Draw();
155  ProfileX_1->SetLineColor(kBlue);
156  ProfileX_2->SetLineColor(kRed);
157  ProfileX_2->Draw("same");
158 
159  TLegend *legend = new TLegend(.12,.62,.64,.90);
160  legend->SetBorderSize(0);
161  legend->SetFillStyle(0);
162  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
163  legend->AddEntry(ProfileX_1, analysis1+" "+viewStr+" view");
164  legend->AddEntry(ProfileX_2, analysis2+" "+viewStr+" view");
165  legend->Draw();
166 
167  gPad->Print("Images/uncorrected_PEcm_vs_w_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_"+viewstr+"view.pdf");
168 
169  //ratios
170  new TCanvas;
171  gPad->SetGrid();
172  const TString hrationame = Form("hratio%s_%s_%s", det.Data(), mcdata.Data(), viewstr.Data());
173  TH1D* hratio = new TH1D(hrationame, "",ProfileX_1->GetXaxis()->GetNbins(), ProfileX_1->GetBinCenter(0), ProfileX_1->GetBinCenter(ProfileX_1->GetXaxis()->GetNbins()));
174  Ratiohist(ProfileX_1, ProfileX_2, hratio, kRed, 24, "W (cm)", "Ratio", 0.85, 1.10);
175 
176  TLegend *legend = new TLegend(.12,.62,.64,.90);
177  legend->SetBorderSize(0);
178  legend->SetFillStyle(0);
179  legend->SetTextSize(0.04);
180  legend->AddEntry((TObject*)0, det+" "+mcdata," ");
181  legend->AddEntry(hratio,analysis2+" to "+analysis1+" "+viewStr+" view");
182  legend->Draw();
183  gPad->Print("Images/uncorrected_PEcm_vs_w_ratio_"+det+"_"+mcdata+"_"+analysis1+"_vs_"+analysis2+"_"+viewstr+"view.pdf");
184  }//end of loop over views
185 } // end of uncorrected_PE()
TH1D * hist2
Definition: plotHisto.C:12
Float_t f2
c2
Definition: demo5.py:33
TH1D * hist1
Definition: plotHisto.C:9
void Ratiohist(TH1 *h1, TH1 *h2, TH1 *hratio, Color_t color, Style_t mstyle, TString labelx, TString labely, Double_t lowy, Double_t highy)
Float_t f1
c1
Definition: demo5.py:24
HistogramAttr1D(TH1 *h, char *title_x, char *title_y, Double_t binlowx, Double_t binhighx, Double_t binlowy, Double_t binhighy, Color_t lcolor)
void uncorrected_PE_tricellvstrajectory ( char *  c1,
TString  det,
TString  mcdata,
TString  analysis1,
TString  analysis2,
TString  hist1,
TString  hist2 
)

Definition at line 187 of file attenuation_calibration_comparisons.C.

References uncorrected_PE().

Referenced by attenuation_calibration_comparisons().

187  {
188  uncorrected_PE(c1, c1, det, mcdata, analysis1, analysis2, hist1, hist2);
189 }//end of uncorrected_PE_tricellvstrajectory()
TH1D * hist2
Definition: plotHisto.C:12
uncorrected_PE(char *c1, char *c2, TString det, TString mcdata, TString analysis1, TString analysis2, TString hist1, TString hist2)
TH1D * hist1
Definition: plotHisto.C:9
c1
Definition: demo5.py:24

Variable Documentation

bool calibrationconstantsplots = false
bool correctedPEplots = true
bool data = true

Definition at line 34 of file attenuation_calibration_comparisons.C.

bool fd = false
bool manytoone = false

Definition at line 37 of file attenuation_calibration_comparisons.C.

bool mc = false

Definition at line 33 of file attenuation_calibration_comparisons.C.

Referenced by ana::AddHistDefSliceTruth(), ana::AddNueHistDefTruth(), ana::AddNumuHistDefSliceTruth(), rwgt::RwgtTest::analyze(), supernova::SimAna::analyze(), supernova::SnovaAna::analyze(), mcchk::RockAna::analyze(), mcchk::NeutrinoAna::analyze(), mcchk::LeptonAna::analyze(), mcchk::DetSimAna::analyze(), slid::LIDTraining::analyze(), attenuation_calibration_comparisons(), caf_numu_reco_minus_true(), caf::Proxy< caf::SRNumuEnergy >::CheckEquals(), caf::Proxy< caf::StandardRecord >::CheckEquals(), data_over_mc_profile(), ana::DataMCAreaNormalizedRatio(), fnex::CorrectedSpectrum::DrawCounts(), fnex::CorrectedSpectrum::DrawDataVsMC(), fnex::CorrectedSpectrum::DrawStacks(), EnergyCont_macro(), fd_plot(), murem::MuonRemove::FillTruthInfo(), jmshower::RecoJMShowerFilter::filter(), fragmentEnergy(), jmshower::RecoJMShower::GetCellDistToTrk(), fnex::CorrectedSpectrum::GetStacksCanvasCopy(), genie::CharmHadronization::Hadronize(), hyperon_macro(), Legend(), ana::CountingExperiment::LoadFrom(), ana::SingleNucAnalysis::LoadFrom(), ana::NumuCC2p2hAnalysis::LoadFrom(), ana::SingleSampleExperiment::LoadFrom(), ana::CrossSectionAnalysis::LoadFrom(), ana::NumuCCIncAnalysis::LoadFrom(), ana::DataMCPair::LoadFrom(), comi::Cana::LoadMCInfo(), MakeEventList(), MakePlot(), earms::ElasticArmsHS::MCVertex(), ana::MichelDecomp::MichelDecomp(), mre_blessed(), caf::Proxy< caf::SRNumuEnergy >::operator=(), caf::Proxy< caf::StandardRecord >::operator=(), ana::OverlayCutFromNuTruthCut(), plot_assess_nd(), plot_assess_pretty_nd(), plot_nd_spectra_2018(), readExfor(), selection_story_plots(), ana::ModularExtrapComponent::SetQuiet(), solve(), and fnex::CorrectedSpectrum::UncorrectedHistMap().

bool nd = true
bool nddatamuoncatcher = true
bool onetoone = true
bool thresholdshadowingplots = false
bool tricellvstrajectory = false
bool uncorrectedPEplots = false