7e92e3af3714081a5b223fbc8c7f926f42aaa342
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCEsdWriterComponent.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   AliHLTTPCEsdWriterComponent.cxx
20     @author Matthias Richter
21     @date   
22     @brief  Writer component to store tracks of the HLT TPC conformal
23             mapping tracker in the AliESD format
24 */
25
26 // see header file for class documentation
27 // or
28 // refer to README to build package
29 // or
30 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
31
32 #include <cassert>
33 #include "AliHLTTPCEsdWriterComponent.h"
34 #include "AliESDEvent.h"
35 #include "AliESDtrack.h"
36 #include "TTree.h"
37 #include "TList.h"
38 #include "AliHLTTPCTrack.h"
39 #include "AliHLTTPCTrackArray.h"
40 #include "AliHLTTPCTrackletDataFormat.h"
41 #include "AliHLTTPCDefinitions.h"
42 #include "AliHLTTPCTransform.h"
43
44 /** ROOT macro for the implementation of ROOT specific class methods */
45 ClassImp(AliHLTTPCEsdWriterComponent)
46
47 AliHLTTPCEsdWriterComponent::AliHLTTPCEsdWriterComponent()
48 {
49   // see header file for class documentation
50   // or
51   // refer to README to build package
52   // or
53   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
54 }
55
56 AliHLTTPCEsdWriterComponent::~AliHLTTPCEsdWriterComponent()
57 {
58   // see header file for class documentation
59 }
60
61 AliHLTTPCEsdWriterComponent::AliWriter::AliWriter()
62   :
63   fTree(NULL),
64   fESD(NULL),
65   fBase(new AliHLTTPCEsdWriterComponent)
66 {
67   // see header file for class documentation
68   // or
69   // refer to README to build package
70   // or
71   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
72 }
73
74 AliHLTTPCEsdWriterComponent::AliWriter::~AliWriter()
75 {
76   // see header file for class documentation
77   if (fBase) delete fBase;
78   fBase=NULL;
79 }
80
81 void AliHLTTPCEsdWriterComponent::AliWriter::GetInputDataTypes(AliHLTComponentDataTypeList& list)
82 {
83   // see header file for class documentation
84   list.push_back(AliHLTTPCDefinitions::fgkTrackSegmentsDataType);
85   list.push_back(AliHLTTPCDefinitions::fgkTracksDataType);
86 }
87
88 int AliHLTTPCEsdWriterComponent::AliWriter::InitWriter()
89 {
90   // see header file for class documentation
91   int iResult=0;
92   fESD = new AliESDEvent;
93   if (fESD) {
94     fESD->CreateStdContent();
95     fTree = new TTree("esdTree", "Tree with HLT ESD objects");
96     if (fTree) {
97       fTree->SetDirectory(0);
98       fESD->WriteToTree(fTree);
99     }
100   }
101   if (fTree==NULL) {
102     iResult=-ENOMEM;
103   }
104   return iResult;
105 }
106
107 int AliHLTTPCEsdWriterComponent::AliWriter::CloseWriter()
108 {
109   // see header file for class documentation
110   int iResult=0;
111   if (fTree) {
112     // the esd structure is written to the user info and is
113     // needed in te ReadFromTree method to read all objects correctly
114     if (fESD) fTree->GetUserInfo()->Add(fESD);
115     WriteObject(kAliHLTVoidEventID, fTree);
116     fTree->GetUserInfo()->Clear();
117     TTree* pTree=fTree;
118     fTree=NULL;
119     delete pTree;
120   } else {
121     HLTWarning("not initialized");
122   }
123
124   if (fESD) {
125     delete fESD;
126   }
127   iResult=AliHLTRootFileWriterComponent::CloseWriter();
128   return iResult;
129 }
130
131 int AliHLTTPCEsdWriterComponent::AliWriter::DumpEvent( const AliHLTComponentEventData& evtData,
132                                             const AliHLTComponentBlockData* blocks, 
133                                             AliHLTComponentTriggerData& /*trigData*/ )
134 {
135   // see header file for class documentation
136   int iResult=0;
137   TTree* pTree=fTree;
138   assert(fBase);
139   if (pTree && fBase) {
140     if (fESD) {
141       AliESDEvent* pESD=fESD;
142
143       iResult=fBase->ProcessBlocks(pTree, pESD, blocks, (int)evtData.fBlockCnt);
144
145     } else {
146       iResult=-ENOMEM;
147     }
148   }
149   return iResult;
150 }
151
152 int AliHLTTPCEsdWriterComponent::AliWriter::ScanArgument(int argc, const char** argv)
153 {
154   // see header file for class documentation
155   int iResult=AliHLTRootFileWriterComponent::ScanArgument(argc, argv);
156   return iResult;
157 }
158
159 int AliHLTTPCEsdWriterComponent::ProcessBlocks(TTree* pTree, AliESDEvent* pESD,
160                                                const AliHLTComponentBlockData* blocks,
161                                                int nBlocks, int* pMinSlice,
162                                                int* pMaxSlice)
163 {
164   // see header file for class documentation
165   int iResult=0;
166   if (pESD && blocks) {
167       const AliHLTComponentBlockData* iter = NULL;
168       AliHLTTPCTrackletData* inPtr=NULL;
169       int bIsTrackSegs=0;
170  
171       for (int ndx=0; ndx<nBlocks && iResult>=0; ndx++) {
172         iter = blocks+ndx;
173         if ( (bIsTrackSegs=(iter->fDataType == AliHLTTPCDefinitions::fgkTrackSegmentsDataType))==1 ||
174              iter->fDataType == AliHLTTPCDefinitions::fgkTracksDataType ) {
175           Int_t minslice=AliHLTTPCDefinitions::GetMinSliceNr(iter->fSpecification);
176           Int_t maxslice=AliHLTTPCDefinitions::GetMaxSliceNr(iter->fSpecification);
177           if (bIsTrackSegs==0) {
178             // slice parameter and data specification ignored, tracks already in global coordinates
179             minslice=-1;
180             maxslice=-1;
181             if (pMinSlice) *pMinSlice=0;
182             if (pMaxSlice) *pMaxSlice=AliHLTTPCTransform::GetNSlice()-1;
183           } else {
184             if (pMinSlice && (*pMinSlice==-1 || *pMinSlice>minslice)) *pMinSlice=minslice;
185             if (pMaxSlice && (*pMaxSlice==-1 || *pMaxSlice<maxslice)) *pMaxSlice=maxslice;
186           }
187           //HLTDebug("dataspec %#x minslice %d", iter->fSpecification, minslice);
188           if (minslice >=-1 && minslice<AliHLTTPCTransform::GetNSlice()) {
189             if (minslice!=maxslice) {
190               HLTWarning("data from multiple sectors in one block: "
191                          "possible mismatch in treatment of local coordinate system");
192             }
193             AliHLTTPCTrackArray tracks;
194             inPtr=(AliHLTTPCTrackletData*)iter->fPtr;
195             HLTDebug("reading block %d (slice %d): %d tracklets", ndx, minslice, inPtr->fTrackletCnt);
196             tracks.FillTracks(inPtr->fTrackletCnt, inPtr->fTracklets, minslice, 0/*don't rotate*/);
197             if ((iResult=Tracks2ESD(&tracks, pESD))>=0) {
198             }
199           } else {
200             HLTError("invalid sector number");
201             iResult=-EBADF;
202           }
203         }
204       }
205       if (iResult>=0 && pTree) {
206         pTree->Fill();
207       }
208
209       pESD->Reset();
210     
211   } else {
212     iResult=-EINVAL;
213   }
214   return iResult;
215 }
216
217 int AliHLTTPCEsdWriterComponent::Tracks2ESD(AliHLTTPCTrackArray* pTracks, AliESDEvent* pESD)
218 {
219   // see header file for class documentation
220   int iResult=0;
221   if (pTracks && pESD) {
222     HLTDebug("converting %d tracks from track array", pTracks->GetNTracks());
223     for (int i=0; i<pTracks->GetNTracks() && iResult>=0; i++) {
224       AliHLTTPCTrack* pTrack=(*pTracks)[i];
225       if (pTrack) {
226         //HLTDebug("convert track %d", i);
227         //pTrack->Print();
228         int iLocal=pTrack->Convert2AliKalmanTrack();
229         if (iLocal>=0) {
230         AliESDtrack iotrack;
231         iotrack.UpdateTrackParams(pTrack,AliESDtrack::kTPCin);
232         iotrack.SetTPCPoints(pTrack->GetPoints());
233         pESD->AddTrack(&iotrack);
234         } else {
235           HLTError("conversion to AliKalmanTrack failed for track %d of %d", i, pTracks->GetNTracks());
236         }
237       } else {
238         HLTError("internal mismatch in array");
239         iResult=-EFAULT;
240       }
241     }
242     
243   } else {
244     iResult=-EINVAL;
245   }
246   return iResult;
247 }
248
249 AliHLTTPCEsdWriterComponent::AliConverter::AliConverter()
250   :
251   fBase(new AliHLTTPCEsdWriterComponent),
252   fWriteTree(1)
253 {
254   // see header file for class documentation
255   // or
256   // refer to README to build package
257   // or
258   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
259 }
260
261 AliHLTTPCEsdWriterComponent::AliConverter::~AliConverter()
262 {
263   // see header file for class documentation
264   if (fBase) delete fBase;
265   fBase=NULL;
266 }
267
268 void AliHLTTPCEsdWriterComponent::AliConverter::GetInputDataTypes(AliHLTComponentDataTypeList& list)
269 {
270   // see header file for class documentation
271   list.push_back(AliHLTTPCDefinitions::fgkTrackSegmentsDataType);
272   list.push_back(AliHLTTPCDefinitions::fgkTracksDataType);
273 }
274
275 AliHLTComponentDataType AliHLTTPCEsdWriterComponent::AliConverter::GetOutputDataType()
276 {
277   // see header file for class documentation
278   return kAliHLTDataTypeESDTree;
279 }
280
281 void AliHLTTPCEsdWriterComponent::AliConverter::GetOutputDataSize(unsigned long& constBase, double& inputMultiplier)
282 {
283   // see header file for class documentation
284   constBase=2000000;
285   inputMultiplier=10.0;
286 }
287
288 int AliHLTTPCEsdWriterComponent::AliConverter::DoInit(int argc, const char** argv)
289 {
290   // see header file for class documentation
291   int iResult=0;
292   TString argument="";
293   int bMissingParam=0;
294   for (int i=0; i<argc && iResult>=0; i++) {
295     argument=argv[i];
296     if (argument.IsNull()) continue;
297
298     // -notree
299     if (argument.CompareTo("-notree")==0) {
300       fWriteTree=0;
301
302       // -tree
303     } else if (argument.CompareTo("-tree")==0) {
304       fWriteTree=1;
305
306     } else {
307       HLTError("unknown argument %s", argument.Data());
308       break;
309     }
310   }
311   if (bMissingParam) {
312     HLTError("missing parameter for argument %s", argument.Data());
313     iResult=-EINVAL;
314   }
315
316   return iResult;
317 }
318
319 int AliHLTTPCEsdWriterComponent::AliConverter::DoDeinit()
320 {
321   // see header file for class documentation
322   return 0;
323 }
324
325 int AliHLTTPCEsdWriterComponent::AliConverter::DoEvent(const AliHLTComponentEventData& evtData, 
326                                                        const AliHLTComponentBlockData* blocks, 
327                                                        AliHLTComponentTriggerData& /*trigData*/,
328                                                        AliHLTUInt8_t* /*outputPtr*/, 
329                                                        AliHLTUInt32_t& /*size*/,
330                                                        AliHLTComponentBlockDataList& /*outputBlocks*/ )
331 {
332   // see header file for class documentation
333   int iResult=0;
334   assert(fBase);
335   AliESDEvent* pESD = new AliESDEvent;
336   if (pESD && fBase) {
337     pESD->CreateStdContent();
338     TTree* pTree = NULL;
339     // TODO: Matthias 06.12.2007
340     // Tried to write the ESD directly instead to a tree, but this did not work
341     // out. Information in the ESD is different, needs investigation.
342     if (fWriteTree)
343       pTree = new TTree("esdTree", "Tree with HLT ESD objects");
344     if (pTree) {
345       pTree->SetDirectory(0);
346       pESD->WriteToTree(pTree);
347     }
348
349     if ((iResult=fBase->ProcessBlocks(pTree, pESD, blocks, (int)evtData.fBlockCnt))>=0) {
350         // TODO: set the specification correctly
351       if (pTree) {
352         // the esd structure is written to the user info and is
353         // needed in te ReadFromTree method to read all objects correctly
354         pTree->GetUserInfo()->Add(pESD);
355         iResult=PushBack(pTree, kAliHLTDataTypeESDTree|kAliHLTDataOriginTPC, 0);
356       } else {
357         iResult=PushBack(pESD, kAliHLTDataTypeESDObject|kAliHLTDataOriginTPC, 0);
358       }
359     }
360     if (pTree) {
361       // clear user info list to prevent objects from being deleted
362       pTree->GetUserInfo()->Clear();
363       delete pTree;
364     }
365
366     delete pESD;
367   }
368   return iResult;
369 }
370