using now the specification provided by the base class
[u/mrichter/AliRoot.git] / HLT / global / AliHLTGlobalHistoComponent.cxx
CommitLineData
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 */
36ClassImp(AliHLTGlobalHistoComponent)
37
38AliHLTGlobalHistoComponent::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
67AliHLTGlobalHistoComponent::~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 76void AliHLTGlobalHistoComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list){
a4c1f5dd 77 // see header file for class documentation
78 list.push_back(kAliHLTAllDataTypes);
79}
80
7a37a034 81AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOutputDataType(){
af747e4b 82 // see header file for class documentation
83 return kAliHLTDataTypeHistogram|kAliHLTDataOriginOut;
84}
85
7a37a034 86TTree* 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 202void AliHLTGlobalHistoComponent::FillHistogramDefinitions(){
a4c1f5dd 203 /// default histogram definitions
204}
205
7a37a034 206int 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
462int 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
474AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOriginDataType() const
475{
476 // get the origin of the output data
4bb46850 477 return kAliHLTDataTypeHistogram|kAliHLTDataOriginHLT;
a4c1f5dd 478}