JFCalibration/BP_DB_Comp.cpp

// to compare the results of backplane pulsing with direct beam

#include "../sls_detector_calibration/jungfrauCommonHeader.h"
#include "../sls_detector_calibration/jungfrauCommonFunctions.h"

#include "TFile.h"
#include "TH2F.h"
#include "TCanvas.h"
#include "TF1.h"
#include "TPaveText.h"

int main(int argc, char* argv[]) {
  jungfrauStyle();
  gROOT->SetBatch(1);
  gStyle->SetOptFit(11);

  if (argc != 2) {
    cout << "Correct usage:" << endl;
    cout << "arg 1: specify module number" << endl;
    exit(1);
  }

  string module_str = argv[1];
  char savename[128];

  TCanvas *mapcanvas = new TCanvas("mapcanvas","",150,10,800,400);
  mapcanvas->SetLeftMargin(0.1);
  mapcanvas->SetRightMargin(0.13);
  mapcanvas->SetTopMargin(0.08);
  mapcanvas->SetBottomMargin(0.15);

  TCanvas *c1 = new TCanvas("c1","");

  // BP data
  sprintf(savename,"data/M%s/BP_ratio_M%s.root", module_str.c_str(), module_str.c_str());
  TFile* saved_file2 = new TFile((const char *)(savename),"READ");
  TH2F* g0overg1_map_BP = (TH2F*)saved_file2->Get("g0overg1_map");

  // DB data
  sprintf(savename,"data/M%s/DB_ratio_M%s.root", module_str.c_str(), module_str.c_str());
  TFile* saved_file = new TFile((const char *)(savename),"READ");
  TH2F* g0overg1_map_DB = (TH2F*)saved_file->Get("g0overg1map");
  TH2F* g0overg1er_map_DB = (TH2F*)saved_file->Get("g0overg1ermap");

  // difference and pull
  TH1F* diff_g0overg1_hist = new TH1F("diff_g0overg1_hist","",100,-10,10);
  TH1F* diff_g0overg1_e = new TH1F("diff_g0overg1_e","",100,-10,10);
  TH1F* diff_g0overg1_d = new TH1F("diff_g0overg1_d","",100,-10,10);
  TH1F* diff_g0overg1_ntd = new TH1F("diff_g0overg1_ntd","",100,-10,10);
  TH1F* diff_g0overg1_q = new TH1F("diff_g0overg1_q","",100,-10,10);
  TH1F* diff_g0overg1_ie = new TH1F("diff_g0overg1_ie","",100,-10,10);
  TH1F* diff_g0overg1_b = new TH1F("diff_g0overg1_b","",100,-10,10);
  TH2F* diff_g0overg1_map = new TH2F("diff_g0overg1_map","",NC,-0.5,NC-0.5,NR,-0.5,NR-0.5);
  TH1F* pull_hist = new TH1F("pull_hist","",100,-20,20);

  // scatter plots
  TH2F* scat = new TH2F("scat","",100,-35,-25,100,-35,-25);
  TH2F* scat_e = new TH2F("scat_e","",100,-35,-25,100,-35,-25);
  TH2F* scat_d = new TH2F("scat_d","",100,-35,-25,100,-35,-25);
  TH2F* scat_ntd = new TH2F("scat_ntd","",100,-35,-25,100,-35,-25);
  TH2F* scat_q = new TH2F("scat_q","",100,-35,-25,100,-35,-25);
  TH2F* scat_ie = new TH2F("scat_ie","",100,-35,-25,100,-35,-25);
  TH2F* scat_b = new TH2F("scat_b","",100,-35,-25,100,-35,-25);

  TH2F* diff_DBer = new TH2F("diff_DBer","",100,-5,5,100,0.2,1);

  double this_BP = 0;
  double this_DB = 0;

  for (int i = 0; i < NCH; i ++) {

    this_BP = g0overg1_map_BP->GetBinContent((i%NC)+1,(i/NC)+1);
    this_DB = g0overg1_map_DB->GetBinContent((i%NC)+1,(i/NC)+1);

    if (this_BP != 0 && this_DB != 0) {
      diff_g0overg1_hist->Fill(this_DB - this_BP);
      diff_g0overg1_map->SetBinContent((i%NC)+1,(i/NC)+1,this_DB - this_BP);
      pull_hist->Fill((this_DB - this_BP) / g0overg1er_map_DB->GetBinContent((i%NC)+1,(i/NC)+1));

      diff_DBer->Fill(this_DB - this_BP, g0overg1er_map_DB->GetBinContent((i%NC)+1,(i/NC)+1), 1);

      scat->Fill(this_DB, this_BP, 1);
      if (isEdge(i)) {
	scat_e->Fill(this_DB, this_BP, 1);
	diff_g0overg1_e->Fill(this_DB - this_BP);
      }
      if (isDouble(i)) {
	scat_d->Fill(this_DB, this_BP, 1);
	diff_g0overg1_d->Fill(this_DB - this_BP);
      }
      if (isNextToDouble(i)) {
	scat_ntd->Fill(this_DB, this_BP, 1);
	diff_g0overg1_ntd->Fill(this_DB - this_BP);
      }
      if (isQuad(i)) {
	scat_q->Fill(this_DB, this_BP, 1);
	diff_g0overg1_q->Fill(this_DB - this_BP);
      }
      if (isInnerEdge(i)) {
	scat_ie->Fill(this_DB, this_BP, 1);
	diff_g0overg1_ie->Fill(this_DB - this_BP);
      }
      if (isBulk(i)) {
	scat_b->Fill(this_DB, this_BP, 1);
	diff_g0overg1_b->Fill(this_DB - this_BP);
      }

    }
  }

  c1->cd();
  diff_g0overg1_hist->GetXaxis()->SetTitle("G0/G1_{DB} - G0/G1_{BP}");
  c1->SetLogy();
  diff_g0overg1_hist->Draw();
  sprintf(savename,"plots/M%s/BP_DB_Comp/diff_g0overg1_hist_log_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));
  c1->SetLogy(0);
  diff_g0overg1_hist->Fit("gaus");
  diff_g0overg1_hist->Draw();
  c1->Update();
  TPaveText *pave = new TPaveText(0.6,0.8,0.94,0.94,"blNDC");
  pave->SetBorderSize(0);
  pave->SetFillStyle(0);
  pave->SetTextSize(0.04);
  pave->SetTextAlign(32);
  TF1* diff_hist_gaus = diff_g0overg1_hist->GetFunction("gaus");
  sprintf(savename,"Mean %0.3f #pm %0.3f", diff_hist_gaus->GetParameter(1), diff_hist_gaus->GetParError(1));
  pave->AddText((const char *)(savename));
  sprintf(savename,"Sigma %0.3f #pm %0.3f", diff_hist_gaus->GetParameter(2), diff_hist_gaus->GetParError(2));
  pave->AddText((const char *)(savename));
  pave->Draw();
  diff_g0overg1_hist->SetStats(kFALSE);
  sprintf(savename,"plots/M%s/BP_DB_Comp/diff_g0overg1_hist_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  mapcanvas->cd();
  diff_g0overg1_map->GetXaxis()->SetTitle("Column");
  diff_g0overg1_map->GetYaxis()->SetTitle("Row");
  diff_g0overg1_map->GetYaxis()->SetTitleOffset(0.7);
  diff_g0overg1_map->Draw("colz");
  diff_g0overg1_map->GetZaxis()->SetRangeUser(-5,5);
  sprintf(savename,"plots/M%s/BP_DB_Comp/diff_g0overg1_map_M%s.png", module_str.c_str(), module_str.c_str());
  mapcanvas->SaveAs((const char *)(savename));

  c1->cd();
  pull_hist->GetXaxis()->SetTitle("G0/G1_{DB} - G0/G1_{BP} / #sigma_{DB}");
  pull_hist->Fit("gaus");
  pull_hist->Draw();
  c1->Update();
  TF1* pull_hist_gaus = pull_hist->GetFunction("gaus");
  pave->Clear();
  sprintf(savename,"Mean %0.3f #pm %0.3f", pull_hist_gaus->GetParameter(1), pull_hist_gaus->GetParError(1));
  pave->AddText((const char *)(savename));
  sprintf(savename,"Sigma %0.3f #pm %0.3f", pull_hist_gaus->GetParameter(2), pull_hist_gaus->GetParError(2));
  pave->AddText((const char *)(savename));
  pave->Draw();
  pull_hist->SetStats(kFALSE);
  sprintf(savename,"plots/M%s/BP_DB_Comp/pull_hist_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  diff_g0overg1_e->SetLineColor(kBlue);
  diff_g0overg1_ie->SetLineColor(kCyan);
  diff_g0overg1_d->SetLineColor(kGreen+2);
  diff_g0overg1_ntd->SetLineColor(kRed);
  diff_g0overg1_q->SetLineColor(kOrange);

  diff_g0overg1_e->Scale(1./diff_g0overg1_e->GetEntries());
  diff_g0overg1_ie->Scale(1./diff_g0overg1_ie->GetEntries());
  diff_g0overg1_d->Scale(1./diff_g0overg1_d->GetEntries());
  diff_g0overg1_ntd->Scale(1./diff_g0overg1_ntd->GetEntries());
  diff_g0overg1_q->Scale(1./diff_g0overg1_q->GetEntries());
  diff_g0overg1_b->Scale(1./diff_g0overg1_b->GetEntries());

  diff_g0overg1_d->GetXaxis()->SetTitle("G0/G1_{DB} - G0/G1_{BP}");
  diff_g0overg1_d->Draw();
  diff_g0overg1_e->Draw("same");
  diff_g0overg1_ie->Draw("same");
  diff_g0overg1_b->Draw("same");
  diff_g0overg1_ntd->Draw("same");
  //diff_g0overg1_q->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/diff_g0overg1_hist_byType_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  c1->SetRightMargin(0.2);

  TF1 *diag = new TF1("diag","[0]*x",-35,-25);
  diag->SetParameter(0,1);

  scat->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat->GetYaxis()->SetTitleOffset(1.3);
  scat->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_e->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_e->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_e->GetYaxis()->SetTitleOffset(1.3);
  scat_e->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_e_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_d->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_d->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_d->GetYaxis()->SetTitleOffset(1.3);
  scat_d->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_d_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_ntd->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_ntd->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_ntd->GetYaxis()->SetTitleOffset(1.3);
  scat_ntd->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_ntd_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_q->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_q->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_q->GetYaxis()->SetTitleOffset(1.3);
  scat_q->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_q_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_ie->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_ie->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_ie->GetYaxis()->SetTitleOffset(1.3);
  scat_ie->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_ie_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  scat_b->GetXaxis()->SetTitle("G0/G1_{DB}");
  scat_b->GetYaxis()->SetTitle("G0/G1_{BP}");
  scat_b->GetYaxis()->SetTitleOffset(1.3);
  scat_b->Draw("colz");
  diag->Draw("same");
  sprintf(savename,"plots/M%s/BP_DB_Comp/scat_b_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

  diff_DBer->GetXaxis()->SetTitle("G0/G1_{DB} - G0/G1_{BP}");
  diff_DBer->GetYaxis()->SetTitle("#sigma_{G0/G1_{DB}}");
  diff_DBer->GetYaxis()->SetTitleOffset(1.3);
  diff_DBer->Draw("colz");
  sprintf(savename,"plots/M%s/BP_DB_Comp/diff_DBer_M%s.png", module_str.c_str(), module_str.c_str());
  c1->SaveAs((const char *)(savename));

}