ac48988f57dd6ecfc783380f2b6904c1b6f25fc3
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowAnalysisWithFittingQDistribution.cxx
1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3 *                                                                        *
4 * Author: The ALICE Off-line Project.                                    *
5 * Contributors are mentioned in the code where appropriate.              *
6 *                                                                        *
7 * Permission to use, copy, modify and distribute this software and its   *
8 * documentation strictly for non-commercial purposes is hereby granted   *
9 * without fee, provided that the above copyright notice appears in all   *
10 * copies and that both the copyright notice and this permission notice   *
11 * appear in the supporting documentation. The authors make no claims     *
12 * about the suitability of this software for any purpose. It is          *
13 * provided "as is" without express or implied warranty.                  * 
14 **************************************************************************/
15
16 /******************************** 
17  * integrated flow estimate by  *
18  *   fitting q-distribution     * 
19  *                              *
20  * author: Ante Bilandzic       * 
21  *          (anteb@nikhef.nl)   *
22  *                              *  
23  *  based on the macro written  *
24  *     by Sergei Voloshin       *
25  *******************************/  
26
27 #define AliFlowAnalysisWithFittingQDistribution_cxx
28
29 #include "Riostream.h"
30 #include "AliFlowCommonConstants.h"
31 #include "AliFlowCommonHist.h"
32 #include "AliFlowCommonHistResults.h"
33 #include "TChain.h"
34 #include "TFile.h"
35 #include "TList.h" 
36 #include "TF1.h"
37 #include "TLegend.h"
38 #include "TParticle.h"
39 #include "TProfile.h"
40 #include "AliFlowEventSimple.h"
41 #include "AliFlowTrackSimple.h"
42 #include "AliFlowAnalysisWithFittingQDistribution.h"
43
44
45 class TH1;
46 class TGraph;
47 class TPave;
48 class TLatex;
49 class TMarker;
50 class TObjArray;
51 class TList;
52 class TSystem;
53 class TROOT;
54 class AliFlowVector;
55 class TVector;
56
57 //================================================================================================================
58
59 ClassImp(AliFlowAnalysisWithFittingQDistribution)
60
61 AliFlowAnalysisWithFittingQDistribution::AliFlowAnalysisWithFittingQDistribution():  
62  fHistList(NULL),
63  fCommonHists(NULL),
64  fCommonHistsResults(NULL),
65  fnBinsPhi(0),
66  fPhiMin(0),
67  fPhiMax(0),
68  fPhiBinWidth(0),
69  fnBinsPt(0),
70  fPtMin(0),
71  fPtMax(0),
72  fPtBinWidth(0),
73  fnBinsEta(0),
74  fEtaMin(0),
75  fEtaMax(0),
76  fEtaBinWidth(0),
77  fHarmonic(2),
78  fAnalysisLabel(NULL),
79  fWeightsList(NULL),
80  fUsePhiWeights(kFALSE),
81  fUsePtWeights(kFALSE),
82  fUseEtaWeights(kFALSE),
83  fUseParticleWeights(NULL),
84  fPhiWeights(NULL),
85  fPtWeights(NULL),
86  fEtaWeights(NULL),
87  // fitting parameters with default values:
88  fFittingParameters(NULL), 
89  fTreshold(5),
90  fvStart(0.05),
91  fvMin(0.0),
92  fvMax(0.25),
93  fSigma2Start(0.75),
94  fSigma2Min(0.5), // (should be kept fixed at 0.5 according to theorists...)
95  fSigma2Max(2.5),
96  fPlotResults(kFALSE),
97  // rest:
98  fLegend(NULL)
99  {
100   // constructor 
101   
102   // base list to hold all output objects:
103   fHistList = new TList();
104   fHistList->SetName("cobjFQD");
105   fHistList->SetOwner(kTRUE);
106   
107   // analysis label;
108   fAnalysisLabel = new TString();
109  
110   // list to hold histograms with phi, pt and eta weights:      
111   fWeightsList = new TList();
112
113   // initialize all arrays:  
114   this->InitializeArrays();
115
116  } // end of constructor
117  
118  
119 //================================================================================================================
120  
121
122 AliFlowAnalysisWithFittingQDistribution::~AliFlowAnalysisWithFittingQDistribution()
123 {
124  // desctructor
125  delete fHistList; 
126 }
127
128
129 //================================================================================================================
130
131
132 void AliFlowAnalysisWithFittingQDistribution::Init()
133 {
134  // access constants and book everything
135  
136  // access constants:
137  this->AccessConstants();
138  
139  // booking:
140  this->BookCommonHistograms();
141  this->BookAndFillWeightsHistograms();
142  this->BookEverythingForDistributions();
143  
144  // store fitting parameters:  
145  this->StoreFittingParameters();
146  
147  // nest lists:
148  fWeightsList->SetName("Weights");
149  fWeightsList->SetOwner(kTRUE);   
150  fHistList->Add(fWeightsList);
151  
152  // set harmonic in common control histograms (to be improved (should I do this somewhere else?)):
153  (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic); 
154  
155 } // end of void AliFlowAnalysisWithFittingQDistribution::Init()
156
157
158 //================================================================================================================
159
160
161 void AliFlowAnalysisWithFittingQDistribution::Make(AliFlowEventSimple* anEvent)
162 {
163  // loop over data
164  
165  // a) fill the common control histograms
166  // b) loop over data and calculate non-weighted and weighted Q-vector and sum of particle weights
167  // c) fill histograms for distribution
168  // d) reset e-b-e quantities
169  
170  // fill the common control histograms
171  fCommonHists->FillControlHistograms(anEvent); 
172  
173  Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
174  Double_t dPt  = 0.; // transverse momentum
175  Double_t dEta = 0.; // pseudorapidity
176
177  Double_t wPhi = 1.; // phi weight
178  Double_t wPt  = 1.; // pt weight
179  Double_t wEta = 1.; // eta weight
180  
181  Int_t nPrim = anEvent->NumberOfTracks();  // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
182                                            // nRP   = # of particles used to determine the reaction plane;
183                                            // nPOI  = # of particles of interest for a detailed flow analysis;
184                                            // rest  = # of particles which are not niether RPs nor POIs.  
185  
186  Int_t n = fHarmonic; // shortcut for the harmonic
187  
188  Double_t dReQ[2] = {0.}; // real part of Q-vector [0=particle weights not used, 1=particle weights used]
189  Double_t dImQ[2] = {0.}; // imaginary part of Q-vector [0=particle weights not used, 1=particle weights used]
190  Double_t dSumOfParticleWeights[2]; // [0=particle weights not used, 1=particle weights used] 
191                                                                                                                                
192  AliFlowTrackSimple *aftsTrack = NULL;                                          
193                                            
194  for(Int_t i=0;i<nPrim;i++) 
195  { 
196   aftsTrack=anEvent->GetTrack(i);
197   if(aftsTrack)
198   {
199    if(!(aftsTrack->InRPSelection())) continue; // consider only tracks which are RPs
200      
201    dPhi = aftsTrack->Phi();
202    dPt  = aftsTrack->Pt();
203    dEta = aftsTrack->Eta();
204   
205    if(fUsePhiWeights && fPhiWeights && fnBinsPhi) // determine phi weight for this particle:
206    {
207     wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
208    }
209    if(fUsePtWeights && fPtWeights && fnBinsPt) // determine pt weight for this particle:
210    {
211     wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth))); 
212    }            
213    if(fUseEtaWeights && fEtaWeights && fEtaBinWidth) // determine eta weight for this particle: 
214    {
215     wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth))); 
216    } 
217     
218    // calculate real and imaginary part of non-weighted and weighted Q-vector and sum of particle weights for this event:
219    for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
220    {
221     // Q-vector:
222     dReQ[pW]+=pow(wPhi*wPt*wEta,pW)*TMath::Cos(n*dPhi); 
223     dImQ[pW]+=pow(wPhi*wPt*wEta,pW)*TMath::Sin(n*dPhi);
224     // sum of particle weights:
225     dSumOfParticleWeights[pW] += pow(wPhi*wPt*wEta,pW); // if pW = 0, this sum gives nRP
226    } 
227    
228   } // end of if(aftsTrack)
229  } // end of for(Int_t i=0;i<nPrim;i++)                                      
230                                            
231  // calculate q = Q\sqrt{sum of particle weights}:
232  // Remark: if particle weights are unit than sum of particle weights = multiplicity
233  Double_t q=0;                                          
234  for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
235  {
236   if(dSumOfParticleWeights[pW])
237   {
238    q = pow(dReQ[pW]*dReQ[pW]+dImQ[pW]*dImQ[pW],0.5)/pow(dSumOfParticleWeights[pW],0.5);
239    // fill histograms:
240    fqDistribution[pW]->Fill(q,1.);
241    fSumOfParticleWeights[pW]->Fill(dSumOfParticleWeights[pW],1.);
242   }
243  } 
244  
245  // reset e-b-e quantities:
246  for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
247  {
248   dReQ[pW] = 0.;
249   dImQ[pW] = 0.;
250   dSumOfParticleWeights[pW] = 0.;
251  }
252
253 }//end of Make()
254
255
256 //================================================================================================================
257
258
259 void AliFlowAnalysisWithFittingQDistribution::GetOutputHistograms(TList *outputListHistos) 
260 {
261  // get pointers to all output histograms (called before Finish()) 
262  
263  if(outputListHistos)   
264  {   
265   // 1.) common control histograms and common histograms for final results:
266   TString commonHistName = "AliFlowCommonHistFQD";
267   commonHistName += fAnalysisLabel->Data();
268   AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(outputListHistos->FindObject(commonHistName.Data()));
269   if(commonHist) this->SetCommonHists(commonHist); 
270   
271   TString commonHistResName = "AliFlowCommonHistResultsFQD";
272   commonHistResName += fAnalysisLabel->Data();
273   AliFlowCommonHistResults *commonHistRes = dynamic_cast<AliFlowCommonHistResults*>
274                                             (outputListHistos->FindObject(commonHistResName.Data()));
275   if(commonHistRes) this->SetCommonHistsResults(commonHistRes); 
276   
277   // 2.) weights: 
278   TList *weightsList = dynamic_cast<TList*>(outputListHistos->FindObject("Weights"));
279   if(weightsList) this->SetWeightsList(weightsList);
280   Bool_t bUsePhiWeights = kFALSE;
281   Bool_t bUsePtWeights = kFALSE;
282   Bool_t bUseEtaWeights = kFALSE;
283   TString fUseParticleWeightsName = "fUseParticleWeightsFQD";
284   fUseParticleWeightsName += fAnalysisLabel->Data();
285   TProfile *useParticleWeights = dynamic_cast<TProfile*>(weightsList->FindObject(fUseParticleWeightsName.Data()));
286   if(useParticleWeights)
287   {
288    this->SetUseParticleWeights(useParticleWeights);  
289    bUsePhiWeights = (Int_t)useParticleWeights->GetBinContent(1);
290    bUsePtWeights = (Int_t)useParticleWeights->GetBinContent(2);
291    bUseEtaWeights = (Int_t)useParticleWeights->GetBinContent(3);
292   }
293   
294   // 3.) distributions and 4.) final results of fitting:
295   TString pWeightsFlag[2] = {"pWeights not used","pWeights used"};
296   TString sigmaFlag[2] = {"#sigma^{2} not fitted","#sigma^{2} fitted"};
297   
298   // q-distribution:
299   TString qDistributionName = "fqDistribution";
300   qDistributionName += fAnalysisLabel->Data();
301   // sum of particle weights:
302   TString sumOfParticleWeightsName = "fSumOfParticleWeightsName"; 
303   sumOfParticleWeightsName += fAnalysisLabel->Data();
304   // final results for integrated flow:
305   TString intFlowName = "fIntFlowFQD";
306   intFlowName += fAnalysisLabel->Data();
307   // sigma^2:
308   TString sigma2Name = "fSigma2";
309   sigma2Name += fAnalysisLabel->Data();
310   // chi^2:
311   TString chi2Name = "fChi2";
312   chi2Name += fAnalysisLabel->Data();
313   
314   TH1D *qDistribution[2] = {NULL};
315   TH1D *sumOfParticleWeights[2] = {NULL};
316   TH1D *intFlow[2][2] = {{NULL}};
317   TH1D *sigma2[2][2] = {{NULL}};
318   TH1D *chi2[2][2] = {{NULL}};
319    
320   for(Int_t pW=0;pW<1+(Int_t)(bUsePhiWeights||bUsePtWeights||bUseEtaWeights);pW++)
321   {
322    // q-distribution:
323    qDistribution[pW] = dynamic_cast<TH1D*>(outputListHistos->FindObject(Form("%s, %s",qDistributionName.Data(),pWeightsFlag[pW].Data())));
324    if(qDistribution[pW])
325    {
326     this->SetqDistribution(qDistribution[pW],pW);
327    } else
328      {
329       cout<<"WARNING: qDistribution[pW] is NULL in AFAWFQD::GOH() !!!!"<<endl;
330       cout<<"pW = "<<pW<<endl;
331      }
332    // sum of particle weights:
333    sumOfParticleWeights[pW] = dynamic_cast<TH1D*>(outputListHistos->FindObject(Form("%s, %s",sumOfParticleWeightsName.Data(),pWeightsFlag[pW].Data())));
334    if(sumOfParticleWeights[pW])
335    {
336     this->SetSumOfParticleWeights(sumOfParticleWeights[pW],pW);
337    } else
338      {
339       cout<<"WARNING: sumOfParticleWeights[pW] is NULL in AFAWFQD::GOH() !!!!"<<endl;
340       cout<<"pW = "<<pW<<endl;
341      }
342    // final results:
343    for(Int_t f=0;f<2;f++)
344    {
345     // final results for integrated flow:
346     intFlow[pW][f] = dynamic_cast<TH1D*>(outputListHistos->FindObject(Form("%s, %s, %s",intFlowName.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data())));
347     if(intFlow[pW][f])
348     {
349      this->SetIntFlow(intFlow[pW][f],pW,f);
350     } else 
351       {
352        cout<<"WARNING: intFlow[pW][f] is NULL in AFAWFQD::GOH() !!!!"<<endl;
353        cout<<"pW = "<<pW<<endl;
354        cout<<"f  = "<<f<<endl;
355       }
356     // sigma^2:
357     sigma2[pW][f] = dynamic_cast<TH1D*>(outputListHistos->FindObject(Form("%s, %s, %s",sigma2Name.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data())));
358     if(sigma2[pW][f])
359     {
360      this->SetSigma2(sigma2[pW][f],pW,f);
361     } else 
362       {
363        cout<<"WARNING: sigma2[pW][f] is NULL in AFAWFQD::GOH() !!!!"<<endl;
364        cout<<"pW = "<<pW<<endl;
365        cout<<"f  = "<<f<<endl;
366       } 
367     // chi^2:
368     chi2[pW][f] = dynamic_cast<TH1D*>(outputListHistos->FindObject(Form("%s, %s, %s",chi2Name.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data())));
369     if(chi2[pW][f])
370     {
371      this->SetChi2(chi2[pW][f],pW,f);
372     } else 
373       {
374        cout<<"WARNING: chi2[pW][f] is NULL in AFAWFQD::GOH() !!!!"<<endl;
375        cout<<"pW = "<<pW<<endl;
376        cout<<"f  = "<<f<<endl;
377       }  
378       
379    } // end of for(Int_t f=0;f<2;f++)
380   } // end of for(Int_t pW=0;pW<1+(Int_t)(bUsePhiWeights||bUsePtWeights||bUseEtaWeights);pW++)
381   
382   // 5.) fitting parameters:
383   // q-distribution:
384   TString fittingParametersName = "fFittingParameters";
385   fittingParametersName += fAnalysisLabel->Data();
386   TProfile *fittingParameters = NULL;
387   fittingParameters = dynamic_cast<TProfile*>(outputListHistos->FindObject(fittingParametersName.Data()));
388   if(fittingParameters)
389   {
390    this->SetFittingParameters(fittingParameters);
391   } else
392     {
393      cout<<"WARNING:fittingParameters is NULL in AFAWFQD::GOH() !!!!"<<endl;
394     }
395   
396  } else // to if(outputListHistos)
397    {
398     cout<<"WARNING: outputListHistos is NULL in AFAWFQD::GOH() !!!!"<<endl;
399     exit(0);
400    } 
401     
402    
403 } // end of void AliFlowAnalysisWithFittingQDistribution::GetOutputHistograms(TList *outputListHistos) 
404
405
406 //================================================================================================================
407
408
409 void AliFlowAnalysisWithFittingQDistribution::Finish(Bool_t doFit)
410 {
411  // calculate the final results
412  
413  // a) acces the constants and fitting paremeters;
414  // b) access the flags;
415  // c) do final fit;
416  // d) fill common hist results;
417  // e) print on the screen the final results.
418  
419  // access the constants and fitting paremeters:
420  this->AccessConstants();
421  this->AccessFittingParameters();
422  
423  // access the flags: 
424  fUsePhiWeights = (Int_t)fUseParticleWeights->GetBinContent(1); 
425  fUsePtWeights = (Int_t)fUseParticleWeights->GetBinContent(2); 
426  fUseEtaWeights = (Int_t)fUseParticleWeights->GetBinContent(3);
427
428  // to be improved (moved somewhere else):
429  if(fPlotResults)
430  {
431   fLegend = new TLegend(0.6,0.55,0.85,0.7); 
432  }
433
434  // do final fit:             
435  if(doFit) 
436  {
437   // particle weights not used:
438   // a) sigma^2 not fitted (fixed to 0.5):
439   this->DoFit(kFALSE,kFALSE);
440   // b) sigma^2 fitted:
441   this->DoFit(kFALSE,kTRUE);
442   // particle weights used:
443   if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)      
444   {
445    // a) sigma^2 not fitted (fixed to 0.5):
446    this->DoFit(kTRUE,kFALSE);  
447    // b) sigma^2 fitted:
448    this->DoFit(kTRUE,kTRUE);  
449   }
450   
451   // fill common hist results (by default fill results obtained with sigma^2 fitted):
452   if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)
453   {
454    this->FillCommonHistResultsIntFlow(kTRUE,kTRUE); 
455   } else
456     {
457      this->FillCommonHistResultsIntFlow(kFALSE,kTRUE);    
458     } 
459   
460   // print final results on the screen:
461   this->PrintFinalResultsForIntegratedFlow();  
462  } // end of if(doFit)
463    
464 } // end of void AliFlowAnalysisWithFittingQDistribution::Finish(Bool_t doFit)
465
466
467 //================================================================================================================
468
469
470 void AliFlowAnalysisWithFittingQDistribution::WriteHistograms(TString* outputFileName)
471 {
472  //store the final results in output .root file
473  TFile *output = new TFile(outputFileName->Data(),"RECREATE");
474  //output->WriteObject(fHistList, "cobjFQD","SingleKey");
475  fHistList->SetName("cobjFQD");
476  fHistList->SetOwner(kTRUE);
477  fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
478  delete output;
479 }
480
481 //================================================================================================================
482
483 void AliFlowAnalysisWithFittingQDistribution::WriteHistograms(TString outputFileName)
484 {
485  //store the final results in output .root file
486  TFile *output = new TFile(outputFileName.Data(),"RECREATE");
487  //output->WriteObject(fHistList, "cobjFQD","SingleKey");
488  fHistList->SetName("cobjFQD");
489  fHistList->SetOwner(kTRUE);
490  fHistList->Write(fHistList->GetName(), TObject::kSingleKey); 
491  delete output;
492 }
493
494
495 //================================================================================================================
496
497
498 void AliFlowAnalysisWithFittingQDistribution::InitializeArrays()
499 {
500  // initialize all arrays
501  
502  for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
503  {
504   fSumOfParticleWeights[pW] = NULL;
505   fqDistribution[pW] = NULL; 
506   for(Int_t f=0;f<2;f++) // sigma^2 not fitted (0) or fitted (1)
507   {
508    fIntFlow[pW][f] = NULL;
509    fSigma2[pW][f] = NULL;
510    fChi2[pW][f] = NULL;
511   }
512  } 
513
514 } // end of void AliFlowAnalysisWithFittingQDistribution::InitializeArrays()
515
516
517 //================================================================================================================
518
519
520 void AliFlowAnalysisWithFittingQDistribution::BookCommonHistograms()
521 {
522  // book common histograms
523  
524  // common control histogram: 
525  TString commonHistName = "AliFlowCommonHistFQD";
526  commonHistName += fAnalysisLabel->Data();
527  fCommonHists = new AliFlowCommonHist(commonHistName.Data());
528  fHistList->Add(fCommonHists);  
529
530  // common histograms for final results:
531  TString commonHistResName = "AliFlowCommonHistResultsFQD";
532  commonHistResName += fAnalysisLabel->Data();
533  fCommonHistsResults = new AliFlowCommonHistResults(commonHistResName.Data());
534  fHistList->Add(fCommonHistsResults); 
535
536 } // end of void AliFlowAnalysisWithFittingQDistribution::BookCommonHistograms(
537
538
539 //================================================================================================================
540
541
542 void AliFlowAnalysisWithFittingQDistribution::BookAndFillWeightsHistograms()
543 {
544  // book and fill histograms which hold phi, pt and eta weights
545
546  if(!fWeightsList)
547  {
548   cout<<"WARNING: fWeightsList is NULL in AFAWFQD::BAFWH() !!!!"<<endl;
549   exit(0);  
550  }
551     
552  TString fUseParticleWeightsName = "fUseParticleWeightsFQD";
553  fUseParticleWeightsName += fAnalysisLabel->Data();
554  fUseParticleWeights = new TProfile(fUseParticleWeightsName.Data(),"0 = particle weight not used, 1 = particle weight used ",3,0,3);
555  fUseParticleWeights->SetLabelSize(0.06);
556  (fUseParticleWeights->GetXaxis())->SetBinLabel(1,"w_{#phi}");
557  (fUseParticleWeights->GetXaxis())->SetBinLabel(2,"w_{p_{T}}");
558  (fUseParticleWeights->GetXaxis())->SetBinLabel(3,"w_{#eta}");
559  fUseParticleWeights->Fill(0.5,(Int_t)fUsePhiWeights);
560  fUseParticleWeights->Fill(1.5,(Int_t)fUsePtWeights);
561  fUseParticleWeights->Fill(2.5,(Int_t)fUseEtaWeights);
562  fWeightsList->Add(fUseParticleWeights); 
563   
564  if(fUsePhiWeights)
565  {
566   if(fWeightsList->FindObject("phi_weights"))
567   {
568    fPhiWeights = dynamic_cast<TH1F*>(fWeightsList->FindObject("phi_weights"));
569    if(fPhiWeights->GetBinWidth(1) != fPhiBinWidth)
570    {
571     cout<<"WARNING: fPhiWeights->GetBinWidth(1) != fPhiBinWidth in AFAWFQD::BAFWH() !!!!        "<<endl;
572     cout<<"         This indicates inconsistent binning in phi histograms throughout the code."<<endl;
573     exit(0);
574    }
575   } else 
576     {
577      cout<<"WARNING: fWeightsList->FindObject(\"phi_weights\") is NULL in AFAWFQD::BAFWH() !!!!"<<endl;
578      exit(0);
579     }
580  } // end of if(fUsePhiWeights)
581  
582  if(fUsePtWeights) 
583  {
584   if(fWeightsList->FindObject("pt_weights"))
585   {
586    fPtWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("pt_weights"));
587    if(fPtWeights->GetBinWidth(1) != fPtBinWidth)
588    {
589     cout<<"WARNING: fPtWeights->GetBinWidth(1) != fPtBinWidth in AFAWFQD::BAFWH() !!!!         "<<endl;
590     cout<<"         This indicates insconsistent binning in pt histograms throughout the code."<<endl;
591     exit(0);
592    }
593   } else 
594     {
595      cout<<"WARNING: fWeightsList->FindObject(\"pt_weights\") is NULL in AFAWFQD::BAFWH() !!!!"<<endl;
596      exit(0);
597     }
598  } // end of if(fUsePtWeights)    
599
600  if(fUseEtaWeights) 
601  {
602   if(fWeightsList->FindObject("eta_weights"))
603   {
604    fEtaWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("eta_weights"));
605    if(fEtaWeights->GetBinWidth(1) != fEtaBinWidth)
606    {
607     cout<<"WARNING: fEtaWeights->GetBinWidth(1) != fEtaBinWidth in AFAWFQD::BAFWH() !!!!        "<<endl;
608     cout<<"         This indicates insconsistent binning in eta histograms throughout the code."<<endl;
609     exit(0);
610    }
611   } else 
612     {
613      cout<<"WARNING: fUseEtaWeights && fWeightsList->FindObject(\"eta_weights\") is NULL in AFAWFQD::BAFWH() !!!!"<<endl;
614      exit(0);
615     }
616  } // end of if(fUseEtaWeights)
617  
618 } // end of AliFlowAnalysisWithFittingQDistribution::BookAndFillWeightsHistograms()
619
620
621 //================================================================================================================================
622
623
624 void AliFlowAnalysisWithFittingQDistribution::AccessConstants()
625 {
626  // access needed common constants from AliFlowCommonConstants
627  
628  fnBinsPhi = AliFlowCommonConstants::GetNbinsPhi();
629  fPhiMin = AliFlowCommonConstants::GetPhiMin();      
630  fPhiMax = AliFlowCommonConstants::GetPhiMax();
631  if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;  
632  fnBinsPt = AliFlowCommonConstants::GetNbinsPt();
633  fPtMin = AliFlowCommonConstants::GetPtMin();        
634  fPtMax = AliFlowCommonConstants::GetPtMax();
635  if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;  
636  fnBinsEta = AliFlowCommonConstants::GetNbinsEta();
637  fEtaMin = AliFlowCommonConstants::GetEtaMin();      
638  fEtaMax = AliFlowCommonConstants::GetEtaMax();
639  if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;  
640  
641 } // end of void AliFlowAnalysisWithFittingQDistribution::AccessConstants()
642
643
644 //================================================================================================================================
645
646
647 void AliFlowAnalysisWithFittingQDistribution::BookEverythingForDistributions()
648 {
649  // book histograms for distributions
650  
651  TString pWeightsFlag[2] = {"pWeights not used","pWeights used"};
652  TString sigmaFlag[2] = {"#sigma^{2} not fitted","#sigma^{2} fitted"};
653  // q-distribution:
654  TString fqDistributionName = "fqDistribution";
655  fqDistributionName += fAnalysisLabel->Data();
656  // sum of particle weights:
657  TString fSumOfParticleWeightsName = "fSumOfParticleWeightsName";
658  fSumOfParticleWeightsName += fAnalysisLabel->Data();
659  // final results for integrated flow:
660  TString fIntFlowName = "fIntFlowFQD";
661  fIntFlowName += fAnalysisLabel->Data();
662  // sigma^2:
663  TString fSigma2Name = "fSigma2";
664  fSigma2Name += fAnalysisLabel->Data();
665  // chi^2:
666  TString fChi2Name = "fChi2";
667  fChi2Name += fAnalysisLabel->Data();
668  
669  for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
670  {
671   // q-distribution:
672   fqDistribution[pW] = new TH1D(Form("%s, %s",fqDistributionName.Data(),pWeightsFlag[pW].Data()),"q-distribution",10000,0,1000);
673   fqDistribution[pW]->SetXTitle("q_{n}=Q_{n}/#sqrt{M}");
674   fqDistribution[pW]->SetYTitle("Counts");
675   fHistList->Add(fqDistribution[pW]);
676   // sum of particle weights: 
677   fSumOfParticleWeights[pW] = new TH1D(Form("%s, %s",fSumOfParticleWeightsName.Data(),pWeightsFlag[pW].Data()),"Sum of particle weights",1000,0,10000);
678   fSumOfParticleWeights[pW]->SetXTitle("#sum_{i=1}^{N} w_{i}");
679   fSumOfParticleWeights[pW]->SetYTitle("Counts");
680   fHistList->Add(fSumOfParticleWeights[pW]);
681   
682   for(Int_t f=0;f<2;f++) // sigma^2 not fitted (0) or fitted (1)
683   {
684    // final results for integrated flow:
685    fIntFlow[pW][f] = new TH1D(Form("%s, %s, %s",fIntFlowName.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data()),"Integrated Flow",1,0,1);
686    fIntFlow[pW][f]->SetLabelSize(0.08);
687    (fIntFlow[pW][f]->GetXaxis())->SetBinLabel(1,"v_{n}");
688    fHistList->Add(fIntFlow[pW][f]);
689    // sigma^2:
690    fSigma2[pW][f] = new TH1D(Form("%s, %s, %s",fSigma2Name.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data()),"#sigma^{2}",1,0,1);
691    fSigma2[pW][f]->SetLabelSize(0.08);
692    (fSigma2[pW][f]->GetXaxis())->SetBinLabel(1,"#sigma^{2}");
693    fHistList->Add(fSigma2[pW][f]);
694    // chi^2:
695    fChi2[pW][f] = new TH1D(Form("%s, %s, %s",fChi2Name.Data(),pWeightsFlag[pW].Data(),sigmaFlag[f].Data()),"#chi^{2} (obtained from Minuit)",1,0,1);
696    fChi2[pW][f]->SetLabelSize(0.08);
697    (fChi2[pW][f]->GetXaxis())->SetLabelOffset(0.01);
698    (fChi2[pW][f]->GetXaxis())->SetBinLabel(1,"#chi^{2}");
699    fHistList->Add(fChi2[pW][f]);
700   } // end of for(Int_t f=0;f<2;f++) // sigma^2 not fixed or fixed
701   
702  } // end of for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++) // pW not used or used
703  
704  // book profile fFittingParameters whiuch will hols fitting parameters:
705  TString fFittingParametersName = "fFittingParameters";
706  fFittingParametersName += fAnalysisLabel->Data(); 
707  fFittingParameters = new TProfile(fFittingParametersName.Data(),"Parameters for fitting q-distribution",8,0,8);
708  fFittingParameters->SetLabelSize(0.05);
709  (fFittingParameters->GetXaxis())->SetBinLabel(1,"Treshold");
710  (fFittingParameters->GetXaxis())->SetBinLabel(2,"starting v_{n}");
711  (fFittingParameters->GetXaxis())->SetBinLabel(3,"min. v_{n}");
712  (fFittingParameters->GetXaxis())->SetBinLabel(4,"max. v_{n}");
713  (fFittingParameters->GetXaxis())->SetBinLabel(5,"starting #sigma^{2}");
714  (fFittingParameters->GetXaxis())->SetBinLabel(6,"min. #sigma^{2}");
715  (fFittingParameters->GetXaxis())->SetBinLabel(7,"max. #sigma^{2}");
716  (fFittingParameters->GetXaxis())->SetBinLabel(8,"plot or not?");
717  fHistList->Add(fFittingParameters);
718  
719 } // end of void AliFlowAnalysisWithFittingQDistribution::BookEverythingForDistributions()
720
721
722 //================================================================================================================================
723
724
725 void AliFlowAnalysisWithFittingQDistribution::DoFit(Bool_t useParticleWeights, Bool_t sigma2Fitted)
726 {
727  // do final fit to q-distribution
728  
729  // shortcuts for flags:
730  Int_t pW = (Int_t)(useParticleWeights);
731  Int_t s2F = (Int_t)(sigma2Fitted);
732  
733  for(Int_t f=0;f<2;f++)
734  {
735   if(!(fqDistribution[pW] && fSumOfParticleWeights[pW] && fIntFlow[pW][f] && fSigma2[pW][f])) 
736   { 
737    cout<<"WARNING: fqDistribution[pW] && fSumOfParticleWeights[pW] && fIntFlow[pW][f] && fSigma2[pw][f] is NULL in AFAWFQD::DoFit() !!!!"<<endl;
738    cout<<"pW = "<<pW<<endl;
739    cout<<"f  = "<<f<<endl;
740    exit(0);
741   }
742  }
743  
744  // average multiplicity and number of events:
745  Double_t AvM = fSumOfParticleWeights[pW]->GetMean(1);
746  //Int_t nEvts = (Int_t)fSumOfParticleWeights[pW]->GetEntries();
747  
748  // for fitting take into account only bins with at least fTreshold entries:
749  Int_t binMin = fqDistribution[pW]->FindFirstBinAbove(fTreshold); // to be improved (add setter for this)  
750  Int_t binMax = fqDistribution[pW]->FindLastBinAbove(fTreshold); // to be improved (add setter for this) 
751  Double_t binWidth = fqDistribution[pW]->GetBinWidth(4); // assuming that all bins have the same width 
752  if(binWidth == 0) 
753  {
754   cout<<"WARNING: binWidth == 0 in AFAWFQD::DoFit()"<<endl;
755   exit(0);
756  }
757  Double_t qmin = (binMin-1)*binWidth; 
758  Double_t qmax = (binMax)*binWidth;
759  Double_t ent = 0.; // number of entries between binMin and binMax:
760  for(Int_t b=binMin;b<=binMax;b++)
761  {
762   ent += fqDistribution[pW]->GetBinContent(b);
763  }
764  Double_t norm = binWidth*ent; // norm (assuming that all bins have the same width)
765
766  // fitting function:
767  TF1 *fittingFun = new TF1("fittingFun","[2]*(x/[1])*exp(-(x*x+[0]*[0])/(2.*[1]))*TMath::BesselI0(x*[0]/[1])",qmin,qmax); 
768  
769  fittingFun->SetParNames("v*sqrt{sum of particle weights}","sigma^2","norm");
770  fittingFun->SetParameters(fvStart*pow(AvM,0.5),fSigma2Start,norm); // to be improved (add setter for starting v)         
771  fittingFun->SetParLimits(0,fvMin*pow(AvM,0.5),fvMax*pow(AvM,0.5)); // to be improved (add setters for vmin and vmax)
772  
773  if(s2F == 0)
774  {
775   fittingFun->FixParameter(1,0.5);
776  } else
777    {
778     fittingFun->SetParLimits(1,fSigma2Min,fSigma2Max); // to be improved (add setters for sigma^2_min and sigma^2_max)          
779    }
780  fittingFun->FixParameter(2,norm);  
781
782  // fitting (do it only if # of entries >50): // to be improved (this is only a pragmatics fix to avoid TMinuit crash)
783  if(ent > 50)
784  {
785   fqDistribution[pW]->Fit("fittingFun","NQ","",qmin,qmax);
786  }
787  // results:
788  Double_t v = 0.; // integrated flow
789  Double_t vError = 0.; // error of integrated flow 
790  Double_t sigma2 = 0.; // sigma^2
791  Double_t sigma2Error = 0.; // error of sigma^2
792  Double_t chi2 = 0; // chi^2 from Minuit
793  
794  if(AvM)
795  { 
796   v = fittingFun->GetParameter(0)/pow(AvM,0.5);
797   vError = fittingFun->GetParError(0)/pow(AvM,0.5);
798   // store the results:
799   fIntFlow[pW][s2F]->SetBinContent(1,v); // s2F = sigma^2 fitted
800   fIntFlow[pW][s2F]->SetBinError(1,vError); // s2F = sigma^2 fitted
801  }
802  
803  if(s2F == 0) // sigma^2 not fitted, but fixed to 0.5
804  {
805   // sigma^2:
806   sigma2 = 0.5;
807   fSigma2[pW][0]->SetBinContent(1,sigma2);  
808   fSigma2[pW][0]->SetBinError(1,0.);
809   // chi^2:
810   chi2 = fittingFun->GetChisquare();
811   fChi2[pW][0]->SetBinContent(1,chi2);  
812   fChi2[pW][0]->SetBinError(1,0.);  
813  } else // sigma^2 fitted
814    {
815     // sigma^2:
816     sigma2 = fittingFun->GetParameter(1);
817     sigma2Error = fittingFun->GetParError(1);
818     fSigma2[pW][1]->SetBinContent(1,sigma2);  
819     fSigma2[pW][1]->SetBinError(1,sigma2Error);    
820     // chi^2:
821     chi2 = fittingFun->GetChisquare();
822     fChi2[pW][1]->SetBinContent(1,chi2);  
823     fChi2[pW][1]->SetBinError(1,0.);  
824    }
825  
826  if(fPlotResults)
827  {
828   // set ranges:
829   Int_t firstNonEmptyBin = fqDistribution[pW]->FindFirstBinAbove(0);
830   Double_t lowRange = fqDistribution[pW]->GetBinLowEdge(firstNonEmptyBin);
831   Int_t lastNonEmptyBin = fqDistribution[pW]->FindLastBinAbove(0);
832   Double_t upperRange = fqDistribution[pW]->GetBinLowEdge(lastNonEmptyBin+10);
833   (fqDistribution[pW]->GetXaxis())->SetRangeUser(lowRange,upperRange);
834   
835   if(s2F == 0)
836   { 
837    // to be improved (perhaps there is a better way to implement this?)
838    fqDistribution[pW]->SetFillColor(16);  
839    fqDistribution[pW]->SetTitle("Fitted q-distribution");
840    fqDistribution[pW]->Draw("");
841    fLegend->AddEntry(fqDistribution[pW],"q-distribution","f");
842    TF1 *fittingFunTemp = (TF1*)fittingFun->Clone("fittingFunTemp");
843    fittingFunTemp->SetLineColor(4); // 4 = blue color
844    fittingFunTemp->Draw("SAME"); 
845    fLegend->AddEntry("fittingFunTemp","#sigma^{2} fixed","l");
846    fLegend->Draw("SAME");      
847   } else
848     {    
849      fittingFun->SetLineColor(2); // 2 = red color   
850      fittingFun->Draw("SAME");
851      fLegend->AddEntry("fittingFun","#sigma^{2} fitted","l");    
852     } 
853  } // end of if(fPlotResults)
854
855 } // end of void AliFlowAnalysisWithFittingQDistribution::DoFit(Bool_t useParticleWeights)
856
857
858 //================================================================================================================================ 
859
860
861 void AliFlowAnalysisWithFittingQDistribution::FillCommonHistResultsIntFlow(Bool_t useParticleWeights, Bool_t sigma2Fitted)
862 {
863  // fill in AliFlowCommonHistResults histograms relevant for 'NONAME' integrated flow (to be improved (name))
864  
865  // shortcuts for the flags:
866  Int_t pW = (Int_t)(useParticleWeights); // 0 = pWeights not useed, 1 = pWeights used
867  Int_t s2F = (Int_t)(sigma2Fitted); // 0 = sigma^2 not fitted (but fixed to 0.5), 1 = sigma^2 fitted
868  
869  if(!fIntFlow[pW][s2F])
870  {
871   cout<<"WARNING: fIntFlow[pW][s2F] is NULL in AFAWFQD::FCHRIF() !!!!"<<endl;
872   cout<<"pW  = "<<pW<<endl;
873   cout<<"s2F = "<<s2F<<endl;
874   exit(0); 
875  }  
876  
877  if(!fSumOfParticleWeights[pW])
878  {
879   cout<<"WARNING: fSumOfParticleWeights[pW] is NULL in AFAWFQD::FCHRIF() !!!!"<<endl;
880   cout<<"pW = "<<pW<<endl;
881   exit(0);
882  }
883  
884  if(!(fCommonHistsResults))
885  {
886   cout<<"WARNING: fCommonHistsResults is NULL in AFAWFQD::FCHRIF() !!!!"<<endl; 
887   exit(0);
888  }
889   
890  // fill integrated flow:
891  Double_t v = fIntFlow[pW][s2F]->GetBinContent(1); 
892  Double_t vError = fIntFlow[pW][s2F]->GetBinError(1);
893  
894  fCommonHistsResults->FillIntegratedFlow(v,vError);   
895  
896  // fill chi (this chi stands for resolution, not to be confused with chi2 used before):
897  Double_t AvM = fSumOfParticleWeights[pW]->GetMean(1);
898  Double_t chi = AvM*pow(v,2); 
899  if(chi>=0)
900  {
901   fCommonHistsResults->FillChi(pow(chi,0.5));   
902   fCommonHistsResults->FillChiRP(pow(chi,0.5));   
903  }
904    
905 } // end of void AliFlowAnalysisWithFittingQDistribution::FillCommonHistResultsIntFlow(Bool_t useParticleWeights, Bool_t sigma2NotFixed) 
906
907
908 //================================================================================================================================ 
909
910
911 void AliFlowAnalysisWithFittingQDistribution::PrintFinalResultsForIntegratedFlow()
912 {
913  // print the final results for integrated flow on the screen
914  
915  // shortcuts: pW  = particle weights 
916  //            s2F = sigma^2 fitted 
917  
918  for(Int_t pW=0;pW<1+(Int_t)(fUsePhiWeights||fUsePtWeights||fUseEtaWeights);pW++)
919  {
920   if(!fSumOfParticleWeights[pW])
921   {
922    cout<<"WARNING: fSumOfParticleWeights[pW] is NULL in AFAWFQD::FCHRIF() !!!!"<<endl;
923    cout<<"pW = "<<pW<<endl;
924    exit(0);
925   }
926   for(Int_t s2F=0;s2F<2;s2F++)
927   {
928    if(!fIntFlow[pW][s2F])
929    {
930     cout<<"WARNING: fIntFlow[pW][s2F] is NULL in AFAWFQD::FCHRIF() !!!!"<<endl;
931     cout<<"pW  = "<<pW<<endl;
932     cout<<"s2F = "<<s2F<<endl;
933     exit(0); 
934    }
935   }  
936  }  
937  
938  if(!(fCommonHistsResults))
939  {
940   cout<<"WARNING: fCommonHistsResults is NULL in AFAWFQD::FCHRIF() !!!!"<<endl; 
941   exit(0);
942  }
943  
944  // shortcut for the harmonic:
945  Int_t n = (Int_t)(fCommonHists->GetHarmonic())->GetBinContent(1); 
946
947  // printing:
948  cout<<" "<<endl;
949  cout<<"***************************************"<<endl;
950  cout<<"***************************************"<<endl;
951  cout<<"      integrated flow by fitting "<<endl;
952  cout<<"           q-distribution:      "<<endl;
953  if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)
954  {
955   cout<<"           (with weights)       "<<endl;
956  } else
957    {
958     cout<<"          (without weights)       "<<endl;
959    }   
960  cout<<endl;
961
962  if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)
963  {
964   cout<<"1.) sigma^2 not fitted: "<<endl;
965   cout<<"  v_"<<n<<"{FQD} = "<<fIntFlow[1][0]->GetBinContent(1)<<" +/- "<<fIntFlow[1][0]->GetBinError(1)<<endl;
966   cout<<"  sigma^2 = 0.5 +/- 0 "<<endl; 
967   cout<<"  chi^2 = "<<fChi2[1][0]->GetBinContent(1)<<" (Minuit)"<<endl; 
968   cout<<" "<<endl;   
969   cout<<"2.) sigma^2 fitted: "<<endl;
970   cout<<"  v_"<<n<<"{FQD} = "<<fIntFlow[1][1]->GetBinContent(1)<<" +/- "<<fIntFlow[1][1]->GetBinError(1)<<endl;
971   cout<<"  sigma^2 = "<<fSigma2[1][1]->GetBinContent(1)<<" +/- "<<fSigma2[1][1]->GetBinError(1)<<endl; 
972   cout<<"  chi^2 = "<<fChi2[1][1]->GetBinContent(1)<<" (Minuit)"<<endl; 
973   cout<<" "<<endl; 
974   cout<<"      nEvts = "<<fSumOfParticleWeights[1]->GetEntries()<<", AvM = "<<fSumOfParticleWeights[1]->GetMean()<<endl;
975   cout<<" "<<endl;
976  } else
977    { 
978     cout<<"1.) sigma^2 not fitted: "<<endl;
979     cout<<endl;
980     cout<<"  v_"<<n<<"{FQD} = "<<fIntFlow[0][0]->GetBinContent(1)<<" +/- "<<fIntFlow[0][0]->GetBinError(1)<<endl;
981     cout<<"  sigma^2 = 0.5 +/- 0 "<<endl; 
982     cout<<"  chi^2 = "<<fChi2[0][0]->GetBinContent(1)<<" (Minuit)"<<endl;  
983     cout<<" "<<endl;   
984     cout<<"2.) sigma^2 fitted: "<<endl;
985     cout<<endl;
986     cout<<"  v_"<<n<<"{FQD} = "<<fIntFlow[0][1]->GetBinContent(1)<<" +/- "<<fIntFlow[0][1]->GetBinError(1)<<endl;
987     cout<<"  sigma^2 = "<<fSigma2[0][1]->GetBinContent(1)<<" +/- "<<fSigma2[0][1]->GetBinError(1)<<endl; 
988     cout<<"  chi^2 = "<<fChi2[0][1]->GetBinContent(1)<<" (Minuit)"<<endl; 
989     cout<<" "<<endl;  
990     cout<<"      nEvts = "<<fSumOfParticleWeights[0]->GetEntries()<<", AvM = "<<fSumOfParticleWeights[0]->GetMean()<<endl;
991     cout<<" "<<endl;
992    }
993     
994  cout<<"***************************************"<<endl;
995  cout<<"***************************************"<<endl; 
996  cout<<endl;
997   
998 } // end of void AliFlowAnalysisWithFittingQDistribution::PrintFinalResultsForIntegratedFlow()
999
1000
1001 //================================================================================================================================ 
1002
1003
1004 void AliFlowAnalysisWithFittingQDistribution::StoreFittingParameters()
1005 {
1006  // store fitting parameters in profile fFittingParameters
1007  
1008  // Binning of fFittingParameters is organized as follows:
1009  // 1st bin: fTreshold
1010  // 2nd bin: fvStart
1011  // 3rd bin: fvMin
1012  // 4th bin: fvMax
1013  // 5th bin: fSigma2Start
1014  // 6th bin: fSigma2Min
1015  // 7th bin: fSigma2Max
1016  
1017  if(!fFittingParameters)
1018  {
1019   cout<<"WARNING: fFittingParameters is NULL in AFAWFQD::SFP() !!!!"<<endl;
1020   exit(0);
1021  }
1022  
1023  fFittingParameters->Reset();
1024  fFittingParameters->Fill(0.5,fTreshold);
1025  fFittingParameters->Fill(1.5,fvStart);
1026  fFittingParameters->Fill(2.5,fvMin);
1027  fFittingParameters->Fill(3.5,fvMax);
1028  fFittingParameters->Fill(4.5,fSigma2Start);
1029  fFittingParameters->Fill(5.5,fSigma2Min);
1030  fFittingParameters->Fill(6.5,fSigma2Max);
1031  fFittingParameters->Fill(7.5,fPlotResults);
1032  
1033 } // end of void AliFlowAnalysisWithFittingQDistribution::StoreFittingParameters()
1034
1035
1036 //================================================================================================================================ 
1037
1038
1039 void AliFlowAnalysisWithFittingQDistribution::AccessFittingParameters()
1040 {
1041  // access fitting parameters:
1042  
1043  if(!fFittingParameters)
1044  {
1045   cout<<"WARNING: fFittingParameters is NULL in AFAWFQD::AFP() !!!!"<<endl;
1046   exit(0);
1047  }
1048  
1049  fTreshold = fFittingParameters->GetBinContent(1);
1050  fvStart = fFittingParameters->GetBinContent(2);
1051  fvMin = fFittingParameters->GetBinContent(3);
1052  fvMax = fFittingParameters->GetBinContent(4);
1053  fSigma2Start = fFittingParameters->GetBinContent(5);
1054  fSigma2Min = fFittingParameters->GetBinContent(6);
1055  fSigma2Max = fFittingParameters->GetBinContent(7);
1056  fPlotResults = (Bool_t) fFittingParameters->GetBinContent(8);
1057  
1058 } // end of void AliFlowAnalysisWithFittingQDistribution::AccessFittingParameters()