update by Gaute:
[u/mrichter/AliRoot.git] / HLT / ITS / AliHLTITSClusterHistoComponent.cxx
1 // $Id$
2 //**************************************************************************
3 //* This file is property of and copyright by the ALICE HLT Project        * 
4 //* ALICE Experiment at CERN, All rights reserved.                         *
5 //*                                                                        *
6 //* Primary Authors: Gaute Ovrebekk <ovrebekk@ift.uib.no>                  *
7 //*                  for The ALICE HLT Project.                            *
8 //*                                                                        *
9 //* Permission to use, copy, modify and distribute this software and its   *
10 //* documentation strictly for non-commercial purposes is hereby granted   *
11 //* without fee, provided that the above copyright notice appears in all   *
12 //* copies and that both the copyright notice and this permission notice   *
13 //* appear in the supporting documentation. The authors make no claims     *
14 //* about the suitability of this software for any purpose. It is          *
15 //* provided "as is" without express or implied warranty.                  *
16 //**************************************************************************
17
18 /** @file   AliHLTITSClusterHistoComponent.cxx
19     @author Gaute Ovrebekk
20     @brief  Component for ploting charge in clusters
21 */
22
23 #if __GNUC__>= 3
24 using namespace std;
25 #endif
26
27 #include "AliHLTITSClusterHistoComponent.h"
28 #include "AliHLTITSClusterDataFormat.h"
29 #include "AliCDBEntry.h"
30 #include "AliCDBManager.h"
31 #include "AliITSRecPoint.h"
32 #include <TFile.h>
33 #include <TString.h>
34 #include "TObjString.h"
35 #include "TObjArray.h"
36
37 //#include <stdlib.h>
38 //#include <cerrno>
39
40 /** ROOT macro for the implementation of ROOT specific class methods */
41 ClassImp(AliHLTITSClusterHistoComponent)
42
43 AliHLTITSClusterHistoComponent::AliHLTITSClusterHistoComponent()
44 :
45 fXY(NULL),                     
46   fXYZ(NULL),                   
47   fCharge(NULL),   
48   fPlotCharge(kFALSE),   
49   fPlotXY(kTRUE),
50   fPlotXYZ(kFALSE) 
51 {
52   // see header file for class documentation
53   // or
54   // refer to README to build package
55   // or
56   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
57
58 }
59
60 AliHLTITSClusterHistoComponent::~AliHLTITSClusterHistoComponent()
61 {
62   // see header file for class documentation
63 }
64
65 // Public functions to implement AliHLTComponent's interface.
66 // These functions are required for the registration process
67
68 const char* AliHLTITSClusterHistoComponent::GetComponentID()
69 {
70   // see header file for class documentation
71   
72   return "ITSClusterHisto";
73 }
74
75 void AliHLTITSClusterHistoComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list)
76 {
77   // see header file for class documentation
78   list.clear();
79   list.push_back( kAliHLTDataTypeTObjArray );
80 }
81
82 AliHLTComponentDataType AliHLTITSClusterHistoComponent::GetOutputDataType()
83 {
84   // see header file for class documentation
85   return kAliHLTDataTypeHistogram;
86
87 }
88
89 void AliHLTITSClusterHistoComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
90 {
91   // see header file for class documentation
92   // XXX TODO: Find more realistic values.
93   constBase = 80000;
94   inputMultiplier = 10;
95 }
96
97 AliHLTComponent* AliHLTITSClusterHistoComponent::Spawn()
98 {
99   // see header file for class documentation
100   return new AliHLTITSClusterHistoComponent;
101 }
102
103 int AliHLTITSClusterHistoComponent::DoInit( int argc, const char** argv )
104 {
105   fPlotCharge=kFALSE;   
106   fPlotXY=kTRUE;
107   fPlotXYZ=kFALSE; 
108      
109   if(fPlotCharge){fCharge = new TH1F("fCharge","Total Charge of clusters",4000,0,4000);}
110   if(fPlotXY){fXY = new TH2F("fXY","Global XY of ITS clusters",1600,-80,80,1600,-80,80);}
111   if(fPlotXYZ){fXYZ = new TH3F("fXYZ","Global XYZ of ITS clusters",1600,-80,80,1600,-80,80,2000,-100,100);}
112   
113   int iResult=0;
114   TString configuration="";
115   TString argument="";
116   for (int i=0; i<argc && iResult>=0; i++) {
117     argument=argv[i];
118     if (!configuration.IsNull()) configuration+=" ";
119     configuration+=argument;
120   }
121   
122   if (!configuration.IsNull()) {
123     iResult=Configure(configuration.Data());
124   }  
125
126   return iResult; 
127 }
128   
129 int AliHLTITSClusterHistoComponent::DoDeinit()
130 {
131   // see header file for class documentation
132   if(fCharge!=NULL) delete fCharge;
133   if(fXY!=NULL) delete fXY;     
134   if(fXYZ!=NULL) delete fXYZ;
135   return 0;
136 }
137
138 int AliHLTITSClusterHistoComponent::DoEvent(const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/)
139 {
140   
141   int TotalSpacePoint = 0;
142   
143   const AliHLTComponentBlockData* iter = NULL;
144   
145   if(!IsDataEvent())
146     return 0;
147   
148   for ( iter = GetFirstInputBlock(kAliHLTDataTypeClusters|kAliHLTDataOriginITSSSD); iter != NULL; iter = GetNextInputBlock() ) {
149     
150     const AliHLTITSClusterData* clusterData = (const AliHLTITSClusterData*) iter->fPtr;
151     Int_t nSpacepoint = (Int_t) clusterData->fSpacePointCnt;
152     TotalSpacePoint += nSpacepoint;
153     AliHLTITSSpacePointData *clusters = (AliHLTITSSpacePointData*) clusterData->fSpacePoints;
154
155     for(int i=0;i<nSpacepoint;i++){
156       Int_t lab[4]={0,0,0,0};
157       Float_t hit[6]={0,0,0,0,0,0};
158       Int_t info[3]={0,0,0};
159       
160       lab[0]=clusters[i].fTracks[0];
161       lab[1]=clusters[i].fTracks[1];
162       lab[2]=clusters[i].fTracks[2];
163       lab[3]=clusters[i].fIndex;
164       hit[0]=clusters[i].fY;
165       hit[1]=clusters[i].fZ;
166       hit[2]=clusters[i].fSigmaY2;
167       hit[3]=clusters[i].fSigmaZ2;
168       hit[4]=clusters[i].fQ;
169       hit[5]=clusters[i].fSigmaYZ;
170       info[0]=clusters[i].fNy;
171       info[1]=clusters[i].fNz;
172       info[2]=clusters[i].fLayer;
173
174       Float_t xyz[3];
175       AliITSRecPoint recpoint(lab,hit,info);
176       recpoint.GetGlobalXYZ(xyz);
177       if(fPlotXY){
178         fXY->Fill(xyz[0],xyz[1]);
179       }
180       if(fPlotXYZ){
181         fXYZ->Fill(xyz[0],xyz[1],xyz[2]);
182       }
183       if(fPlotCharge){
184         fCharge->Fill(recpoint.GetQ());
185       }
186     }
187   }
188   
189   if(fPlotCharge){
190     AliHLTUInt32_t fSpecification = 0x0;
191     PushBack( (TObject*) fCharge,kAliHLTDataTypeHistogram,fSpecification);
192   }
193   if(fPlotXY){
194     AliHLTUInt32_t fSpecification = 0x0;
195     PushBack( (TObject*) fXY,kAliHLTDataTypeHistogram,fSpecification);
196   }
197   if(fPlotXYZ){
198     AliHLTUInt32_t fSpecification = 0x0;
199     PushBack( (TObject*) fXYZ,kAliHLTDataTypeHistogram,fSpecification);
200   }
201   
202   HLTInfo("ITSClusterHisto found %d Total Spacepoints", TotalSpacePoint);
203   
204   return 0;
205 }
206
207 int AliHLTITSClusterHistoComponent::Configure(const char* arguments)
208 {
209   
210   int iResult=0;
211   
212   if (!arguments) return iResult;
213   
214   TString allArgs=arguments;
215   TString argument;
216   
217   TObjArray* pTokens=allArgs.Tokenize(" ");
218   
219   if (pTokens) {
220     for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
221       argument=((TObjString*)pTokens->At(i))->GetString();
222       if (argument.IsNull()) continue;
223       
224       if (argument.CompareTo("-plot-all")==0) {
225         HLTInfo("Ploting all historgams");
226         fPlotXY = kTRUE;
227         fPlotXYZ = kTRUE;
228         fPlotCharge = kTRUE;
229         continue;
230       }
231       
232       else if (argument.CompareTo("-plot-xy")==0) {
233         HLTInfo("Ploting Global XY");
234         fPlotXY = kTRUE;
235         continue;
236       }
237
238       else if (argument.CompareTo("-plot-xyz")==0) {
239         HLTInfo("Ploting Global XYZ");
240         //fPlotXYZ = kTRUE;
241         continue;
242       }
243       else if (argument.CompareTo("-plot-charge")==0) {
244         HLTInfo("Ploting charge of clusters");
245         fPlotCharge = kTRUE;
246         continue;
247       }
248
249       else {
250         HLTError("unknown argument %s", argument.Data());
251         iResult=-EINVAL;
252         break;
253       }
254     }
255     delete pTokens;
256   }
257   
258   if(!fCharge && fPlotCharge){fCharge = new TH1F("fCharge","Total Charge of clusters",4000,0,4000);}
259   if(!fXY && fPlotXY){fXY = new TH2F("fXY","Global XY of ITS clusters",1600,-80,80,1600,-80,80);}
260   if(!fXYZ && fPlotXYZ){fXYZ = new TH3F("fXYZ","Global XYZ of ITS clusters",1600,-80,80,1600,-80,80,2000,-100,100);}
261   
262   return iResult;
263 }
264
265 int AliHLTITSClusterHistoComponent::Reconfigure(const char* cdbEntry, const char* chainId)
266 {
267   // see header file for class documentation
268   int iResult=0;
269   
270   const char* path="HLT/ConfigITS/HistoComponent";
271   const char* defaultNotify="";
272   if (cdbEntry) {
273     path=cdbEntry;
274     defaultNotify=" (default)";
275   }
276   if (path) {
277     HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
278     AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path);
279     if (pEntry) {
280       TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
281       if (pString) {
282         HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
283         iResult=Configure(pString->GetString().Data());
284       } else {
285         HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
286       }
287     } else {
288       HLTError("can not fetch object \"%s\" from CDB", path);
289     }
290   }
291   
292   return iResult;
293 }