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