PWGLF/FORWARD/analysis2/tests/TestFitELoss.C

   1 /**
   2  * Test script to fit the energy loss spectra
   3  *
   4  * @deprecated
   5  * This is a simple test script
   6  *
   7  * @ingroup pwg2_forward_scripts_tests
   8  */
   9 #ifndef __CINT__
  10 # include "AliForwardUtil.h"
  11 # include <TFile.h>
  12 # include <TList.h>
  13 # include <TH1.h>
  14 # include <TF1.h>
  15 # include <TFitResult.h>
  16 # include <TMath.h>
  17 # include <TStyle.h>
  18 # include <TArrow.h>
  19 # include <TCanvas.h>
  20 #else
  21 class TCanvas;
  22 class TFile;
  23 class TH1;
  24 class TList;
  25 class TF1;
  26 #endif
  27
  28 //____________________________________________________________________
  29 /**
  30  *
  31  *
  32  * @param ef
  33  * @param d
  34  * @param r
  35  * @param etabin
  36  *
  37  * @return
  38  *
  39  * @deprecated
  40  * This is a simple test script
  41  *
  42  * @ingroup pwg2_forward_scripts_tests
  43  */
  44 TH1* GetEDist(TList* ef, UShort_t d, Char_t r, UShort_t etabin)
  45 {
  46   TList* dL = static_cast<TList*>(ef->FindObject(Form("FMD%d%c",d,r)));
  47   if (!dL) {
  48     Error("GetEDist", "Couldn't find list FMD%d%c",d,r);
  49     ef->ls();
  50     return 0;
  51   }
  52   if (etabin > 999) {
  53     TH1* hist = static_cast<TH1*>(dL->FindObject(Form("FMD%d%c_edist",d,r)));
  54     if (hist) {
  55       Error("GetEDist", "Couldn't find EDists histogram for FMD%d%c",d,r);
  56       return 0;
  57     }
  58   }
  59
  60   TList* edL = static_cast<TList*>(dL->FindObject("EDists"));
  61   if (!edL) {
  62     Error("GetEDist", "Couldn't find list EDists for FMD%d%c",d,r);
  63     return 0;
  64   }
  65
  66   TH1* hist = static_cast<TH1*>(edL->FindObject(Form("FMD%d%c_etabin%03d",
  67                                                      d, r, etabin)));
  68   if (!hist) {
  69     Error("GetEDist", "Couldn't find histogra FMD%d%c_etabin%03d",
  70           d,r, etabin);
  71     return 0;
  72   }
  73
  74   return hist;
  75 }
  76
  77 //____________________________________________________________________
  78 /**
  79  *
  80  *
  81  * @param ef
  82  * @param d
  83  * @param r
  84  * @param eta
  85  *
  86  * @return
  87  *
  88  * @deprecated
  89  * This is a simple test script
  90  *
  91  * @ingroup pwg2_forward_scripts_tests
  92 */
  93 TH1* GetEDist(TList* ef, UShort_t d, Char_t r, Float_t eta)
  94 {
  95   if (!ef) {
  96     Error("GetEDist", "EF not set");
  97     return 0;
  98   }
  99   TAxis* etaAxis = static_cast<TAxis*>(ef->FindObject("etaAxis"));
 100   if (!etaAxis) {
 101     Error("GetEDist", "Couldn't find eta axis in list");
 102     return 0;
 103   }
 104
 105   UShort_t bin = etaAxis->FindBin(eta);
 106   if (bin <= 0 || bin > etaAxis->GetNbins()) {
 107     Error("GetEDist", "eta=%f out of range [%f,%f] - getting ring histo",
 108           eta, etaAxis->GetXmin(), etaAxis->GetXmax());
 109     return GetEDist(ef, d, r, UShort_t(1000));
 110   }
 111
 112   return GetEDist(ef, d, r, bin);
 113 }
 114
 115 //____________________________________________________________________
 116 /**
 117  *
 118  *
 119  * @param file
 120  *
 121  * @return
 122  *
 123  * @deprecated
 124  * This is a simple test script
 125  *
 126  * @ingroup pwg2_forward_scripts_tests
 127  */
 128 TList* GetEF(TFile* file)
 129 {
 130   TList* forward = static_cast<TList*>(file->Get("PWG2forwardDnDeta/Forward"));
 131   if (!forward) {
 132     Error("GetEF", "Failed to get list PWG2forwardDnDeta/Forward from file");
 133     return 0;
 134   }
 135   TList* ef = static_cast<TList*>(forward->FindObject("fmdEnergyFitter"));
 136   if (!ef) {
 137     Error("GetEF", "Failed to get energy fitter list");
 138     return 0;
 139   }
 140
 141   return ef;
 142 }
 143
 144 //____________________________________________________________________
 145 TList* ef = 0;
 146
 147 //____________________________________________________________________
 148 /**
 149  *
 150  *
 151  *
 152  * @return
 153  *
 154  * @deprecated
 155  * This is a simple test script
 156  *
 157  * @ingroup pwg2_forward_scripts_tests
 158  */
 159 TList*  CheckEF()
 160 {
 161   if (ef) return ef;
 162
 163   TFile* file = TFile::Open("AnalysisResults.root", "READ");
 164   if (!file) {
 165     Error("Fit1", "Failed to open file");
 166     return 0;
 167   }
 168   return GetEF(file);
 169 }
 170
 171 //____________________________________________________________________
 172 TCanvas* c = 0;
 173
 174 //____________________________________________________________________
 175 /**
 176  *
 177  *
 178  *
 179  * @return
 180  *
 181  * @deprecated
 182  * This is a simple test script
 183  *
 184  * @ingroup pwg2_forward_scripts_tests
 185  */
 186 TCanvas* CheckC()
 187 {
 188   if (c) return c;
 189
 190   gStyle->SetOptFit(1111);
 191   gStyle->SetCanvasColor(0);
 192   gStyle->SetCanvasBorderSize(0);
 193   gStyle->SetCanvasBorderMode(0);
 194   gStyle->SetCanvasDefW(800);
 195   gStyle->SetCanvasDefH(800);
 196   gStyle->SetPadTopMargin(0.05);
 197   gStyle->SetPadRightMargin(0.05);
 198   gStyle->SetTitleBorderSize(0);
 199   gStyle->SetTitleFillColor(0);
 200   gStyle->SetTitleStyle(0);
 201   gStyle->SetStatBorderSize(1);
 202   gStyle->SetStatColor(0);
 203   gStyle->SetStatStyle(0);
 204   gStyle->SetStatX(0.95);
 205   gStyle->SetStatY(0.95);
 206   gStyle->SetStatW(0.15);
 207   gStyle->SetStatH(0.15);
 208
 209   c = new TCanvas("c", "c");
 210   c->SetLogy();
 211
 212   return c;
 213 }
 214
 215 //____________________________________________________________________
 216 /**
 217  *
 218  *
 219  * @param f
 220  *
 221  * @deprecated
 222  * This is a simple test script
 223  *
 224  * @ingroup pwg2_forward_scripts_tests
 225  */
 226 void PrintFit(TF1* f)
 227 {
 228   Int_t    nu     = f->GetNDF();
 229   Double_t chi2   = f->GetChisquare();
 230   Double_t chi2nu = (nu > 0 ? chi2/nu : 0);
 231   Printf("%s: chi^2/nu=%f/%d=%f [%f,%f]",
 232          f->GetName(), chi2, nu, chi2nu,
 233          f->GetXmin(), f->GetXmax());
 234   for (Int_t i = 0; i < f->GetNpar(); i++) {
 235     Double_t v = f->GetParameter(i);
 236     Double_t e = f->GetParError(i);
 237     Double_t r = 100*(v == 0 ? 1 : e / v);
 238     Printf("%32s = %14.7f +/- %14.7f (%5.1f%%)",f->GetParName(i),v,e,r);
 239   }
 240 }
 241
 242 //____________________________________________________________________
 243 /**
 244  *
 245  *
 246  * @param n
 247  * @param d
 248  * @param r
 249  * @param eta
 250  *
 251  * @deprecated
 252  * This is a simple test script
 253  *
 254  * @ingroup pwg2_forward_scripts_tests
 255  */
 256 void TestFitELoss(Int_t n, UShort_t d, Char_t r, Float_t eta)
 257 {
 258   TList* ef1 = CheckEF();
 259   TCanvas* c1 = CheckC();
 260   if (!ef1) return;
 261   if (!c1)  return;
 262
 263   TH1* dist = GetEDist(ef1, d, r, eta);
 264   if (!dist) return;
 265
 266   AliForwardUtil::ELossFitter f(0.4, 10, 4);
 267   TF1* landau1 = new TF1(*f.Fit1Particle(dist, 0));
 268   landau1->SetName("Landau1");
 269   PrintFit(landau1);
 270   TF1* landaun = new TF1(*f.FitNParticle(dist, n, 0));
 271   landau1->SetName(Form("Landau%d", n));
 272   PrintFit(landaun);
 273   landau1->SetRange(0,10);
 274   landaun->SetRange(0,10);
 275   landau1->SetLineWidth(4);
 276   landaun->SetLineWidth(4);
 277   landau1->SetLineColor(kBlack);
 278   landaun->SetLineColor(kBlack);
 279
 280   dist->GetListOfFunctions()->Clear();
 281   dist->GetListOfFunctions()->Add(landau1);
 282   dist->GetListOfFunctions()->Add(landaun);
 283   dist->GetListOfFunctions()->ls();
 284   dist->Draw();
 285   landau1->Draw("same");
 286   landaun->Draw("same");
 287
 288   Double_t mp = landaun->GetParameter(1);
 289   Double_t xi = landaun->GetParameter(2);
 290   for (Int_t i  = 1; i <= n; i++) {
 291     Double_t x  = i * (mp + xi * TMath::Log(i));
 292     Double_t y  = landaun->Eval(x);
 293     Double_t y1 = y < 0.05 ? 1 : 0.01;
 294     TArrow* a = new TArrow(x,y1,x,y,0.03,"|>");
 295     Info("FitSteer", "Delta_{p,%d}=%f", i, x);
 296     a->SetLineWidth(2);
 297     a->SetAngle(30);
 298     a->Draw();
 299   }
 300
 301   c1->cd();
 302 }
 303 //
 304 // EOF
 305 //