PWG2/FORWARD/analysis2/corrs/CompareSecMaps.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 CompareSecMaps(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::kSecondaryMap, fn1, fn2, o1, o2);
  30   if (!o1 || !o2) return;
  31   AliFMDCorrSecondaryMap* obj1 = static_cast<AliFMDCorrSecondaryMap*>(o1);
  32   AliFMDCorrSecondaryMap* obj2 = static_cast<AliFMDCorrSecondaryMap*>(o2);
  33   UShort_t nVtx = obj1->GetVertexAxis().GetNbins();
  34
  35   // --- Make canvas -------------------------------------------------
  36   Canvas* c = new Canvas("secMapComparison", "Ratio of secondary maps", n1, n2);
  37   c->Open();
  38
  39   // --- Loop over the data ------------------------------------------
  40   for (UShort_t d = 1; d <= 3; d++) {
  41     UShort_t nR = (d == 1 ? 1 : 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(d, r, v);
  53         TH2* h2 = obj2->GetCorrection(d, r, v);
  54
  55         if (!h1) {
  56           Error("CompareSecMaps",
  57                 "Map for FMD%d%c, vtxbin %3d not found in first",
  58                 d, r, v);
  59           continue;
  60         }
  61         if (!h2) {
  62           Error("CompareSecMaps",
  63                 "Map 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         hists.AddAt(ratio, v-1);
  84       }
  85       c->Print(d, r);
  86
  87       // --- Make 1D profiles ----------------------------------------
  88       c->Clear(nVtx, d, r);
  89       for (UShort_t v=1; v <= nVtx; v++) {
  90         c->cd(v);
  91         TH2* hist = static_cast<TH2*>(hists.At(v-1));
  92         TH1* prof = hist->ProjectionX();
  93         prof->Scale(1. / nS);
  94         prof->SetStats(0);
  95         prof->SetMinimum(0.8);
  96         prof->SetMaximum(1.2);
  97
  98         prof->Draw();
  99         prof->Fit("pol0","Q");
 100
 101         TF1* f = prof->GetFunction("pol0");
 102
 103         TLatex* l = new TLatex(0.5, 0.4, Form("A = %f #pm %f",
 104                                               f->GetParameter(0),
 105                                               f->GetParError(0)));
 106         l->SetTextAlign(22);
 107         l->SetNDC();
 108         l->Draw();
 109         l->DrawLatex(0.5, 0.3, Form("#chi^2/NDF = %f / %d = %f",
 110                                     f->GetChisquare(),
 111                                     f->GetNDF(),
 112                                     f->GetChisquare() / f->GetNDF()));
 113         Double_t dist = TMath::Abs(1 - f->GetParameter(0));
 114         l->DrawLatex(0.5, 0.35, Form("|1 - A| = %f %s #deltaA",
 115                                      dist, dist <= f->GetParError(0) ?
 116                                      "#leq" : ">"));
 117
 118         TLine* l1 = new TLine(-4, 1, 6, 1);
 119         l1->SetLineColor(kRed);
 120         l1->SetLineStyle(2);
 121         l1->Draw();
 122       }
 123
 124       c->Print(d, r, "profiles");
 125     }
 126   }
 127
 128   // --- Close stuff -------------------------------------------------
 129   c->Close();
 130   // file1->Close();
 131   // file2->Close();
 132 }
 133
 134
 135 //____________________________________________________________________
 136 //
 137 // EOF
 138 //