]>
Commit | Line | Data |
---|---|---|
ab886e13 | 1 | |
2 | /************************************************************************** | |
3 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
4 | * * | |
5 | * Author: The ALICE Off-line Project. * | |
6 | * Contributors are mentioned in the code where appropriate. * | |
7 | * * | |
8 | * Permission to use, copy, modify and distribute this software and its * | |
9 | * documentation strictly for non-commercial purposes is hereby granted * | |
10 | * without fee, provided that the above copyright notice appears in all * | |
11 | * copies and that both the copyright notice and this permission notice * | |
12 | * appear in the supporting documentation. The authors make no claims * | |
13 | * about the suitability of this software for any purpose. It is * | |
14 | * provided "as is" without express or implied warranty. * | |
15 | **************************************************************************/ | |
16 | ||
ba06aaec | 17 | //------------------------------------------------------------------------------ |
18 | // Implementation of AliPerfAnalyzeInvPt class. It analyzes the output of | |
ab886e13 | 19 | // AliPerformancePtCalib.cxx and AliPerformancePtCalibMC.cxx: via |
20 | // AliPerformancePtCalib::Analyse() or AliPerformancePtCalibMC::Analyse() | |
21 | // | |
22 | // Projection of charge/pt vs theta and vs phi resp. histoprams will be fitted with either | |
23 | // polynomial or gaussian fit function to extract minimum position of charge/pt. | |
ba06aaec | 24 | // Fit options and theta, phi bins can be set by user. |
25 | // Attention: use the Set* functions of AliPerformancePtCalib.h and AliPerformancePtCalibMC.h | |
26 | // when running AliPerformancePtCalib*::Analyse() | |
27 | // The result of the analysis (histograms/graphs) are stored in the folder which is | |
28 | // a data member of AliPerformancePtCalib*. | |
29 | // | |
30 | // Author: S.Schuchmann 11/13/2009 | |
31 | //------------------------------------------------------------------------------ | |
32 | ||
33 | /* | |
34 | ||
ab886e13 | 35 | // after running the performance task of AliPerformancePtCalib* or AliPerformancePtCalibMC* |
36 | // read the file, and get component: | |
ba06aaec | 37 | |
38 | TFile f("Output.root"); | |
ab886e13 | 39 | AliPerformancePtCalib *compObj = (AliPerformancePtCalib*)coutput->FindObject("AliPerformancePtCalib"); |
40 | //or | |
41 | AliPerformancePtCalibMC *compObj = (AliPerformancePtCalibMC*)coutput->FindObject("AliPerformancePtCalibMC"); | |
42 | ||
43 | // set phi and theta bins for fitting and analyse comparison data | |
44 | compObj->SetProjBinsTheta(thetaBins,nThetaBins,minPhi, maxPhi); | |
45 | compObj->SetProjBinsPhi(phiBins,nPhiBins,minTheta,maxTheta); | |
46 | compObj->SetMakeFitOption(kFALSE,exclRange,fitRange); | |
47 | compObj->SetDoRebin(rebin); | |
ba06aaec | 48 | compObj->Analyse(); |
ab886e13 | 49 | //for details see functions of class AliPerformancePtCalib* |
ba06aaec | 50 | |
51 | // the output histograms/graphs will be stored in the folder "folderRes" | |
52 | compObj->GetAnalysisFolder()->ls("*"); | |
53 | ||
54 | // user can save whole comparison object (or only folder with anlysed histograms) | |
55 | // in the seperate output file (e.g.) | |
56 | TFile fout("Analysed_InvPt.root","recreate"); | |
57 | compObj->Write(); // compObj->GetAnalysisFolder()->Write(); | |
58 | fout.Close(); | |
59 | ||
60 | */ | |
61 | ||
acb9d358 | 62 | #include "TNamed.h" |
ba06aaec | 63 | #include "TF1.h" |
acb9d358 | 64 | #include "TH1F.h" |
65 | #include "TH2F.h" | |
ba06aaec | 66 | #include "TH1D.h" |
ab886e13 | 67 | #include "TMath.h" |
acb9d358 | 68 | #include "TGraphErrors.h" |
ba06aaec | 69 | #include "TCanvas.h" |
70 | #include "TObjArray.h" | |
acb9d358 | 71 | |
acb9d358 | 72 | #include "AliPerfAnalyzeInvPt.h" |
73 | ||
74 | ClassImp(AliPerfAnalyzeInvPt) | |
75 | ||
76 | ||
77 | // fit functions | |
78 | //____________________________________________________________________________________________________________________________________________ | |
ba06aaec | 79 | Double_t AliPerfAnalyzeInvPt::Polynomial(Double_t *x, const Double_t *par) |
acb9d358 | 80 | { |
ba06aaec | 81 | // fit function for fitting of 1/pt with a polynomial of 4th order, rejecting points between +/-par[4] |
acb9d358 | 82 | |
6bcacaf1 | 83 | |
acb9d358 | 84 | if (x[0] > -par[4] && x[0] < par[4]) { |
85 | TF1::RejectPoint(); | |
86 | return 0; | |
87 | } | |
6bcacaf1 | 88 | |
acb9d358 | 89 | return par[2]+par[0]*pow((x[0]-par[1]),2)+par[3]*pow((x[0]-par[1]),4); |
6bcacaf1 | 90 | |
acb9d358 | 91 | } |
92 | //____________________________________________________________________________________________________________________________________________ | |
ba06aaec | 93 | Double_t AliPerfAnalyzeInvPt::PolynomialRejP(Double_t *x, const Double_t *par) |
acb9d358 | 94 | { |
ba06aaec | 95 | // fit function for fitting of 1/pt with a polynomial of 4th order to improve result (fit range and rejection zone is adjusted to first guess of minimum position) |
6bcacaf1 | 96 | |
acb9d358 | 97 | Double_t pos = par[5]; |
98 | Double_t neg = - par[5]; | |
99 | pos += par[4]; | |
100 | neg += par[4]; | |
6bcacaf1 | 101 | |
acb9d358 | 102 | if (x[0] > neg && x[0] < pos) { |
103 | TF1::RejectPoint(); | |
104 | return 0; | |
105 | } | |
6bcacaf1 | 106 | |
acb9d358 | 107 | return par[2]+par[0]*pow((x[0]-par[1]),2)+par[3]*pow((x[0]-par[1]),4); |
6bcacaf1 | 108 | |
acb9d358 | 109 | } |
110 | //____________________________________________________________________________________________________________________________________________ | |
ba06aaec | 111 | Double_t AliPerfAnalyzeInvPt::InvGauss(Double_t *x, const Double_t *par) |
acb9d358 | 112 | { |
ba06aaec | 113 | // fit function for fitting of 1/pt with gaussian, rejecting points between +/-par[4] |
acb9d358 | 114 | if (x[0] > -par[6] && x[0] < par[6]) { |
115 | TF1::RejectPoint(); | |
116 | return 0; | |
117 | } | |
118 | ||
ba06aaec | 119 | return par[3]+par[0]*(TMath::Exp(-0.5*(TMath::Power((x[0]-par[1])/par[2], 2.0)))+par[4]*pow((x[0]-par[1]),2)+par[5]*pow((x[0]-par[1]),4)) ; |
acb9d358 | 120 | } |
121 | //____________________________________________________________________________________________________________________________________________ | |
ba06aaec | 122 | Double_t AliPerfAnalyzeInvPt::InvGaussRejP(Double_t *x, const Double_t *par) |
acb9d358 | 123 | { |
ba06aaec | 124 | // fit function for fitting of 1/pt with gaussian to improve result (fit range and rejection zone is adjusted to first guess of minimum position) |
acb9d358 | 125 | Double_t pos = par[7];//0.12; |
126 | Double_t neg = - par[7];//0.12; | |
127 | pos += par[6]; | |
128 | neg += par[6]; | |
129 | ||
130 | if (x[0] > neg && x[0] < pos) { | |
131 | TF1::RejectPoint(); | |
132 | return 0; | |
133 | } | |
134 | ||
135 | ||
136 | return par[3]+par[0]*TMath::Exp(-0.5*(TMath::Power((x[0]-par[1])/par[2], 2.0)))+par[4]*pow((x[0]-par[1]),2)+par[5]*pow((x[0]-par[1]),4) ; | |
137 | } | |
138 | ||
139 | ||
140 | //_____________________________________________________________________________________________________________________________________________ | |
141 | AliPerfAnalyzeInvPt::AliPerfAnalyzeInvPt(): | |
142 | TNamed("AliPerfAnalyzeInvPt","AliPerfAnalyzeInvPt"), | |
143 | fNThetaBins(0), | |
144 | fNPhiBins(0), | |
145 | fRange(0), | |
146 | fExclRange(0), | |
147 | fFitGaus(0) , | |
6bcacaf1 | 148 | fDoRebin(0), |
149 | fRebin(0), | |
ab886e13 | 150 | //histograms,graphs and functions |
acb9d358 | 151 | fHistH2InvPtTheta(0), |
152 | fHistH2InvPtPhi(0), | |
153 | fGrMinPosTheta(0), | |
154 | fGrMinPosPhi(0), | |
155 | fFitMinPos(0), | |
156 | fFitMinPosRejP(0), | |
157 | fFitInvGauss(0), | |
158 | fFitInvGaussRejP(0) | |
159 | { | |
160 | // Default constructor | |
161 | ||
ab886e13 | 162 | fFitGaus = kFALSE;// flag for gaussian fit |
163 | fNThetaBins = 0; // theta bins for projections | |
164 | fNPhiBins = 0; //phi bins for projections | |
165 | fRange = 0; //fit range | |
166 | fExclRange = 0; //fit exclusion range | |
167 | fDoRebin = kFALSE; // flag for rebin | |
168 | fRebin = 0; // bins for rebin | |
acb9d358 | 169 | |
acb9d358 | 170 | for(Int_t i=0;i<100;i++){ |
171 | ||
172 | fHistFitTheta[i] = NULL; | |
173 | fHistFitPhi[i] = NULL; | |
174 | } | |
ba06aaec | 175 | |
176 | ||
acb9d358 | 177 | } |
178 | //_____________________________________________________________________________________________________________________________________________ | |
c6d45264 | 179 | AliPerfAnalyzeInvPt::AliPerfAnalyzeInvPt(const Char_t* name="AliAnalyzeInvPt",const Char_t* title="AliAnalyzeInvPt"): |
ba06aaec | 180 | TNamed(name, title), |
acb9d358 | 181 | fNThetaBins(0), |
182 | fNPhiBins(0), | |
183 | fRange(0), | |
184 | fExclRange(0), | |
185 | fFitGaus(0) , | |
6bcacaf1 | 186 | fDoRebin(0), |
187 | fRebin(0), | |
ab886e13 | 188 | //histograms,graphs and functions |
acb9d358 | 189 | fHistH2InvPtTheta(0), |
190 | fHistH2InvPtPhi(0), | |
191 | fGrMinPosTheta(0), | |
192 | fGrMinPosPhi(0), | |
193 | fFitMinPos(0), | |
194 | fFitMinPosRejP(0), | |
195 | fFitInvGauss(0), | |
196 | fFitInvGaussRejP(0) | |
ba06aaec | 197 | { |
ab886e13 | 198 | fFitGaus = kFALSE;// flag for gaussian fit |
199 | fNThetaBins = 0; // theta bins for projections | |
200 | fNPhiBins = 0; //phi bins for projections | |
201 | fRange = 0; //fit range | |
202 | fExclRange = 0; //fit exclusion range | |
203 | fDoRebin = kFALSE; // flag for rebin | |
204 | fRebin = 0; // bins for rebin | |
ba06aaec | 205 | |
ba06aaec | 206 | for(Int_t i=0;i<100;i++){ |
acb9d358 | 207 | |
ba06aaec | 208 | fHistFitTheta[i] = NULL; |
209 | fHistFitPhi[i] = NULL; | |
acb9d358 | 210 | } |
ba06aaec | 211 | |
212 | } | |
acb9d358 | 213 | |
214 | ||
ba06aaec | 215 | //______________________________________________________________________________________________________________________________________ |
216 | void AliPerfAnalyzeInvPt::InitGraphs(Double_t *binsXTheta,Double_t *fitParamTheta,Double_t *errFitParamTheta,Double_t *binsXPhi,Double_t *fitParamPhi,Double_t *errFitParamPhi){ | |
acb9d358 | 217 | |
ba06aaec | 218 | // initialize graphs for storing fit parameter |
acb9d358 | 219 | |
ba06aaec | 220 | fGrMinPosTheta = new TGraphErrors(fNThetaBins,binsXTheta,fitParamTheta,0, errFitParamTheta); |
221 | fGrMinPosTheta->SetMarkerStyle(20); | |
222 | fGrMinPosTheta->SetMarkerColor(2); | |
223 | fGrMinPosTheta->SetLineColor(2); | |
ab886e13 | 224 | fGrMinPosTheta->GetYaxis()->SetTitle("min pos (GeV/c)^{-1}"); |
ba06aaec | 225 | fGrMinPosTheta->GetXaxis()->SetTitle("#theta bin no."); |
226 | fGrMinPosTheta->GetYaxis()->SetTitleOffset(1.2); | |
ab886e13 | 227 | fGrMinPosTheta->SetTitle("#theta bins"); |
ba06aaec | 228 | |
229 | fGrMinPosPhi = new TGraphErrors(fNPhiBins,binsXPhi,fitParamPhi,0,errFitParamPhi); | |
230 | fGrMinPosPhi->SetMarkerStyle(20); | |
231 | fGrMinPosPhi->SetMarkerColor(4); | |
232 | fGrMinPosPhi->SetLineColor(4); | |
ab886e13 | 233 | fGrMinPosPhi->GetYaxis()->SetTitle("min pos (GeV/c)^{-1}"); |
ba06aaec | 234 | fGrMinPosPhi->GetXaxis()->SetTitle("#phi bin no."); |
235 | fGrMinPosPhi->GetYaxis()->SetTitleOffset(1.2); | |
ab886e13 | 236 | fGrMinPosPhi->SetTitle("#phi bins"); |
acb9d358 | 237 | } |
238 | ||
239 | //______________________________________________________________________________________________________________________________________ | |
240 | ||
241 | void AliPerfAnalyzeInvPt::InitFitFcn(){ | |
ba06aaec | 242 | // fit functions |
243 | fFitMinPos = new TF1("fFitMinPos", Polynomial,-4.0,4.0,5); | |
244 | fFitMinPos->SetLineColor(4); | |
245 | fFitMinPos->SetLineWidth(1); | |
246 | fFitMinPos->SetParameter(0,1.0); | |
247 | fFitMinPos->SetParameter(1,0.0); | |
248 | fFitMinPos->SetParameter(2,0.0); | |
249 | ||
250 | fFitMinPosRejP = new TF1("fFitMinPosRejP",PolynomialRejP,-4.0,4.0,6); | |
251 | fFitMinPosRejP->SetLineColor(2); | |
252 | fFitMinPosRejP->SetLineWidth(1); | |
253 | fFitMinPosRejP->SetParameter(0,1.0); | |
254 | fFitMinPosRejP->SetParameter(1,0.0); | |
255 | fFitMinPosRejP->SetParameter(2,0.0); | |
acb9d358 | 256 | |
257 | ||
ba06aaec | 258 | fFitInvGauss = new TF1("fFitInvGauss", InvGauss,-4.0,4.0,6); |
259 | fFitInvGauss->SetLineColor(4); | |
260 | fFitInvGauss->SetLineWidth(1); | |
261 | fFitInvGauss->SetParameter(2,1.0); | |
acb9d358 | 262 | |
ba06aaec | 263 | fFitInvGaussRejP = new TF1("fFitInvGaussRejP", InvGaussRejP,-4.0,4.0,7); |
264 | fFitInvGaussRejP->SetLineColor(2); | |
265 | fFitInvGaussRejP->SetLineWidth(1); | |
266 | fFitInvGaussRejP->SetParameter(2,1.0); | |
acb9d358 | 267 | |
ba06aaec | 268 | } |
acb9d358 | 269 | //______________________________________________________________________________________________________________________________________ |
ba06aaec | 270 | void AliPerfAnalyzeInvPt::StartAnalysis(const TH2F *histThetaInvPt, const TH2F *histPhiInvPt, TObjArray* aFolderObj){ |
ab886e13 | 271 | //start analysis: fitting charge/pt spectra |
acb9d358 | 272 | |
273 | ||
274 | ||
275 | if(!histThetaInvPt) { | |
276 | Printf("warning: no 1/pt histogram to analyse in theta bins!"); | |
277 | } | |
278 | if(!histPhiInvPt) { | |
279 | Printf("warning: no 1/pt histogram to analyse in phit bins!"); | |
280 | } | |
281 | ||
ab886e13 | 282 | Double_t thetaBins[9] = {0.77,0.97,1.17,1.37,1.57,1.77,1.97,2.17,2.37}; // default theta bins |
283 | Double_t phiBins[13] = {0.0,0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5.0,5.5,6.5}; // default phi bins | |
acb9d358 | 284 | |
285 | Int_t nThetaBins = 9; | |
286 | Int_t nPhiBins = 13; | |
287 | Double_t range = 1.0; // fit range | |
288 | Double_t exclRange = 0.13; // range of point rejection in fit | |
289 | Bool_t fitGaus=kFALSE; // fit with Gaussian or with polynomial (kFALSE) | |
290 | ||
ba06aaec | 291 | if(!fNThetaBins){ |
acb9d358 | 292 | fNThetaBins = nThetaBins; |
293 | for(Int_t k = 0;k<fNThetaBins;k++){ | |
294 | fThetaBins[k] = thetaBins[k]; | |
295 | } | |
296 | Printf("warning: for theta bins no user intput! bins are set to default values."); | |
297 | } | |
298 | ||
ba06aaec | 299 | if(!fNPhiBins){ |
acb9d358 | 300 | fNPhiBins = nPhiBins; |
301 | for(Int_t k = 0;k<fNPhiBins;k++){ | |
302 | fPhiBins[k] = phiBins[k]; | |
303 | } | |
304 | ||
305 | Printf("warning: for phi bins no user intput! bins are set to default values."); | |
306 | ||
307 | } | |
308 | ||
ba06aaec | 309 | if(!fRange){ |
acb9d358 | 310 | fRange = range; |
311 | fExclRange = exclRange; | |
312 | fFitGaus = fitGaus; | |
313 | Printf("warning: no fit range is set. default fitting conditions are used."); | |
314 | } | |
315 | ||
316 | ||
ab886e13 | 317 | //fit parameter arrays |
acb9d358 | 318 | Double_t fitParamTheta[100],fitParamPhi[100]; |
319 | Double_t errFitParamTheta[100], errFitParamPhi[100]; | |
320 | Double_t binsXTheta[100],binsXPhi[100]; | |
ab886e13 | 321 | |
322 | //copy histograms | |
ba06aaec | 323 | fHistH2InvPtTheta = (TH2F*)histThetaInvPt->Clone("invPtVsTheta"); |
324 | fHistH2InvPtPhi = (TH2F*)histPhiInvPt->Clone("invPtVsPhi"); | |
ab886e13 | 325 | |
326 | // initialize fit functions | |
acb9d358 | 327 | InitFitFcn(); |
ab886e13 | 328 | |
329 | // canvas for 2D histograms (input) | |
330 | TCanvas *thephiCan =new TCanvas("thephiCan","theta and phi vs invPt",800,400); | |
331 | thephiCan->Divide(2,1); | |
332 | thephiCan->cd(1); | |
333 | fHistH2InvPtTheta->SetTitle("#theta vs charge/pt for selected #phi range"); | |
334 | fHistH2InvPtTheta->GetXaxis()->SetTitle("charge/p_{t} (GeV/c)^{-1}"); | |
335 | fHistH2InvPtTheta->GetYaxis()->SetTitle("#theta (rad)"); | |
336 | fHistH2InvPtTheta->DrawCopy("col1z"); | |
337 | ||
338 | thephiCan->cd(2); | |
339 | fHistH2InvPtPhi->SetTitle("#phi vs charge/pt for selected #theta range"); | |
340 | fHistH2InvPtPhi->GetXaxis()->SetTitle("charge/p_{t} (GeV/c)^{-1}"); | |
341 | fHistH2InvPtPhi->GetYaxis()->SetTitle("#phi (rad)"); | |
342 | fHistH2InvPtPhi->DrawCopy("col1z"); | |
343 | ||
344 | // canvas for 1D histograms (projections) | |
acb9d358 | 345 | TCanvas *theCan =new TCanvas("theCan","invPt theta bins",1200,900); |
346 | theCan->Divide((fNThetaBins+2)/2,2); | |
347 | TCanvas *phiCan =new TCanvas("phiCan","invPt phi bins",1200,900); | |
348 | phiCan->Divide((fNPhiBins+2)/2,2); | |
ab886e13 | 349 | |
acb9d358 | 350 | Int_t countPad = 1; |
351 | ||
ab886e13 | 352 | // analyse charge/pt in bins of theta |
acb9d358 | 353 | for(Int_t i=0;i<fNThetaBins;i++){ |
ab886e13 | 354 | |
355 | //set name for each projection | |
acb9d358 | 356 | TString name = "fit_theta_"; |
357 | name +=i; | |
ab886e13 | 358 | fHistH2InvPtTheta->SetName(name.Data()); |
359 | ||
360 | //find bins for projection | |
acb9d358 | 361 | Int_t firstBin= fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i]); |
362 | if(i>0) firstBin +=1; | |
363 | Int_t lastBin = fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i+1]); | |
364 | if( i == fNThetaBins-1) { | |
365 | firstBin= fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[0]); | |
366 | lastBin = fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i]); | |
367 | } | |
ab886e13 | 368 | |
369 | //make projection | |
370 | fHistFitTheta[i] = (TH1D*)fHistH2InvPtTheta->ProjectionX("_px",firstBin,lastBin,"e"); | |
acb9d358 | 371 | |
372 | Char_t titleTheta[50]; | |
ab886e13 | 373 | if(TMath::Abs(i-(fNThetaBins-1))<0.5) sprintf(titleTheta,"charge/pt (GeV/c) integrated over #theta"); |
374 | else sprintf(titleTheta,"charge/pt (GeV/c) for #theta range: %1.3f - %1.3f",fThetaBins[i],fThetaBins[i+1]); | |
acb9d358 | 375 | |
376 | fHistFitTheta[i]->SetTitle(titleTheta); | |
377 | ||
378 | Double_t invPtMinPos = 0; | |
379 | Double_t invPtMinPosErr = 0; | |
380 | Double_t invPtMinPosImpr = 0; | |
381 | Double_t invPtMinPosErrImpr = 0; | |
382 | ||
6bcacaf1 | 383 | if(fDoRebin) fHistFitTheta[i]->Rebin(fRebin); |
ab886e13 | 384 | |
ba06aaec | 385 | //start fitting |
386 | if(!fFitGaus){ | |
ab886e13 | 387 | Printf("making polynomial fit in charge/pt in theta bins"); |
acb9d358 | 388 | theCan->cd(countPad); |
389 | MakeFit(fHistFitTheta[i],fFitMinPos, invPtMinPos,invPtMinPosErr, fExclRange,fRange); | |
390 | MakeFitBetter(fHistFitTheta[i],fFitMinPosRejP, invPtMinPosImpr, invPtMinPosErrImpr, invPtMinPos,fExclRange,fRange); | |
391 | ||
ab886e13 | 392 | //plot projection and fit |
acb9d358 | 393 | fHistFitTheta[i]->DrawCopy(); |
394 | fFitMinPos->DrawCopy("L,same"); | |
395 | fFitMinPosRejP->DrawCopy("L,same"); | |
396 | } | |
397 | else{ | |
ab886e13 | 398 | Printf("making gauss fit in charge/pt in theta bins"); |
acb9d358 | 399 | theCan->cd(countPad); |
400 | MakeFitInvGauss(fHistFitTheta[i],fFitInvGauss, invPtMinPos, invPtMinPosErr,fExclRange,fRange); | |
401 | MakeFitInvGaussBetter(fHistFitTheta[i],fFitInvGaussRejP, invPtMinPosImpr, invPtMinPosErrImpr,invPtMinPos,fExclRange,fRange); | |
ab886e13 | 402 | |
403 | //plot projection and fit | |
acb9d358 | 404 | fHistFitTheta[i]->DrawCopy(); |
405 | fFitInvGauss->DrawCopy("L,same"); | |
406 | fFitInvGaussRejP->DrawCopy("L,same"); | |
407 | } | |
ab886e13 | 408 | //add objects for analysis folder |
acb9d358 | 409 | aFolderObj->Add(fHistFitTheta[i]); |
ab886e13 | 410 | |
411 | //store fit parameter | |
acb9d358 | 412 | fitParamTheta[i] = invPtMinPosImpr; |
413 | errFitParamTheta[i] = invPtMinPosErrImpr; | |
ab886e13 | 414 | |
415 | //count bins and pad number for displaying | |
acb9d358 | 416 | binsXTheta[i] = i+1.0; |
417 | countPad++; | |
418 | } | |
419 | ||
420 | ||
421 | countPad = 1; | |
422 | ||
423 | ||
ab886e13 | 424 | // analyse charge/pt in bins of phi |
acb9d358 | 425 | |
426 | for(Int_t i=0;i<fNPhiBins;i++){ | |
ab886e13 | 427 | |
428 | //set name for each projection | |
acb9d358 | 429 | TString name = "fit_phi_"; |
430 | name +=i; | |
acb9d358 | 431 | fHistH2InvPtPhi->SetName(name.Data()); |
ab886e13 | 432 | |
433 | //find bins for projection | |
acb9d358 | 434 | Int_t firstBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i]); |
435 | if(i>0) firstBin +=1; | |
436 | Int_t lastBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i+1]); | |
b90a2249 | 437 | if(TMath::Abs(i-(fNPhiBins-1))<0.5){ |
acb9d358 | 438 | firstBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[0]); |
439 | lastBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i]); | |
440 | } | |
ab886e13 | 441 | |
442 | //make projection | |
443 | fHistFitPhi[i] = (TH1D*) fHistH2InvPtPhi->ProjectionX("_px",firstBin,lastBin,"e"); | |
acb9d358 | 444 | |
445 | Char_t titlePhi[50]; | |
ab886e13 | 446 | if(TMath::Abs(i-(fNPhiBins-1))<0.5) sprintf(titlePhi,"charge/pt (GeV/c) integrated over #phi"); |
447 | else sprintf(titlePhi,"charge/pt (GeV/c) for #phi range: %1.3f - %1.3f",fPhiBins[i],fPhiBins[i+1]); | |
acb9d358 | 448 | |
449 | fHistFitPhi[i]->SetTitle(titlePhi); | |
450 | ||
451 | Double_t invPtMinPos = 0; | |
452 | Double_t invPtMinPosErr = 0; | |
453 | Double_t invPtMinPosImpr = 0; | |
454 | Double_t invPtMinPosErrImpr = 0; | |
ab886e13 | 455 | if(fDoRebin) fHistFitPhi[i]->Rebin(fRebin); |
456 | ||
457 | //start fitting | |
ba06aaec | 458 | if(!fFitGaus){ |
ab886e13 | 459 | Printf("making polynomial fit in charge/pt in phi bins"); |
acb9d358 | 460 | phiCan->cd(countPad); |
461 | MakeFit(fHistFitPhi[i],fFitMinPos, invPtMinPos, invPtMinPosErr,fExclRange,fRange); | |
462 | MakeFitBetter(fHistFitPhi[i],fFitMinPosRejP, invPtMinPosImpr, invPtMinPosErrImpr,invPtMinPos,fExclRange,fRange); | |
ab886e13 | 463 | |
464 | //plot projection and fit | |
acb9d358 | 465 | fHistFitPhi[i]->DrawCopy(); |
466 | fFitMinPos->DrawCopy("L,same"); | |
467 | fFitMinPosRejP->DrawCopy("L,same"); | |
468 | ||
469 | } | |
470 | else { | |
ab886e13 | 471 | Printf("making gauss fit in charge/pt in phi bins"); |
acb9d358 | 472 | phiCan->cd(countPad); |
473 | MakeFitInvGauss(fHistFitPhi[i],fFitInvGauss, invPtMinPos, invPtMinPosErr, exclRange,fRange); | |
474 | MakeFitInvGaussBetter(fHistFitPhi[i],fFitInvGaussRejP, invPtMinPosImpr, invPtMinPosErrImpr, invPtMinPos, fExclRange,fRange); | |
ab886e13 | 475 | |
476 | //plot projection and fit | |
acb9d358 | 477 | fHistFitPhi[i]->DrawCopy(); |
478 | fFitInvGauss->DrawCopy("L,same"); | |
479 | fFitInvGaussRejP->DrawCopy("L,same"); | |
480 | } | |
ab886e13 | 481 | //add objects for analysis folder |
acb9d358 | 482 | aFolderObj->Add(fHistFitPhi[i]); |
ab886e13 | 483 | |
484 | //store fit parameter | |
acb9d358 | 485 | fitParamPhi[i] = invPtMinPosImpr; |
486 | errFitParamPhi[i] = invPtMinPosErrImpr; | |
ab886e13 | 487 | |
488 | //count bins and pad number for displaying | |
acb9d358 | 489 | binsXPhi[i] = i+1.0; |
490 | countPad++; | |
491 | } | |
ab886e13 | 492 | |
493 | //initialize graphs for displaying the fit parameter | |
ba06aaec | 494 | InitGraphs(binsXTheta,fitParamTheta,errFitParamTheta,binsXPhi,fitParamPhi,errFitParamPhi); |
acb9d358 | 495 | |
ba06aaec | 496 | //plot fit values = minimum positions of charge/pt |
acb9d358 | 497 | TCanvas *canFitVal = new TCanvas("canFitVal","min pos histos",800,400); |
498 | canFitVal->Divide(2,1); | |
499 | ||
500 | canFitVal->cd(1); | |
6bcacaf1 | 501 | fGrMinPosTheta->Draw("AP"); |
acb9d358 | 502 | canFitVal->cd(2); |
6bcacaf1 | 503 | fGrMinPosPhi->Draw("AP"); |
acb9d358 | 504 | |
ab886e13 | 505 | Printf("AliPerfAnalyzeInvPt::StartAnalysis: ***** NOTE: last bin is always fit result of integral over all angle BINS which have been set by user! angle ranges of 2D projections are kept.*****"); |
ba06aaec | 506 | |
507 | //add objects to folder | |
acb9d358 | 508 | aFolderObj->Add(fGrMinPosTheta); |
509 | aFolderObj->Add(fGrMinPosPhi); | |
510 | aFolderObj->Add(fHistH2InvPtTheta); | |
511 | aFolderObj->Add(fHistH2InvPtPhi); | |
512 | ||
513 | ||
514 | } | |
515 | ||
516 | ||
517 | //____________________________________________________________________________________________________________________________________________ | |
518 | ||
ab886e13 | 519 | void AliPerfAnalyzeInvPt::MakeFit(TH1D *hproy, TF1 * fitpb, Double_t &mean, Double_t &errMean, Double_t &excl,Double_t &range) |
acb9d358 | 520 | { |
ab886e13 | 521 | // fit charge/pt and extract minimum position |
ba06aaec | 522 | |
acb9d358 | 523 | fitpb->SetRange(-range,range); |
524 | fitpb->SetParLimits(1,-0.05,0.05); | |
525 | fitpb->SetParameter(0,1.0); | |
526 | fitpb->FixParameter(4,excl); | |
527 | fitpb->FixParameter(2,0.0); | |
528 | ||
ba06aaec | 529 | hproy->Fit(fitpb,"RM"); |
acb9d358 | 530 | mean = fitpb->GetParameter(1); |
531 | errMean = fitpb->GetParError(1); | |
ba06aaec | 532 | if(!mean) errMean = 0.0; |
acb9d358 | 533 | } |
534 | //____________________________________________________________________________________________________________________________________________ | |
ab886e13 | 535 | void AliPerfAnalyzeInvPt::MakeFitBetter(TH1D *hproy, TF1 * fitpb2, Double_t &mean, Double_t &errMean, Double_t &f, Double_t &excl, Double_t &range) |
acb9d358 | 536 | { |
ba06aaec | 537 | // adjust fit range to minimum position of AliPerfAnalyzeInvPt::MakeFit and fit 1/pt and extract new minimum position |
538 | ||
acb9d358 | 539 | fitpb2->FixParameter(5,excl); |
540 | fitpb2->FixParameter(4,f); | |
541 | // fitpb2->FixParameter(2,0.0); | |
ba06aaec | 542 | fitpb2->SetRange(-range+f,range+f); |
acb9d358 | 543 | fitpb2->SetParLimits(1,-0.05,0.05); |
544 | fitpb2->SetParameter(0,1.0); | |
ba06aaec | 545 | hproy->Fit(fitpb2,"RM"); |
acb9d358 | 546 | mean = fitpb2->GetParameter(1); |
547 | errMean = fitpb2->GetParError(1); | |
ba06aaec | 548 | if(!mean) errMean = 0.0; |
acb9d358 | 549 | |
550 | } | |
551 | //____________________________________________________________________________________________________________________________________________ | |
ab886e13 | 552 | void AliPerfAnalyzeInvPt::MakeFitInvGauss(TH1D *hproy, TF1 * fitpb, Double_t &mean, Double_t &errMean, Double_t &excl,Double_t &range) |
acb9d358 | 553 | { |
ab886e13 | 554 | // fit charge/pt and extract minimum position |
ba06aaec | 555 | |
acb9d358 | 556 | fitpb->FixParameter(6,excl); |
557 | fitpb->SetRange(-range,range); | |
558 | fitpb->SetParameter(0,-1.0); | |
559 | fitpb->SetParameter(2,1.0); | |
560 | fitpb->SetParameter(3,25000.0); | |
561 | fitpb->SetParameter(4,-1.0); | |
562 | fitpb->SetParLimits(1,-0.02,0.02); | |
ba06aaec | 563 | hproy->Fit(fitpb,"RM"); |
acb9d358 | 564 | mean = fitpb->GetParameter(1); |
565 | errMean = fitpb->GetParError(1); | |
ba06aaec | 566 | if(!mean) errMean = 0.0; |
acb9d358 | 567 | |
568 | } | |
569 | //____________________________________________________________________________________________________________________________________________ | |
ab886e13 | 570 | void AliPerfAnalyzeInvPt::MakeFitInvGaussBetter(TH1D *hproy, TF1 * fitpb2, Double_t &mean, Double_t &errMean, Double_t &f, Double_t &excl, Double_t &range) |
acb9d358 | 571 | { |
ba06aaec | 572 | // adjust fit range to minimum position of AliPerfAnalyzeInvPt::MakeFitInvGauss and fit 1/pt and extract new minimum position |
573 | ||
acb9d358 | 574 | fitpb2->FixParameter(7,excl); |
575 | fitpb2->FixParameter(6,f); | |
576 | fitpb2->SetRange(-range+f,range+f); | |
577 | fitpb2->SetParameter(0,-1.0); | |
578 | fitpb2->SetParameter(2,1.0); | |
579 | fitpb2->SetParameter(3,25000.0); | |
580 | fitpb2->SetParameter(4,-1.0); | |
581 | fitpb2->SetParLimits(1,-0.02,0.02); | |
ba06aaec | 582 | hproy->Fit(fitpb2,"RM"); |
acb9d358 | 583 | mean = fitpb2->GetParameter(1); |
584 | errMean = fitpb2->GetParError(1); | |
ba06aaec | 585 | if(!mean) errMean = 0.0; |
acb9d358 | 586 | |
587 | } | |
588 | ||
589 | ||
590 | //____________________________________________________________________________________________________________________________________________ | |
acb9d358 | 591 | |
ba06aaec | 592 | // set variables for analysis |
593 | void AliPerfAnalyzeInvPt::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins){ | |
acb9d358 | 594 | |
ba06aaec | 595 | // set theta bins for Analyse() |
596 | //set phi bins as array and set number of this array which is equal to number of bins analysed | |
597 | //the last analysed bin will always be the projection from first to last bin in the array | |
598 | if(nphBins){ | |
599 | fNPhiBins = nphBins; | |
acb9d358 | 600 | |
ba06aaec | 601 | for(Int_t k = 0;k<fNPhiBins;k++){ |
602 | fPhiBins[k] = phiBinArray[k]; | |
603 | } | |
604 | Printf("AliPerfAnalyzeInvPt::SetProjBinsPhi:number of bins in phi set to %i",fNPhiBins); | |
acb9d358 | 605 | } |
ba06aaec | 606 | else Printf("AliPerfAnalyzeInvPt::SetProjBinsPhi: number of bins in theta NOT set. Default values are used."); |
acb9d358 | 607 | } |
608 | //____________________________________________________________________________________________________________________________________________ | |
ba06aaec | 609 | void AliPerfAnalyzeInvPt::SetProjBinsTheta(const Double_t *thetaBinArray, const Int_t nthBins){ |
610 | // set theta bins for Analyse() | |
611 | //set theta bins as array and set number of this array which is equal to number of bins analysed | |
612 | //the last analysed bin will always be the projection from first to last bin in the array | |
613 | if(nthBins){ | |
614 | fNThetaBins = nthBins; | |
615 | for(Int_t k = 0;k<fNThetaBins;k++){ | |
616 | fThetaBins[k] = thetaBinArray[k]; | |
617 | } | |
618 | Printf("AliPerfAnalyzeInvPt::SetProjBinsTheta: number of bins in theta set to %i",fNThetaBins); | |
acb9d358 | 619 | } |
ba06aaec | 620 | else Printf("AliPerfAnalyzeInvPt::SetProjBinsTheta: number of bins in theta NOT set. Default values are used."); |
acb9d358 | 621 | } |
622 | //____________________________________________________________________________________________________________________________________________ | |
623 | void AliPerfAnalyzeInvPt::SetMakeFitOption(const Bool_t setGausFit, const Double_t exclusionR,const Double_t fitR ){ | |
624 | ||
ba06aaec | 625 | //set the fit options: |
626 | //for usage of gaussian function instead of polynomial (default) set setGausFit=kTRUE | |
627 | //set the range of rejection of points around 0 via exclusionR | |
628 | //set the fit range around 0 with fitR | |
629 | if(fitR){ | |
630 | fFitGaus = setGausFit; | |
631 | fExclRange = exclusionR; | |
632 | fRange = fitR; | |
acb9d358 | 633 | |
ab886e13 | 634 | if(fFitGaus) Printf("AliPerfAnalyzeInvPt::SetMakeFitOption: set MakeGausFit with fit range %2.3f and exclusion range in fabs(1/pt): %2.3f",fRange,fExclRange); |
635 | else Printf("AliPerfAnalyzeInvPt::SetMakeFitOption: set standard polynomial fit with fit range %2.3f and exclusion range in fabs(1/pt): %2.3f",fRange,fExclRange); | |
ba06aaec | 636 | } |
637 | else Printf(" AliPerfAnalyzeInvPt::SetMakeFitOption: no user input. Set standard polynomial fit with fit range 1.0 and exclusion range 0.13."); | |
acb9d358 | 638 | } |