]>
Commit | Line | Data |
---|---|---|
1 | #include "AlidNdEtaCorrectionSelector.h" | |
2 | ||
3 | #include <TStyle.h> | |
4 | #include <TSystem.h> | |
5 | #include <TCanvas.h> | |
6 | #include <TParticle.h> | |
7 | #include <TParticlePDG.h> | |
8 | ||
9 | #include <AliLog.h> | |
10 | #include <AliGenEventHeader.h> | |
11 | #include <AliTracker.h> | |
12 | ||
13 | #include "../esdTrackCuts/AliESDtrackCuts.h" | |
14 | #include "dNdEtaCorrection.h" | |
15 | ||
16 | ClassImp(AlidNdEtaCorrectionSelector) | |
17 | ||
18 | AlidNdEtaCorrectionSelector::AlidNdEtaCorrectionSelector(TTree *) : | |
19 | AliSelector(), | |
20 | fEsdTrackCuts(0), | |
21 | fdNdEtaCorrection(0), | |
22 | fdNdEtaCorrectionFinal(0) | |
23 | { | |
24 | // | |
25 | // Constructor. Initialization of pointers | |
26 | // | |
27 | } | |
28 | ||
29 | AlidNdEtaCorrectionSelector::~AlidNdEtaCorrectionSelector() | |
30 | { | |
31 | // | |
32 | // Destructor | |
33 | // | |
34 | ||
35 | // histograms are in the output list and deleted when the output | |
36 | // list is deleted by the TSelector dtor | |
37 | } | |
38 | ||
39 | void AlidNdEtaCorrectionSelector::Begin(TTree * tree) | |
40 | { | |
41 | // The Begin() function is called at the start of the query. | |
42 | // When running with PROOF Begin() is only called on the client. | |
43 | // The tree argument is deprecated (on PROOF 0 is passed). | |
44 | ||
45 | AliSelector::Begin(tree); | |
46 | } | |
47 | ||
48 | void AlidNdEtaCorrectionSelector::SlaveBegin(TTree * tree) | |
49 | { | |
50 | // The SlaveBegin() function is called after the Begin() function. | |
51 | // When running with PROOF SlaveBegin() is called on each slave server. | |
52 | // The tree argument is deprecated (on PROOF 0 is passed). | |
53 | ||
54 | AliSelector::SlaveBegin(tree); | |
55 | ||
56 | fdNdEtaCorrection = new dNdEtaCorrection(); | |
57 | ||
58 | if (fChain) | |
59 | fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fChain->GetUserInfo()->FindObject("AliESDtrackCuts")); | |
60 | ||
61 | if (!fEsdTrackCuts) | |
62 | AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from user info"); | |
63 | } | |
64 | ||
65 | Bool_t AlidNdEtaCorrectionSelector::Notify() | |
66 | { | |
67 | // The Notify() function is called when a new file is opened. This | |
68 | // can be either for a new TTree in a TChain or when when a new TTree | |
69 | // is started when using PROOF. Typically here the branch pointers | |
70 | // will be retrieved. It is normaly not necessary to make changes | |
71 | // to the generated code, but the routine can be extended by the | |
72 | // user if needed. | |
73 | ||
74 | if (AliSelector::Notify() == kFALSE) | |
75 | return kFALSE; | |
76 | ||
77 | return kTRUE; | |
78 | } | |
79 | ||
80 | Bool_t AlidNdEtaCorrectionSelector::IsPrimary(const TParticle* aParticle, Int_t aTotalPrimaries) | |
81 | { | |
82 | // | |
83 | // Returns if the given particle is a primary particle | |
84 | // This function or a equivalent should be available in some common place of AliRoot | |
85 | // | |
86 | ||
87 | // if the particle has a daughter primary, we do not want to count it | |
88 | if (aParticle->GetFirstDaughter() != -1 && aParticle->GetFirstDaughter() < aTotalPrimaries) | |
89 | return kFALSE; | |
90 | ||
91 | Int_t pdgCode = TMath::Abs(aParticle->GetPdgCode()); | |
92 | ||
93 | // skip quarks and gluon | |
94 | if (pdgCode > 10 && pdgCode != 21) | |
95 | return kTRUE; | |
96 | ||
97 | return kFALSE; | |
98 | } | |
99 | ||
100 | Bool_t AlidNdEtaCorrectionSelector::Process(Long64_t entry) | |
101 | { | |
102 | // The Process() function is called for each entry in the tree (or possibly | |
103 | // keyed object in the case of PROOF) to be processed. The entry argument | |
104 | // specifies which entry in the currently loaded tree is to be processed. | |
105 | // It can be passed to either TTree::GetEntry() or TBranch::GetEntry() | |
106 | // to read either all or the required parts of the data. When processing | |
107 | // keyed objects with PROOF, the object is already loaded and is available | |
108 | // via the fObject pointer. | |
109 | // | |
110 | // This function should contain the "body" of the analysis. It can contain | |
111 | // simple or elaborate selection criteria, run algorithms on the data | |
112 | // of the event and typically fill histograms. | |
113 | ||
114 | // WARNING when a selector is used with a TChain, you must use | |
115 | // the pointer to the current TTree to call GetEntry(entry). | |
116 | // The entry is always the local entry number in the current tree. | |
117 | // Assuming that fChain is the pointer to the TChain being processed, | |
118 | // use fChain->GetTree()->GetEntry(entry). | |
119 | ||
120 | if (AliSelector::Process(entry) == kFALSE) | |
121 | return kFALSE; | |
122 | ||
123 | // check prerequesites | |
124 | if (!fESD) | |
125 | { | |
126 | AliDebug(AliLog::kError, "ESD branch not available"); | |
127 | return kFALSE; | |
128 | } | |
129 | ||
130 | if (!fHeader) | |
131 | { | |
132 | AliDebug(AliLog::kError, "Header branch not available"); | |
133 | return kFALSE; | |
134 | } | |
135 | ||
136 | if (!fEsdTrackCuts) | |
137 | { | |
138 | AliDebug(AliLog::kError, "fESDTrackCuts not available"); | |
139 | return kFALSE; | |
140 | } | |
141 | ||
142 | // ######################################################## | |
143 | // get the EDS vertex | |
144 | const AliESDVertex* vtxESD = fESD->GetVertex(); | |
145 | ||
146 | // the vertex should be reconstructed | |
147 | if (strcmp(vtxESD->GetName(),"default")==0) | |
148 | return kTRUE; | |
149 | ||
150 | Double_t vtx_res[3]; | |
151 | vtx_res[0] = vtxESD->GetXRes(); | |
152 | vtx_res[1] = vtxESD->GetYRes(); | |
153 | vtx_res[2] = vtxESD->GetZRes(); | |
154 | ||
155 | // the resolution should be reasonable??? | |
156 | if (vtx_res[2]==0 || vtx_res[2]>0.1) | |
157 | return kTRUE; | |
158 | ||
159 | // ######################################################## | |
160 | // get the MC vertex | |
161 | AliGenEventHeader* genHeader = fHeader->GenEventHeader(); | |
162 | ||
163 | TArrayF vtxMC(3); | |
164 | genHeader->PrimaryVertex(vtxMC); | |
165 | ||
166 | // ######################################################## | |
167 | // loop over mc particles | |
168 | TTree* particleTree = GetKinematics(); | |
169 | TParticle* particle = 0; | |
170 | particleTree->SetBranchAddress("Particles", &particle); | |
171 | ||
172 | Int_t nPrim = fHeader->GetNprimary(); | |
173 | Int_t nTotal = fHeader->GetNtrack(); | |
174 | ||
175 | for (Int_t i_mc = nTotal - nPrim; i_mc < nTotal; ++i_mc) | |
176 | { | |
177 | particleTree->GetEntry(i_mc); | |
178 | ||
179 | if (!particle) | |
180 | continue; | |
181 | ||
182 | if (strcmp(particle->GetName(),"XXX") == 0) | |
183 | { | |
184 | AliDebug(AliLog::kWarning, Form("WARNING: There is a particle named XXX (%d).", i_mc)); | |
185 | continue; | |
186 | } | |
187 | ||
188 | TParticlePDG* pdgPart = particle->GetPDG(); | |
189 | ||
190 | if (strcmp(pdgPart->ParticleClass(),"Unknown") == 0) | |
191 | { | |
192 | AliDebug(AliLog::kError, Form("WARNING: There is a particle with an unknown particle class (%d pdg code %d).", i_mc, particle->GetPdgCode())); | |
193 | continue; | |
194 | } | |
195 | ||
196 | if (IsPrimary(particle, nPrim) == kFALSE) | |
197 | continue; | |
198 | ||
199 | if (pdgPart->Charge() == 0) | |
200 | continue; | |
201 | ||
202 | fdNdEtaCorrection->FillGene(vtxMC[2], particle->Eta()); | |
203 | ||
204 | }// end of mc particle | |
205 | ||
206 | // ######################################################## | |
207 | // loop over esd tracks | |
208 | Int_t nTracks = fESD->GetNumberOfTracks(); | |
209 | for (Int_t t=0; t<nTracks; t++) | |
210 | { | |
211 | AliESDtrack* esdTrack = fESD->GetTrack(t); | |
212 | ||
213 | // cut the esd track? | |
214 | if (!fEsdTrackCuts->AcceptTrack(esdTrack)) | |
215 | continue; | |
216 | ||
217 | // using the eta of the mc particle | |
218 | Int_t label = TMath::Abs(esdTrack->GetLabel()); | |
219 | if (label == 0) | |
220 | { | |
221 | AliDebug(AliLog::kWarning, Form("WARNING: cannot find corresponding mc part for track %d.", t)); | |
222 | continue; | |
223 | } | |
224 | particleTree->GetEntry(nTotal - nPrim + label); | |
225 | ||
226 | fdNdEtaCorrection->FillMeas(vtxMC[2], particle->Eta()); | |
227 | ||
228 | } // end of track loop | |
229 | ||
230 | return kTRUE; | |
231 | } | |
232 | ||
233 | void AlidNdEtaCorrectionSelector::SlaveTerminate() | |
234 | { | |
235 | // The SlaveTerminate() function is called after all entries or objects | |
236 | // have been processed. When running with PROOF SlaveTerminate() is called | |
237 | // on each slave server. | |
238 | ||
239 | AliSelector::SlaveTerminate(); | |
240 | ||
241 | // Add the histograms to the output on each slave server | |
242 | if (!fOutput) | |
243 | { | |
244 | AliDebug(AliLog::kError, "ERROR: Output list not initialized"); | |
245 | return; | |
246 | } | |
247 | ||
248 | fOutput->Add(fdNdEtaCorrection->GetGeneratedHistogram()); | |
249 | fOutput->Add(fdNdEtaCorrection->GetMeasuredHistogram()); | |
250 | } | |
251 | ||
252 | void AlidNdEtaCorrectionSelector::Terminate() | |
253 | { | |
254 | // The Terminate() function is the last function to be called during | |
255 | // a query. It always runs on the client, it can be used to present | |
256 | // the results graphically or save the results to file. | |
257 | ||
258 | AliSelector::Terminate(); | |
259 | ||
260 | fdNdEtaCorrectionFinal = new dNdEtaCorrection(); | |
261 | TH2F* measuredHistogram = dynamic_cast<TH2F*> (fOutput->FindObject("etaVsVtx_meas")); | |
262 | TH2F* generatedHistogram = dynamic_cast<TH2F*> (fOutput->FindObject("etaVsVtx_gene")); | |
263 | if (!measuredHistogram || !generatedHistogram) | |
264 | { | |
265 | AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p", (void*) generatedHistogram, (void*) measuredHistogram)); | |
266 | return; | |
267 | } | |
268 | fdNdEtaCorrectionFinal->SetGeneratedHistogram(generatedHistogram); | |
269 | fdNdEtaCorrectionFinal->SetMeasuredHistogram(measuredHistogram); | |
270 | ||
271 | fdNdEtaCorrectionFinal->Finish(); | |
272 | ||
273 | TFile* fout = new TFile("correction_map.root","RECREATE"); | |
274 | ||
275 | fEsdTrackCuts->SaveHistograms("esd_track_cuts"); | |
276 | fdNdEtaCorrectionFinal->SaveHistograms(); | |
277 | ||
278 | fout->Write(); | |
279 | fout->Close(); | |
280 | ||
281 | fdNdEtaCorrectionFinal->DrawHistograms(); | |
282 | } |