ALIROOT-5433 Transition to CDHv3 in HLT
[u/mrichter/AliRoot.git] / HLT / TPCLib / HWCFemulator / AliHLTTPCHWCFConsistencyControlComponent.cxx
1 //****************************************************************************
2 //* This file is property of and copyright by the ALICE HLT Project          * 
3 //* ALICE Experiment at CERN, All rights reserved.                           *
4 //*                                                                          *
5 //* Primary Authors: Sergey Gorbunov, Torsten Alt                            *
6 //* Developers:      Sergey Gorbunov <sergey.gorbunov@fias.uni-frankfurt.de> *
7 //*                  Torsten Alt <talt@cern.ch>                              *
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   AliHLTTPCHWCFConsistencyControlComponent.cxx
20 //  @author Sergey Gorbunov <sergey.gorbunov@fias.uni-frankfurt.de>
21 //  @author Torsten Alt <talt@cern.ch> 
22 //  @brief  Comparison of TPC clusters produced by FPGA clusterfinder and by FPGA Emulator 
23 //  @brief  ( see AliHLTTPCHWCFEmulator class )
24 //  @note
25
26
27 #include "AliHLTTPCHWCFConsistencyControlComponent.h"
28
29 #include "AliHLTTPCDefinitions.h"
30 #include "AliHLTTPCHWCFDataTypes.h"
31 #include "AliHLTTPCClusterMCData.h"
32
33 #include "AliGRPObject.h"
34 #include "AliCDBEntry.h"
35 #include "AliCDBManager.h"
36 #include "AliHLTCDHWrapper.h"
37 #include <cstdlib>
38 #include <cerrno>
39 #include "TString.h"
40 #include "TObjString.h"
41 #include "TObjArray.h"
42 #include "TH1F.h"
43
44 #include <sys/time.h>
45 #include "TFile.h"
46
47 using namespace std;
48
49 AliHLTTPCHWCFConsistencyControlComponent::AliHLTTPCHWCFConsistencyControlComponent()
50   :
51   AliHLTProcessor(),
52   fNDismatch(0),
53   fNBlocks(0),
54   fBenchmark("TPCHWConsistencyControl"), 
55   fHistHeaderAll(0),
56   fHistHeaderGood(0),
57   fHistClusterAll(0),
58   fHistClusterGood(0),
59   fProfHeader(0),
60   fProfCluster(0)
61 {
62   // see header file for class documentation
63   // or
64   // refer to README to build package
65   // or
66   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
67 }
68
69
70 AliHLTTPCHWCFConsistencyControlComponent::AliHLTTPCHWCFConsistencyControlComponent(const AliHLTTPCHWCFConsistencyControlComponent&)
71   :
72   AliHLTProcessor(),
73   fNDismatch(0),
74   fNBlocks(0),
75   fBenchmark("TPCHWConsistencyControl"),
76   fHistHeaderAll(0),
77   fHistHeaderGood(0),
78   fHistClusterAll(0),
79   fHistClusterGood(0),
80   fProfHeader(0),
81   fProfCluster(0)
82 {
83   // dummy
84 }
85
86 AliHLTTPCHWCFConsistencyControlComponent& AliHLTTPCHWCFConsistencyControlComponent::operator=(const AliHLTTPCHWCFConsistencyControlComponent&)
87 {
88   // dummy
89   return *this;
90 }
91
92 AliHLTTPCHWCFConsistencyControlComponent::~AliHLTTPCHWCFConsistencyControlComponent()
93 {
94   // see header file for class documentation
95   DoDeinit();
96 }
97
98 // Public functions to implement AliHLTComponent's interface.
99 // These functions are required for the registration process
100
101 const char* AliHLTTPCHWCFConsistencyControlComponent::GetComponentID()
102 {
103   // see header file for class documentation
104   return "TPCHWCFConsistenyControl";
105 }
106
107 void AliHLTTPCHWCFConsistencyControlComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
108 {
109   // see header file for class documentation
110   list.clear();
111   list.push_back( AliHLTTPCDefinitions::fgkHWClustersDataType| kAliHLTDataOriginTPC );
112 }
113
114 AliHLTComponentDataType AliHLTTPCHWCFConsistencyControlComponent::GetOutputDataType()
115 {
116   // see header file for class documentation
117   return kAliHLTDataTypeHistogram  | kAliHLTDataOriginTPC;
118 }
119
120
121 void AliHLTTPCHWCFConsistencyControlComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
122 {
123   // see header file for class documentation
124   // XXX TODO: Find more realistic values.  
125   constBase = 1000;
126   inputMultiplier = 0;
127 }
128
129
130 AliHLTComponent* AliHLTTPCHWCFConsistencyControlComponent::Spawn()
131 {
132   // see header file for class documentation
133   return new AliHLTTPCHWCFConsistencyControlComponent();
134 }
135
136 void AliHLTTPCHWCFConsistencyControlComponent::GetOCDBObjectDescription( TMap* const /*targetMap*/){
137 // Get a list of OCDB object description needed for the particular component
138   /*
139   if (!targetMap) return;
140   
141   // OCDB entries for component arguments
142   targetMap->Add(new TObjString("HLT/ConfigTPC/TPCHWClusterFinder"), new TObjString("component arguments, empty at the moment"));
143   */
144 }
145
146
147 int AliHLTTPCHWCFConsistencyControlComponent::DoInit( int argc, const char** argv )
148 {
149   // see header file for class documentation
150
151   fHistHeaderAll = new TH1F("hHWCFHeaderAll", "fHistHeaderAll",6,0.,6.);
152   fHistHeaderGood = new TH1F("hHWCFHeaderGood", "fHistHeaderGood",6,0.,6.);
153   fHistClusterAll = new TH1F("hHWCFClusterAll", "fHistClusterAll",8,0.,8.);
154   fHistClusterGood = new TH1F("hHWCFClusterGood", "fHistClusterGood",8,0.,8.);
155   fProfHeader = new TH1F("pHWCFHeader", "HWCF: Consistency of header data", 6, 0., 6.);
156   fProfCluster = new TH1F("pHWCFClusters", "HWCF: Consisteny of cluster data", 8, 0., 8.);
157
158   TString arguments = "";
159   for ( int i = 0; i < argc; i++ ) {
160     if ( !arguments.IsNull() ) arguments += " ";
161     arguments += argv[i];
162   }
163
164   return Configure( NULL, NULL, arguments.Data()  );
165 }
166
167 int AliHLTTPCHWCFConsistencyControlComponent::Reconfigure( const char* cdbEntry, const char* chainId )
168 {
169   // Reconfigure the component from OCDB
170
171   return Configure( cdbEntry, chainId, NULL );
172 }
173
174 int AliHLTTPCHWCFConsistencyControlComponent::ScanConfigurationArgument(int argc, const char** argv)
175 {
176   // see header file for class documentation
177   TString arguments = "";
178   for ( int i = 0; i < argc; i++ ) {
179     if ( !arguments.IsNull() ) arguments += " ";
180     arguments += argv[i];
181   }
182   return ReadConfigurationString(arguments);
183 }
184
185
186
187 void AliHLTTPCHWCFConsistencyControlComponent::SetDefaultConfiguration()
188 {
189   // Set default configuration for the FPGA ClusterFinder Emulator component
190   // Some parameters can be later overwritten from the OCDB
191
192   fNDismatch = 0;
193   fNBlocks = 0;
194   fBenchmark.Reset();
195   fBenchmark.SetTimer(0,"total");
196   fBenchmark.SetTimer(1,"reco");    
197 }
198
199 int AliHLTTPCHWCFConsistencyControlComponent::ReadConfigurationString(  const char* arguments )
200 {
201   // Set configuration parameters for the FPGA ClusterFinder Emulator component
202   // from the string
203
204   int iResult = 0;
205   if ( !arguments ) return iResult;
206
207   TString allArgs = arguments;
208   TString argument;
209   int bMissingParam = 0;
210
211   TObjArray* pTokens = allArgs.Tokenize( " " );
212
213   int nArgs =  pTokens ? pTokens->GetEntries() : 0;
214
215   for ( int i = 0; i < nArgs; i++ ) {
216     argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
217     if ( argument.IsNull() ) continue;
218
219     //if ( argument.CompareTo( "-deconvolute-time" ) == 0 ) {
220     //if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
221     //fDoDeconvTime  = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
222     //HLTInfo( "Time deconvolution is set to: %d", fDoDeconvTime );
223     //continue;
224     //}
225
226     HLTError( "Unknown option \"%s\"", argument.Data() );
227     iResult = -EINVAL;
228   }
229   delete pTokens;
230
231   if ( bMissingParam ) {
232     HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
233     iResult = -EINVAL;
234   }
235
236   return iResult;
237 }
238
239
240 int AliHLTTPCHWCFConsistencyControlComponent::ReadCDBEntry( const char* /*cdbEntry*/, const char* /*chainId*/ )
241 {
242   // Read configuration from OCDB
243   return 0;
244   /*
245   const char* defaultNotify = "";
246
247   if ( !cdbEntry ){
248     cdbEntry = "HLT/ConfigTPC/TPCHWClusterFinder";
249     defaultNotify = " (default)";
250     chainId = 0;
251   }
252
253   HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
254   AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
255   
256   if ( !pEntry ) {
257     HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
258     return -EINVAL;
259   }
260
261   TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
262
263   if ( !pString ) {
264     HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
265     return -EINVAL;
266   }
267
268   HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
269
270   return  ReadConfigurationString( pString->GetString().Data() );
271   */
272 }
273
274
275 int AliHLTTPCHWCFConsistencyControlComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
276 {
277   // Configure the component
278   // There are few levels of configuration,
279   // parameters which are set on one step can be overwritten on the next step
280
281   //* read hard-coded values
282
283   SetDefaultConfiguration();
284
285   //* read the default CDB entry
286
287   int iResult1 = ReadCDBEntry( NULL, chainId );
288
289   //* read the actual CDB entry if required
290
291   int iResult2 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
292
293   //* read extra parameters from input (if they are)
294
295   int iResult3 = 0;
296
297   if ( commandLine && commandLine[0] != '\0' ) {
298     HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
299     iResult3 = ReadConfigurationString( commandLine );
300   }
301
302   return iResult1 ? iResult1 : ( iResult2 ? iResult2 : iResult3 );
303 }
304
305
306 int AliHLTTPCHWCFConsistencyControlComponent::DoDeinit()
307 {
308   // see header file for class documentation 
309   if( fHistHeaderAll ) delete fHistHeaderAll;
310   if( fHistHeaderGood ) delete fHistHeaderGood;
311   if( fHistClusterAll ) delete fHistClusterAll;
312   if( fHistClusterGood ) delete fHistClusterGood;
313   if( fProfHeader ) delete fProfHeader;
314   if( fProfCluster ) delete fProfCluster;
315
316   fHistHeaderAll = 0;
317   fHistHeaderGood = 0;
318   fHistClusterAll = 0;
319   fHistClusterGood = 0;
320   fProfHeader = 0;
321   fProfCluster = 0;
322   return 0;
323 }
324
325
326 int AliHLTTPCHWCFConsistencyControlComponent::DoEvent( const AliHLTComponentEventData& evtData, 
327                                                        const AliHLTComponentBlockData* blocks, 
328                                                        AliHLTComponentTriggerData& /*trigData*/, AliHLTUInt8_t* /*outputPtr*/, 
329                                                        AliHLTUInt32_t& size, 
330                                                        vector<AliHLTComponentBlockData>& /*outputBlocks*/ )
331 {
332   // see header file for class documentation
333
334   int iResult=0;
335   //AliHLTUInt32_t maxSize = size;
336   size = 0;
337   
338   if(!IsDataEvent()){
339     return 0;
340   }
341
342   fBenchmark.StartNewEvent();
343   fBenchmark.Start(0);
344   bool *checkedBlocks = new bool[evtData.fBlockCnt];
345   for( unsigned long i=0; i<evtData.fBlockCnt; i++ ) checkedBlocks[i] = 0;
346   
347   for ( unsigned long ndx1 = 0; ndx1 < evtData.fBlockCnt; ndx1++ ){
348     const AliHLTComponentBlockData* iter1 = blocks+ndx1;
349     fBenchmark.AddInput(iter1->fSize);
350     
351     if (  iter1->fDataType != (AliHLTTPCDefinitions::fgkHWClustersDataType | kAliHLTDataOriginTPC) 
352           ) continue;
353  
354     if( checkedBlocks[ndx1] ) continue;// block already checked
355
356     checkedBlocks[ndx1] = 1;
357
358     int slice1 = AliHLTTPCDefinitions::GetMinSliceNr( *iter1 );
359     int patch1 = AliHLTTPCDefinitions::GetMinPatchNr( *iter1 );
360
361     AliHLTCDHWrapper cdh1(iter1->fPtr);
362    
363     bool header1OK = 1;
364     if( iter1->fSize>0 ){
365       header1OK = ( iter1->fSize % sizeof(AliHLTUInt32_t) == 0 ) 
366         && (iter1->fPtr != NULL )
367         && ( iter1->fSize >= cdh1.GetHeaderSize() ) ;
368     }
369     
370     int  nMatchedBlocks=0;
371         
372     bool sameHLTHd = 1;     
373     bool sameCDH = 1;     
374     bool sameRCU = 1;      
375     bool sameSize = 1;
376    
377     bool sameFlag = 1;
378     bool sameQMax = 1;
379     bool sameCharge = 1;
380     bool sameRow = 1;
381     bool sameFloat[4] = {1,1,1,1};
382     
383     for ( unsigned long ndx2 = ndx1+1; ndx2 < evtData.fBlockCnt; ndx2++ ){      
384       const AliHLTComponentBlockData* iter2 = blocks+ndx2;
385       if (  iter2->fDataType != (AliHLTTPCDefinitions::fgkHWClustersDataType | kAliHLTDataOriginTPC)
386             ) continue;
387             
388       int slice2 = AliHLTTPCDefinitions::GetMinSliceNr( *iter2 );
389       int patch2 = AliHLTTPCDefinitions::GetMinPatchNr( *iter2 );
390       if( slice1!=slice2 || patch1!=patch2 ) continue;
391       
392       if( checkedBlocks[ndx2] ) continue;
393       checkedBlocks[ndx2] = 1;
394
395       nMatchedBlocks++;
396       
397       AliHLTCDHWrapper cdh2(iter1->fPtr);
398
399       bool header2OK = 1;
400       if( iter2->fSize>0 ){
401         header2OK = ( iter2->fSize % sizeof(AliHLTUInt32_t) == 0 ) 
402           && (iter2->fPtr != NULL )
403           && ( iter2->fSize >= cdh2.GetHeaderSize() ) ;
404       }
405  
406       fHistHeaderAll->Fill(2);
407       sameHLTHd = ( header1OK == header2OK );
408       if( sameHLTHd ) fHistHeaderGood->Fill(2);
409       if( !header1OK || !header2OK ) continue;
410       
411       int nWordsHeader =  cdh2.GetHeaderSize()/sizeof(AliHLTUInt32_t);
412       int nWords1 = iter1->fSize/sizeof(AliHLTUInt32_t);
413       int nWords2 = iter2->fSize/sizeof(AliHLTUInt32_t);
414       const AliHLTUInt32_t *p1 = (const AliHLTUInt32_t *) iter1->fPtr;
415       const AliHLTUInt32_t *p2 = (const AliHLTUInt32_t *) iter2->fPtr;
416
417       // compare CDH headers
418       
419        for(AliHLTInt32_t i=0; i<nWordsHeader; i++ ){
420         if( p1[i] != p2[i] ) sameCDH = 0;
421       }
422
423       fHistHeaderAll->Fill(3);
424       if( sameCDH ) fHistHeaderGood->Fill(3);
425
426       // find rcu headers
427       
428       int startRCU1 = nWordsHeader;
429       int startRCU2 = nWordsHeader;
430
431       for( ; startRCU1<nWords1; startRCU1+=6 ){
432         if( p1[startRCU1]>>30 == 0x2 ) break;     
433       }
434
435       for( ; startRCU2<nWords2; startRCU2+=6 ){
436         if( p2[startRCU2]>>30 == 0x2 ) break;     
437       } 
438       
439
440       // compare RCU headers
441
442       if( startRCU1 < nWords1 || startRCU2 < nWords2 ){
443         if( nWords1 - startRCU1 != nWords2 - startRCU2 ) sameRCU = 0;
444         else{
445           for( AliHLTInt32_t i1=startRCU1, i2=startRCU2; (i1<nWords1) && (i2<nWords2); i1++,i2++ ){
446             if( p1[i1]!=p2[i2] ) sameRCU = 0;
447           }
448         }
449         fHistHeaderAll->Fill(4);
450         if( sameRCU ) fHistHeaderGood->Fill(4);
451       }
452       
453       sameSize = ( startRCU1 == startRCU2 );
454
455       fHistHeaderAll->Fill(5);
456       if( sameSize ) fHistHeaderGood->Fill(5);
457
458       // compare clusters
459
460       if( startRCU1 == startRCU2 ){
461         for( AliHLTInt32_t i=nWordsHeader; i<startRCU1; i+=6){
462           AliHLTUInt32_t header1 = p1[i];
463           AliHLTUInt32_t charge1 = p1[i+1];
464           AliHLTUInt32_t header2 = p2[i];
465           AliHLTUInt32_t charge2 = p2[i+1];
466
467
468           AliHLTUInt32_t flag1 = header1 >> 30;
469           AliHLTUInt32_t qmax1 = ( header1 & 0xFFFFFF ) >> 6;
470           AliHLTUInt32_t row1 = (header1>>24) & 0x3f;
471           AliHLTUInt32_t flag2 = header2 >> 30;
472           AliHLTUInt32_t qmax2 = ( header2 & 0xFFFFFF ) >> 6;
473           AliHLTUInt32_t row2 = (header2>>24) & 0x3f;
474           
475           fHistClusterAll->Fill(0);
476           if( flag1 == flag2 ) fHistClusterGood->Fill(0);
477           else sameFlag = 0;
478           
479           if( flag1!=0x3 || flag2!=0x3 ) continue;
480           
481           // compare cluster qMax
482           
483           fHistClusterAll->Fill(1);
484           if( qmax1 == qmax2 ) fHistClusterGood->Fill(1);
485           else sameQMax = 0;
486           
487           // compare cluster row index
488           
489           fHistClusterAll->Fill(2);
490           if( row1 == row2 ) fHistClusterGood->Fill(2);
491           else sameRow = 0;
492           
493           // compare cluster charge
494           
495           fHistClusterAll->Fill(3);
496           if( charge1 == charge2 ) fHistClusterGood->Fill(3);
497           else sameCharge = 0;
498
499           // compare floating point data
500
501           for( int j=0; j<4; j++ ){
502             AliHLTFloat64_t f1     = *((AliHLTFloat32_t*)&p1[i+2+j]);
503             AliHLTFloat64_t f2     = *((AliHLTFloat32_t*)&p2[i+2+j]);
504             double w = (fabs(f1) + fabs(f2))/2;
505             if( w>1.e-20 ){
506               fHistClusterAll->Fill(4+j);
507               if( fabs(f1 - f2) < 1.e-6*w ) fHistClusterGood->Fill(4+j);
508               else sameFloat[j] = 0;
509             }
510           }
511         }
512       }
513     }
514         
515     fNBlocks++;
516       
517     bool err = 0;
518     TString warn;
519     warn.Form("HWCF consistency check for slice %d patch %d. Values  not matched:", slice1, patch1 );
520     
521     if( nMatchedBlocks!=1 ){ 
522       err=1; 
523       TString x;
524       x.Form(" NMatchedBlocks(%d)",nMatchedBlocks);
525       warn+=x;
526     }
527     if( !sameHLTHd ){ err=1; warn+=", HLT header";}
528     if( !sameCDH ){ err=1; warn+=", CDH header";}
529     if( !sameRCU ){ err=1; warn+=", RCU header";}
530     if( !sameSize ){ err=1; warn+=", N Clusters";}
531     if( !sameFlag ){ err=1; warn+=", Cluster Header";}
532     if( !sameQMax ){ err=1; warn+=", Cluster Qmax";}
533     if( !sameCharge ){ err=1; warn+=", Cluster Charge";}
534     if( !sameRow ){ err=1; warn+=", Cluster Row";}
535     if( !sameFloat[0] ){ err=1; warn+=", Cluster Pad value";}
536     if( !sameFloat[1] ){ err=1; warn+=", Cluster Time value";}
537     if( !sameFloat[2] ){ err=1; warn+=", Cluster PadErr value";}
538     if( !sameFloat[3] ){ err=1; warn+=", Cluster TimeErr value";}
539       
540     if( err ){
541       fNDismatch++;             
542       HLTWarning(warn.Data());
543     } else {
544       //warn+=" NO ";
545       //HLTWarning(warn.Data());     
546     }
547        
548     fHistHeaderAll->Fill(0);
549     if( nMatchedBlocks >=1 ) fHistHeaderGood->Fill(0);
550     fHistHeaderAll->Fill(1);
551     if( nMatchedBlocks <=1 ) fHistHeaderGood->Fill(1);
552
553   } // first block
554     
555   
556   delete[] checkedBlocks;
557    
558   HLTInfo("HWCF consistency check: %.10f %s of %ld data blocks are OK",
559              (double)(fNBlocks-fNDismatch)/(double)fNBlocks*100.,"%",fNBlocks);
560   
561   //if( fNDismatch>0 ){
562   //HLTWarning("HWCF inconsistency: %ld of %ld data blocks are not OK",
563   //fNDismatch,fNBlocks);      
564   //}  
565
566   fProfHeader->Divide(fHistHeaderGood, fHistHeaderAll,1,1.,"b");
567   fProfCluster->Divide(fHistClusterGood, fHistClusterAll, 1, 1, "b");
568
569   PushBack( (TObject*) fProfHeader, kAliHLTDataTypeHistogram|kAliHLTDataOriginTPC,0);
570   fBenchmark.AddOutput(GetLastObjectSize());
571   PushBack( (TObject*) fProfCluster, kAliHLTDataTypeHistogram|kAliHLTDataOriginTPC,0);
572   fBenchmark.AddOutput(GetLastObjectSize());
573   
574
575   fBenchmark.Stop(0);  
576   HLTInfo(fBenchmark.GetStatistics());
577   return iResult;
578 }
579
580