PWG2/FORWARD/analysis2/corrs/CompareVtxBias.C

   1 /**
   2  * Make ratio of two specific maps
   3  *
   4  * @param d        Detector
   5  * @param r        Ring
   6  * @param v        Vertex bin (1 based)
   7  * @param first    First correction
   8  * @param second   Second correction
   9  *
  10  * @return Ratio of the two, or null
  11  */
  12
  13 //____________________________________________________________________
  14 void
  15 CompareVtxBias(const char* fn1,   const char* fn2,
  16                const char* n1=0,  const char* n2=0,
  17                bool load=true)
  18 {
  19
  20   // --- Load Utilities ----------------------------------------------
  21   if (load) {
  22     gROOT->Macro("$ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/LoadLibs.C");
  23     gROOT->LoadMacro("$ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/CompareCorrs.C");
  24   }
  25
  26   // --- Get Objects -------------------------------------------------
  27   TObject* o1 = 0;
  28   TObject* o2 = 0;
  29   GetObjects(AliForwardCorrectionManager::kVertexBias, fn1, fn2, o1, o2);
  30   if (!o1 || !o2) return;
  31   AliFMDCorrVertexBias* obj1 = static_cast<AliFMDCorrVertexBias*>(o1);
  32   AliFMDCorrVertexBias* obj2 = static_cast<AliFMDCorrVertexBias*>(o2);
  33   UShort_t nVtx = obj1->GetVertexAxis().GetNbins();
  34
  35   // --- Make canvas -------------------------------------------------
  36   Canvas* c = new Canvas("vtxBiasComparison", "Ratio of vertex bias", n1, n2);
  37   c->Open();
  38
  39   // --- Loop over the data ------------------------------------------
  40   UShort_t d  = 0;
  41   UShort_t nR = 2;
  42   for (UShort_t q = 0; q < nR; q++) {
  43     Char_t   r  = (q == 0 ? 'I' : 'O');
  44     UShort_t nS = (q == 0 ?  20 :  40);
  45
  46     // --- Make 2D ratios ------------------------------------------
  47     c->Clear(nVtx, d, r);
  48     TList hists;
  49     for (UShort_t v=1; v <= nVtx; v++) {
  50       TVirtualPad* p = c->cd(v);
  51
  52       TH2* h1 = obj1->GetCorrection(r, v);
  53       TH2* h2 = obj2->GetCorrection(r, v);
  54
  55       if (!h1) {
  56         Error("CompareVtxBias",
  57               "Bias for FMD%d%c, vtxbin %3d not found in first",
  58               d, r, v);
  59         continue;
  60       }
  61       if (!h2) {
  62         Error("CompareVtxBias",
  63               "Bias for FMD%d%c, vtxbin %3d not found in second",
  64               d, r, v);
  65         continue;
  66       }
  67
  68       Double_t vl    = obj1->GetVertexAxis().GetBinLowEdge(v);
  69       Double_t vh    = obj1->GetVertexAxis().GetBinUpEdge(v);
  70       TH2*     ratio =
  71         static_cast<TH2*>(h1->Clone(Form("tmpFMD%d%c_%3d",d,r,v)));
  72       ratio->SetName(Form("FMD%d%c_vtx%03d_ratio", d, r, v));
  73       ratio->SetTitle(Form("%+5.1f<v_{z}<%-+5.1f", vl, vh));
  74       ratio->Divide(h2);
  75       ratio->SetStats(0);
  76       ratio->SetDirectory(0);
  77       ratio->SetZTitle("ratio");
  78
  79       if (ratio->GetMaximum()-ratio->GetMinimum() > 10)
  80         p->SetLogz();
  81
  82       ratio->Draw("colz");
  83       // h2->Draw("colz");
  84       hists.AddAt(ratio, v-1);
  85     }
  86     c->Print(d, r);
  87
  88     // --- Make 1D profiles ----------------------------------------
  89     c->Clear(nVtx, d, r);
  90     for (UShort_t v=1; v <= nVtx; v++) {
  91       c->cd(v);
  92       TH2* hist = static_cast<TH2*>(hists.At(v-1));
  93       TH1* prof = hist->ProjectionX();
  94       prof->Scale(1. / nS);
  95       prof->SetStats(0);
  96       prof->SetMinimum(0.8);
  97       prof->SetMaximum(1.2);
  98
  99       prof->Draw();
 100       prof->Fit("pol0","Q");
 101
 102       TF1* f = prof->GetFunction("pol0");
 103
 104       TLatex* l = new TLatex(0.5, 0.4, Form("A = %f #pm %f",
 105                                             f->GetParameter(0),
 106                                               f->GetParError(0)));
 107       l->SetTextAlign(22);
 108       l->SetNDC();
 109       l->Draw();
 110       l->DrawLatex(0.5, 0.3, Form("#chi^2/NDF = %f / %d = %f",
 111                                   f->GetChisquare(),
 112                                   f->GetNDF(),
 113                                   f->GetChisquare() / f->GetNDF()));
 114       Double_t dist = TMath::Abs(1 - f->GetParameter(0));
 115       l->DrawLatex(0.5, 0.35, Form("|1 - A| = %f %s #deltaA",
 116                                    dist, dist <= f->GetParError(0) ?
 117                                    "#leq" : ">"));
 118
 119       TLine* l1 = new TLine(-4, 1, 6, 1);
 120       l1->SetLineColor(kRed);
 121       l1->SetLineStyle(2);
 122       l1->Draw();
 123     }
 124
 125     c->Print(d, r, "profiles");
 126   }
 127
 128   // --- Close stuff -------------------------------------------------
 129   c->Close();
 130   // file1->Close();
 131   // file2->Close();
 132 }
 133
 134
 135 //____________________________________________________________________
 136 //
 137 // EOF
 138 //