]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/TPCLib/AliHLTTPCClusterFinderComponent.cxx
Adding a cut on ADC to speed up things (Dieggo)
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCClusterFinderComponent.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  *                  Timm Steinbeck <timm@kip.uni-heidelberg.de>           *
9  *                  Jochen Thaeder <thaeder@kip.uni-heidelberg.de>        *
10  *                  for The ALICE HLT Project.                            *
11  *                                                                        *
12  * Permission to use, copy, modify and distribute this software and its   *
13  * documentation strictly for non-commercial purposes is hereby granted   *
14  * without fee, provided that the above copyright notice appears in all   *
15  * copies and that both the copyright notice and this permission notice   *
16  * appear in the supporting documentation. The authors make no claims     *
17  * about the suitability of this software for any purpose. It is          *
18  * provided "as is" without express or implied warranty.                  *
19  **************************************************************************/
20
21 /** @file   AliHLTTPCClusterFinderComponent.cxx
22     @author Timm Steinbeck, Matthias Richter, Jochen Thaeder, Kenneth Aamodt
23     @date   
24     @brief  The TPC cluster finder processing component
25 */
26
27 // see header file for class documentation                                   //
28 // or                                                                        //
29 // refer to README to build package                                          //
30 // or                                                                        //
31 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt                          //
32
33 #if __GNUC__>= 3
34 using namespace std;
35 #endif
36 #include "AliHLTTPCClusterFinderComponent.h"
37 #include "AliHLTTPCDigitReaderPacked.h"
38 #include "AliHLTTPCDigitReaderUnpacked.h"
39 #include "AliHLTTPCDigitReaderRaw.h"
40 #include "AliHLTTPCDigitReaderDecoder.h"
41 #include "AliHLTTPCClusterFinder.h"
42 #include "AliHLTTPCSpacePointData.h"
43 #include "AliHLTTPCClusterDataFormat.h"
44 #include "AliHLTTPCTransform.h"
45 #include "AliHLTTPCClusters.h"
46 #include "AliHLTTPCDefinitions.h"
47 #include "AliCDBEntry.h"
48 #include "AliCDBManager.h"
49
50 #include <cstdlib>
51 #include <cerrno>
52 #include "TString.h"
53 #include "TObjString.h"
54 #include <sys/time.h>
55
56 // this is a global object used for automatic component registration, do not use this
57 // use fPackedSwitch = true for packed inputtype "gkDDLPackedRawDataType"
58 // use fPackedSwitch = false for unpacked inputtype "gkUnpackedRawDataType"
59 AliHLTTPCClusterFinderComponent gAliHLTTPCClusterFinderComponentPacked(AliHLTTPCClusterFinderComponent::kClusterFinderPacked);
60 AliHLTTPCClusterFinderComponent gAliHLTTPCClusterFinderComponentUnpacked(AliHLTTPCClusterFinderComponent::kClusterFinderUnpacked);
61 AliHLTTPCClusterFinderComponent gAliHLTTPCClusterFinderComponentDecoder(AliHLTTPCClusterFinderComponent::kClusterFinderDecoder);
62
63 /** ROOT macro for the implementation of ROOT specific class methods */
64 ClassImp(AliHLTTPCClusterFinderComponent)
65
66 AliHLTTPCClusterFinderComponent::AliHLTTPCClusterFinderComponent(int mode)
67   :
68   fClusterFinder(NULL),
69   fReader(NULL),
70   fClusterDeconv(true),
71   fXYClusterError(-1),
72   fZClusterError(-1),
73   fModeSwitch(mode),
74   fUnsorted(0),
75   fPatch(0),
76   fPadArray(NULL),
77   fGetActivePads(0)
78 {
79   // see header file for class documentation
80   // or
81   // refer to README to build package
82   // or
83   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
84 }
85
86 AliHLTTPCClusterFinderComponent::~AliHLTTPCClusterFinderComponent()
87 {
88   // see header file for class documentation
89 }
90
91 // Public functions to implement AliHLTComponent's interface.
92 // These functions are required for the registration process
93
94 const char* AliHLTTPCClusterFinderComponent::GetComponentID()
95 {
96   // see header file for class documentation
97   switch(fModeSwitch){
98   case kClusterFinderPacked:
99     return "TPCClusterFinderPacked";
100     break;
101   case kClusterFinderUnpacked:
102     return "TPCClusterFinderUnpacked";
103     break;
104   case kClusterFinderDecoder:
105     return "TPCClusterFinderDecoder";
106     break;
107   }
108   HLTFatal("unknown digit reader type");
109   return "";
110 }
111
112 void AliHLTTPCClusterFinderComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
113 {
114   // see header file for class documentation
115   list.clear(); 
116   switch(fModeSwitch){
117   case kClusterFinderPacked:
118     list.push_back( kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC );
119     break;
120   case kClusterFinderUnpacked:
121     list.push_back( AliHLTTPCDefinitions::fgkUnpackedRawDataType );
122     break;
123   case kClusterFinderDecoder:
124     list.push_back( kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC );
125     break;
126   }
127 }
128
129 AliHLTComponentDataType AliHLTTPCClusterFinderComponent::GetOutputDataType()
130 {
131   // see header file for class documentation
132   return kAliHLTMultipleDataType;
133 }
134
135 int AliHLTTPCClusterFinderComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
136
137 {
138   // see header file for class documentation
139   tgtList.clear();
140   tgtList.push_back(AliHLTTPCDefinitions::fgkClustersDataType);
141   tgtList.push_back(AliHLTTPCDefinitions::fgkActivePadsDataType);
142   return tgtList.size();
143 }
144
145 void AliHLTTPCClusterFinderComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
146 {
147   // see header file for class documentation
148   // XXX TODO: Find more realistic values.  
149   constBase = 0;
150   switch(fModeSwitch){
151   case 0:
152     inputMultiplier = (6 * 0.4);
153     break;
154   case 1:
155     inputMultiplier = 0.4;
156     break;
157   case 2:
158     inputMultiplier = (6 * 0.4);
159     break;
160   }
161 }
162
163 AliHLTComponent* AliHLTTPCClusterFinderComponent::Spawn()
164 {
165   // see header file for class documentation
166   return new AliHLTTPCClusterFinderComponent(fModeSwitch);
167 }
168         
169 int AliHLTTPCClusterFinderComponent::DoInit( int argc, const char** argv )
170 {
171   // see header file for class documentation
172   if ( fClusterFinder )
173     return EINPROGRESS;
174
175   fClusterFinder = new AliHLTTPCClusterFinder();
176
177   Int_t rawreadermode =  -1;
178   Int_t sigthresh = -1;
179   Double_t sigmathresh= -1;
180   Float_t occulimit = 1.0;
181   Int_t oldRCUFormat=0;
182   // Data Format version numbers:
183   // 0: RCU Data format as delivered during TPC commissioning, pads/padrows are sorted, RCU trailer is one 32 bit word.
184   // 1: As 0, but pads/padrows are delivered "as is", without sorting
185   // 2: As 0, but RCU trailer is 3 32 bit words.
186   // 3: As 1, but RCU trailer is 3 32 bit words.
187   // -1: use offline raw reader
188
189   Int_t i = 0;
190   Char_t* cpErr;
191
192   while ( i < argc ) {      
193
194     // -- raw reader mode option
195     if ( !strcmp( argv[i], "rawreadermode" ) ) {
196       if ( argc <= i+1 ) {
197         Logging( kHLTLogError, "HLT::TPCClusterFinder::DoInit", "Missing rawreadermode", "Raw Reader Mode not specified" );
198         return ENOTSUP;
199       }
200
201       // Decodes the rawreader mode: either number or string and returns the rawreadermode
202       // -1 on failure, -2 for offline
203       rawreadermode = AliHLTTPCDigitReaderRaw::DecodeMode( argv[i+1] );
204
205       if (rawreadermode == -1 ) {
206         Logging( kHLTLogError, "HLT::TPCClusterFinder::DoInit", "Missing rawreadermode", "Cannot convert rawreadermode specifier '%s'.", argv[i+1] );
207         return EINVAL;
208       }
209
210       i += 2;
211       continue;
212     }
213
214     // -- pp run option
215     if ( !strcmp( argv[i], "pp-run" ) ) {
216       fClusterDeconv = false;
217       i++;
218       continue;
219     }
220
221     // -- zero suppression threshold
222     if ( !strcmp( argv[i], "adc-threshold" ) ) {
223       sigthresh = strtoul( argv[i+1], &cpErr ,0);
224       if ( *cpErr ) {
225         HLTError("Cannot convert threshold specifier '%s'.", argv[i+1]);
226         return EINVAL;
227       }
228       i+=2;
229       continue;
230     }
231
232     // -- pad occupancy limit
233     if ( !strcmp( argv[i], "occupancy-limit" ) ) {
234       occulimit = strtod( argv[i+1], &cpErr);
235       if ( *cpErr ) {
236         HLTError("Cannot convert occupancy specifier '%s'.", argv[i+1]);
237         return EINVAL;
238       }
239       i+=2;
240       continue;
241     }
242
243     // -- number of timebins (default 1024)
244     if ( !strcmp( argv[i], "timebins" ) ) {
245       TString parameter(argv[i+1]);
246       parameter.Remove(TString::kLeading, ' '); // remove all blanks
247       if (parameter.IsDigit()) {
248         AliHLTTPCTransform::SetNTimeBins(parameter.Atoi());
249         HLTInfo("number of timebins set to %d, zbin=%f", AliHLTTPCTransform::GetNTimeBins(), AliHLTTPCTransform::GetZWidth());
250       } else {
251         HLTError("Cannot timebin specifier '%s'.", argv[i+1]);
252         return EINVAL;
253       }
254       i+=2;
255       continue;
256     }
257
258     // -- checking for rcu format
259     if ( !strcmp( argv[i], "oldrcuformat" ) ) {
260       oldRCUFormat = strtoul( argv[i+1], &cpErr ,0);
261       if ( *cpErr ){
262         HLTError("Cannot convert oldrcuformat specifier '%s'. Should  be 0(off) or 1(on), must be integer", argv[i+1]);
263         return EINVAL;
264       }
265       i+=2;
266       continue;
267     }
268       
269     // -- checking for unsorted clusterfinding
270     if ( !strcmp( argv[i], "unsorted" ) ) {
271       fUnsorted = strtoul( argv[i+1], &cpErr ,0);
272       if ( *cpErr ){
273         HLTError("Cannot convert unsorted specifier '%s'. Should  be 0(off) or 1(on), must be integer", argv[i+1]);
274         return EINVAL;
275       }
276       i+=2;
277       continue;
278     }
279       
280     // -- checking for active pads, used in 2007 December run
281     if ( !strcmp( argv[i], "activepads" ) ) {
282       fGetActivePads = strtoul( argv[i+1], &cpErr ,0);
283       if ( *cpErr ){
284         HLTError("Cannot convert activepads specifier '%s'. Should  be 0(off) or 1(on), must be integer", argv[i+1]);
285         return EINVAL;
286       }
287       i+=2;
288       continue;
289     }
290
291     // -- checking for nsigma-threshold, used in 2007 December run in ZeroSuppression
292     if ( !strcmp( argv[i], "nsigma-threshold" ) ) {
293        sigmathresh = strtoul( argv[i+1], &cpErr ,0);
294       if ( *cpErr ){
295         HLTError("Cannot convert nsigma-threshold specifier '%s'. Must be integer", argv[i+1]);
296         return EINVAL;
297       }
298       i+=2;
299       continue;
300     }
301
302     Logging(kHLTLogError, "HLT::TPCClusterFinder::DoInit", "Unknown Option", "Unknown option '%s'", argv[i] );
303     return EINVAL;
304
305   }
306
307   // Choose reader
308   if (fModeSwitch==kClusterFinderPacked) { 
309     if (rawreadermode == -2) {
310       HLTDebug("using AliHLTTPCDigitReaderPacked");
311       fReader = new AliHLTTPCDigitReaderPacked();
312       if(oldRCUFormat==1){
313         fReader->SetOldRCUFormat(kTRUE);
314       }
315       else if(oldRCUFormat!=0){
316         HLTWarning("Wrong oldrcuformat specifier %d; oldrcuformat set to default(kFALSE)",oldRCUFormat);
317       }
318       if(fUnsorted==1){
319         fReader->SetUnsorted(kTRUE);
320       }
321       fClusterFinder->SetReader(fReader);
322     } 
323     else {
324 #if defined(HAVE_TPC_MAPPING)
325       HLTDebug("using AliHLTTPCDigitReaderRaw mode %d", rawreadermode);
326       fReader = new AliHLTTPCDigitReaderRaw(rawreadermode);
327       fClusterFinder->SetReader(fReader);
328 #else //! defined(HAVE_TPC_MAPPING)
329       HLTFatal("DigitReaderRaw not available - check your build");
330       return -ENODEV;
331 #endif //defined(HAVE_TPC_MAPPING)
332     }
333   }
334   else if(fModeSwitch==kClusterFinderUnpacked){
335     HLTDebug("using AliHLTTPCDigitReaderUnpacked");
336     fReader = new AliHLTTPCDigitReaderUnpacked();
337     fClusterFinder->SetReader(fReader);
338   }
339   else if(fModeSwitch==kClusterFinderDecoder){
340     fReader = new AliHLTTPCDigitReaderDecoder();
341     fClusterFinder->SetReader(fReader);
342   }
343   else{
344     HLTFatal("No mode set for clusterfindercomponent");
345   }
346   // if pp-run use occupancy limit else set to 1. ==> use all 
347   if ( !fClusterDeconv )
348     fClusterFinder->SetOccupancyLimit(occulimit);
349   else 
350     fClusterFinder->SetOccupancyLimit(1.0);
351       
352   // Variables to setup the Clusterfinder
353   // TODO: this sounds strange and has to be verified; is the cluster finder not working when
354   // fClusterDeconv = false ?
355   fClusterDeconv = true;
356   fXYClusterError = -1;
357   fZClusterError = -1;
358
359  
360   fClusterFinder->SetDeconv( fClusterDeconv );
361   fClusterFinder->SetXYError( fXYClusterError );
362   fClusterFinder->SetZError( fZClusterError );
363   if ( (fXYClusterError>0) && (fZClusterError>0) )
364     fClusterFinder->SetCalcErr( false );
365   fClusterFinder->SetSignalThreshold(sigthresh);
366   fClusterFinder->SetNSigmaThreshold(sigmathresh);
367
368   return 0;
369 }
370
371 int AliHLTTPCClusterFinderComponent::DoDeinit()
372 {
373   // see header file for class documentation
374
375   if ( fClusterFinder )
376     delete fClusterFinder;
377   fClusterFinder = NULL;
378  
379   if ( fReader )
380     delete fReader;
381   fReader = NULL;
382     
383   return 0;
384 }
385
386 int AliHLTTPCClusterFinderComponent::DoEvent( const AliHLTComponentEventData& evtData, 
387                                               const AliHLTComponentBlockData* blocks, 
388                                               AliHLTComponentTriggerData& /*trigData*/, AliHLTUInt8_t* outputPtr, 
389                                               AliHLTUInt32_t& size, 
390                                               vector<AliHLTComponentBlockData>& outputBlocks )
391 {
392   // see header file for class documentation
393
394   //  == init iter (pointer to datablock)
395   const AliHLTComponentBlockData* iter = NULL;
396   unsigned long ndx;
397
398   //  == OUTdatatype pointer
399   AliHLTTPCClusterData* outPtr;
400
401   AliHLTUInt8_t* outBPtr;
402   UInt_t offset, mysize, nSize, tSize = 0;
403
404   outBPtr = outputPtr;
405   outPtr = (AliHLTTPCClusterData*)outBPtr;
406
407   Int_t slice, patch, row[2];
408   unsigned long maxPoints, realPoints = 0;
409
410   for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ )
411     {
412       iter = blocks+ndx;
413       mysize = 0;
414       offset = tSize;
415
416
417       if (fModeSwitch==0 || fModeSwitch==2) {
418         HLTDebug("Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s",
419                  evtData.fEventID, evtData.fEventID, 
420                  DataType2Text( iter->fDataType).c_str(), 
421                  DataType2Text(kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC).c_str());
422
423         if (iter->fDataType == AliHLTTPCDefinitions::fgkDDLPackedRawDataType &&
424             GetEventCount()<2) {
425           HLTWarning("data type %s is depricated, use %s (kAliHLTDataTypeDDLRaw)!",
426                      DataType2Text(AliHLTTPCDefinitions::fgkDDLPackedRawDataType).c_str(),
427                      DataType2Text(kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC).c_str());
428           }
429
430         if ( iter->fDataType != (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC) &&
431              iter->fDataType != AliHLTTPCDefinitions::fgkDDLPackedRawDataType ) continue;
432
433       }
434       else if(fModeSwitch==1){
435         HLTDebug("Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s",
436                  evtData.fEventID, evtData.fEventID, 
437                  DataType2Text( iter->fDataType).c_str(), 
438                  DataType2Text(AliHLTTPCDefinitions::fgkUnpackedRawDataType).c_str());
439
440         if ( iter->fDataType != AliHLTTPCDefinitions::fgkUnpackedRawDataType ) continue;
441
442       }
443         
444       slice = AliHLTTPCDefinitions::GetMinSliceNr( *iter );
445       patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
446       row[0] = AliHLTTPCTransform::GetFirstRow( patch );
447       row[1] = AliHLTTPCTransform::GetLastRow( patch );
448
449
450       if(fUnsorted){
451         fClusterFinder->SetUnsorted(fUnsorted);
452         fClusterFinder->SetPatch(patch);
453       }
454       /*      if(fUnsorted){
455               if(fPadArray==NULL){
456               fClusterFinder->SetUnsorted(fUnsorted);
457               //fPadArray = new AliHLTTPCPadArray(patch);
458               //fPadArray->InitializeVector(fModeSwitch);
459               }
460               else if(fPadArray->GetPatch()!=patch||fPadArray->GetPatch()==-1){
461               if (GetEventCount()<3) {
462               HLTWarning("pad array not initialized for data of specification 0x%08x, block skipped", iter->fSpecification);
463               } else if ((GetEventCount()%5000)==0) { // assuming 0.5 to 1kHz this gives a message rate of 0.1 to 0.5 Hz
464               HLTWarning("reminder: pad array not initialized for data of specification 0x%08x", iter->fSpecification);
465               }
466               continue;
467               }
468               }
469       */
470       
471
472       outPtr = (AliHLTTPCClusterData*)outBPtr;
473
474       maxPoints = (size-tSize-sizeof(AliHLTTPCClusterData))/sizeof(AliHLTTPCSpacePointData);
475
476       fClusterFinder->InitSlice( slice, patch, row[0], row[1], maxPoints );
477       fClusterFinder->SetOutputArray( (AliHLTTPCSpacePointData*)outPtr->fSpacePoints );
478         
479       if(fUnsorted){
480         fClusterFinder->ReadDataUnsorted(iter->fPtr, iter->fSize, fModeSwitch);
481
482         fClusterFinder->FindClusters();
483       }
484       else{
485         fClusterFinder->Read(iter->fPtr, iter->fSize );
486         fClusterFinder->ProcessDigits();
487       }
488       realPoints = fClusterFinder->GetNumberOfClusters();
489         
490       outPtr->fSpacePointCnt = realPoints;
491       nSize = sizeof(AliHLTTPCSpacePointData)*realPoints;
492       mysize += nSize+sizeof(AliHLTTPCClusterData);
493
494       Logging( kHLTLogDebug, "HLT::TPCClusterFinder::DoEvent", "Spacepoints", 
495                "Number of spacepoints: %lu Slice/Patch/RowMin/RowMax: %d/%d/%d/%d.",
496                realPoints, slice, patch, row[0], row[1] );
497       AliHLTComponentBlockData bd;
498       FillBlockData( bd );
499       bd.fOffset = offset;
500       bd.fSize = mysize;
501       bd.fSpecification = iter->fSpecification;
502       bd.fDataType = AliHLTTPCDefinitions::fgkClustersDataType;
503       //AliHLTSubEventDescriptor::FillBlockAttributes( bd.fAttributes );
504       outputBlocks.push_back( bd );
505         
506       tSize += mysize;
507       outBPtr += mysize;
508       outPtr = (AliHLTTPCClusterData*)outBPtr;
509         
510       /*      if(fGetActivePads){
511               AliHLTTPCPadArray::AliHLTTPCActivePads* outPtrActive;
512               UInt_t activePadsSize, activePadsN = 0;
513               outPtrActive = (AliHLTTPCPadArray::AliHLTTPCActivePads*)outBPtr;
514               offset=tSize;
515               Int_t maxActivePads = (size-tSize)/sizeof(AliHLTTPCPadArray::AliHLTTPCActivePads);
516               activePadsSize= fClusterFinder->GetActivePads((AliHLTTPCPadArray::AliHLTTPCActivePads*)outPtrActive,maxActivePads)*sizeof(AliHLTTPCPadArray::AliHLTTPCActivePads);
517         
518               AliHLTComponentBlockData bdActive;
519               FillBlockData( bdActive );
520               bdActive.fOffset = offset;
521               bdActive.fSize = activePadsSize;
522               bdActive.fSpecification = iter->fSpecification;
523               bdActive.fDataType = AliHLTTPCDefinitions::fgkActivePadsDataType;
524               outputBlocks.push_back( bdActive );
525         
526               tSize+=activePadsSize;
527               outBPtr += activePadsSize;
528               outPtrActive = (AliHLTTPCPadArray::AliHLTTPCActivePads*)outBPtr;
529               }
530       */
531  
532
533       if ( tSize > size )
534         {
535           Logging( kHLTLogFatal, "HLT::TPCClusterFinder::DoEvent", "Too much data", 
536                    "Data written over allowed buffer. Amount written: %lu, allowed amount: %lu.",
537                    tSize, size );
538           return EMSGSIZE;
539         }
540     }
541     
542   size = tSize;
543
544   return 0;
545 }
546
547 int AliHLTTPCClusterFinderComponent::Reconfigure(const char* cdbEntry, const char* chainId)
548 {
549   // see header file for class documentation
550   const char* path="HLT/ConfigTPC";
551   if (cdbEntry) path=cdbEntry;
552   if (path) {
553     HLTInfo("reconfigure from entry %s, chain id %s", path, (chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
554     AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path/*,GetRunNo()*/);
555     if (pEntry) {
556       TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
557       if (pString) {
558         HLTInfo("received configuration object: %s", pString->GetString().Data());
559       } else {
560         HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
561       }
562     } else {
563       HLTError("can not fetch object \"%s\" from CDB", path);
564     }
565   }
566   return 0;
567 }