eb67f03ef831fb214680d1d468f87ae0dd112dc5
[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 "AliTRDclusterizerHLT.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"),
69   fClusterizer(NULL),
70   fRecoParam(NULL),
71   fCDB(NULL),
72   fMemReader(NULL),
73   fReconstructor(NULL),
74   fGeometryFileName(""),
75   fUseHLTClusters(kFALSE)
76 {
77   // Default constructor
78
79   fGeometryFileName = getenv("ALICE_ROOT");
80   fGeometryFileName += "/HLT/TRD/geometry.root";
81 }
82
83 AliHLTTRDClusterizerComponent::~AliHLTTRDClusterizerComponent()
84 {
85   // Destructor
86   // Work is Done in DoDeInit()
87 }
88
89 /**
90  * Convert AliTRDcluster to AliHLTTRDCluster 
91  * Add HLTCluster to the output, defined by pointer
92  * Fill block desctiptors 
93  * Return size of the added to ouput objects
94  */
95 //============================================================================
96 UInt_t AliHLTTRDClusterizerComponent::AddToOutput(TClonesArray* inClusterArray, AliHLTUInt8_t* outBlockPtr)
97 {
98   AliTRDcluster* cluster = 0;
99   AliHLTUInt32_t addedSize = 0;
100   //  == OUTdatatype pointer
101   AliHLTTRDCluster * outPtr = (AliHLTTRDCluster*)outBlockPtr;
102   
103   if (inClusterArray){
104     Int_t nbEntries  = inClusterArray->GetEntries();
105     for (Int_t iCluster = 0; iCluster<nbEntries; iCluster++){
106       //cout << "Geting cluster #" << iCluster << endl;
107       UInt_t blockSize=0;
108       
109       cluster = dynamic_cast<AliTRDcluster*>(inClusterArray->At(iCluster));
110
111       AliHLTTRDCluster *hltCluster = new (outPtr) AliHLTTRDCluster(cluster);
112       //cout << Form("cluster added at 0x%x (%i)\n",outPtr,outPtr);
113
114       blockSize = sizeof(*hltCluster);
115
116       addedSize += blockSize;
117       outBlockPtr += blockSize;
118       outPtr = (AliHLTTRDCluster*)outBlockPtr;
119     }
120   }
121   return addedSize;
122   
123 }
124
125 const char* AliHLTTRDClusterizerComponent::GetComponentID()
126 {
127   // Return the component ID const char *
128   return "TRDClusterizer"; // The ID of this component
129 }
130
131 void AliHLTTRDClusterizerComponent::GetInputDataTypes( vector<AliHLTComponent_DataType>& list)
132 {
133   // Get the list of input data
134   list.clear(); // We do not have any requirements for our input data type(s).
135   list.push_back( AliHLTTRDDefinitions::fgkDDLRawDataType );
136 }
137
138 AliHLTComponent_DataType AliHLTTRDClusterizerComponent::GetOutputDataType()
139 {
140   // Get the output data type
141   return AliHLTTRDDefinitions::fgkClusterDataType;
142 }
143
144 void AliHLTTRDClusterizerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
145 {
146   // Get the output data size
147   constBase = 0;
148   inputMultiplier = ((double)fOutputPercentage)/100.0;
149 }
150
151 AliHLTComponent* AliHLTTRDClusterizerComponent::Spawn()
152 {
153   // Spawn function, return new instance of this class
154   return new AliHLTTRDClusterizerComponent;
155 };
156
157 int AliHLTTRDClusterizerComponent::DoInit( int argc, const char** argv )
158 {
159   // perform initialization. We check whether our relative output size is specified in the arguments.
160   fOutputPercentage = 100;
161   Int_t iRawDataVersion = 2;
162   int i = 0;
163   char* cpErr;
164   Bool_t bWriteClusters = kFALSE;
165   
166
167   Int_t iRecoParamType = -1; // default will be the low flux
168
169   // the data type will become obsolete as soon as the formats are established
170   Int_t iRecoDataType = -1; // default will be simulation
171   
172   while ( i < argc )
173     {
174       HLTDebug("argv[%d] == %s", i, argv[i] );
175       if ( !strcmp( argv[i], "output_percentage" ) )
176         {
177           if ( i+1>=argc )
178             {
179               HLTError("Missing output_percentage parameter");
180               return ENOTSUP;
181             }
182           HLTDebug("argv[%d+1] == %s", i, argv[i+1] );
183           fOutputPercentage = strtoul( argv[i+1], &cpErr, 0 );
184           if ( *cpErr )
185             {
186               HLTError("Cannot convert output_percentage parameter '%s'", argv[i+1] );
187               return EINVAL;
188             }
189           HLTInfo("Output percentage set to %lu %%", fOutputPercentage );
190           i += 2;
191         }
192       else if ( strcmp( argv[i], "-cdb" ) == 0)
193         {
194           if ( i+1 >= argc )
195             {
196               HLTError("Missing -cdb argument");
197               return ENOTSUP;         
198             }
199           fStrorageDBpath = argv[i+1];
200           HLTInfo("DB storage is %s", fStrorageDBpath.c_str() );          
201           i += 2;
202           
203         }      
204
205       else if ( strcmp( argv[i], "-lowflux" ) == 0)
206         {
207           iRecoParamType = 0;     
208           HLTDebug("Low flux reco selected.");
209           i++;
210           
211         }      
212
213       else if ( strcmp( argv[i], "-highflux" ) == 0)
214         {
215           iRecoParamType = 1;     
216           HLTDebug("High flux reco selected.");
217           i++;
218           
219         }      
220
221       else if ( strcmp( argv[i], "-cosmics" ) == 0)
222         {
223           iRecoParamType = 2;     
224           HLTDebug("Cosmic test reco selected.");
225           i++;
226           
227         }      
228
229       // raw data type - sim or experiment
230       else if ( strcmp( argv[i], "-simulation" ) == 0)
231         {
232           iRecoDataType = 0;
233           i++;
234           
235         }
236
237       else if ( strcmp( argv[i], "-experiment" ) == 0)
238         {
239           iRecoDataType = 1;
240           i++;
241           
242         }
243
244       else if ( strcmp( argv[i], "-rawver" ) == 0)
245         {
246           if ( i+1 >= argc )
247             {
248               HLTError("Missing -rawver argument");
249               return ENOTSUP;         
250             }
251           iRawDataVersion = atoi( argv[i+1] );
252           HLTInfo("Raw data version is %d", iRawDataVersion );    
253           i += 2;
254           
255         }      
256
257       else if ( strcmp( argv[i], "-geometry" ) == 0)
258         {
259           if ( i+1 >= argc )
260             {
261               HLTError("Missing -geometry argument");
262               return ENOTSUP;         
263             }
264           fGeometryFileName = argv[i+1];
265           HLTInfo("GeomFile storage is %s", fGeometryFileName.c_str() );          
266           i += 2;
267           
268         }      
269       else if ( strcmp( argv[i], "-writeClusters" ) == 0)
270         {
271           bWriteClusters = kTRUE;
272           i++;
273         }
274       else if ( strcmp( argv[i], "-useHLTClusters" ) == 0)
275         {
276           fUseHLTClusters = kTRUE;
277           i++;
278           HLTInfo("Using AliHLTCluster to pass data further in the chain");
279         }
280       else{
281         HLTError("Unknown option '%s'", argv[i] );
282         return EINVAL;
283       }
284       
285     }
286
287   // THE "REAL" INIT COMES HERE
288
289   if (iRecoParamType < 0 || iRecoParamType > 2)
290     {
291       HLTWarning("No reco param selected. Use -lowflux or -highflux flag. Defaulting to low flux.");
292       iRecoParamType = 0;
293     }
294
295   if (iRecoParamType == 0)
296     {
297       fRecoParam = AliTRDrecoParam::GetLowFluxParam();
298       HLTDebug("Low flux params init.");
299     }
300
301   if (iRecoParamType == 1)
302     {
303       fRecoParam = AliTRDrecoParam::GetHighFluxParam();
304       HLTDebug("High flux params init.");
305     }
306
307   if (iRecoParamType == 2)
308     {
309       fRecoParam = AliTRDrecoParam::GetCosmicTestParam();
310       HLTDebug("Cosmic Test params init.");
311     }
312
313   if (fRecoParam == 0)
314     {
315       HLTError("No reco params initialized. Sniffing big trouble!");
316       return -1;
317     }
318
319   fReconstructor = new AliTRDReconstructor();
320   fReconstructor->SetRecoParam(fRecoParam);
321   fReconstructor->SetStreamLevel(0, AliTRDReconstructor::kClusterizer); // default value
322   if (bWriteClusters)
323     {
324       HLTInfo("Writing clusters. I.e. output is a TTree with clusters");
325       fReconstructor->SetOption("cw,sl_cf_0");
326     }
327   else
328     {
329       HLTInfo("Not writing clusters. I.e. output is a TClonesArray of clusters");
330       fReconstructor->SetOption("!cw,sl_cf_0");
331     }
332   
333   
334   // init the raw data type to be used...
335   // the switch here will become obsolete as soon as the data structures is fixed 
336   // both: in sim and reality
337   if (iRecoDataType < 0 || iRecoDataType > 1)
338     {
339       HLTWarning("No data type selected. Use -simulation or -experiment flag. Defaulting to simulation.");
340       iRecoDataType = 0;
341     }
342
343   if (iRecoDataType == 0)
344     {
345       AliTRDrawStreamBase::SetRawStreamVersion(AliTRDrawStreamBase::kTRDsimStream);
346       HLTDebug("Data type expected is SIMULATION!");
347     }
348
349   if (iRecoDataType == 1)
350     {
351       AliTRDrawStreamBase::SetRawStreamVersion(AliTRDrawStreamBase::kTRDrealStream);
352       HLTDebug("Data type expected is EXPERIMENT!");
353     }
354
355   // the DATA BASE STUFF
356   fCDB = AliCDBManager::Instance();
357   if (!fCDB)
358     {
359       HLTError("Could not get CDB instance", "fCDB 0x%x", fCDB);
360     }
361   else
362     {
363       fCDB->SetRun(0); // THIS HAS TO BE RETRIEVED !!!
364       fCDB->SetDefaultStorage(fStrorageDBpath.c_str());
365       HLTDebug("CDB instance; fCDB 0x%x", fCDB);
366     }
367
368   if((AliGeomManager::GetGeometry()) == NULL){
369     
370     if ( TFile::Open(fGeometryFileName.c_str())) {
371       AliGeomManager::LoadGeometry(fGeometryFileName.c_str());
372     }
373     else {
374       HLTError("Cannot load geometry from file %s",fGeometryFileName.c_str());
375       return EINVAL;
376     }
377   }
378   else
379     HLTInfo("Geometry Already Loaded");
380   
381   fMemReader = new AliRawReaderMemory;
382
383   fClusterizer = new AliTRDclusterizerHLT("TRDCclusterizer", "TRDCclusterizer");
384   fClusterizer->SetReconstructor(fReconstructor);
385
386   fClusterizer->SetRawVersion(iRawDataVersion);
387   fClusterizer->InitClusterTree();
388   return 0;
389 }
390
391 int AliHLTTRDClusterizerComponent::DoDeinit()
392 {
393   // Deinitialization of the component
394   delete fMemReader;
395   fMemReader = 0;
396   delete fClusterizer;
397   fClusterizer = 0;
398   
399   fReconstructor->SetClusters(0x0);
400   delete fReconstructor;
401   fReconstructor = 0x0;
402   return 0;
403
404
405   if (fCDB)
406     {
407       HLTDebug("destroy fCDB");
408       fCDB->Destroy();
409       fCDB = 0;
410     }
411
412   if (fRecoParam)
413     {
414       HLTDebug("Deleting fRecoParam");
415       delete fRecoParam;
416       fRecoParam = 0;
417     }
418 }
419
420 int AliHLTTRDClusterizerComponent::DoEvent( const AliHLTComponentEventData& evtData, 
421                                             const AliHLTComponentBlockData* blocks, 
422                                             AliHLTComponent_TriggerData& /*trigData*/, 
423                                             AliHLTUInt8_t* outputPtr, 
424                                             AliHLTUInt32_t& size, 
425                                             vector<AliHLTComponent_BlockData>& outputBlocks )
426 {
427   // Process an event
428   HLTDebug( "NofBlocks %lu", evtData.fBlockCnt );
429   // Process an event
430   AliHLTUInt32_t totalSize = 0, offset = 0;
431   Bool_t bWriteClusters = fReconstructor->IsWritingClusters();
432
433   //implement a usage of the following
434   //   AliHLTUInt32_t triggerDataStructSize = trigData.fStructSize;
435   //   AliHLTUInt32_t triggerDataSize = trigData.fDataSize;
436   //   void *triggerData = trigData.fData;
437   //HLTDebug( "Trigger data received. Struct size %d Data size %d Data location 0x%x", trigData.fStructSize, trigData.fDataSize, (UInt_t*)trigData.fData);
438
439   // Loop over all input blocks in the event
440   AliHLTComponentDataType expectedDataType = (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTRD);
441   for ( unsigned long i = 0; i < evtData.fBlockCnt; i++ )
442     {
443       if (evtData.fEventID == 1)
444         CALLGRIND_START_INSTRUMENTATION();
445       
446       const AliHLTComponentBlockData &block = blocks[i];
447       offset = totalSize;
448       // lets not use the internal TRD data types here : AliHLTTRDDefinitions::fgkDDLRawDataType
449       // which is depreciated - we use HLT global defs instead
450       //      if ( block.fDataType != (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTRD) )
451       AliHLTComponentDataType inputDataType = block.fDataType;
452       if ( inputDataType != expectedDataType)
453         {
454           HLTDebug( "Block # %i/%i; Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s; Skipping",
455                     i, evtData.fBlockCnt,
456                     evtData.fEventID, evtData.fEventID, 
457                     DataType2Text(inputDataType).c_str(), 
458                     DataType2Text(expectedDataType).c_str());
459           continue;
460         }
461       else 
462         {
463         HLTDebug("We get the right data type: Block # %i/%i; Event 0x%08LX (%Lu) Received datatype: %s",
464                     i, evtData.fBlockCnt,
465                     evtData.fEventID, evtData.fEventID, 
466                     DataType2Text(inputDataType).c_str());
467         }
468       
469       //      fMemReader->Reset();
470       fMemReader->SetMemory((UChar_t*) block.fPtr, block.fSize);
471
472       AliHLTUInt32_t spec = block.fSpecification;
473       
474       Int_t id = 1024;
475       
476       for ( Int_t ii = 0; ii < 18 ; ii++ ) {
477         if ( spec & 0x00000001 ) {
478           id += ii;
479           break;
480         }
481         spec = spec >> 1 ;
482       }
483
484       fMemReader->SetEquipmentID( id ); 
485       
486       fClusterizer->ResetTree();  
487       Bool_t iclustered = fClusterizer->Raw2ClustersChamber(fMemReader);
488       if (iclustered == kTRUE)
489         {
490           HLTDebug( "Clustered successfully");
491         }
492       else
493         {
494           HLTError("Clustering ERROR");
495           return -1;
496         }
497
498       // put the tree into output
499       //fcTree->Print();
500       if (bWriteClusters)
501         {
502           TTree *fcTree = fClusterizer->GetClusterTree();
503           if (fcTree)
504             {
505               Int_t nbEntries = fcTree->GetEntriesFast();
506               if (nbEntries > 0){
507                 HLTDebug("fcTree: Entries - %i; Size - %i",nbEntries,sizeof(fcTree));
508                 PushBack(fcTree, AliHLTTRDDefinitions::fgkClusterDataType, spec);
509               }
510               
511             }
512           fClusterizer->RecPoints()->Delete();
513           
514         }
515       else
516         {
517           TClonesArray *clustersArray = fClusterizer->RecPoints();
518           fClusterizer->SetClustersOwner(kFALSE);
519           if (clustersArray){
520             Int_t nbEntries = clustersArray->GetEntriesFast();
521             if (nbEntries > 0){
522               if (fUseHLTClusters){
523                 // Using low-level interface 
524                 // with interface classes
525                 AliHLTUInt32_t addedSize = AddToOutput(clustersArray, outputPtr);
526                 totalSize += addedSize;
527                 if ( totalSize > size )
528                   {
529                     HLTError("Too much data; Data written over allowed buffer. Amount written: %lu, allowed amount: %lu.",
530                              totalSize, size );
531                     return EMSGSIZE;
532                   }
533                 
534                 // Fill block 
535                 AliHLTComponentBlockData bd;
536                 FillBlockData( bd );
537                 bd.fOffset = offset;
538                 bd.fSize = addedSize;
539                 //bd.fSpecification = spec;
540                 bd.fSpecification = gkAliEventTypeData;
541                 bd.fDataType = AliHLTTRDDefinitions::fgkClusterDataType;
542                 outputBlocks.push_back( bd );
543                 HLTDebug( "Block ; size %i; dataType %s; spec 0x%x ",
544                           bd.fSize, DataType2Text(bd.fDataType).c_str(), spec);
545
546               }
547               else{
548                 // Using high-level interface
549                 // pass TClonesArray with streamers
550                 clustersArray->BypassStreamer(kFALSE);
551                 HLTDebug("clustersArray: Entries - %i; Size - %i",clustersArray->GetEntriesFast(),sizeof(clustersArray));
552                 PushBack(clustersArray, AliHLTTRDDefinitions::fgkClusterDataType, spec);
553                 //PrintObject(clustersArray);
554               }
555               
556               clustersArray->Delete();
557             }
558             delete clustersArray;
559           }
560           else 
561             HLTWarning("Array of clusters is empty!");
562         }
563     }//  for ( unsigned long i = 0; i < evtData.fBlockCnt; i++ )
564   fReconstructor->SetClusters(0x0);
565
566   size = totalSize;
567   HLTDebug("Event is done. size written to the output is %i", size);
568   return 0;
569 }
570
571
572 void AliHLTTRDClusterizerComponent::PrintObject( TClonesArray* inClustersArray)
573 {
574   AliTRDcluster* cluster=0x0;
575   
576   for (Int_t i=0; i < inClustersArray->GetEntriesFast(); i++){
577     cluster = dynamic_cast<AliTRDcluster*>(inClustersArray->At(i));
578     HLTDebug("cluster[%i]",i);
579     HLTDebug("  PadCol = %i; PadRow = %i; PadTime = %i", cluster->GetPadCol(), cluster->GetPadRow(), cluster->GetPadTime());
580     HLTDebug("  Detector = %i, Amplitude = %f, Center = %f", cluster->GetDetector(), cluster->GetQ(), cluster->GetCenter());
581     HLTDebug("  LocalTimeBin =  %i; NPads = %i; maskedPosition: %s, status: %s", cluster->GetLocalTimeBin(), cluster->GetNPads(),cluster->GetPadMaskedPosition(),cluster->GetPadMaskedPosition());
582   }
583   
584 }
585