Renames and new scripts
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / scripts / DrawCorrELoss.C
CommitLineData
1c762251 1/**
2 * @file
3 *
4 * Scripts to draw energy loss fits from correction object file
5 *
6 * @ingroup pwg2_forward_analysis_scripts
7 */
8/**
9 * Clear canvas
10 *
11 * @param c Canvas to clear
12 *
13 * @ingroup pwg2_forward_analysis_scripts
14 */
cbd6464b 15void
16ClearCanvas(TCanvas* c)
17{
18 c->SetLeftMargin(.1);
19 c->SetRightMargin(.05);
20 c->SetBottomMargin(.1);
21 c->SetTopMargin(.05);
22 c->Clear();
23}
24
1c762251 25/**
26 * Draw energy loss fits to a multi-page PDF.
27 *
28 * @par Input:
29 * The input file is expected to contain a AliFMDCorrELossFit object
0f84fefb 30 * named @i elossfits in the top level directory.
1c762251 31 *
32 * @para Output:
33 * A multi-page PDF. Note, that the PDF generated by ROOT in this way
34 * is broken (cannot be read by Acrobat Reader on Windows and MacOSX)
35 * and one should pass it through a filter to correct these problems.
36 *
37 * @param fname File name
38 * @param option Drawing options
39 *
40 * @ingroup pwg2_forward_analysis_scripts
41 */
cbd6464b 42void
0f84fefb 43DrawCorrELoss(const char* fname, const char* option="err")
cbd6464b 44{
1c762251 45 //__________________________________________________________________
46 // Load libraries and object
cbd6464b 47 gROOT->Macro("$ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/LoadLibs.C");
48
49 TFile* file = TFile::Open(fname, "READ");
50 if (!file) {
51 Error("DrawELossFits", "Failed to open %s", fname);
52 return;
53 }
54 TString pname(fname);
55 pname.ReplaceAll(".root", ".pdf");
56
57 AliFMDCorrELossFit* fits =
58 static_cast<AliFMDCorrELossFit*>(file->Get("elossfits"));
59 if (!fits) {
60 Error("DrawELossFits", "Object 'elossfits' not found in %s", fname);
61 return;
62 }
63
1c762251 64 //__________________________________________________________________
65 // Create a canvas
cbd6464b 66 TCanvas* c = new TCanvas("c", "c", 800 / TMath::Sqrt(2), 800);
67 c->SetFillColor(0);
68 c->SetBorderSize(0);
69 c->SetBorderMode(0);
70 c->Print(Form("%s[", pname.Data()));
71
72 gStyle->SetOptStat(0);
73 gStyle->SetTitleColor(0);
74 gStyle->SetTitleStyle(0);
75 gStyle->SetTitleBorderSize(0);
76 gStyle->SetTitleX(.7);
77 gStyle->SetTitleY(1);
78 gStyle->SetTitleW(.3);
79 gStyle->SetTitleH(.1);
80 gStyle->SetFrameFillColor(kWhite);
81 gStyle->SetFrameBorderSize(1);
82 gStyle->SetFrameBorderMode(1);
83
84 ClearCanvas(c);
1c762251 85 //__________________________________________________________________
86 // Create a title page
cbd6464b 87 TLatex* ll = new TLatex(.5,.8, fname);
88 ll->SetTextAlign(22);
89 ll->SetTextSize(0.05);
90 ll->SetNDC();
91 ll->Draw();
92
93 TLatex* l = new TLatex(.5,.8, fname);
94 l->SetNDC();
95 l->SetTextSize(0.03);
96 l->SetTextFont(132);
97 l->SetTextAlign(12);
98 l->DrawLatex(0.2, 0.70, "1^{st} page is a summary of fit parameters");
99 l->DrawLatex(0.2, 0.67, "2^{nd} page is a summary of relative errors");
100 l->DrawLatex(0.2, 0.64, "Subsequent pages shows the fitted functions");
101 l->DrawLatex(0.3, 0.60, "Black line is the full fitted function");
102 l->DrawLatex(0.3, 0.57, "Coloured lines are the individual N-mip components");
103 l->DrawLatex(0.3, 0.54, "Full drawn lines correspond to used components");
104 l->DrawLatex(0.3, 0.51, "Dashed lines correspond to ignored components");
105 l->DrawLatex(0.2, 0.47, "Each component has the form");
106 l->DrawLatex(0.3, 0.42, "f_{n}(x; #Delta, #xi, #sigma') = "
107 "#int_{-#infty}^{+#infty}d#Delta' "
108 "landau(x; #Delta, #xi)gaus(x; #Delta', #sigma')");
109 l->DrawLatex(0.2, 0.37, "The full function is given by");
110 l->DrawLatex(0.3, 0.32, "f_{N}(x; #Delta, #xi, #sigma', #bf{a}) = "
111 "#sum_{i=1}^{N} a_{i} "
112 "f_{i}(x; #Delta_{i}, #xi_{i}, #sigma_{i}')");
113 l->DrawLatex(0.3, 0.26, "#Delta_{i} = i (#Delta_{1} + #xi_{1} log(i))");
114 l->DrawLatex(0.3, 0.23, "#xi_{i} = i #xi_{1}");
115 l->DrawLatex(0.3, 0.20, "#sigma_{i} = #sqrt{i} #sigma_{1}");
116 l->DrawLatex(0.3, 0.17, "#sigma_{n} #dot{=} 0");
117 l->DrawLatex(0.3, 0.14, "#sigma'^{2} = #sigma^{2}_{n} + #sigma^{2}");
118 l->DrawLatex(0.3, 0.11, "a_{1} = 1");
119 c->Print(pname.Data(), "Title:Title page");
120
121 ClearCanvas(c);
122
1c762251 123 //__________________________________________________________________
124 // Draw overview page
cbd6464b 125 fits->Draw(option);
126 c->Print(pname.Data(), "Title:Fit overview");
127
128 ClearCanvas(c);
129
1c762251 130 //__________________________________________________________________
131 // Draw relative parameter errors
cbd6464b 132 fits->Draw("rel");
133 c->Print(pname.Data(), "Title:Relative parameter errors");
134
1c762251 135 //__________________________________________________________________
136 // Draw all fits individually
cbd6464b 137 Int_t nPad = 6;
138 for (UShort_t d=1; d<=3; d++) {
139 UShort_t nQ = (d == 1 ? 1 : 2);
140 for (UShort_t q = 0; q < nQ; q++) {
141 Char_t r = (q == 0 ? 'I' : 'O');
142
143 TObjArray* ra = fits->GetRingArray(d, r);
144 if (!ra) continue;
145
146 AliFMDCorrELossFit::ELossFit* fit = 0;
147 TIter next(ra);
148 Int_t i = 0;
149 while ((fit = static_cast<AliFMDCorrELossFit::ELossFit*>(next()))) {
150 if ((i % nPad) == 0) {
151 ClearCanvas(c);
152 c->Divide(2,nPad/2,0,0);
153 }
154 TVirtualPad* p = c->cd((i % nPad) + 1);
155 p->SetLogy();
156 p->SetGridx();
157 p->SetGridy();
158 p->SetFillColor(kWhite);
159 fit->Draw("comp");
160
161 if ((i % nPad) == (nPad-1))
162 c->Print(pname.Data(),
163 Form("Title:FMD%d%c page %d", d, r, (i/nPad)+1));
164 i++;
165 }
166 if (i % nPad != 0)
167 c->Print(pname.Data(),
168 Form("Title:FMD%d%c page %d", d, r, (i/nPad)+1));
169 }
170 }
171
1c762251 172 //__________________________________________________________________
173 // Close output file
cbd6464b 174 c->Print(Form("%s]", pname.Data()));
175}
1c762251 176//
177// EOF
178//