HLT TRD update by Konstantin
[u/mrichter/AliRoot.git] / HLT / TRD / AliHLTTRDClusterizerComponent.cxx
1 // $Id$
2
3 /**************************************************************************
4  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5  *                                                                        *
6  * Authors: Matthias Richter <Matthias.Richter@ift.uib.no>                *
7  *          Timm Steinbeck <timm@kip.uni-heidelberg.de>                   *
8  *          for The ALICE Off-line 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   AliHLTTRDClusterizerComponent.cxx
20     @author Timm Steinbeck, Matthias Richter
21     @date   
22     @brief  A TRDClusterizer processing component for the HLT. */
23
24 // see header file for class documentation                                   //
25 // or                                                                        //
26 // refer to README to build package                                          //
27 // or                                                                        //
28 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt                          //
29
30 #if __GNUC__ >= 3
31 using namespace std;
32 #endif
33
34 #include "TTree.h"
35 #include "TFile.h"
36 #include "TBranch.h"
37
38 #include "AliHLTTRDClusterizerComponent.h"
39 #include "AliHLTTRDDefinitions.h"
40 #include "AliHLTTRDCluster.h"
41
42 #include "AliGeomManager.h"
43 #include "AliTRDReconstructor.h"
44 #include "AliCDBManager.h"
45 #include "AliHLTTRDClusterizer.h"
46 #include "AliTRDrecoParam.h"
47 #include "AliTRDrawStreamBase.h"
48 #include "AliTRDcluster.h"
49
50 #include "AliRawReaderMemory.h"
51
52 #ifdef HAVE_VALGRIND_CALLGRIND_H
53 #include <valgrind/callgrind.h>
54 #else
55 #define CALLGRIND_START_INSTRUMENTATION() do { } while (0)
56 #define CALLGRIND_STOP_INSTRUMENTATION() do { } while (0)
57 #endif
58
59 #include <cstdlib>
60 #include <cerrno>
61 #include <string>
62
63 ClassImp(AliHLTTRDClusterizerComponent);
64    
65 AliHLTTRDClusterizerComponent::AliHLTTRDClusterizerComponent():
66   AliHLTProcessor(),
67   fOutputPercentage(100), // By default we copy to the output exactly what we got as input
68   fStrorageDBpath("local://$ALICE_ROOT/OCDB"),
69   fClusterizer(NULL),
70   fRecoParam(NULL),
71   fCDB(NULL),
72   fMemReader(NULL),
73   fReconstructor(NULL),
74   fGeometryFileName("")
75 {
76   // Default constructor
77
78   fGeometryFileName = getenv("ALICE_ROOT");
79   fGeometryFileName += "/HLT/TRD/geometry.root";
80 }
81
82 AliHLTTRDClusterizerComponent::~AliHLTTRDClusterizerComponent()
83 {
84   // Destructor
85   // Work is Done in DoDeInit()
86 }
87
88
89 const char* AliHLTTRDClusterizerComponent::GetComponentID()
90 {
91   // Return the component ID const char *
92   return "TRDClusterizer"; // The ID of this component
93 }
94
95 void AliHLTTRDClusterizerComponent::GetInputDataTypes( vector<AliHLTComponent_DataType>& list)
96 {
97   // Get the list of input data
98   list.clear(); // We do not have any requirements for our input data type(s).
99   list.push_back( AliHLTTRDDefinitions::fgkDDLRawDataType );
100 }
101
102 AliHLTComponent_DataType AliHLTTRDClusterizerComponent::GetOutputDataType()
103 {
104   // Get the output data type
105   return AliHLTTRDDefinitions::fgkClusterDataType;
106 }
107
108 void AliHLTTRDClusterizerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
109 {
110   // Get the output data size
111   constBase = 0;
112   inputMultiplier = ((double)fOutputPercentage)/100.0;
113 }
114
115 AliHLTComponent* AliHLTTRDClusterizerComponent::Spawn()
116 {
117   // Spawn function, return new instance of this class
118   return new AliHLTTRDClusterizerComponent;
119 };
120
121 int AliHLTTRDClusterizerComponent::DoInit( int argc, const char** argv )
122 {
123   // perform initialization. We check whether our relative output size is specified in the arguments.
124   fOutputPercentage = 100;
125   Int_t iRawDataVersion = 2;
126   int i = 0;
127   char* cpErr;
128
129   Int_t iRecoParamType = -1; // default will be the low flux
130
131   // the data type will become obsolete as soon as the formats are established
132   Int_t iRecoDataType = -1; // default will be simulation
133   
134   while ( i < argc )
135     {
136       HLTDebug("argv[%d] == %s", i, argv[i] );
137       if ( !strcmp( argv[i], "output_percentage" ) )
138         {
139           if ( i+1>=argc )
140             {
141               HLTError("Missing output_percentage parameter");
142               return ENOTSUP;
143             }
144           HLTDebug("argv[%d+1] == %s", i, argv[i+1] );
145           fOutputPercentage = strtoul( argv[i+1], &cpErr, 0 );
146           if ( *cpErr )
147             {
148               HLTError("Cannot convert output_percentage parameter '%s'", argv[i+1] );
149               return EINVAL;
150             }
151           HLTInfo("Output percentage set to %lu %%", fOutputPercentage );
152           i += 2;
153         }
154       else if ( strcmp( argv[i], "-cdb" ) == 0)
155         {
156           if ( i+1 >= argc )
157             {
158               HLTError("Missing -cdb argument");
159               return ENOTSUP;         
160             }
161           fStrorageDBpath = argv[i+1];
162           HLTInfo("DB storage is %s", fStrorageDBpath.c_str() );          
163           i += 2;
164           
165         }      
166
167       else if ( strcmp( argv[i], "-lowflux" ) == 0)
168         {
169           iRecoParamType = 0;     
170           HLTDebug("Low flux reco selected.");
171           i++;
172           
173         }      
174
175       else if ( strcmp( argv[i], "-highflux" ) == 0)
176         {
177           iRecoParamType = 1;     
178           HLTDebug("High flux reco selected.");
179           i++;
180           
181         }      
182
183       else if ( strcmp( argv[i], "-cosmics" ) == 0)
184         {
185           iRecoParamType = 2;     
186           HLTDebug("Cosmic test reco selected.");
187           i++;
188           
189         }      
190
191       // raw data type - sim or experiment
192       else if ( strcmp( argv[i], "-simulation" ) == 0)
193         {
194           iRecoDataType = 0;
195           i++;
196           
197         }
198
199       else if ( strcmp( argv[i], "-experiment" ) == 0)
200         {
201           iRecoDataType = 1;
202           i++;
203           
204         }
205
206       else if ( strcmp( argv[i], "-rawver" ) == 0)
207         {
208           if ( i+1 >= argc )
209             {
210               HLTError("Missing -rawver argument");
211               return ENOTSUP;         
212             }
213           iRawDataVersion = atoi( argv[i+1] );
214           HLTInfo("Raw data version is %d", iRawDataVersion );    
215           i += 2;
216           
217         }      
218
219       else if ( strcmp( argv[i], "-geometry" ) == 0)
220         {
221           if ( i+1 >= argc )
222             {
223               HLTError("Missing -geometry argument");
224               return ENOTSUP;         
225             }
226           fGeometryFileName = argv[i+1];
227           HLTInfo("GeomFile storage is %s", fGeometryFileName.c_str() );          
228           i += 2;
229         }      
230       else{
231         HLTError("Unknown option '%s'", argv[i] );
232         return EINVAL;
233       }
234       
235     }
236
237   // THE "REAL" INIT COMES HERE
238
239   if (iRecoParamType < 0 || iRecoParamType > 2)
240     {
241       HLTWarning("No reco param selected. Use -lowflux or -highflux flag. Defaulting to low flux.");
242       iRecoParamType = 0;
243     }
244
245   if (iRecoParamType == 0)
246     {
247       fRecoParam = AliTRDrecoParam::GetLowFluxParam();
248       HLTDebug("Low flux params init.");
249     }
250
251   if (iRecoParamType == 1)
252     {
253       fRecoParam = AliTRDrecoParam::GetHighFluxParam();
254       HLTDebug("High flux params init.");
255     }
256
257   if (iRecoParamType == 2)
258     {
259       fRecoParam = AliTRDrecoParam::GetCosmicTestParam();
260       HLTDebug("Cosmic Test params init.");
261     }
262
263   if (fRecoParam == 0)
264     {
265       HLTError("No reco params initialized. Sniffing big trouble!");
266       return -1;
267     }
268
269   fReconstructor = new AliTRDReconstructor();
270   fReconstructor->SetRecoParam(fRecoParam);
271   fReconstructor->SetStreamLevel(0, AliTRDReconstructor::kClusterizer); // default value
272   HLTInfo("Not writing clusters. I.e. output is a TClonesArray of clusters");
273   fReconstructor->SetOption("hlt,!cw,sl_cf_0");
274   
275   
276   // init the raw data type to be used...
277   // the switch here will become obsolete as soon as the data structures is fixed 
278   // both: in sim and reality
279   if (iRecoDataType < 0 || iRecoDataType > 1)
280     {
281       HLTWarning("No data type selected. Use -simulation or -experiment flag. Defaulting to simulation.");
282       iRecoDataType = 0;
283     }
284
285   if (iRecoDataType == 0)
286     {
287       AliTRDrawStreamBase::SetRawStreamVersion(AliTRDrawStreamBase::kTRDsimStream);
288       HLTDebug("Data type expected is SIMULATION!");
289     }
290
291   if (iRecoDataType == 1)
292     {
293       AliTRDrawStreamBase::SetRawStreamVersion(AliTRDrawStreamBase::kTRDrealStream);
294       HLTDebug("Data type expected is EXPERIMENT!");
295     }
296
297   // the DATA BASE STUFF
298   fCDB = AliCDBManager::Instance();
299   if (!fCDB)
300     {
301       HLTError("Could not get CDB instance", "fCDB 0x%x", fCDB);
302     }
303   else
304     {
305       fCDB->SetRun(0); // THIS HAS TO BE RETRIEVED !!!
306       fCDB->SetDefaultStorage(fStrorageDBpath.c_str());
307       HLTDebug("CDB instance; fCDB 0x%x", fCDB);
308     }
309
310   if((AliGeomManager::GetGeometry()) == NULL){
311     
312     if ( TFile::Open(fGeometryFileName.c_str())) {
313       AliGeomManager::LoadGeometry(fGeometryFileName.c_str());
314     }
315     else {
316       HLTError("Cannot load geometry from file %s",fGeometryFileName.c_str());
317       return EINVAL;
318     }
319   }
320   else
321     HLTInfo("Geometry Already Loaded");
322   
323   fMemReader = new AliRawReaderMemory;
324
325   fClusterizer = new AliHLTTRDClusterizer("TRDCclusterizer", "TRDCclusterizer");
326   fClusterizer->SetReconstructor(fReconstructor);
327   fClusterizer->SetAddLabels(kFALSE);
328   fClusterizer->SetRawVersion(iRawDataVersion);
329   return 0;
330 }
331
332 int AliHLTTRDClusterizerComponent::DoDeinit()
333 {
334   // Deinitialization of the component
335   delete fMemReader;
336   fMemReader = 0;
337   delete fClusterizer;
338   fClusterizer = 0;
339   
340   fReconstructor->SetClusters(0x0);
341   delete fReconstructor;
342   fReconstructor = 0x0;
343   return 0;
344
345
346   if (fCDB)
347     {
348       HLTDebug("destroy fCDB");
349       fCDB->Destroy();
350       fCDB = 0;
351     }
352
353   if (fRecoParam)
354     {
355       HLTDebug("Deleting fRecoParam");
356       delete fRecoParam;
357       fRecoParam = 0;
358     }
359 }
360
361 int AliHLTTRDClusterizerComponent::DoEvent( const AliHLTComponentEventData& evtData, 
362                                             const AliHLTComponentBlockData* blocks, 
363                                             AliHLTComponent_TriggerData& /*trigData*/, 
364                                             AliHLTUInt8_t* outputPtr, 
365                                             AliHLTUInt32_t& size, 
366                                             vector<AliHLTComponent_BlockData>& outputBlocks )
367 {
368   // Process an event
369   HLTDebug( "NofBlocks %lu", evtData.fBlockCnt );
370   // Process an event
371   AliHLTUInt32_t totalSize = 0, offset = 0;
372
373   //implement a usage of the following
374   //   AliHLTUInt32_t triggerDataStructSize = trigData.fStructSize;
375   //   AliHLTUInt32_t triggerDataSize = trigData.fDataSize;
376   //   void *triggerData = trigData.fData;
377   //HLTDebug( "Trigger data received. Struct size %d Data size %d Data location 0x%x", trigData.fStructSize, trigData.fDataSize, (UInt_t*)trigData.fData);
378
379   // Loop over all input blocks in the event
380   AliHLTComponentDataType expectedDataType = (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTRD);
381   for ( unsigned long i = 0; i < evtData.fBlockCnt; i++ )
382     {
383       if (evtData.fEventID == 1)
384         CALLGRIND_START_INSTRUMENTATION();
385       
386       const AliHLTComponentBlockData &block = blocks[i];
387       offset = totalSize;
388       // lets not use the internal TRD data types here : AliHLTTRDDefinitions::fgkDDLRawDataType
389       // which is depreciated - we use HLT global defs instead
390       //      if ( block.fDataType != (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTRD) )
391       AliHLTComponentDataType inputDataType = block.fDataType;
392       if ( inputDataType != expectedDataType)
393         {
394           HLTDebug( "Block # %i/%i; Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s; Skipping",
395                     i, evtData.fBlockCnt,
396                     evtData.fEventID, evtData.fEventID, 
397                     DataType2Text(inputDataType).c_str(), 
398                     DataType2Text(expectedDataType).c_str());
399           continue;
400         }
401       else 
402         {
403         HLTDebug("We get the right data type: Block # %i/%i; Event 0x%08LX (%Lu) Received datatype: %s",
404                     i, evtData.fBlockCnt,
405                     evtData.fEventID, evtData.fEventID, 
406                     DataType2Text(inputDataType).c_str());
407         }
408       
409       //      fMemReader->Reset();
410       fMemReader->SetMemory((UChar_t*) block.fPtr, block.fSize);
411
412       AliHLTUInt32_t spec = block.fSpecification;
413       
414       Int_t id = 1024;
415       
416       for ( Int_t ii = 0; ii < 18 ; ii++ ) {
417         if ( spec & 0x00000001 ) {
418           id += ii;
419           break;
420         }
421         spec = spec >> 1 ;
422       }
423
424       fMemReader->SetEquipmentID( id ); 
425       
426       fClusterizer->SetMemBlock(outputPtr);
427       Bool_t iclustered = fClusterizer->Raw2ClustersChamber(fMemReader);
428       if (iclustered == kTRUE)
429         {
430           HLTDebug( "Clustered successfully");
431         }
432       else
433         {
434           HLTError("Clustering ERROR");
435           return -1;
436         }
437
438       // put the tree into output
439       //fcTree->Print();
440       
441       AliHLTUInt32_t addedSize = fClusterizer->GetAddedSize();
442         if (addedSize > 0){
443           // Using low-level interface 
444           // with interface classes
445           totalSize += addedSize;
446           if ( totalSize > size )
447             {
448               HLTError("Too much data; Data written over allowed buffer. Amount written: %lu, allowed amount: %lu.",
449                        totalSize, size );
450               return EMSGSIZE;
451             }
452                 
453           // Fill block 
454           AliHLTComponentBlockData bd;
455           FillBlockData( bd );
456           bd.fOffset = offset;
457           bd.fSize = addedSize;
458           //bd.fSpecification = spec;
459           bd.fSpecification = gkAliEventTypeData;
460           bd.fDataType = AliHLTTRDDefinitions::fgkClusterDataType;
461           outputBlocks.push_back( bd );
462           HLTDebug( "Block ; size %i; dataType %s; spec 0x%x ",
463                     bd.fSize, DataType2Text(bd.fDataType).c_str(), spec);
464               
465         }
466         else 
467           HLTWarning("Array of clusters is empty!");
468     }
469   fReconstructor->SetClusters(0x0);
470
471   size = totalSize;
472   HLTDebug("Event is done. size written to the output is %i", size);
473   return 0;
474 }
475
476
477 void AliHLTTRDClusterizerComponent::PrintObject( TClonesArray* inClustersArray)
478 {
479   AliTRDcluster* cluster=0x0;
480   
481   for (Int_t i=0; i < inClustersArray->GetEntriesFast(); i++){
482     cluster = dynamic_cast<AliTRDcluster*>(inClustersArray->At(i));
483     HLTDebug("cluster[%i]",i);
484     HLTDebug("  PadCol = %i; PadRow = %i; PadTime = %i", cluster->GetPadCol(), cluster->GetPadRow(), cluster->GetPadTime());
485     HLTDebug("  Detector = %i, Amplitude = %f, Center = %f", cluster->GetDetector(), cluster->GetQ(), cluster->GetCenter());
486     HLTDebug("  LocalTimeBin =  %i; NPads = %i; maskedPosition: %s, status: %s", cluster->GetLocalTimeBin(), cluster->GetNPads(),cluster->GetPadMaskedPosition(),cluster->GetPadMaskedPosition());
487   }
488   
489 }
490