Added some more scripts
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / scripts / MakeELossFit.C
CommitLineData
ab0f914c 1#ifndef __CINT__
2#include "AliForwardCorrectionManager.h"
3#include "AliFMDCorrELossFit.h"
4#include <TAxis.h>
5#include <TFile.h>
6#include <TList.h>
7#include <TError.h>
8#include <TH1.h>
9#include <TF1.h>
10#include <TClonesArray.h>
11#else
12class TAxis;
13class AliFMDCirrELossFit;
14class TH1;
15
16#endif
17
18/**
19 * Class to make correction object and save to file
20 *
21 * Run like
22 *
23 * @verbatim
24 * gROOT->Macro("$ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/LoadLibs.C");
25 * gROOT->LoadMacro("$ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/Compile.C");
26 * Compile("$ALICE_ROOT/PWG2/FORWARD/analysis2/MakeELossFit.C");
27 * MakeELossFit mef(sys, cms, field, mc, "AnalysisResults.root");
28 * mef.Run();
29 * @endverbatim
30 * where @e sys, the collision system, is one of
31 * - 1: pp
32 * - 2: PbPb
33 *
34 * @e cms is the center of mass energy per nuclean in GeV (e.g., 2760
35 * for first PbPb data), @e field is (signed) L3 magnetic in kG, and
36 * @e mc is wether this correction applies to MC data or not.
37 *
38 * The class generates a file like
39 * @verbatim
40 * elossfits<sys>_<cms>GeV_<field>kG_<realmc>.root
41 * @endverbatim
42 * in the working directory. This file can be moved to
43 * @verbatim
44 * $(ALICE_ROOT)/PWG2/FORWARD/corrections/ELossFits
45 * @endverbatim
46 * and stored in SVN for later use.
47 *
48 * @ingroup pwg2_forward_analysis_scripts
49 */
50
51class MakeELossFit
52{
53protected:
54public:
55 TList* fFitter;
56 TAxis* fAxis;
57 TClonesArray fFits;
58 UShort_t fSys;
59 UShort_t fCMS;
60 Short_t fField;
61 Bool_t fMC;
62
63 //__________________________________________________________________
64 MakeELossFit(UShort_t sys,
65 UShort_t cms,
66 Short_t field,
67 Bool_t mc,
68 const char* filename="AnalysisResults.root")
69 : fFitter(0),
70 fAxis(0),
71 fFits("AliFMDCorrELossFit::ELossFit"),
72 fSys(sys),
73 fCMS(cms),
74 fField(field),
75 fMC(mc)
76 {
77 TFile* file = TFile::Open(filename, "READ");
78 if (!file) {
79 Error("MakeELossFit", "Failed to open file %s", filename);
80 return;
81 }
82 TList* forward = static_cast<TList*>(file->Get("Forward"));
83 // static_cast<TList*>(file->Get("PWG2forwardDnDeta/Forward"));
84 if (!forward) {
85 Error("MakeELossFit", "Couldn't get forward list from %s", filename);
86 return;
87 }
88
89 fFitter = static_cast<TList*>(forward->FindObject("fmdEnergyFitter"));
90 if (!fFitter) {
91 Error("MakeELossFit", "Couldn't get fitter folder");
92 return;
93 }
94
95 fAxis = static_cast<TAxis*>(fFitter->FindObject("etaAxis"));
96 if (!fAxis) {
97 Error("MakeELossFit", "Couldn't get eta axis");
98 return;
99 }
100 file->Close();
101
102#if 0
103 AliForwardCorrectionManager& mgr = AliForwardCorrectionManager::Instance();
104 mgr.Init(sys, cms, field, mc, 0, true);
105#endif
106 }
107 //__________________________________________________________________
108 TList* GetRing(UShort_t d, Char_t r) const
109 {
110 TList* rL = static_cast<TList*>(fFitter->FindObject(Form("FMD%d%c",d,r)));
111 if (!rL) {
112 Warning("DrawFits", "List FMD%d%c not found", d, r);
113 return 0;
114 }
115 return rL;
116 }
117 //__________________________________________________________________
118 TList* GetEDists(UShort_t d, Char_t r) const
119 {
120 TList* rList = GetRing(d, r);
121 if (!rList)
122 return 0;
123 // rList->ls();
124
125 TList* edists = static_cast<TList*>(rList->FindObject("EDists"));
126 if (!edists) {
127 Error("PrintFits", "Couldn't get EDists list for FMD%d%c", d,r);
128 return 0;
129 }
130 return edists;
131 }
132 //__________________________________________________________________
133 TH1* GetEDist(UShort_t d, Char_t r, UShort_t etabin)
134 {
135 TList* eList = GetEDists(d, r);
136 if (!eList) {
137 Warning("GetEDist", "No list for FMD%d%c", d, r);
138 return 0;
139 }
140 // eList->ls();
141
142 TString cmp(Form("FMD%d%c_etabin%03d", d, r, etabin));
143
144 TIter next(eList);
145 TObject* o = 0;
146 while ((o = next())) {
147 if (!cmp.CompareTo(o->GetName())) {
148 return static_cast<TH1*>(o);
149 }
150 }
151 return 0;
152 }
153#ifndef __CINT__
154 //__________________________________________________________________
155 AliFMDCorrELossFit::ELossFit* FindBestFit(TH1* dist)
156 {
157 TList* funcs = dist->GetListOfFunctions();
158 TIter next(funcs);
159 TF1* func = 0;
160 fFits.Clear();
161 Int_t i = 0;
162 Info("FindBestFit", "%s", dist->GetName());
163 while ((func = static_cast<TF1*>(next()))) {
164 AliFMDCorrELossFit::ELossFit* fit =
165 new(fFits[i++]) AliFMDCorrELossFit::ELossFit(0,*func);
166 fit->CalculateQuality(10, .20, 1e-7);
167 }
168 fFits.Sort(false);
169 // fFits.Print();
170 return static_cast<AliFMDCorrELossFit::ELossFit*>(fFits.At(0));
171 }
172#endif
173 //__________________________________________________________________
174 Bool_t Run()
175 {
176 if (!fFitter || !fAxis) {
177 Error("Run", "Missing objects");
178 return kFALSE;
179 }
180 AliFMDCorrELossFit* obj = new AliFMDCorrELossFit;
181 obj->SetEtaAxis(*fAxis);
182
183 for (UShort_t d = 1; d <= 3; d++) {
184 Info("Run", "detector is FMD%d", d);
185 UShort_t nQ = (d == 1 ? 1 : 2);
186 for (UShort_t q = 0; q < nQ; q++) {
187 Char_t r = (q == 0 ? 'I' : 'O');
188 Int_t nBin = fAxis->GetNbins();
189 Info("Run", " ring is FMD%d%c - %d bins", d, r, nBin);
190
191 Bool_t oneSeen = kFALSE;
192 for (UShort_t b = 1; b <= nBin; b++) {
193 TH1* h = GetEDist(d, r, b);
194 if (oneSeen && !h) break;
195 if (!h) continue;
196 if (!oneSeen) oneSeen = true;
197
198 AliFMDCorrELossFit::ELossFit* best = FindBestFit(h);
199 best->Print("");
200 best->fDet = d;
201 best->fRing = r;
202 best->fBin = b;
203 // Double_t eta = fAxis->GetBinCenter(b);
204 Info("Run", " bin=%3d->%6.4f", b, eta);
205 obj->SetFit(d, r, b, new AliFMDCorrELossFit::ELossFit(*best));
206 }
207 }
208 }
209 AliForwardCorrectionManager& mgr = AliForwardCorrectionManager::Instance();
210 TString fname(mgr.GetFileName(AliForwardCorrectionManager::kELossFits,
211 fSys, fCMS, fField, fMC));
212 TFile* output = TFile::Open(fname.Data(), "RECREATE");
213 if (!output) {
214 Warning("Run", "Failed to open output file %s", fname.Data());
215 return kFALSE;
216 }
217 obj->Write(mgr.GetObjectName(AliForwardCorrectionManager::kELossFits));
218 output->Write();
219 output->Close();
220 Info("Run", "File %s created. It should be copied to %s and stored in SVN",
221 fname.Data(),mgr.GetFileDir(AliForwardCorrectionManager::kELossFits));
222
223 return kTRUE;
224 }
225};
226//
227// EOF
228//