]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWG2/SPECTRA/LambdaK0PbPb/MultYields2QA.C
e3497883797bb5a1bf877e2d6038672fd8bfb9e7
[u/mrichter/AliRoot.git] / PWG2 / SPECTRA / LambdaK0PbPb / MultYields2QA.C
1 #include "AliMassFitControl.h"
2
3 // Two functions now.
4 // MultYields3 takes a 3D histogram and calls the 2D function MultYields2
5 // MultYields2 could be called directly if we only have a 2D histogram
6
7 void MultYields3(TH3F *PtMassMult, Int_t particleMode, Int_t MultBin, Char_t* label){
8
9   // Leave open possibility to choose different values depending on MultBin - LSB
10   Float_t MultLo[1] = {0};
11   Float_t MultHi[1] = {300};
12
13   //Make 2D projections from the 3D histogram
14   PtMassMult->GetZaxis()->SetRange(MultLo[MultBin],MultHi[MultBin]);
15   TH2F* hParMass = (TH2F*)PtMassMult->Project3D("XY");// FIX:MF
16   hParMass->SetTitle("PtMass");
17   
18   //hParMass->Draw(); // Drawing invokes default c1 canvas making it inaccessible in MultYields2 - TO FIX
19   MultYields2(hParMass, particleMode,MultBin,label);
20   
21 }
22
23 void MultYields2QA(TH2F *hParMass, Int_t particleMode, Int_t ihist,Int_t Nev = 1, Int_t MultBin, Char_t* label){
24
25   hParMass->Draw();
26   ////////////////////////// minhist is needed by AliMassFitControl.h (in case the minimum on the parameter's axis is arbitrary value); 
27   Double_t minhist = 0;
28   minhist = hParMass->GetBinLowEdge(1);
29   cout<<" value of the first bin of parameter : " <<minhist << endl;
30
31   // Do .L on command line first
32   //gROOT->LoadMacro("macros/PtMassAna2.C");
33   /* Modifications to produce a single uncorrected spectrum in a particular mult. bin for K0 or Lambda
34      Old ratio code preserved in MultYieldsRatio.C
35      Dec 2003 */
36   Char_t* part; //for name of particle
37   Float_t BR; //branching ratio - check them
38   if (particleMode==0){
39     part = "K0";
40     BR=0.686;
41   } else if (particleMode==1){
42     part = "Lambda";
43     BR=0.639;
44   } else if (particleMode==2){
45     part = "AntiLambda";
46     BR=0.639;
47   } else if (particleMode==3){
48     part = "Lambda+antiLambda";
49     BR=0.639;
50   } else if (particleMode==4){
51     part = "Xi";
52     BR=1.0; // Should be Lam ->p pi Br
53   } else if (particleMode ==6){
54     part = "Omega";
55     BR=1.0; // Should be Lam->p pi  * Om -> K Lam
56   }
57   
58
59
60   TString title[1]={"MinimumBias"};
61   //Make 2D projections from the 3D histogram
62   //Minbias (i.e. everything)
63   
64   TObjArray *controllerArray = new TObjArray(40,1); // 2nd arg means can count from 1!
65
66   //Here probably need switch-case, depending on mult bin and particle
67   // LoPt, HiPt, polynomial order, rebinning factor
68
69   /////  LAMBDA and LAMBDA+ANTILAMBDA (Combination)
70
71
72
73   if(particleMode == 1 || particleMode == 3){ // Lambda or Lambda+Anti-Lambda
74     if(ihist ==0){
75       /*     controllerArray->AddLast(new AliMassFitControl(0.2,0.3, 2,2, 1.095,1.17)); //1
76              controllerArray->AddLast(new AliMassFitControl(0.3,0.4, 2,2, 1.095,1.17)); //2
77              controllerArray->AddLast(new AliMassFitControl(0.4,0.5, 2,2, 1.095,1.17)); //3
78              controllerArray->AddLast(new AliMassFitControl(0.5,0.6, 2,2, 1.095,1.17)); //4
79              controllerArray->AddLast(new AliMassFitControl(0.6,0.7, 2,2, 1.095,1.17)); //5
80              controllerArray->AddLast(new AliMassFitControl(0.7,0.8, 2,2, 1.1,1.17)); //6
81       */      controllerArray->AddLast(new AliMassFitControl(0.8,0.9, 2,2, 1.1,1.17)); //7
82       controllerArray->AddLast(new AliMassFitControl(0.9,1.0, 2,2, 1.1,1.18)); //8
83       controllerArray->AddLast(new AliMassFitControl(1.0,1.1, 2,2, 1.1,1.17)); //9
84       controllerArray->AddLast(new AliMassFitControl(1.1,1.2, 2,2, 1.1,1.14)); //10
85       controllerArray->AddLast(new AliMassFitControl(1.2,1.3, 2,2, 1.1,1.14)); //11 
86       controllerArray->AddLast(new AliMassFitControl(1.3,1.4, 2,2, 1.1,1.17)); //12
87       controllerArray->AddLast(new AliMassFitControl(1.4,1.5, 2,2, 1.095,1.17)); //13
88       controllerArray->AddLast(new AliMassFitControl(1.5,1.6, 2,2, 1.095,1.17)); //14
89       controllerArray->AddLast(new AliMassFitControl(1.6,1.7, 2,2, 1.095,1.17)); //15
90       controllerArray->AddLast(new AliMassFitControl(1.7,1.8, 2,2, 1.095,1.18)); //16
91       controllerArray->AddLast(new AliMassFitControl(1.8,1.9, 2,2, 1.095,1.17)); //17
92       controllerArray->AddLast(new AliMassFitControl(1.9,2.0, 2,2, 1.095,1.17)); //18
93       controllerArray->AddLast(new AliMassFitControl(2.0,2.2, 2,2, 1.095,1.17)); //19
94       controllerArray->AddLast(new AliMassFitControl(2.2,2.4, 2,2, 1.095,1.17)); //20
95       controllerArray->AddLast(new AliMassFitControl(2.4,2.6, 2,2, 1.095,1.17)); //21
96       controllerArray->AddLast(new AliMassFitControl(2.6,2.8, 2,2, 1.095,1.17)); //22
97       controllerArray->AddLast(new AliMassFitControl(2.8,3.0, 2,2, 1.095,1.17)); //23
98       controllerArray->AddLast(new AliMassFitControl(3.0,3.2, 2,2, 1.095,1.17)); //24
99       controllerArray->AddLast(new AliMassFitControl(3.2,3.4, 2,2, 1.095,1.17)); //25
100       controllerArray->AddLast(new AliMassFitControl(3.4,3.6, 2,2, 1.095,1.17)); //26
101       controllerArray->AddLast(new AliMassFitControl(3.6,3.8, 2,2, 1.095,1.17)); //27
102       controllerArray->AddLast(new AliMassFitControl(3.8,4.0, 2,2, 1.095,1.167)); //28
103       controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 2,2, 1.095,1.17)); //29
104       controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 2,2, 1.095,1.17)); //30
105       controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 2,2, 1.095,1.17)); //31  bin05
106       controllerArray->AddLast(new AliMassFitControl(5.5,6.5, 2,2, 1.095,1.17)); //32
107       controllerArray->AddLast(new AliMassFitControl(6.5,8.0, 2,2, 1.095,1.17)); //33
108       //controllerArray->AddLast(new AliMassFitControl(8.0,12.0, 2,2, 1.095,1.17));//34
109       
110
111     }
112     if(ihist == 1){
113       cout << "histogram : " <<1<<endl;
114       for(int i = 0; i<49; i++)
115         //     for(int i = 0; i<33; i++) //for 05 centrality
116         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
117         controllerArray->AddLast(new AliMassFitControl(minhist,0.1+i*0.2,0.1+(i+1.0)*0.2, 2,1, 1.095,1.15));
118     }
119     if(ihist == 2){
120       cout << "histogram : " <<2<<endl;
121       for(int i = 0; i<49; i++)
122         //      for(int i = 0; i<33; i++) //for 05 centrality
123         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
124         controllerArray->AddLast(new AliMassFitControl(minhist,0.1+i*0.2,0.1+(i+1.0)*0.2, 2,1, 1.095,1.17));
125     }
126     if(ihist == 3){
127       cout << "histogram : " <<3<<endl;
128       for(int i = 0; i<=40; i++)
129         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
130         controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 1.095,1.17));
131     }
132     if(ihist == 4){
133       cout << "histogram : " <<4<<endl;
134       for(int i = 0; i<=40; i++)
135         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
136         controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 1.095,1.17));
137     }
138     if(ihist == 5){
139       cout << "histogram : " <<5<<endl;
140       for(int i = 0; i<=33; i++)
141         //      if(i==11 || i==15 || i==16)
142         //        controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 0.43,0.57));
143         //        else
144         //controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 0.44,0.56));
145         controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 1.095,1.17));   //00
146     }
147     if(ihist == 6){
148       cout << "histogram : " <<6<<endl;
149       for(int i = 1; i<40; i++){
150         //   if(i==23) 
151         controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 1.095,1.17));
152         /*        if(  i==17 || i==18 || i==23) 
153                   controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 0.446,0.55));
154                   else
155                   controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 0.44,0.56));
156         */        }
157     }
158   }
159   /// ANTI LAMBDA ---->
160   else if (particleMode == 2){ // Anti-Lambdas
161     controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 0,4, 1.095,1.17));
162     controllerArray->AddLast(new AliMassFitControl(5.5,6.0, 0,6, 1.085,1.17));
163   } // end if anti-Lambda
164   else if (particleMode == 0){ // K0s case
165     if(ihist == 0){
166       /*      controllerArray->AddLast(new AliMassFitControl(0.2,0.3, 2,2, 0.428,0.56)); //1
167               controllerArray->AddLast(new AliMassFitControl(0.3,0.4, 2,2, 0.428,0.56)); //2
168               controllerArray->AddLast(new AliMassFitControl(0.4,0.5, 2,2, 0.428,0.56)); //3
169               controllerArray->AddLast(new AliMassFitControl(0.5,0.6, 2,2, 0.428,0.46)); //4
170               controllerArray->AddLast(new AliMassFitControl(0.6,0.7, 2,2, 0.428,0.56)); //5
171               controllerArray->AddLast(new AliMassFitControl(0.7,0.8, 2,2, 0.438,0.56)); //6
172       */      controllerArray->AddLast(new AliMassFitControl(0.8,0.9, 2,2, 0.438,0.56)); //7
173       controllerArray->AddLast(new AliMassFitControl(0.9,1.0, 2,2, 0.438,0.56)); //8
174       controllerArray->AddLast(new AliMassFitControl(1.0,1.1, 2,2, 0.438,0.56)); //9
175       controllerArray->AddLast(new AliMassFitControl(1.1,1.2, 2,2, 0.438,0.56)); //10
176       controllerArray->AddLast(new AliMassFitControl(1.2,1.3, 2,2, 0.438,0.56)); //11 
177       controllerArray->AddLast(new AliMassFitControl(1.3,1.4, 2,2, 0.438,0.56)); //12
178       controllerArray->AddLast(new AliMassFitControl(1.4,1.5, 2,2, 0.438,0.56)); //13
179       controllerArray->AddLast(new AliMassFitControl(1.5,1.6, 2,2, 0.438,0.56)); //14
180       controllerArray->AddLast(new AliMassFitControl(1.6,1.7, 2,2, 0.438,0.56)); //15
181       controllerArray->AddLast(new AliMassFitControl(1.7,1.8, 2,2, 0.438,0.56)); //16
182       controllerArray->AddLast(new AliMassFitControl(1.8,1.9, 2,2, 0.438,0.56)); //17
183       controllerArray->AddLast(new AliMassFitControl(1.9,2.0, 2,2, 0.438,0.56)); //18
184       controllerArray->AddLast(new AliMassFitControl(2.0,2.2, 2,2, 0.4385,0.535)); //19
185       controllerArray->AddLast(new AliMassFitControl(2.2,2.4, 2,2, 0.438,0.56)); //20
186       controllerArray->AddLast(new AliMassFitControl(2.4,2.6, 2,2, 0.438,0.56)); //21
187       controllerArray->AddLast(new AliMassFitControl(2.6,2.8, 2,2, 0.428,0.56)); //22
188       controllerArray->AddLast(new AliMassFitControl(2.8,3.0, 2,2, 0.428,0.56)); //23
189       controllerArray->AddLast(new AliMassFitControl(3.0,3.2, 2,2, 0.428,0.56)); //24
190       controllerArray->AddLast(new AliMassFitControl(3.2,3.4, 2,2, 0.428,0.56)); //25
191       controllerArray->AddLast(new AliMassFitControl(3.4,3.6, 2,2, 0.428,0.56)); //26
192       controllerArray->AddLast(new AliMassFitControl(3.6,3.8, 2,2, 0.428,0.56)); //27
193       controllerArray->AddLast(new AliMassFitControl(3.8,4.0, 2,2, 0.428,0.56)); //28
194       controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 2,2, 0.428,0.56)); //29
195       controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 2,2, 0.428,0.56)); //30
196       controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 2,2, 0.428,0.56)); //31
197       controllerArray->AddLast(new AliMassFitControl(5.5,6.5, 2,2, 0.428,0.56)); //32
198       controllerArray->AddLast(new AliMassFitControl(6.5,8.0, 2,2, 0.428,0.56)); //33
199       controllerArray->AddLast(new AliMassFitControl(8.0,12.0, 2,2, 0.428,0.56));//34
200       
201     }
202     if(ihist == 1){
203       cout << "histogram : " <<1<<endl;
204       for(int i = 0; i<49; i++)
205         //     for(int i = 0; i<33; i++) //for 05 centrality
206         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
207         controllerArray->AddLast(new AliMassFitControl(minhist,0.1+i*0.2,0.1+(i+1.0)*0.2, 2,1, 0.445,0.56));
208     }
209     if(ihist == 2){
210       cout << "histogram : " <<2<<endl;
211       for(int i = 0; i<42; i++){
212         //      for(int i = 0; i<33; i++) //for 05 centrality
213         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
214         if(i==2)controllerArray->AddLast(new AliMassFitControl(minhist,0.1+i*0.2,0.1+(i+1.0)*0.2, 1,1, 0.445,0.56));
215         else
216           controllerArray->AddLast(new AliMassFitControl(minhist,0.1+i*0.2,0.1+(i+1.0)*0.2, 2,1, 0.445,0.56));
217       }
218     }
219     if(ihist == 3){
220       cout << "histogram : " <<3<<endl;
221       //for(int i = 0; i<=28; i++) //bin05
222       for(int i = 0; i<=40; i++)
223         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
224         controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.445,0.56));
225     }
226     if(ihist == 4){
227       cout << "histogram : " <<4<<endl;
228       for(int i = 0; i<=40; i++)
229         //      controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.44,0.56));
230         controllerArray->AddLast(new AliMassFitControl(minhist,i,i+1.0, 2,1, 0.445,0.56));
231     }
232     if(ihist == 5){
233       cout << "histogram : " <<5<<endl;
234       for(int i = 0; i<=33; i++)
235         //      if(i==11 || i==15 || i==16)
236         //        controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 0.43,0.57));
237         //        else
238         //controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 0.44,0.56));
239         controllerArray->AddLast(new AliMassFitControl(minhist,0+i*0.03,(i+1)*0.03, 2,1, 0.43,0.56));   //00
240     }
241     if(ihist == 6){
242       cout << "histogram : " <<6<<endl;
243       for(int i = 1; i<40; i++){
244         //   if(i==23) 
245         controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 0.445,0.55));
246         /*        if(  i==17 || i==18 || i==23) 
247                   controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 0.446,0.55));
248                   else
249                   controllerArray->AddLast(new AliMassFitControl(minhist,0.998 +i*0.00005,0.998 +(i+1)*0.00005, 2,1, 0.44,0.56));
250         */        }
251     }
252   }  else if (particleMode == 4) { //Xi case
253     //controllerArray->AddLast(new AliMassFitControl(0.5,0.7, 1,1, 1.28,1.45)); //signal not visible with 
254     controllerArray->AddLast(new AliMassFitControl(0.7,0.9, 1,1, 1.285,1.45));
255     controllerArray->AddLast(new AliMassFitControl(0.9,1.1, 1,1, 1.285,1.45));
256     controllerArray->AddLast(new AliMassFitControl(1.1,1.3, 1,1, 1.27,1.45));
257     controllerArray->AddLast(new AliMassFitControl(1.3,1.5, 1,1, 1.27,1.45));
258     controllerArray->AddLast(new AliMassFitControl(1.5,1.7, 1,1, 1.27,1.45));
259     controllerArray->AddLast(new AliMassFitControl(1.7,2.0, 1,1, 1.27,1.45));
260     controllerArray->AddLast(new AliMassFitControl(2.0,2.5, 1,1, 1.27,1.44));
261     controllerArray->AddLast(new AliMassFitControl(2.5,3.0, 1,1, 1.27,1.45));
262     controllerArray->AddLast(new AliMassFitControl(3.0,3.5, 1,1, 1.27,1.44));
263     controllerArray->AddLast(new AliMassFitControl(3.5,4.0, 1,1, 1.27,1.45));
264     controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 2,1, 1.27,1.44));
265     controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 2,1, 1.27,1.45));
266   }
267   controllerArray->Sort();
268
269   // Make the proper label to pass in for use in labelling the diagnostics
270   // saved canvas files and histos
271   Char_t fulllabel[80];
272   //sprintf(fulllabel,"%s%s",title[MultBin].Data(),label);
273   sprintf(fulllabel,"%s",label);
274
275   //Slice up projection into mass histograms to extract yield
276   TH1F* hYield = (TH1F*)PtMassAna2(hParMass,particleMode,ihist,controllerArray->GetEntries(),controllerArray,fulllabel);
277
278
279   // CORRECTIONS Nev -comming from "FitSpectrum.C" (events with Zvertex < |10cm|) 
280   hYield->Scale(1.0/Nev); //Divide by the number of events
281
282   //hYield->Scale(Veff[MultBin]); //Multiply by the vertex efficiency effectively increaing number of events 
283   //(since Veff<1) therefore decreases yield
284   //hYield->Scale(1.0/BR);  //Divide by branching ratio (again increases yield since BR<1)
285   //hYield->Scale(1.0/(2*TMath::Pi())); // Always plot 1/2pi ...
286
287   Char_t yieldTitle[80];
288   //sprintf(yieldTitle,"Uncorrected %s yield: %s",part,title[MultBin].Data());
289   sprintf(yieldTitle,"Uncorrected %s yield",part);
290   hYield->SetTitle(yieldTitle);
291   if(ihist == 0)hYield->SetXTitle("p_{t} / [GeV/c]");
292   if(ihist == 1)hYield->SetXTitle("DCA / [cm]");
293   if(ihist == 2)hYield->SetXTitle("DCA / [cm]");
294   if(ihist == 3)hYield->SetXTitle("Radius / [cm]");
295   if(ihist == 4)hYield->SetXTitle("Decay Lenght / [cm]");
296   if(ihist == 5)hYield->SetXTitle("V0 Daughters / [cm]");
297   if(ihist == 6)hYield->SetXTitle("Cos of pointing angle");
298   hYield->SetYTitle("1/Nev.dN/dp_{t}");
299
300   // Create plots
301
302   Char_t fileNameBase[80];
303   if(ihist == 0)sprintf(fileNameBase,"Masses%s",part);
304   else sprintf(fileNameBase,"Masses%s",label);
305   Char_t fileNamePng[80];
306   sprintf(fileNamePng,"%s.png",fileNameBase);
307   Char_t fileNameEps[80];
308   sprintf(fileNameEps,"%s.eps",fileNameBase);
309   Char_t fileNamePdf[80];
310   sprintf(fileNamePdf,"%s.pdf",fileNameBase);
311
312   c1->SaveAs(fileNamePng);
313   c1->SaveAs(fileNameEps);
314   c1->SaveAs(fileNamePdf);
315
316   //c1->Clear();
317
318   //c1->SetLogy();
319   TCanvas *cYield = new TCanvas("Yield","Corrected Yield",600,400);
320   cYield->cd();
321   //cYield->SetLogy();
322   hYield->SetStats(kFALSE);
323   hYield->Draw();
324   // cYield->SetLogy();
325   cYield->Update();
326
327   //  hRC_MB->SetMarkerStyle(20);
328   //  hRC_MB->SetMarkerColor(4);
329   //  hRC_MB->Scale(NBinMB/NBin3);
330   Char_t fnametext[80];
331   if(ihist == 0)sprintf(fnametext,"Yield%s",part);
332   else sprintf(fnametext,"Yield%s",label);
333   Char_t fnamePng[80];
334   sprintf(fnamePng,"%s.png",fnametext);
335   c1->SaveAs(fnamePng);
336   Char_t fnameEps[80];
337   sprintf(fnameEps,"%s.eps",fnametext);
338   c1->SaveAs(fnameEps);
339
340   // This section for yield scaled by number of binary collisions.
341   // Could add array of values and do scaling according to 'MultBin' index
342   TH1F* hScYield = hYield->Clone("ScYield");
343   Char_t scYieldTitle[80];
344   //sprintf(scYieldTitle,"<N_{bin}> Scaled %s",hYield->GetTitle());
345   //hScYield->SetTitle(scYieldTitle);
346   //SCALING for scaled yield only  divide by mean Nbin (scaled yield is therefore smaller)
347   //hScYield->Scale(1/NBin[MultBin]);
348
349   Char_t fnameRoot[80];
350   sprintf(fnameRoot,"%s.root",fnametext);
351   TFile *YieldFile = new TFile(fnameRoot,"RECREATE");
352   hYield->Write();
353   hScYield->Write();
354   YieldFile->Close();
355
356   controllerArray->Delete();
357 }