a4c1f5dd |
1 | // $Id$ |
2 | |
3 | //************************************************************************** |
4 | //* This file is property of and copyright by the ALICE HLT Project * |
5 | //* ALICE Experiment at CERN, All rights reserved. * |
6 | //* * |
7 | //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> * |
8 | //* for The ALICE HLT Project. * |
9 | //* * |
10 | //* Permission to use, copy, modify and distribute this software and its * |
11 | //* documentation strictly for non-commercial purposes is hereby granted * |
12 | //* without fee, provided that the above copyright notice appears in all * |
13 | //* copies and that both the copyright notice and this permission notice * |
14 | //* appear in the supporting documentation. The authors make no claims * |
15 | //* about the suitability of this software for any purpose. It is * |
16 | //* provided "as is" without express or implied warranty. * |
17 | //************************************************************************** |
18 | |
19 | /// @file AliHLTGlobalHistoComponent.cxx |
20 | /// @author Matthias Richter |
21 | /// @date 2010-09-16 |
22 | /// @brief A histogramming component for global ESD properties based |
23 | /// on the AliHLTTTreeProcessor |
24 | |
25 | #include "AliHLTGlobalHistoComponent.h" |
26 | #include "AliESDEvent.h" |
7a37a034 |
27 | #include "AliESDv0.h" |
28 | #include "AliKFParticle.h" |
29 | #include "AliKFVertex.h" |
30 | |
a4c1f5dd |
31 | #include "TTree.h" |
32 | #include "TString.h" |
a4c1f5dd |
33 | #include <cassert> |
34 | |
35 | /** ROOT macro for the implementation of ROOT specific class methods */ |
36 | ClassImp(AliHLTGlobalHistoComponent) |
37 | |
38 | AliHLTGlobalHistoComponent::AliHLTGlobalHistoComponent() |
39 | : AliHLTTTreeProcessor() |
40 | , fEvent(0) |
41 | , fNofTracks(0) |
1b89b93a |
42 | , fNofV0s(0) |
1c1af02a |
43 | , fNofUPCpairs(0) |
1b89b93a |
44 | , fNofContributors(0) |
54da0d34 |
45 | , fVertexX(-99) |
46 | , fVertexY(-99) |
47 | , fVertexZ(-99) |
1b89b93a |
48 | , fVertexStatus(kFALSE) |
a4c1f5dd |
49 | , fTrackVariables() |
e7f3baeb |
50 | , fTrackVariablesInt() |
7a37a034 |
51 | , fV0Variables() |
1c1af02a |
52 | , fUPCVariables() |
7a37a034 |
53 | , fNEvents(0) |
54 | , fNGammas(0) |
55 | , fNKShorts(0) |
56 | , fNLambdas(0) |
57 | , fNPi0s(0) |
a4c1f5dd |
58 | { |
59 | // see header file for class documentation |
60 | // or |
61 | // refer to README to build package |
62 | // or |
63 | // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt |
64 | |
65 | } |
66 | |
67 | AliHLTGlobalHistoComponent::~AliHLTGlobalHistoComponent() |
68 | { |
69 | // see header file for class documentation |
70 | fTrackVariables.Reset(); |
e7f3baeb |
71 | fTrackVariablesInt.Reset(); |
7a37a034 |
72 | fV0Variables.Reset(); |
1c1af02a |
73 | fUPCVariables.Reset(); |
a4c1f5dd |
74 | } |
75 | |
7a37a034 |
76 | void AliHLTGlobalHistoComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list){ |
a4c1f5dd |
77 | // see header file for class documentation |
78 | list.push_back(kAliHLTAllDataTypes); |
79 | } |
80 | |
7a37a034 |
81 | AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOutputDataType(){ |
af747e4b |
82 | // see header file for class documentation |
83 | return kAliHLTDataTypeHistogram|kAliHLTDataOriginOut; |
84 | } |
85 | |
7a37a034 |
86 | TTree* AliHLTGlobalHistoComponent::CreateTree(int /*argc*/, const char** /*argv*/){ |
a4c1f5dd |
87 | // create the tree and branches |
88 | int iResult=0; |
89 | TTree* pTree=new TTree("ESDproperties", "HLT ESD properties"); |
90 | if (!pTree) return NULL; |
91 | |
92 | const char* trackVariableNames = { |
e7f3baeb |
93 | // Note the black at the end of each name! |
a4c1f5dd |
94 | "Track_pt " |
95 | "Track_phi " |
96 | "Track_eta " |
54da0d34 |
97 | "Track_p " |
98 | "Track_theta " |
99 | "Track_Nclusters " |
97491309 |
100 | "Track_status " |
54da0d34 |
101 | "Track_charge " |
97491309 |
102 | "Track_DCAr " |
103 | "Track_DCAz " |
104 | "Track_dEdx " |
a4c1f5dd |
105 | }; |
54da0d34 |
106 | |
e7f3baeb |
107 | const char* trackIntVariableNames = { |
108 | // Note the black at the end of each name! |
109 | "Track_status " |
110 | }; |
111 | |
7a37a034 |
112 | const char* V0VariableNames = { |
113 | // Note the black at the end of each name! |
114 | "V0_AP " |
115 | "V0_pt " |
116 | "clust1 " |
117 | "clust2 " |
118 | "dev1 " |
119 | "dev2 " |
120 | "devPrim " |
121 | "length " |
122 | "sigmaLength " |
123 | "r " |
1c1af02a |
124 | }; |
125 | |
126 | const char* UPCVariableNames = { |
127 | // Note the black at the end of each name! |
128 | "nClusters_1 " |
129 | "nClusters_2 " |
130 | "polarity_1 " |
131 | "polarity_2 " |
49da4a19 |
132 | "px_1 " |
133 | "py_1 " |
134 | "px_2 " |
135 | "py_2 " |
1c1af02a |
136 | }; |
7a37a034 |
137 | |
af747e4b |
138 | int maxTrackCount = 20000; // FIXME: make configurable |
7a37a034 |
139 | int maxV0Count = 100000; |
1c1af02a |
140 | int maxUPCCount = 1; |
af747e4b |
141 | |
a4c1f5dd |
142 | if ((iResult=fTrackVariables.Init(maxTrackCount, trackVariableNames))<0) { |
e7f3baeb |
143 | HLTError("failed to initialize internal structure for track properties (float)"); |
144 | } |
145 | if ((iResult=fTrackVariablesInt.Init(maxTrackCount, trackIntVariableNames))<0) { |
146 | HLTError("failed to initialize internal structure for track properties (int)"); |
a4c1f5dd |
147 | } |
7a37a034 |
148 | if ((iResult=fV0Variables.Init(maxV0Count, V0VariableNames))<0) { |
149 | HLTError("failed to initialize internal structure for V0 properties (float)"); |
150 | } |
1c1af02a |
151 | if ((iResult=fUPCVariables.Init(maxUPCCount, UPCVariableNames))<0) { |
152 | HLTError("failed to initialize internal structure for UPC properties (float)"); |
153 | } |
a4c1f5dd |
154 | |
155 | if (iResult>=0) { |
1b89b93a |
156 | pTree->Branch("event", &fEvent, "event/I"); |
157 | pTree->Branch("trackcount", &fNofTracks, "trackcount/I"); |
158 | pTree->Branch("vertexX", &fVertexX, "vertexX/F"); |
159 | pTree->Branch("vertexY", &fVertexY, "vertexY/F"); |
160 | pTree->Branch("vertexZ", &fVertexZ, "vertexZ/F"); |
7a37a034 |
161 | pTree->Branch("V0", &fNofV0s, "V0/I"); |
1c1af02a |
162 | pTree->Branch("UPC", &fNofUPCpairs, "UPC/I"); |
1b89b93a |
163 | pTree->Branch("nContributors",&fNofContributors, "nContributors/I"); |
164 | pTree->Branch("vertexStatus", &fVertexStatus, "vertexStatus/I"); |
e7f3baeb |
165 | |
166 | int i=0; |
167 | // FIXME: this is a bit ugly since type 'f' and 'i' are specified |
168 | // explicitely. Would be better to use a function like |
169 | // AliHLTGlobalHistoVariables::GetType but could not get this working |
170 | for (i=0; i<fTrackVariables.Variables(); i++) { |
a4c1f5dd |
171 | TString specifier=fTrackVariables.GetKey(i); |
172 | float* pArray=fTrackVariables.GetArray(specifier); |
173 | specifier+="[trackcount]/f"; |
174 | pTree->Branch(fTrackVariables.GetKey(i), pArray, specifier.Data()); |
175 | } |
e7f3baeb |
176 | for (i=0; i<fTrackVariablesInt.Variables(); i++) { |
177 | TString specifier=fTrackVariablesInt.GetKey(i); |
178 | int* pArray=fTrackVariablesInt.GetArray(specifier); |
179 | specifier+="[trackcount]/i"; |
180 | pTree->Branch(fTrackVariablesInt.GetKey(i), pArray, specifier.Data()); |
7a37a034 |
181 | } |
182 | for (i=0; i<fV0Variables.Variables(); i++) { |
183 | TString specifier=fV0Variables.GetKey(i); |
184 | float* pArray=fV0Variables.GetArray(specifier); |
185 | specifier+="[V0]/f"; |
186 | pTree->Branch(fV0Variables.GetKey(i), pArray, specifier.Data()); |
e7f3baeb |
187 | } |
1c1af02a |
188 | for (i=0; i<fUPCVariables.Variables(); i++) { |
189 | TString specifier=fUPCVariables.GetKey(i); |
190 | float* pArray=fUPCVariables.GetArray(specifier); |
191 | specifier+="[UPC]/f"; |
192 | pTree->Branch(fUPCVariables.GetKey(i), pArray, specifier.Data()); |
193 | } |
a4c1f5dd |
194 | } else { |
195 | delete pTree; |
196 | pTree=NULL; |
197 | } |
7a37a034 |
198 | |
a4c1f5dd |
199 | return pTree; |
200 | } |
201 | |
7a37a034 |
202 | void AliHLTGlobalHistoComponent::FillHistogramDefinitions(){ |
a4c1f5dd |
203 | /// default histogram definitions |
204 | } |
205 | |
7a37a034 |
206 | int AliHLTGlobalHistoComponent::FillTree(TTree* pTree, const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/ ){ |
207 | |
a4c1f5dd |
208 | /// fill the tree from the ESD |
209 | int iResult=0; |
210 | if (!IsDataEvent()) return 0; |
211 | |
212 | ResetVariables(); |
213 | |
214 | // fetch ESD from input stream |
215 | const TObject* obj = GetFirstInputObject(kAliHLTAllDataTypes, "AliESDEvent"); |
216 | AliESDEvent* esd = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj)); |
217 | esd->GetStdContent(); |
218 | |
219 | // fill track variables |
1b89b93a |
220 | fNofTracks = esd->GetNumberOfTracks(); |
221 | fVertexX = esd->GetPrimaryVertexTracks()->GetX(); |
222 | fVertexY = esd->GetPrimaryVertexTracks()->GetY(); |
223 | fVertexZ = esd->GetPrimaryVertexTracks()->GetZ(); |
224 | fNofV0s = esd->GetNumberOfV0s(); |
225 | fNofContributors = esd->GetPrimaryVertexTracks()->GetNContributors(); |
226 | fVertexStatus = esd->GetPrimaryVertexTracks()->GetStatus(); |
1c1af02a |
227 | fNofUPCpairs = 1; |
54da0d34 |
228 | |
a4c1f5dd |
229 | for (int i=0; i<fNofTracks; i++) { |
230 | AliESDtrack *esdTrack = esd->GetTrack(i); |
231 | if (!esdTrack) continue; |
232 | |
97491309 |
233 | Float_t DCAr, DCAz = -99; |
234 | esdTrack->GetImpactParametersTPC(DCAr, DCAz); |
235 | |
54da0d34 |
236 | fTrackVariables.Fill("Track_pt" , esdTrack->Pt() ); |
237 | fTrackVariables.Fill("Track_phi" , esdTrack->Phi()*TMath::RadToDeg() ); |
238 | fTrackVariables.Fill("Track_eta" , esdTrack->Theta() ); |
239 | fTrackVariables.Fill("Track_p" , esdTrack->P() ); |
240 | fTrackVariables.Fill("Track_theta" , esdTrack->Theta()*TMath::RadToDeg() ); |
241 | fTrackVariables.Fill("Track_Nclusters" , esdTrack->GetTPCNcls() ); |
97491309 |
242 | fTrackVariables.Fill("Track_status" , esdTrack->GetStatus() ); |
54da0d34 |
243 | fTrackVariables.Fill("Track_charge" , esdTrack->Charge() ); |
97491309 |
244 | fTrackVariables.Fill("Track_DCAr" , DCAr ); |
245 | fTrackVariables.Fill("Track_DCAz" , DCAz ); |
246 | fTrackVariables.Fill("Track_dEdx" , esdTrack->GetTPCsignal() ); |
1c1af02a |
247 | fTrackVariablesInt.Fill("Track_status" , esdTrack->GetStatus() ); |
7a37a034 |
248 | } |
1b89b93a |
249 | //HLTInfo("added parameters for %d tracks", fNofTracks); |
1c1af02a |
250 | |
251 | if(fNofTracks==2){ |
252 | AliESDtrack *esdTrack1 = esd->GetTrack(0); |
253 | if(!esdTrack1) return 0; |
254 | AliESDtrack *esdTrack2 = esd->GetTrack(1); |
255 | if(!esdTrack2) return 0; |
256 | |
257 | if(esdTrack1->Charge()*esdTrack2->Charge()<0){ |
258 | |
259 | fUPCVariables.Fill("nClusters_1", esdTrack1->GetTPCNcls() ); |
260 | fUPCVariables.Fill("nClusters_2", esdTrack2->GetTPCNcls() ); |
261 | fUPCVariables.Fill("polarity_1", esdTrack1->Charge() ); |
262 | fUPCVariables.Fill("polarity_2", esdTrack2->Charge() ); |
49da4a19 |
263 | fUPCVariables.Fill("px_1", esdTrack1->Px() ); |
264 | fUPCVariables.Fill("py_1", esdTrack1->Py() ); |
265 | fUPCVariables.Fill("px_2", esdTrack2->Px() ); |
266 | fUPCVariables.Fill("py_2", esdTrack2->Py() ); |
1c1af02a |
267 | } |
268 | } |
269 | |
7a37a034 |
270 | AliKFParticle::SetField( esd->GetMagneticField() ); |
271 | |
8dabf5a9 |
272 | //const double kKsMass = 0.49767; |
273 | //const double kLambdaMass = 1.11568; |
7a37a034 |
274 | //const double kPi0Mass = 0.13498; |
275 | |
1c1af02a |
276 | //std::vector<AliKFParticle> vGammas; |
7a37a034 |
277 | |
278 | |
279 | for (int i=0; i<fNofV0s; i++) { |
280 | AliESDv0 *esdV0 = esd->GetV0(i); |
281 | if (!esdV0) continue; |
282 | |
283 | AliESDtrack *t1 = esd->GetTrack( esd->GetV0(i)->GetNindex()); |
284 | AliESDtrack *t2 = esd->GetTrack( esd->GetV0(i)->GetPindex()); |
285 | |
286 | AliKFParticle kf1( *t1, 11 ); |
287 | AliKFParticle kf2( *t2, 11 ); |
288 | |
289 | AliKFVertex primVtx( *esd->GetPrimaryVertexTracks() ); |
290 | double dev1 = kf1.GetDeviationFromVertex( primVtx ); |
291 | double dev2 = kf2.GetDeviationFromVertex( primVtx ); |
292 | |
293 | AliKFParticle v0( kf1, kf2 ); |
294 | double devPrim = v0.GetDeviationFromVertex( primVtx ); |
295 | primVtx+=v0; |
296 | v0.SetProductionVertex( primVtx ); |
297 | |
298 | Double_t length, sigmaLength; |
299 | if( v0.GetDecayLength( length, sigmaLength ) ) continue; |
300 | |
301 | double dx = v0.GetX()-primVtx.GetX(); |
302 | double dy = v0.GetY()-primVtx.GetY(); |
303 | double r = sqrt(dx*dx + dy*dy); |
304 | |
305 | |
306 | // Armenteros-Podolanski plot |
307 | |
308 | double pt=0, ap=0; |
309 | //{ |
310 | AliKFParticle kf01 = kf1, kf02 = kf2; |
311 | kf01.SetProductionVertex(v0); |
312 | kf02.SetProductionVertex(v0); |
313 | kf01.TransportToProductionVertex(); |
314 | kf02.TransportToProductionVertex(); |
315 | double px1=kf01.GetPx(), py1=kf01.GetPy(), pz1=kf01.GetPz(); |
316 | double px2=kf02.GetPx(), py2=kf02.GetPy(), pz2=kf02.GetPz(); |
317 | double px = px1+px2, py = py1+py2, pz = pz1+pz2; |
318 | double p = sqrt(px*px+py*py+pz*pz); |
319 | double l1 = (px*px1 + py*py1 + pz*pz1)/p; |
320 | double l2 = (px*px2 + py*py2 + pz*pz2)/p; |
321 | pt = sqrt(px1*px1+py1*py1+pz1*pz1 - l1*l1); |
322 | ap = (l2-l1)/(l1+l2); |
323 | //} |
324 | |
325 | // if( |
326 | // t1->GetTPCNcls()>=fAPCuts[0] |
327 | // && t2->GetTPCNcls()>=fAPCuts[0] |
328 | // && dev1>=fAPCuts[1] |
329 | // && dev2>=fAPCuts[1] |
330 | // && devPrim <= fAPCuts[2] |
331 | // && length >= fAPCuts[3]*sigmaLength |
332 | // && length >= fAPCuts[4] |
333 | // && r <= fAPCuts[5] |
334 | // ) |
335 | // //{ |
336 | // //if( fAP ) fAP->Fill( ap, pt ); |
337 | // //} |
338 | // |
339 | // // Gamma finder |
340 | // |
341 | // bool isGamma = 0; |
342 | // |
343 | // if( |
344 | // t1->GetTPCNcls()>=fGammaCuts[0] |
345 | // && t2->GetTPCNcls()>=fGammaCuts[0] |
346 | // && dev1>=fGammaCuts[1] |
347 | // && dev2>=fGammaCuts[1] |
348 | // && devPrim <= fGammaCuts[2] |
349 | // && length >= fGammaCuts[3]*sigmaLength |
350 | // && length >= fGammaCuts[4] |
351 | // && r <= fGammaCuts[5] |
352 | // ){ |
353 | // double mass, error; |
354 | // v0.GetMass(mass,error); |
355 | // //if( fGamma ) fGamma->Fill( mass ); |
356 | // |
357 | // if( TMath::Abs(mass)<=fGammaCuts[6]*error || TMath::Abs(mass)<=fGammaCuts[7] ){ |
358 | // AliKFParticle gamma = v0; |
359 | // gamma.SetMassConstraint(0); |
360 | // // if( fGammaXY |
361 | // // && t1->GetTPCNcls()>=60 |
362 | // // && t2->GetTPCNcls()>=60 |
363 | // // ) fGammaXY->Fill(gamma.GetX(), gamma.GetY()); |
364 | // isGamma = 1; |
365 | // fNGammas++; |
366 | // vGammas.push_back( gamma ); |
367 | // } |
368 | // } |
369 | // |
370 | // if( isGamma ) continue; |
371 | // |
372 | // |
373 | // // KShort finder |
374 | // |
375 | // bool isKs = 0; |
376 | // |
377 | // if( |
378 | // t1->GetTPCNcls()>=fKsCuts[0] |
379 | // && t2->GetTPCNcls()>=fKsCuts[0] |
380 | // && dev1>=fKsCuts[1] |
381 | // && dev2>=fKsCuts[1] |
382 | // && devPrim <= fKsCuts[2] |
383 | // && length >= fKsCuts[3]*sigmaLength |
384 | // && length >= fKsCuts[4] |
385 | // && r <= fKsCuts[5] |
386 | // ){ |
387 | // |
388 | // AliKFParticle piN( *t1, 211 ); |
389 | // AliKFParticle piP( *t2, 211 ); |
390 | // |
391 | // AliKFParticle Ks( piN, piP ); |
392 | // Ks.SetProductionVertex( primVtx ); |
393 | // |
394 | // double mass, error; |
395 | // Ks.GetMass( mass, error); |
396 | // //if( fKShort ) fKShort->Fill( mass ); |
397 | // if( TMath::Abs( mass - kKsMass )<=fKsCuts[6]*error || TMath::Abs( mass - kKsMass )<=fKsCuts[7] ){ |
398 | // isKs = 1; |
399 | // fNKShorts++; |
400 | // } |
401 | // } |
402 | // |
403 | // if( isKs ) continue; |
404 | // |
405 | // // Lambda finder |
406 | // //printf("QQQQQQQQQQQQQQQQq :%f\n",fLambdaCuts[0]); |
407 | // if( |
408 | // t1->GetTPCNcls()>=fLambdaCuts[0] |
409 | // && t2->GetTPCNcls()>=fLambdaCuts[0] |
410 | // && dev1>=fLambdaCuts[1] |
411 | // && dev2>=fLambdaCuts[1] |
412 | // && devPrim <= fLambdaCuts[2] |
413 | // && length >= fLambdaCuts[3]*sigmaLength |
414 | // && length >= fLambdaCuts[4] |
415 | // && r <= fLambdaCuts[5] |
416 | // && TMath::Abs( ap )>.4 |
417 | // ){ |
418 | // |
419 | // AliKFParticle kP, kpi; |
420 | // if( ap<0 ){ |
421 | // kP = AliKFParticle( *t2, 2212 ); |
422 | // kpi = AliKFParticle( *t1, 211 ); |
423 | // } else { |
424 | // kP = AliKFParticle( *t1, 2212 ); |
425 | // kpi = AliKFParticle( *t2, 211 ); |
426 | // } |
427 | // |
428 | // AliKFParticle lambda = AliKFParticle( kP, kpi ); |
429 | // lambda.SetProductionVertex( primVtx ); |
430 | // //double mass, error; |
431 | // lambda.GetMass( Lmass, Lerror); |
432 | // //if( fLambda ) fLambda->Fill( mass ); |
433 | // if( TMath::Abs( Lmass - kLambdaMass )<=fLambdaCuts[6]*Lerror || TMath::Abs( Lmass - kLambdaMass )<=fLambdaCuts[7] ){ |
434 | // fNLambdas++; |
435 | // } |
436 | // } |
437 | |
438 | |
439 | fV0Variables.Fill("V0_AP", ap); |
440 | fV0Variables.Fill("V0_pt", pt); |
441 | fV0Variables.Fill("clust1", t1->GetTPCNcls()); |
442 | fV0Variables.Fill("clust2", t2->GetTPCNcls()); |
443 | fV0Variables.Fill("dev1", dev1); |
444 | fV0Variables.Fill("dev2", dev2); |
445 | fV0Variables.Fill("devPrim", devPrim); |
446 | fV0Variables.Fill("length", length); |
447 | fV0Variables.Fill("sigmaLength", sigmaLength); |
448 | fV0Variables.Fill("r", r); |
449 | |
450 | } // end of loop over V0s |
451 | |
a4c1f5dd |
452 | if (iResult<0) { |
453 | // fill an empty event |
454 | ResetVariables(); |
455 | } |
456 | fEvent++; |
457 | |
458 | pTree->Fill(); |
459 | return iResult; |
460 | } |
461 | |
462 | int AliHLTGlobalHistoComponent::ResetVariables() |
463 | { |
464 | /// reset all filling variables |
465 | fNofTracks=0; |
4868ab3f |
466 | fNofV0s=0; |
a4c1f5dd |
467 | fTrackVariables.ResetCount(); |
e7f3baeb |
468 | fTrackVariablesInt.ResetCount(); |
c3176951 |
469 | fV0Variables.ResetCount(); |
1c1af02a |
470 | fUPCVariables.ResetCount(); |
a4c1f5dd |
471 | return 0; |
472 | } |
473 | |
474 | AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOriginDataType() const |
475 | { |
476 | // get the origin of the output data |
4bb46850 |
477 | return kAliHLTDataTypeHistogram|kAliHLTDataOriginHLT; |
a4c1f5dd |
478 | } |