d6bee46ef7dc7589634a23eff6d62242f1c69b03
[u/mrichter/AliRoot.git] / HLT / global / AliHLTGlobalHistoComponent.cxx
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"
27 #include "AliESDv0.h"
28 #include "AliKFParticle.h"
29 #include "AliKFVertex.h"
30
31 #include "TTree.h"
32 #include "TString.h"
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)
42   , fNofV0s(0)
43   , fNofContributors(0)
44   , fVertexX(-99)
45   , fVertexY(-99)
46   , fVertexZ(-99)
47   , fVertexStatus(kFALSE)
48   , fMaxTrackCount(20000)
49   , fMaxV0Count(1000)
50   , fFillV0(kFALSE)
51   , fTrackVariables()
52   , fTrackVariablesInt()
53   , fV0Variables()
54 {
55   // see header file for class documentation
56   // or
57   // refer to README to build package
58   // or
59   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
60 }
61
62 AliHLTGlobalHistoComponent::~AliHLTGlobalHistoComponent(){
63   // see header file for class documentation
64   fTrackVariables.Reset();
65   fTrackVariablesInt.Reset();
66   fV0Variables.Reset();
67 }
68
69 void AliHLTGlobalHistoComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list){
70   // see header file for class documentation
71   list.push_back(kAliHLTAllDataTypes);
72 }
73
74 AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOutputDataType(){
75   // see header file for class documentation
76   return kAliHLTDataTypeHistogram|kAliHLTDataOriginOut;
77 }
78
79 TTree* AliHLTGlobalHistoComponent::CreateTree(int /*argc*/, const char** /*argv*/){
80 // create the tree and branches
81  
82   int iResult=0;
83   TTree* pTree = new TTree("ESDproperties", "HLT ESD properties");
84   if (!pTree) return NULL;
85
86   const char* trackVariableNames = {
87     // Note the black at the end of each name!
88     "Track_pt "
89     "Track_phi "
90     "Track_eta "
91     "Track_p "
92     "Track_theta "
93     "Track_Nclusters "
94     "Track_status "
95     "Track_charge "
96     "Track_DCAr "
97     "Track_DCAz "
98     "Track_dEdx "
99   };
100   
101   const char* trackIntVariableNames = {
102     // Note the black at the end of each name!
103     "Track_status "
104   };
105
106   const char* V0VariableNames = {
107     // Note the black at the end of each name!
108     "V0_AP "
109     "V0_pt "  
110     "clust1 "
111     "clust2 "
112     "dev1 "
113     "dev2 "
114     "devPrim "
115     "length "
116     "sigmaLength "
117     "r "
118   };
119      
120   if((iResult=fTrackVariables.Init(fMaxTrackCount, trackVariableNames))<0){
121     HLTError("failed to initialize internal structure for track properties (float)");
122   }
123   if((iResult=fTrackVariablesInt.Init(fMaxTrackCount, trackIntVariableNames))<0){
124     HLTError("failed to initialize internal structure for track properties (int)");
125   }
126   if((iResult=fV0Variables.Init(fMaxV0Count, V0VariableNames))<0){
127       HLTError("failed to initialize internal structure for V0 properties (float)");
128   }
129   
130   if(iResult>=0){     
131      pTree->Branch("event",        &fEvent,           "event/I");
132      pTree->Branch("trackcount",   &fNofTracks,       "trackcount/I");
133      pTree->Branch("vertexX",      &fVertexX,         "vertexX/F");
134      pTree->Branch("vertexY",      &fVertexY,         "vertexY/F");
135      pTree->Branch("vertexZ",      &fVertexZ,         "vertexZ/F");
136      pTree->Branch("nContributors",&fNofContributors, "nContributors/I");
137      pTree->Branch("vertexStatus", &fVertexStatus,    "vertexStatus/I");
138      if(fFillV0==kTRUE) pTree->Branch("V0", &fNofV0s, "V0/I");
139
140      int i=0;
141      // FIXME: this is a bit ugly since type 'f' and 'i' are specified
142      // explicitely. Would be better to use a function like
143      // AliHLTGlobalHistoVariables::GetType but could not get this working
144     
145      for(i=0; i<fTrackVariables.Variables(); i++){
146          TString specifier=fTrackVariables.GetKey(i);
147          float* pArray=fTrackVariables.GetArray(specifier);
148          specifier+="[trackcount]/f";
149          pTree->Branch(fTrackVariables.GetKey(i), pArray, specifier.Data());
150      }
151      for(i=0; i<fTrackVariablesInt.Variables(); i++){
152          TString specifier=fTrackVariablesInt.GetKey(i);
153          int* pArray=fTrackVariablesInt.GetArray(specifier);
154          specifier+="[trackcount]/i";
155          pTree->Branch(fTrackVariablesInt.GetKey(i), pArray, specifier.Data());
156      }    
157      if(fFillV0==kTRUE){
158         for(i=0; i<fV0Variables.Variables(); i++){
159             TString specifier=fV0Variables.GetKey(i);
160             float* pArray=fV0Variables.GetArray(specifier);
161             specifier+="[V0]/f";
162             pTree->Branch(fV0Variables.GetKey(i), pArray, specifier.Data());
163         }
164      }
165   } else {
166     delete pTree;
167     pTree = NULL;
168   }
169   
170   return pTree;
171 }
172
173 void AliHLTGlobalHistoComponent::FillHistogramDefinitions(){
174   /// default histogram definitions
175 }
176
177 int AliHLTGlobalHistoComponent::FillTree(TTree* pTree, const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/ ){
178
179   /// fill the tree from the ESD
180   int iResult=0;
181   if (!IsDataEvent()) return 0;
182
183   ResetVariables();
184
185   // fetch ESD from input stream
186   const TObject* obj = GetFirstInputObject(kAliHLTAllDataTypes, "AliESDEvent");
187   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj));
188   esd->GetStdContent();
189
190   // fill track variables
191   fNofTracks       = esd->GetNumberOfTracks();
192   fVertexX         = esd->GetPrimaryVertexTracks()->GetX();
193   fVertexY         = esd->GetPrimaryVertexTracks()->GetY();
194   fVertexZ         = esd->GetPrimaryVertexTracks()->GetZ();
195   fNofContributors = esd->GetPrimaryVertexTracks()->GetNContributors();
196   fVertexStatus    = esd->GetPrimaryVertexTracks()->GetStatus();
197   fNofV0s = esd->GetNumberOfV0s();
198   
199   if(fNofV0s > fMaxV0Count){
200      HLTWarning("Found V0s are %d, while component argument is %d, the respective TTree branch is not filled properly. Need to reconfigure.\n", fNofV0s, fMaxV0Count);
201   }
202   
203   for(int i=0; i<fNofTracks; i++){    
204       AliESDtrack *esdTrack = esd->GetTrack(i);
205       if (!esdTrack) continue;
206       
207       Float_t DCAr, DCAz = -99;
208       esdTrack->GetImpactParametersTPC(DCAr, DCAz);
209       
210       fTrackVariables.Fill("Track_pt"        , esdTrack->Pt()                      );
211       fTrackVariables.Fill("Track_phi"       , esdTrack->Phi()*TMath::RadToDeg()   );
212       fTrackVariables.Fill("Track_eta"       , esdTrack->Theta()                   );
213       fTrackVariables.Fill("Track_p"         , esdTrack->P()                       );
214       fTrackVariables.Fill("Track_theta"     , esdTrack->Theta()*TMath::RadToDeg() );
215       fTrackVariables.Fill("Track_Nclusters" , esdTrack->GetTPCNcls()              );
216       fTrackVariables.Fill("Track_status"    , esdTrack->GetStatus()               );
217       fTrackVariables.Fill("Track_charge"    , esdTrack->Charge()                  );
218       fTrackVariables.Fill("Track_DCAr"      , DCAr                                );
219       fTrackVariables.Fill("Track_DCAz"      , DCAz                                );   
220       fTrackVariables.Fill("Track_dEdx"      , esdTrack->GetTPCsignal()            );   
221       fTrackVariablesInt.Fill("Track_status" , esdTrack->GetStatus()               );      
222   }
223   
224   if(fFillV0==kTRUE){
225      for(int i=0; i<fNofV0s; i++){
226      
227      
228          AliESDv0 *esdV0 = esd->GetV0(i);
229          if(!esdV0) continue;
230         
231             
232         //fV0Variables.Fill("V0_AP", ap);
233         //fV0Variables.Fill("V0_pt", pt); 
234         //fV0Variables.Fill("clust1", t1->GetTPCNcls()); 
235         //fV0Variables.Fill("clust2", t2->GetTPCNcls()); 
236         //fV0Variables.Fill("dev1", dev1); 
237         //fV0Variables.Fill("dev2", dev2); 
238         //fV0Variables.Fill("devPrim", devPrim); 
239         //fV0Variables.Fill("length", length); 
240         //fV0Variables.Fill("sigmaLength", sigmaLength); 
241         //fV0Variables.Fill("r", r); 
242           
243      } // end of loop over V0s
244   }
245   
246   if(iResult<0){
247     // fill an empty event
248     ResetVariables();
249   }
250   
251   fEvent++;
252   pTree->Fill();
253   return iResult;
254 }
255
256 int AliHLTGlobalHistoComponent::ResetVariables(){
257 /// reset all filling variables
258   fNofTracks=0;
259   fNofV0s=0;
260   fTrackVariables.ResetCount();
261   fTrackVariablesInt.ResetCount();
262   fV0Variables.ResetCount();
263   return 0;
264 }
265
266 AliHLTComponentDataType AliHLTGlobalHistoComponent::GetOriginDataType() const{
267 // get the origin of the output data
268   return kAliHLTDataTypeHistogram|kAliHLTDataOriginHLT;
269 }
270
271 int AliHLTGlobalHistoComponent::ScanConfigurationArgument(int argc, const char** argv){
272 /// inherited from AliHLTComponent, scan argument
273
274   if (argv==NULL || argc<1) return 0;
275
276   int i=0;
277   TString argument = argv[i];
278
279   // -max-track-count
280   if(argument.CompareTo("-max-track-count")==0){    
281      if (++i>=argc) return -EPROTO; // missing parameter
282      argument = argv[i];
283      fMaxTrackCount = argument.Atoi();    
284      HLTInfo("got %s with parameter %s", argument.Data(), argv[i]);
285      return ++i; 
286   }
287  
288   // -max-V0-count 
289   if(argument.CompareTo("-max-V0-count")==0){    
290      if (++i>=argc) return -EPROTO; // missing parameter
291      argument = argv[i];
292      fMaxV0Count = argument.Atoi();    
293      HLTInfo("got %s with parameter %s", argument.Data(), argv[i]);
294      return ++i;
295   }
296   
297   // -fill-V0
298   if(argument.CompareTo("-fill-V0")==0){
299      fFillV0 = kTRUE;
300      HLTInfo("got %s", argument.Data());
301      return ++i;
302   }
303   // no recognized argument, forward to base class
304   return AliHLTTTreeProcessor::ScanConfigurationArgument(argc, argv);
305 }
306
307 int AliHLTGlobalHistoComponent::Reconfigure(const char* cdbEntry, const char* /*chainId*/){  
308 // see header file for class documentation
309
310   TString cdbPath;
311   if (cdbEntry) {
312     cdbPath=cdbEntry;
313   } else {
314     cdbPath="HLT/ConfigHLT/";
315     cdbPath+=GetComponentID();
316   }
317
318   return ConfigureFromCDBTObjString(cdbPath.Data());
319 }
320
321