memory consumption measures only memory used by preprocessor
[u/mrichter/AliRoot.git] / PWG3 / vertexingOld / AliD0toKpiPlots.C
CommitLineData
9af36e3a 1void AliD0toKpiPlots(const Char_t *inName="AliD0toKpi.root",
2 const Char_t *outName="D0histograms.root") {
3 //--------------------------------------------------------------------------
4 // This macro histograms many variables of D0->Kpi candidates
5 //
6 // Andrea Dainese, andrea.dainese@lnl.infn.it
7 //--------------------------------------------------------------------------
8
c75b58db 9 gSystem->Load("libANALYSIS");
9af36e3a 10 gSystem->Load("libAOD.so");
11 gSystem->Load("libPWG3base.so");
3d1d2038 12 gSystem->Load("libPWG3vertexingOld.so");
9af36e3a 13
14 // set of cuts
15 Double_t D0Cuts[9] = {0.1, // mass [GeV]
16 1000000., // dca [micron]
17 1.1, // cosThetaStar
18 0., // pT K [GeV/c]
19 0., // pT Pi [GeV/c]
20 100000., // d0K upper [micron]
21 100000., // d0Pi upper [micron]
22 10000000000., // d0d0 [micron^2]
23 -1.1}; // cosThetaPointing
24
25 // number of events (for normalization)
26 Bool_t normalize = kFALSE;
27 Double_t events = 1.;
28
29
30
31 // define histograms
32 TH1F *hptK = new TH1F("hptK","\"K\" p_{t} distribution",50,0,10);
33 hptK->SetXTitle("p_{t} [GeV]");
34
35 TH1F *hptPi = new TH1F("hptPi","\"#pi\" p_{t} distribution",50,0,10);
36 hptPi->SetXTitle("p_{t} [GeV]");
37
38 TH1F *hDCA = new TH1F("hDCA","DCA",50,0,1000);
39 hDCA->SetXTitle("dca [#mu m]");
40
41 TH1F *hptD0 = new TH1F("hptD0","D^{0} p_{t} distribution",40,0,40);
42 hptD0->SetXTitle("p_{t} [GeV]");
43
44 TH1F *hyD0 = new TH1F("hyD0","D^{0} rapidity distribution",50,-2,2);
45 hyD0->SetXTitle("y");
46
47 TH1F *hCPtaD0 = new TH1F("hCPtaD0","cosine of pointing angle distribution",100,-1,1);
48 hCPtaD0->SetXTitle("cos #theta_{point}");
49
50 TH1F *hCPtaXY = new TH1F("hCPtaXY","cosine of pointing angle in (x,y) plane",100,-1,1);
51 hCPtaXY->SetXTitle("cos #theta_{point}");
52
53 TH1F *hCts = new TH1F("hCts","cosine of decay angle",50,-1.2,1.2);
54 hCts->SetXTitle("cos #theta^{*}");
55
56 TH2F *hCtsVsPtK = new TH2F("hCtsVsPtK","cosine of decay angle VS \"K\" p_{t}",50,0,5,50,-1,1);
57 hCtsVsPtK->SetYTitle("cos #theta^{*}");
58 hCtsVsPtK->SetXTitle("p_{t} [GeV]");
59
60 TH1F *hd0d0 = new TH1F("hd0d0","Product of the impact parameters",100,-100000,100000);
61 hd0d0->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]");
62
63 TH1F *hd0K = new TH1F("hd0K","Impact parameter of \"K\"",100,-5000,5000);
64 hd0K->SetXTitle("d_{0}^{K} [#mu m]");
65
66 TH1F *hd0Pi = new TH1F("hd0Pi","Impact parameter of \"#pi\"",100,-5000,5000);
67 hd0Pi->SetXTitle("d_{0}^{#pi} [#mu m]");
68
69 TH2F *hCPtaVsd0d0 = new TH2F("hCPtaVsd0d0","cos #theta_{point} vs d_{0}^{K} #times d_{0}^{#pi}",100,-100000,100000,100,-1,1);
70 hCPtaVsd0d0->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]");
71 hCPtaVsd0d0->SetYTitle("cos #theta_{point}");
72
73 TH2F *hCPtaVsd0d0zoom = new TH2F("hCPtaVsd0d0zoom","cos #theta_{point} vs d_{0}^{K} #times d_{0}^{#pi}",100,-100000,0,100,.9,1);
74 hCPtaVsd0d0zoom->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]");
75 hCPtaVsd0d0zoom->SetYTitle("cos #theta_{point}");
76
77 TH2F *hd0d0VSptD0 = new TH2F("hd0d0VSptD0","d_{0}^{K} #times d_{0}^{#pi} VS D^{0} p_{t}",50,0,25,100,-120000,120000);
78 hd0d0VSptD0->SetYTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]");
79 hd0d0VSptD0->SetXTitle("D^{0} p_{t} [GeV]");
80
81 TH1F *hMass = new TH1F("hMass","Invariant mass distribution",50,1.765,1.965);
82 hMass->SetXTitle("M[K,#pi] [GeV]");
83
84 TH2F *hArm = new TH2F("hArm","Armenteros plot",50,-2,2,50,0,1);
85 hArm->SetXTitle("#alpha");
86 hArm->SetYTitle("q_{t}");
87
88 // open input file and get tree
89 TFile *inFile = TFile::Open(inName);
90
91 TTree *treeD0 = (TTree*)inFile->Get("TreeD0");
92 AliD0toKpi *D = 0;
93 treeD0->SetBranchAddress("D0toKpi",&D);
94 Int_t entries = (Int_t)treeD0->GetEntries();
95
96 printf("+++\n+++ Number of D0 in tree: %d\n+++\n",entries);
97
98 Double_t MD0,MD0bar,ctsD0,ctsD0bar,ctsPiD0,ctsPiD0bar;
99 Double_t WgtD0,WgtD0bar;
100 Double_t sampleABC=0.;
101 Int_t okD0=0,okD0bar=0;
102 Int_t nSel = 0;
103 Int_t ptbin;
104
105 // loop on D0
106 for(Int_t i=0; i<entries; i++) {
107 if(i%10000==0) printf(" candidate %d of %d\n",i,entries);
108
109 // get event from tree
110 treeD0->GetEvent(i);
111 //--- select the PID strategy & compute weights
112 // D->ApplyPID("TOFparam_PbPb");
113 // D->ComputeWgts();
114 // get weights for the three samples A+B+C
115 // D->GetWgts(WgtD0,WgtD0bar,"ABC");
116 WgtD0 = 1.; WgtD0bar = 1.;
117
118 // normalize to 1 event
119 if(normalize) { WgtD0 /= events; WgtD0bar /= events; }
120
121 // check if candidate passes selection (as D0 or D0bar)
122 D->Select(D0Cuts,okD0,okD0bar);
123
124 // set weights to 0 if the candidate doesn't pass selection
125 if(!okD0) WgtD0=0.;
126 if(!okD0bar) WgtD0bar=0.;
127 if(okD0 || okD0bar) nSel++;
128
129 // count selected candidates
130 sampleABC += WgtD0 + WgtD0bar;
131
132 // inv mass and cosThetaStar
133 D->InvMass(MD0,MD0bar);
134 D->CosThetaStar(ctsD0,ctsD0bar);
135
136 // fill histograms
137 hptK->Fill(D->PtChild(1),WgtD0);
138 hptK->Fill(D->PtChild(0),WgtD0bar);
139 hptPi->Fill(D->PtChild(0),WgtD0);
140 hptPi->Fill(D->PtChild(1),WgtD0bar);
141 hd0K->Fill(D->Getd0Child(1),WgtD0);
142 hd0K->Fill(D->Getd0Child(0),WgtD0bar);
143 hd0Pi->Fill(D->Getd0Child(0),WgtD0);
144 hd0Pi->Fill(D->Getd0Child(1),WgtD0bar);
145 hMass->Fill(MD0,WgtD0);
146 hMass->Fill(MD0bar,WgtD0bar);
147 hCts->Fill(ctsD0,WgtD0);
148 hCts->Fill(ctsD0bar,WgtD0bar);
149 hCtsVsPtK->Fill(D->PtChild(1),ctsD0,WgtD0);
150 hCtsVsPtK->Fill(D->PtChild(0),ctsD0bar,WgtD0bar);
151 hDCA->Fill(D->GetDCA(),WgtD0+WgtD0bar);
152 hptD0->Fill(D->Pt(),WgtD0+WgtD0bar);
153 hyD0->Fill(D->Rapidity(),WgtD0+WgtD0bar);
154 hd0d0->Fill(D->ProdImpParams(),WgtD0+WgtD0bar);
155 hCPtaD0->Fill(D->CosPointing(),WgtD0+WgtD0bar);
156 hCPtaXY->Fill(D->CosPointingXY(),WgtD0+WgtD0bar);
157 hCPtaVsd0d0->Fill(D->ProdImpParams(),D->CosPointing(),WgtD0+WgtD0bar);
158 hd0d0VSptD0->Fill(D->Pt(),D->ProdImpParams(),WgtD0+WgtD0bar);
159 hCPtaVsd0d0zoom->Fill(D->ProdImpParams(),D->CosPointing(),WgtD0+WgtD0bar);
160 hArm->Fill(D->Alpha(),D->Qt(),WgtD0+WgtD0bar);
161
162
163 } // end loop on D0 candidates
164
165 inFile->Close();
166
167 printf("\n\n --- Total number of candidates passing selection: %d\n\n --- Sum of weights sample A+B+C: %f\n\n",nSel,sampleABC);
168
169 // draw histograms
170 TCanvas *c1 = new TCanvas("c1","pt K & pi",0,0,700,700);
171 c1->SetLogy();
172 hptK->Draw();
173 hptPi->Draw("same");
174
175 TCanvas *c2 = new TCanvas("c2","pt D0",0,0,700,700);
176 c2->SetLogy();
177 hptD0->Draw();
178
179 TCanvas *c3 = new TCanvas("c3","rapidity D0",0,0,700,700);
180 hyD0->Draw();
181
182 TCanvas *c4 = new TCanvas("c4","pointing angle",0,0,700,700);
183 hCPtaD0->Draw();
184
185 TCanvas *c5 = new TCanvas("c5","d0 x d0",0,0,700,700);
186 c5->SetLogy();
187 hd0d0->Draw();
188
189 TCanvas *c6 = new TCanvas("c6","pointing angle VS d0d0",0,0,700,700);
190 c6->SetLogz();
191 hCPtaVsd0d0->Draw("box");
192
193 TCanvas *c7 = new TCanvas("c7","mass",0,0,700,700);
194 hMass->Draw();
195
196 TCanvas *c8 = new TCanvas("c8","armenteros",0,0,700,700);
197 hArm->Draw("box");
198
199 TCanvas *c9 = new TCanvas("c9","decay angle",0,0,700,700);
200 hCts->Draw();
201
202 TCanvas *c10 = new TCanvas("c10","dca",0,0,700,700);
203 c10->SetLogy();
204 hDCA->Draw();
205
206 TCanvas *c11 = new TCanvas("c11","d0 K & pi",0,0,700,700);
207 c11->SetLogy();
208 hd0K->Draw();
209 hd0Pi->Draw("same");
210
211 // write all histograms to file
212 TFile *outFile = new TFile(outName,"recreate");
213 hMass->Write();
214 hDCA->Write();
215 hCts->Write();
216 hCtsVsPtK->Write();
217 hArm->Write();
218 hCPtaVsd0d0->Write();
219 hd0d0VSptD0->Write();
220 hCPtaVsd0d0zoom->Write();
221 hd0d0->Write();
222 hCPtaD0->Write();
223 hCPtaXY->Write();
224 hptK->Write();
225 hptPi->Write();
226 hptD0->Write();
227 hyD0->Write();
228 hd0K->Write();
229 hd0Pi->Write();
230 outFile->Close();
231
232 return;
233}
234
235
236
237