[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TOF / pPb502 / macros / TZEROcalib.C

#include "CommonDefs.C"

Int_t ncentbins = 5;
Int_t nvertexbins = 20;

TZEROcalib(const Char_t *filename = "data.root", const Char_t *calibfilename = NULL, Int_t evMax = kMaxInt)
{

  /* include path for ACLic */
  gSystem->AddIncludePath("-I$ALICE_ROOT/include");
  gSystem->AddIncludePath("-I$ALICE_ROOT/TOF");
  /* load libraries */
  gSystem->Load("libANALYSIS");
  gSystem->Load("libANALYSISalice");
  /* build analysis task class */
  gROOT->LoadMacro("AliAnalysisParticle.cxx+g");
  gROOT->LoadMacro("AliAnalysisEvent.cxx+g");
  gROOT->LoadMacro("AliAnalysisTrack.cxx+g");

  /* open file, get tree and connect */
  TFile *filein = TFile::Open(filename);
  TTree *treein = (TTree *)filein->Get("aodTree");
  printf("got \"aodTree\": %d entries\n", treein->GetEntries());
  AliAnalysisEvent *analysisEvent = new AliAnalysisEvent();
  TClonesArray *analysisTrackArray = new TClonesArray("AliAnalysisTrack");
  AliAnalysisTrack *analysisTrack = NULL;
  treein->SetBranchAddress("AnalysisEvent", &analysisEvent);
  treein->SetBranchAddress("AnalysisTrack", &analysisTrackArray);

  /* open calibfile */
  TH1 *hCentrality_TZEROA_mean = NULL;
  TH1 *hCentrality_TZEROA_sigma = NULL;
  TH1 *hCentrality_TZEROC_mean = NULL;
  TH1 *hCentrality_TZEROC_sigma = NULL;
  TH1 *hCentrality_TOF_mean = NULL;
  TH1 *hCentrality_TOF_TZEROA_mean = NULL;
  TH1 *hCentrality_TOF_TZEROC_mean = NULL;
  TH1 *hCentrality_TOF_TZEROTOF_mean = NULL;
  if (calibfilename) {
    TFile *calibfile = TFile::Open(calibfilename);
    hCentrality_TZEROA_mean = (TH1 *)calibfile->Get("hCentrality_TZEROA_mean");
    hCentrality_TZEROA_sigma = (TH1 *)calibfile->Get("hCentrality_TZEROA_sigma");
    hCentrality_TZEROC_mean = (TH1 *)calibfile->Get("hCentrality_TZEROC_calib");
    hCentrality_TZEROC_sigma = (TH1 *)calibfile->Get("hCentrality_TZEROC_sigma");
    hCentrality_TOF_mean = (TH1 *)calibfile->Get("hCentrality_TOF_mean");
    hCentrality_TOF_TZEROA_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROA_mean");
    hCentrality_TOF_TZEROC_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROC_mean");
    hCentrality_TOF_TZEROTOF_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROTOF_mean");
  }

  /* histos */
  TH1F *hCentrality = new TH1F("hCentrality", "", NcentralityBins, centralityBin);
  TH1F *hCentrality_AC = new TH1F("hCentrality_AC", "", NcentralityBins, centralityBin);
  TH1F *hCentrality_A = new TH1F("hCentrality_A", "", NcentralityBins, centralityBin);
  TH1F *hCentrality_C = new TH1F("hCentrality_C", "", NcentralityBins, centralityBin);
  TH1F *hCentrality_NONE = new TH1F("hCentrality_NONE", "", NcentralityBins, centralityBin);

  TH2F *hCentrality_TZEROA = new TH2F("hCentrality_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROC = new TH2F("hCentrality_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROOR = new TH2F("hCentrality_TZEROOR", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROAND = new TH2F("hCentrality_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROTOF = new TH2F("hCentrality_TZEROTOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZERODIFF = new TH2F("hCentrality_TZERODIFF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);

  TH2F *hCentrality_TZEROTOF_TZEROA = new TH2F("hCentrality_TZEROTOF_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROTOF_TZEROC = new TH2F("hCentrality_TZEROTOF_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROTOF_TZEROAND = new TH2F("hCentrality_TZEROTOF_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TZEROA_TZEROC = new TH2F("hCentrality_TZEROA_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);

  TH2F *hVertex_TZEROA = new TH2F("hVertex_TZEROA", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TZEROC = new TH2F("hVertex_TZEROC", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TZEROOR = new TH2F("hVertex_TZEROOR", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TZEROAND = new TH2F("hVertex_TZEROAND", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TZEROTOF = new TH2F("hVertex_TZEROTOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TZERODIFF = new TH2F("hVertex_TZERODIFF", "", nvertexbins, -10., 10., 200, -2440., 2440.);


  TH2F *hCentrality_TOF = new TH2F("hCentrality_TOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TOF_TZEROA = new TH2F("hCentrality_TOF_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TOF_TZEROC = new TH2F("hCentrality_TOF_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TOF_TZEROAND = new TH2F("hCentrality_TOF_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_TOF_TZEROTOF = new TH2F("hCentrality_TOF_TZEROTOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
  TH2F *hCentrality_Resolution = new TH2F("hCentrality_Resolution", "", NcentralityBins, centralityBin, 20, 0., 100.);

  TH2F *hResolution_TZERODIFF = new TH2F("hResolution_TZERODIFF", "", 20, 0., 100., 200, -2440., 2440.);
  TH2F *hResolution_TOF_TZEROA = new TH2F("hResolution_TOF_TZEROA", "", 20, 0., 100., 200, -2440., 2440.);
  TH2F *hResolution_TOF_TZEROC = new TH2F("hResolution_TOF_TZEROC", "", 20, 0., 100., 200, -2440., 2440.);
  TH2F *hResolution_TOF_TZEROAND = new TH2F("hResolution_TOF_TZEROAND", "", 20, 0., 100., 200, -2440., 2440.);
  TH2F *hResolution_TOF_TZEROTOF = new TH2F("hResolution_TOF_TZEROTOF", "", 20, 0., 100., 200, -2440., 2440.);

  TH2F *hVertex_TOF = new TH2F("hVertex_TOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TOF_TZEROA = new TH2F("hVertex_TOF_TZEROA", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TOF_TZEROC = new TH2F("hVertex_TOF_TZEROC", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TOF_TZEROOR = new TH2F("hVertex_TOF_TZEROOR", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TOF_TZEROAND = new TH2F("hVertex_TOF_TZEROAND", "", nvertexbins, -10., 10., 200, -2440., 2440.);
  TH2F *hVertex_TOF_TZEROTOF = new TH2F("hVertex_TOF_TZEROTOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);

  TH2F *hDiffDiff = new TH2F("hDiffDiff", "", 200, -2440., 2440., 200, -2440., 2440.);

  /* loop over events */
  for (Int_t iev = 0; iev < treein->GetEntries() && iev < evMax; iev++) {
    /* get event */
    treein->GetEvent(iev);
    if (iev % 100000 == 0) printf("iev = %d\n", iev);
    /* check event */
    if (!analysisEvent->AcceptEvent(acceptEventType)) continue;
    
    /*** ACCEPTED EVENT ***/

    /* get vertex position */
    Double_t vertexz = analysisEvent->GetVertexZ();

    /* get centrality */
    Double_t cent = analysisEvent->GetCentralityPercentile(centralityEstimator);
    hCentrality->Fill(cent);

    Int_t icent;
    for (icent = 0; icent < NcentralityBins; icent++)
      if (cent < centralityBin[icent + 1])
	break;

    Double_t TZEROA_mean = hCentrality_TZEROA_mean ? hCentrality_TZEROA_mean->GetBinContent(icent + 1) : 0.;
    Double_t TZEROA_sigma = hCentrality_TZEROA_sigma ? hCentrality_TZEROA_sigma->GetBinContent(icent + 1) : 1000.;
    Double_t TZEROC_mean  = hCentrality_TZEROC_mean ? hCentrality_TZEROC_mean->GetBinContent(icent + 1) : 0.;
    Double_t TZEROC_sigma = hCentrality_TZEROC_sigma ? hCentrality_TZEROC_sigma->GetBinContent(icent + 1) : 1000.;

    Double_t TOF_mean = hCentrality_TOF_mean ? hCentrality_TOF_mean->GetBinContent(icent + 1) : 0.;
    Double_t TOF_TZEROA_mean = hCentrality_TOF_TZEROA_mean ? hCentrality_TOF_TZEROA_mean->GetBinContent(icent + 1) : 0.;
    Double_t TOF_TZEROC_mean = hCentrality_TOF_TZEROC_mean ? hCentrality_TOF_TZEROC_mean->GetBinContent(icent + 1) : 0.;
    Double_t TOF_TZEROTOF_mean = hCentrality_TOF_TZEROTOF_mean ? hCentrality_TOF_TZEROTOF_mean->GetBinContent(icent + 1) : 0.;
    
    
    /* TZERO */
    Double_t TZEROA = analysisEvent->GetTimeZeroT0(1) - TZEROA_shift;
    Bool_t hasTZEROA = TMath::Abs(TZEROA - TZEROA_mean) < 3. * TZEROA_sigma;
    Double_t TZEROC = analysisEvent->GetTimeZeroT0(2) - TZEROC_shift;
    Bool_t hasTZEROC = TMath::Abs(TZEROC - TZEROC_mean) < 3. * TZEROC_sigma;
    /* vertex correction */
    //    TZEROA += -TZEROvertexCorr * vertexz;
    //    TZEROC += TZEROvertexCorr * vertexz;
    /* alignment to TOF */
    TZEROA += TOF_TZEROA_mean - TOF_mean;
    TZEROC += TOF_TZEROC_mean - TOF_mean;
    /* TZEROAND */
    Double_t TZEROAND = (TZEROA + TZEROC) * 0.5;
    Bool_t hasTZEROAND = hasTZEROA && hasTZEROC;
    /* TZEROTOF */
    Double_t TZEROTOF = analysisEvent->GetTimeZeroTOF()[9];
    Bool_t hasTZEROTOF = analysisEvent->GetTimeZeroTOFSigma()[9] < 150.;
    TZEROTOF += TOF_TZEROTOF_mean - TOF_mean;

    if (hasTZEROA) {
      hCentrality_TZEROA->Fill(cent, TZEROA);
      hVertex_TZEROA->Fill(vertexz, TZEROA);
    }
    if (hasTZEROC) {
      hCentrality_TZEROC->Fill(cent, TZEROC);
      hVertex_TZEROC->Fill(vertexz, TZEROC);
    }
    if (hasTZEROAND) {
      hCentrality_TZEROAND->Fill(cent, TZEROAND);
      hVertex_TZEROAND->Fill(vertexz, TZEROAND);
    }
    if (hasTZEROTOF) {
      hCentrality_TZEROTOF->Fill(cent, TZEROTOF);
      hVertex_TZEROTOF->Fill(vertexz, TZEROTOF);
    }

    if (hasTZEROA && hasTZEROC && hasTZEROAND && hasTZEROTOF) {
      hCentrality_TZEROTOF_TZEROA->Fill(cent, TZEROTOF - TZEROA);
      hCentrality_TZEROTOF_TZEROC->Fill(cent, TZEROTOF - TZEROC);
      hCentrality_TZEROTOF_TZEROAND->Fill(cent, TZEROTOF - TZEROAND);
      hCentrality_TZEROA_TZEROC->Fill(cent, TZEROA - TZEROC);
    }

    /* loop over tracks */
    for (Int_t itrk = 0; itrk < analysisTrackArray->GetEntries(); itrk++) {
      /* get track */
      AliAnalysisTrack *analysisTrack = (AliAnalysisTrack *)analysisTrackArray->At(itrk);
      if (!analysisTrack) continue;
      /* check accepted track */
      if (!analysisTrack->AcceptTrack()) continue;
      /* check momentum */
      if (analysisTrack->GetP() < 0.9 || analysisTrack->GetP() > 1.1) continue;
      /* check TOF pid */
      if (!analysisTrack->HasTOFPID()) continue;
      
      TZEROTOF = analysisEvent->GetTimeZeroTOF(analysisTrack->GetP());
      TZEROTOF += TOF_TZEROTOF_mean - TOF_mean;
      hasTZEROTOF = analysisEvent->GetTimeZeroTOFSigma(analysisTrack->GetP()) < 150.;

      Double_t time = analysisTrack->GetTOFTime() - TOF_mean;
      Double_t texp = analysisTrack->GetTOFExpTime(2);
      Double_t deltat = time - texp;

      
      hCentrality_TOF->Fill(cent, deltat);
      hVertex_TOF->Fill(vertexz, deltat);
      if (hasTZEROA) {
	hCentrality_TOF_TZEROA->Fill(cent, deltat - TZEROA);
	hVertex_TOF_TZEROA->Fill(vertexz, deltat - TZEROA);
      }
      if (hasTZEROC) {
	hCentrality_TOF_TZEROC->Fill(cent, deltat - TZEROC);
	hVertex_TOF_TZEROC->Fill(vertexz, deltat - TZEROC);
      }
      if (hasTZEROAND) {
	hCentrality_TOF_TZEROAND->Fill(cent, deltat - TZEROAND);
	hVertex_TOF_TZEROAND->Fill(vertexz, deltat - TZEROAND);
      }
      if (hasTZEROTOF) {
	hCentrality_TOF_TZEROTOF->Fill(cent, deltat - TZEROTOF);
	hVertex_TOF_TZEROTOF->Fill(vertexz, deltat - TZEROTOF);
      }

    }

  } /* end of loop over events */
  
  /* output */
  if (!calibfilename) TFile *fileout = TFile::Open(Form("TZEROcalib.%s", filename), "RECREATE");
  else TFile *fileout = TFile::Open(Form("TZEROcheck.%s", filename), "RECREATE");

  hCentrality->Write();
  hCentrality_TZEROA->Write();
  hCentrality_TZEROC->Write();
  hCentrality_TZEROAND->Write();
  hCentrality_TZEROTOF->Write();
  hCentrality_TZERODIFF->Write();

  hVertex_TZEROA->Write();
  hVertex_TZEROC->Write();
  hVertex_TZEROAND->Write();
  hVertex_TZEROTOF->Write();
  hVertex_TZERODIFF->Write();

  hCentrality_TOF->Write();
  hCentrality_TOF_TZEROA->Write();
  hCentrality_TOF_TZEROC->Write();
  hCentrality_TOF_TZEROAND->Write();
  hCentrality_TOF_TZEROTOF->Write();
  hCentrality_Resolution->Write();

  hCentrality_TZEROTOF_TZEROA->Write();
  hCentrality_TZEROTOF_TZEROC->Write();
  hCentrality_TZEROTOF_TZEROAND->Write();
  hCentrality_TZEROA_TZEROC->Write();

  hVertex_TOF->Write();
  hVertex_TOF_TZEROA->Write();
  hVertex_TOF_TZEROC->Write();
  hVertex_TOF_TZEROAND->Write();
  hVertex_TOF_TZEROTOF->Write();

  hResolution_TZERODIFF->Write();
  hResolution_TOF_TZEROA->Write();
  hResolution_TOF_TZEROC->Write();
  hResolution_TOF_TZEROAND->Write();

  hDiffDiff->Write();

  fileout->Close();

  if (!calibfilename) {
    TZEROcalibration(Form("TZEROcalib.%s", filename));
    TZEROcalib(filename, "TZEROcalibration.root");
    TZEROresolution(Form("TZEROcheck.%s", filename));
  }

}

TZEROcalibration(const Char_t *filename = "TZEROcalib.data.root", Float_t min = 2., Float_t max = 1.)
{

  gROOT->LoadMacro("FitPeak.C");

  TFile *fin = TFile::Open(filename);
  TH2 *hCentrality_TOF = (TH2 *)fin->Get("hCentrality_TOF");
  TH2 *hCentrality_TZEROA = (TH2 *)fin->Get("hCentrality_TZEROA");
  TH2 *hCentrality_TZEROC = (TH2 *)fin->Get("hCentrality_TZEROC");

  TH2 *hCentrality_TOF_TZEROTOF = (TH2 *)fin->Get("hCentrality_TOF_TZEROTOF");
  TH2 *hCentrality_TOF_TZEROA = (TH2 *)fin->Get("hCentrality_TOF_TZEROA");
  TH2 *hCentrality_TOF_TZEROC = (TH2 *)fin->Get("hCentrality_TOF_TZEROC");
  TH2 *hCentrality_TOF_TZEROTOF = (TH2 *)fin->Get("hCentrality_TOF_TZEROTOF");

  TFile *fout = TFile::Open("TZEROcalibration.root", "RECREATE");
  FitPeak(hCentrality_TOF, 200., 2., 1., "hCentrality_TOF")->Write();
  FitPeak(hCentrality_TZEROA, 200., 2., 1., "hCentrality_TZEROA")->Write();
  FitPeak(hCentrality_TZEROC, 200., 2., 1., "hCentrality_TZEROC")->Write();

  FitPeak(hCentrality_TOF_TZEROA, 200., 2., 1., "hCentrality_TOF_TZEROA")->Write();
  FitPeak(hCentrality_TOF_TZEROC, 200., 2., 1., "hCentrality_TOF_TZEROC")->Write();
  FitPeak(hCentrality_TOF_TZEROTOF, 200., 2., 1., "hCentrality_TOF_TZEROTOF")->Write();

  fout->Close();
}

TZEROresolution(const Char_t *filename, Float_t min = 2., Float_t max = 2., Char_t *what = "Centrality")
{

  TFile *fin = TFile::Open(filename);
  hdiff = (TH2 *)fin->Get(Form("h%s_TZEROA_TZEROC", what));
  hA = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROA", what));
  hC = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROC", what));
  hAND = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROAND", what));

  gROOT->LoadMacro("FitPeak.C");

  
  TH1 *hdiff_sigma = (TH1 *)FitPeak(hdiff, 100., min, max, "sigma1")->At(1);
  
  TH1 *hA_sigma = (TH1 *)FitPeak(hA, 100., min, max, "sigma2")->At(1);
  
  TH1 *hC_sigma = (TH1 *)FitPeak(hC, 100., min, max, "sigma3")->At(1);
  
  TH1 *hAND_sigma = (TH1 *)FitPeak(hAND, 100., min, max, "sigma4")->At(1);

  TH1 *hTZEROTOF_reso = (TH1 *)hdiff_sigma->Clone("hTZEROTOF_reso");
  hTZEROTOF_reso->Reset();
  hTZEROTOF_reso->SetMarkerStyle(20);
  hTZEROTOF_reso->SetMarkerColor(2);
  hTZEROTOF_reso->SetLineColor(2);
  hTZEROTOF_reso->SetTitle("T0-TOF");
  TH1 *hTZEROA_reso = (TH1 *)hdiff_sigma->Clone("hTZEROA_reso");
  hTZEROA_reso->Reset();
  hTZEROA_reso->SetMarkerStyle(21);
  hTZEROA_reso->SetMarkerColor(8);
  hTZEROA_reso->SetLineColor(8);
  hTZEROA_reso->SetTitle("T0-TZEROA");
  TH1 *hTZEROC_reso = (TH1 *)hdiff_sigma->Clone("hTZEROC_reso");
  hTZEROC_reso->Reset();
  hTZEROC_reso->SetMarkerStyle(22);
  hTZEROC_reso->SetMarkerColor(4);
  hTZEROC_reso->SetLineColor(4);
  hTZEROC_reso->SetTitle("T0-TZEROC");
  TH1 *hTZEROAND_reso = (TH1 *)hdiff_sigma->Clone("hTZEROAND_reso");
  hTZEROAND_reso->Reset();
  hTZEROAND_reso->SetMarkerStyle(23);
  hTZEROAND_reso->SetMarkerColor(kYellow+1);
  hTZEROAND_reso->SetLineColor(kYellow+1);
  hTZEROAND_reso->SetTitle("T0-TZEROAND");

  for (Int_t i = 0; i < hTZEROTOF_reso->GetNbinsX(); i++) {

    Double_t diff = hdiff_sigma->GetBinContent(i+1);
    Double_t diff_e = hdiff_sigma->GetBinError(i+1);
    Double_t diff2 = diff*diff;
    Double_t diff2_e = 2.*diff*diff_e;
    Double_t a = hA_sigma->GetBinContent(i+1);
    Double_t a_e = hA_sigma->GetBinError(i+1);
    Double_t a2 = a*a;
    Double_t a2_e = 2.*a*a_e;
    Double_t c = hC_sigma->GetBinContent(i+1);
    Double_t c_e = hC_sigma->GetBinError(i+1);
    Double_t c2 = c*c;
    Double_t c2_e = 2.*c*c_e;
    Double_t ac = hAND_sigma->GetBinContent(i+1);
    Double_t ac_e = hAND_sigma->GetBinError(i+1);
    Double_t ac2 = ac*ac;
    Double_t ac2_e = 2.*ac*ac_e;

    if (diff == 0.) continue;
    if (a2 + c2 - diff2 < 0.) continue;

    printf("a=%f, c=%f, diff=%f, ac=%f\n", a, c, diff, ac);

    Double_t tofreso2 = (a2 + c2 - diff2) * 0.5;
    Double_t tofreso2_e = (a2_e + c2_e + diff2_e) * 0.5;
    Double_t tofreso = TMath::Sqrt(tofreso2);
    Double_t tofreso_e = tofreso2_e / (2. * tofreso);
    hTZEROTOF_reso->SetBinContent(i+1, tofreso);
    hTZEROTOF_reso->SetBinError(i+1, tofreso_e);

    if (a2 - tofreso*tofreso < 0.) continue;
    if (c2 - tofreso*tofreso < 0.) continue;
    if (ac2 - tofreso*tofreso < 0.) continue;

    Double_t tzeroareso2 = a2 - tofreso2;
    Double_t tzeroareso2_e = a2_e + tofreso2_e;
    Double_t tzeroareso = TMath::Sqrt(tzeroareso2);
    Double_t tzeroareso_e = tzeroareso2_e / (2. * tzeroareso);
    hTZEROA_reso->SetBinContent(i+1, tzeroareso);
    hTZEROA_reso->SetBinError(i+1, tzeroareso_e);

    Double_t tzerocreso2 = c2 - tofreso2;
    Double_t tzerocreso2_e = c2_e + tofreso2_e;
    Double_t tzerocreso = TMath::Sqrt(tzerocreso2);
    Double_t tzerocreso_e = tzerocreso2_e / (2. * tzerocreso);
    hTZEROC_reso->SetBinContent(i+1, tzerocreso);
    hTZEROC_reso->SetBinError(i+1, tzerocreso_e);

    Double_t tzeroacreso2 = ac2 - tofreso2;
    Double_t tzeroacreso2_e = ac2_e + tofreso2_e;
    Double_t tzeroacreso = TMath::Sqrt(tzeroacreso2);
    Double_t tzeroacreso_e = tzeroacreso2_e / (2. * tzeroacreso);
    hTZEROAND_reso->SetBinContent(i+1, tzeroacreso);
    hTZEROAND_reso->SetBinError(i+1, tzeroacreso_e);

    //   Double_t tzerocreso = TMath::Sqrt(c2 - tofreso*tofreso);
    //    hTZEROC_reso->SetBinContent(i+1, tzerocreso);

    //    Double_t tzeroacreso = TMath::Sqrt(ac2 - tofreso*tofreso);
    //    hTZEROAND_reso->SetBinContent(i+1, tzeroacreso);

    printf("TOF=%f, TZEROA=%f, TZEROC=%f, TZEROAND=%f\n", tofreso, tzeroareso, tzerocreso, tzeroacreso);
  }

  hTZEROTOF_reso->Draw();
  hTZEROA_reso->Draw("same");
  hTZEROC_reso->Draw("same");
  hTZEROAND_reso->Draw("same");

  TFile *fout = TFile::Open("TZEROresolution.root", "RECREATE");
  hTZEROTOF_reso->Write();
  hTZEROA_reso->Write();
  hTZEROC_reso->Write();
  hTZEROAND_reso->Write();
  fout->Close();
}
Commit	Line	Data
59e49925	1	#include "CommonDefs.C"
	2
	3	Int_t ncentbins = 5;
	4	Int_t nvertexbins = 20;
	5
	6	TZEROcalib(const Char_t filename = "data.root", const Char_t calibfilename = NULL, Int_t evMax = kMaxInt)
	7	{
	8
	9	/* include path for ACLic */
	10	gSystem->AddIncludePath("-I$ALICE_ROOT/include");
	11	gSystem->AddIncludePath("-I$ALICE_ROOT/TOF");
	12	/* load libraries */
	13	gSystem->Load("libANALYSIS");
	14	gSystem->Load("libANALYSISalice");
	15	/* build analysis task class */
	16	gROOT->LoadMacro("AliAnalysisParticle.cxx+g");
	17	gROOT->LoadMacro("AliAnalysisEvent.cxx+g");
	18	gROOT->LoadMacro("AliAnalysisTrack.cxx+g");
	19
	20	/* open file, get tree and connect */
	21	TFile *filein = TFile::Open(filename);
	22	TTree treein = (TTree )filein->Get("aodTree");
	23	printf("got \"aodTree\": %d entries\n", treein->GetEntries());
	24	AliAnalysisEvent *analysisEvent = new AliAnalysisEvent();
	25	TClonesArray *analysisTrackArray = new TClonesArray("AliAnalysisTrack");
	26	AliAnalysisTrack *analysisTrack = NULL;
	27	treein->SetBranchAddress("AnalysisEvent", &analysisEvent);
	28	treein->SetBranchAddress("AnalysisTrack", &analysisTrackArray);
	29
	30	/* open calibfile */
	31	TH1 *hCentrality_TZEROA_mean = NULL;
	32	TH1 *hCentrality_TZEROA_sigma = NULL;
	33	TH1 *hCentrality_TZEROC_mean = NULL;
	34	TH1 *hCentrality_TZEROC_sigma = NULL;
	35	TH1 *hCentrality_TOF_mean = NULL;
	36	TH1 *hCentrality_TOF_TZEROA_mean = NULL;
	37	TH1 *hCentrality_TOF_TZEROC_mean = NULL;
	38	TH1 *hCentrality_TOF_TZEROTOF_mean = NULL;
	39	if (calibfilename) {
	40	TFile *calibfile = TFile::Open(calibfilename);
	41	hCentrality_TZEROA_mean = (TH1 *)calibfile->Get("hCentrality_TZEROA_mean");
	42	hCentrality_TZEROA_sigma = (TH1 *)calibfile->Get("hCentrality_TZEROA_sigma");
	43	hCentrality_TZEROC_mean = (TH1 *)calibfile->Get("hCentrality_TZEROC_calib");
	44	hCentrality_TZEROC_sigma = (TH1 *)calibfile->Get("hCentrality_TZEROC_sigma");
	45	hCentrality_TOF_mean = (TH1 *)calibfile->Get("hCentrality_TOF_mean");
	46	hCentrality_TOF_TZEROA_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROA_mean");
	47	hCentrality_TOF_TZEROC_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROC_mean");
	48	hCentrality_TOF_TZEROTOF_mean = (TH1 *)calibfile->Get("hCentrality_TOF_TZEROTOF_mean");
	49	}
	50
	51	/* histos */
	52	TH1F *hCentrality = new TH1F("hCentrality", "", NcentralityBins, centralityBin);
	53	TH1F *hCentrality_AC = new TH1F("hCentrality_AC", "", NcentralityBins, centralityBin);
	54	TH1F *hCentrality_A = new TH1F("hCentrality_A", "", NcentralityBins, centralityBin);
	55	TH1F *hCentrality_C = new TH1F("hCentrality_C", "", NcentralityBins, centralityBin);
	56	TH1F *hCentrality_NONE = new TH1F("hCentrality_NONE", "", NcentralityBins, centralityBin);
	57
	58	TH2F *hCentrality_TZEROA = new TH2F("hCentrality_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	59	TH2F *hCentrality_TZEROC = new TH2F("hCentrality_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	60	TH2F *hCentrality_TZEROOR = new TH2F("hCentrality_TZEROOR", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	61	TH2F *hCentrality_TZEROAND = new TH2F("hCentrality_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	62	TH2F *hCentrality_TZEROTOF = new TH2F("hCentrality_TZEROTOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	63	TH2F *hCentrality_TZERODIFF = new TH2F("hCentrality_TZERODIFF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
	64
65	TH2F *hCentrality_TZEROTOF_TZEROA = new TH2F("hCentrality_TZEROTOF_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
66	TH2F *hCentrality_TZEROTOF_TZEROC = new TH2F("hCentrality_TZEROTOF_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
67	TH2F *hCentrality_TZEROTOF_TZEROAND = new TH2F("hCentrality_TZEROTOF_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
68	TH2F *hCentrality_TZEROA_TZEROC = new TH2F("hCentrality_TZEROA_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
69
70	TH2F *hVertex_TZEROA = new TH2F("hVertex_TZEROA", "", nvertexbins, -10., 10., 200, -2440., 2440.);
71	TH2F *hVertex_TZEROC = new TH2F("hVertex_TZEROC", "", nvertexbins, -10., 10., 200, -2440., 2440.);
72	TH2F *hVertex_TZEROOR = new TH2F("hVertex_TZEROOR", "", nvertexbins, -10., 10., 200, -2440., 2440.);
73	TH2F *hVertex_TZEROAND = new TH2F("hVertex_TZEROAND", "", nvertexbins, -10., 10., 200, -2440., 2440.);
74	TH2F *hVertex_TZEROTOF = new TH2F("hVertex_TZEROTOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
75	TH2F *hVertex_TZERODIFF = new TH2F("hVertex_TZERODIFF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
76
77
78	TH2F *hCentrality_TOF = new TH2F("hCentrality_TOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
79	TH2F *hCentrality_TOF_TZEROA = new TH2F("hCentrality_TOF_TZEROA", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
80	TH2F *hCentrality_TOF_TZEROC = new TH2F("hCentrality_TOF_TZEROC", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
81	TH2F *hCentrality_TOF_TZEROAND = new TH2F("hCentrality_TOF_TZEROAND", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
82	TH2F *hCentrality_TOF_TZEROTOF = new TH2F("hCentrality_TOF_TZEROTOF", "", NcentralityBins, centralityBin, 200, -2440., 2440.);
83	TH2F *hCentrality_Resolution = new TH2F("hCentrality_Resolution", "", NcentralityBins, centralityBin, 20, 0., 100.);
84
85	TH2F *hResolution_TZERODIFF = new TH2F("hResolution_TZERODIFF", "", 20, 0., 100., 200, -2440., 2440.);
86	TH2F *hResolution_TOF_TZEROA = new TH2F("hResolution_TOF_TZEROA", "", 20, 0., 100., 200, -2440., 2440.);
87	TH2F *hResolution_TOF_TZEROC = new TH2F("hResolution_TOF_TZEROC", "", 20, 0., 100., 200, -2440., 2440.);
88	TH2F *hResolution_TOF_TZEROAND = new TH2F("hResolution_TOF_TZEROAND", "", 20, 0., 100., 200, -2440., 2440.);
89	TH2F *hResolution_TOF_TZEROTOF = new TH2F("hResolution_TOF_TZEROTOF", "", 20, 0., 100., 200, -2440., 2440.);
90
91	TH2F *hVertex_TOF = new TH2F("hVertex_TOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
92	TH2F *hVertex_TOF_TZEROA = new TH2F("hVertex_TOF_TZEROA", "", nvertexbins, -10., 10., 200, -2440., 2440.);
93	TH2F *hVertex_TOF_TZEROC = new TH2F("hVertex_TOF_TZEROC", "", nvertexbins, -10., 10., 200, -2440., 2440.);
94	TH2F *hVertex_TOF_TZEROOR = new TH2F("hVertex_TOF_TZEROOR", "", nvertexbins, -10., 10., 200, -2440., 2440.);
95	TH2F *hVertex_TOF_TZEROAND = new TH2F("hVertex_TOF_TZEROAND", "", nvertexbins, -10., 10., 200, -2440., 2440.);
96	TH2F *hVertex_TOF_TZEROTOF = new TH2F("hVertex_TOF_TZEROTOF", "", nvertexbins, -10., 10., 200, -2440., 2440.);
97
98	TH2F *hDiffDiff = new TH2F("hDiffDiff", "", 200, -2440., 2440., 200, -2440., 2440.);
99
100	/* loop over events */
101	for (Int_t iev = 0; iev < treein->GetEntries() && iev < evMax; iev++) {
102	/* get event */
103	treein->GetEvent(iev);
104	if (iev % 100000 == 0) printf("iev = %d\n", iev);
105	/* check event */
106	if (!analysisEvent->AcceptEvent(acceptEventType)) continue;
107
108	/* ACCEPTED EVENT */
109
110	/* get vertex position */
111	Double_t vertexz = analysisEvent->GetVertexZ();
112
113	/* get centrality */
114	Double_t cent = analysisEvent->GetCentralityPercentile(centralityEstimator);
115	hCentrality->Fill(cent);
116
117	Int_t icent;
118	for (icent = 0; icent < NcentralityBins; icent++)
119	if (cent < centralityBin[icent + 1])
120	break;
121
122	Double_t TZEROA_mean = hCentrality_TZEROA_mean ? hCentrality_TZEROA_mean->GetBinContent(icent + 1) : 0.;
123	Double_t TZEROA_sigma = hCentrality_TZEROA_sigma ? hCentrality_TZEROA_sigma->GetBinContent(icent + 1) : 1000.;
124	Double_t TZEROC_mean = hCentrality_TZEROC_mean ? hCentrality_TZEROC_mean->GetBinContent(icent + 1) : 0.;
125	Double_t TZEROC_sigma = hCentrality_TZEROC_sigma ? hCentrality_TZEROC_sigma->GetBinContent(icent + 1) : 1000.;
126
127	Double_t TOF_mean = hCentrality_TOF_mean ? hCentrality_TOF_mean->GetBinContent(icent + 1) : 0.;
128	Double_t TOF_TZEROA_mean = hCentrality_TOF_TZEROA_mean ? hCentrality_TOF_TZEROA_mean->GetBinContent(icent + 1) : 0.;
129	Double_t TOF_TZEROC_mean = hCentrality_TOF_TZEROC_mean ? hCentrality_TOF_TZEROC_mean->GetBinContent(icent + 1) : 0.;
130	Double_t TOF_TZEROTOF_mean = hCentrality_TOF_TZEROTOF_mean ? hCentrality_TOF_TZEROTOF_mean->GetBinContent(icent + 1) : 0.;
131
132
133	/* TZERO */
134	Double_t TZEROA = analysisEvent->GetTimeZeroT0(1) - TZEROA_shift;
135	Bool_t hasTZEROA = TMath::Abs(TZEROA - TZEROA_mean) < 3. * TZEROA_sigma;
136	Double_t TZEROC = analysisEvent->GetTimeZeroT0(2) - TZEROC_shift;
137	Bool_t hasTZEROC = TMath::Abs(TZEROC - TZEROC_mean) < 3. * TZEROC_sigma;
138	/* vertex correction */
139	// TZEROA += -TZEROvertexCorr * vertexz;
140	// TZEROC += TZEROvertexCorr * vertexz;
141	/* alignment to TOF */
142	TZEROA += TOF_TZEROA_mean - TOF_mean;
143	TZEROC += TOF_TZEROC_mean - TOF_mean;
144	/* TZEROAND */
145	Double_t TZEROAND = (TZEROA + TZEROC) * 0.5;
146	Bool_t hasTZEROAND = hasTZEROA && hasTZEROC;
147	/* TZEROTOF */
148	Double_t TZEROTOF = analysisEvent->GetTimeZeroTOF()[9];
149	Bool_t hasTZEROTOF = analysisEvent->GetTimeZeroTOFSigma()[9] < 150.;
150	TZEROTOF += TOF_TZEROTOF_mean - TOF_mean;
151
152	if (hasTZEROA) {
153	hCentrality_TZEROA->Fill(cent, TZEROA);
154	hVertex_TZEROA->Fill(vertexz, TZEROA);
155	}
156	if (hasTZEROC) {
157	hCentrality_TZEROC->Fill(cent, TZEROC);
158	hVertex_TZEROC->Fill(vertexz, TZEROC);
159	}
160	if (hasTZEROAND) {
161	hCentrality_TZEROAND->Fill(cent, TZEROAND);
162	hVertex_TZEROAND->Fill(vertexz, TZEROAND);
163	}
164	if (hasTZEROTOF) {
165	hCentrality_TZEROTOF->Fill(cent, TZEROTOF);
166	hVertex_TZEROTOF->Fill(vertexz, TZEROTOF);
167	}
168
169	if (hasTZEROA && hasTZEROC && hasTZEROAND && hasTZEROTOF) {
170	hCentrality_TZEROTOF_TZEROA->Fill(cent, TZEROTOF - TZEROA);
171	hCentrality_TZEROTOF_TZEROC->Fill(cent, TZEROTOF - TZEROC);
172	hCentrality_TZEROTOF_TZEROAND->Fill(cent, TZEROTOF - TZEROAND);
173	hCentrality_TZEROA_TZEROC->Fill(cent, TZEROA - TZEROC);
174	}
175
176	/* loop over tracks */
177	for (Int_t itrk = 0; itrk < analysisTrackArray->GetEntries(); itrk++) {
178	/* get track */
179	AliAnalysisTrack analysisTrack = (AliAnalysisTrack )analysisTrackArray->At(itrk);
180	if (!analysisTrack) continue;
181	/* check accepted track */
182	if (!analysisTrack->AcceptTrack()) continue;
183	/* check momentum */
184	if (analysisTrack->GetP() < 0.9 \|\| analysisTrack->GetP() > 1.1) continue;
185	/* check TOF pid */
186	if (!analysisTrack->HasTOFPID()) continue;
187
188	TZEROTOF = analysisEvent->GetTimeZeroTOF(analysisTrack->GetP());
189	TZEROTOF += TOF_TZEROTOF_mean - TOF_mean;
190	hasTZEROTOF = analysisEvent->GetTimeZeroTOFSigma(analysisTrack->GetP()) < 150.;
191
192	Double_t time = analysisTrack->GetTOFTime() - TOF_mean;
193	Double_t texp = analysisTrack->GetTOFExpTime(2);
194	Double_t deltat = time - texp;
195
196
197	hCentrality_TOF->Fill(cent, deltat);
198	hVertex_TOF->Fill(vertexz, deltat);
199	if (hasTZEROA) {
200	hCentrality_TOF_TZEROA->Fill(cent, deltat - TZEROA);
201	hVertex_TOF_TZEROA->Fill(vertexz, deltat - TZEROA);
202	}
203	if (hasTZEROC) {
204	hCentrality_TOF_TZEROC->Fill(cent, deltat - TZEROC);
205	hVertex_TOF_TZEROC->Fill(vertexz, deltat - TZEROC);
206	}
207	if (hasTZEROAND) {
208	hCentrality_TOF_TZEROAND->Fill(cent, deltat - TZEROAND);
209	hVertex_TOF_TZEROAND->Fill(vertexz, deltat - TZEROAND);
210	}
211	if (hasTZEROTOF) {
212	hCentrality_TOF_TZEROTOF->Fill(cent, deltat - TZEROTOF);
213	hVertex_TOF_TZEROTOF->Fill(vertexz, deltat - TZEROTOF);
214	}
215
216	}
217
218	} /* end of loop over events */
219
220	/* output */
221	if (!calibfilename) TFile *fileout = TFile::Open(Form("TZEROcalib.%s", filename), "RECREATE");
222	else TFile *fileout = TFile::Open(Form("TZEROcheck.%s", filename), "RECREATE");
223
224	hCentrality->Write();
225	hCentrality_TZEROA->Write();
226	hCentrality_TZEROC->Write();
227	hCentrality_TZEROAND->Write();
228	hCentrality_TZEROTOF->Write();
229	hCentrality_TZERODIFF->Write();
230
231	hVertex_TZEROA->Write();
232	hVertex_TZEROC->Write();
233	hVertex_TZEROAND->Write();
234	hVertex_TZEROTOF->Write();
235	hVertex_TZERODIFF->Write();
236
237	hCentrality_TOF->Write();
238	hCentrality_TOF_TZEROA->Write();
239	hCentrality_TOF_TZEROC->Write();
240	hCentrality_TOF_TZEROAND->Write();
241	hCentrality_TOF_TZEROTOF->Write();
242	hCentrality_Resolution->Write();
243
244	hCentrality_TZEROTOF_TZEROA->Write();
245	hCentrality_TZEROTOF_TZEROC->Write();
246	hCentrality_TZEROTOF_TZEROAND->Write();
247	hCentrality_TZEROA_TZEROC->Write();
248
249	hVertex_TOF->Write();
250	hVertex_TOF_TZEROA->Write();
251	hVertex_TOF_TZEROC->Write();
252	hVertex_TOF_TZEROAND->Write();
253	hVertex_TOF_TZEROTOF->Write();
254
255	hResolution_TZERODIFF->Write();
256	hResolution_TOF_TZEROA->Write();
257	hResolution_TOF_TZEROC->Write();
258	hResolution_TOF_TZEROAND->Write();
259
260	hDiffDiff->Write();
261
262	fileout->Close();
263
264	if (!calibfilename) {
265	TZEROcalibration(Form("TZEROcalib.%s", filename));
266	TZEROcalib(filename, "TZEROcalibration.root");
267	TZEROresolution(Form("TZEROcheck.%s", filename));
268	}
269
270	}
271
272	TZEROcalibration(const Char_t *filename = "TZEROcalib.data.root", Float_t min = 2., Float_t max = 1.)
273	{
274
275	gROOT->LoadMacro("FitPeak.C");
276
277	TFile *fin = TFile::Open(filename);
278	TH2 hCentrality_TOF = (TH2 )fin->Get("hCentrality_TOF");
279	TH2 hCentrality_TZEROA = (TH2 )fin->Get("hCentrality_TZEROA");
280	TH2 hCentrality_TZEROC = (TH2 )fin->Get("hCentrality_TZEROC");
281
282	TH2 hCentrality_TOF_TZEROTOF = (TH2 )fin->Get("hCentrality_TOF_TZEROTOF");
283	TH2 hCentrality_TOF_TZEROA = (TH2 )fin->Get("hCentrality_TOF_TZEROA");
284	TH2 hCentrality_TOF_TZEROC = (TH2 )fin->Get("hCentrality_TOF_TZEROC");
285	TH2 hCentrality_TOF_TZEROTOF = (TH2 )fin->Get("hCentrality_TOF_TZEROTOF");
286
287	TFile *fout = TFile::Open("TZEROcalibration.root", "RECREATE");
288	FitPeak(hCentrality_TOF, 200., 2., 1., "hCentrality_TOF")->Write();
289	FitPeak(hCentrality_TZEROA, 200., 2., 1., "hCentrality_TZEROA")->Write();
290	FitPeak(hCentrality_TZEROC, 200., 2., 1., "hCentrality_TZEROC")->Write();
291
292	FitPeak(hCentrality_TOF_TZEROA, 200., 2., 1., "hCentrality_TOF_TZEROA")->Write();
293	FitPeak(hCentrality_TOF_TZEROC, 200., 2., 1., "hCentrality_TOF_TZEROC")->Write();
294	FitPeak(hCentrality_TOF_TZEROTOF, 200., 2., 1., "hCentrality_TOF_TZEROTOF")->Write();
295
296	fout->Close();
297	}
298
299	TZEROresolution(const Char_t filename, Float_t min = 2., Float_t max = 2., Char_t what = "Centrality")
300	{
301
302	TFile *fin = TFile::Open(filename);
303	hdiff = (TH2 *)fin->Get(Form("h%s_TZEROA_TZEROC", what));
304	hA = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROA", what));
305	hC = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROC", what));
306	hAND = (TH2 *)fin->Get(Form("h%s_TZEROTOF_TZEROAND", what));
307
308	gROOT->LoadMacro("FitPeak.C");
309
310
311	TH1 hdiff_sigma = (TH1 )FitPeak(hdiff, 100., min, max, "sigma1")->At(1);
312
313	TH1 hA_sigma = (TH1 )FitPeak(hA, 100., min, max, "sigma2")->At(1);
314
315	TH1 hC_sigma = (TH1 )FitPeak(hC, 100., min, max, "sigma3")->At(1);
316
317	TH1 hAND_sigma = (TH1 )FitPeak(hAND, 100., min, max, "sigma4")->At(1);
318
319	TH1 hTZEROTOF_reso = (TH1 )hdiff_sigma->Clone("hTZEROTOF_reso");
320	hTZEROTOF_reso->Reset();
321	hTZEROTOF_reso->SetMarkerStyle(20);
322	hTZEROTOF_reso->SetMarkerColor(2);
323	hTZEROTOF_reso->SetLineColor(2);
324	hTZEROTOF_reso->SetTitle("T0-TOF");
325	TH1 hTZEROA_reso = (TH1 )hdiff_sigma->Clone("hTZEROA_reso");
326	hTZEROA_reso->Reset();
327	hTZEROA_reso->SetMarkerStyle(21);
328	hTZEROA_reso->SetMarkerColor(8);
329	hTZEROA_reso->SetLineColor(8);
330	hTZEROA_reso->SetTitle("T0-TZEROA");
331	TH1 hTZEROC_reso = (TH1 )hdiff_sigma->Clone("hTZEROC_reso");
332	hTZEROC_reso->Reset();
333	hTZEROC_reso->SetMarkerStyle(22);
334	hTZEROC_reso->SetMarkerColor(4);
335	hTZEROC_reso->SetLineColor(4);
336	hTZEROC_reso->SetTitle("T0-TZEROC");
337	TH1 hTZEROAND_reso = (TH1 )hdiff_sigma->Clone("hTZEROAND_reso");
338	hTZEROAND_reso->Reset();
339	hTZEROAND_reso->SetMarkerStyle(23);
340	hTZEROAND_reso->SetMarkerColor(kYellow+1);
341	hTZEROAND_reso->SetLineColor(kYellow+1);
342	hTZEROAND_reso->SetTitle("T0-TZEROAND");
343
344	for (Int_t i = 0; i < hTZEROTOF_reso->GetNbinsX(); i++) {
345
346	Double_t diff = hdiff_sigma->GetBinContent(i+1);
347	Double_t diff_e = hdiff_sigma->GetBinError(i+1);
348	Double_t diff2 = diff*diff;
349	Double_t diff2_e = 2.diffdiff_e;
350	Double_t a = hA_sigma->GetBinContent(i+1);
351	Double_t a_e = hA_sigma->GetBinError(i+1);
352	Double_t a2 = a*a;
353	Double_t a2_e = 2.aa_e;
354	Double_t c = hC_sigma->GetBinContent(i+1);
355	Double_t c_e = hC_sigma->GetBinError(i+1);
356	Double_t c2 = c*c;
357	Double_t c2_e = 2.cc_e;
358	Double_t ac = hAND_sigma->GetBinContent(i+1);
359	Double_t ac_e = hAND_sigma->GetBinError(i+1);
360	Double_t ac2 = ac*ac;
361	Double_t ac2_e = 2.acac_e;
362
363	if (diff == 0.) continue;
364	if (a2 + c2 - diff2 < 0.) continue;
365
366	printf("a=%f, c=%f, diff=%f, ac=%f\n", a, c, diff, ac);
367
368	Double_t tofreso2 = (a2 + c2 - diff2) * 0.5;
369	Double_t tofreso2_e = (a2_e + c2_e + diff2_e) * 0.5;
370	Double_t tofreso = TMath::Sqrt(tofreso2);
371	Double_t tofreso_e = tofreso2_e / (2. * tofreso);
372	hTZEROTOF_reso->SetBinContent(i+1, tofreso);
373	hTZEROTOF_reso->SetBinError(i+1, tofreso_e);
374
375	if (a2 - tofreso*tofreso < 0.) continue;
376	if (c2 - tofreso*tofreso < 0.) continue;
377	if (ac2 - tofreso*tofreso < 0.) continue;
378
379	Double_t tzeroareso2 = a2 - tofreso2;
380	Double_t tzeroareso2_e = a2_e + tofreso2_e;
381	Double_t tzeroareso = TMath::Sqrt(tzeroareso2);
382	Double_t tzeroareso_e = tzeroareso2_e / (2. * tzeroareso);
383	hTZEROA_reso->SetBinContent(i+1, tzeroareso);
384	hTZEROA_reso->SetBinError(i+1, tzeroareso_e);
385
386	Double_t tzerocreso2 = c2 - tofreso2;
387	Double_t tzerocreso2_e = c2_e + tofreso2_e;
388	Double_t tzerocreso = TMath::Sqrt(tzerocreso2);
389	Double_t tzerocreso_e = tzerocreso2_e / (2. * tzerocreso);
390	hTZEROC_reso->SetBinContent(i+1, tzerocreso);
391	hTZEROC_reso->SetBinError(i+1, tzerocreso_e);
392
393	Double_t tzeroacreso2 = ac2 - tofreso2;
394	Double_t tzeroacreso2_e = ac2_e + tofreso2_e;
395	Double_t tzeroacreso = TMath::Sqrt(tzeroacreso2);
396	Double_t tzeroacreso_e = tzeroacreso2_e / (2. * tzeroacreso);
397	hTZEROAND_reso->SetBinContent(i+1, tzeroacreso);
398	hTZEROAND_reso->SetBinError(i+1, tzeroacreso_e);
399
400	// Double_t tzerocreso = TMath::Sqrt(c2 - tofreso*tofreso);
401	// hTZEROC_reso->SetBinContent(i+1, tzerocreso);
402
403	// Double_t tzeroacreso = TMath::Sqrt(ac2 - tofreso*tofreso);
404	// hTZEROAND_reso->SetBinContent(i+1, tzeroacreso);
405
406	printf("TOF=%f, TZEROA=%f, TZEROC=%f, TZEROAND=%f\n", tofreso, tzeroareso, tzerocreso, tzeroacreso);
407	}
408
409	hTZEROTOF_reso->Draw();
410	hTZEROA_reso->Draw("same");
411	hTZEROC_reso->Draw("same");
412	hTZEROAND_reso->Draw("same");
413
414	TFile *fout = TFile::Open("TZEROresolution.root", "RECREATE");
415	hTZEROTOF_reso->Write();
416	hTZEROA_reso->Write();
417	hTZEROC_reso->Write();
418	hTZEROAND_reso->Write();
419	fout->Close();
420	}