2 //**************************************************************************
3 //* This file is property of and copyright by the ALICE HLT Project *
4 //* ALICE Experiment at CERN, All rights reserved. *
6 //* Primary Authors: Gaute Øvrebekk <st05886@alf.uib.no> *
7 //* for The ALICE HLT Project. *
9 //* Permission to use, copy, modify and distribute this software and its *
10 //* documentation strictly for non-commercial purposes is hereby granted *
11 //* without fee, provided that the above copyright notice appears in all *
12 //* copies and that both the copyright notice and this permission notice *
13 //* appear in the supporting documentation. The authors make no claims *
14 //* about the suitability of this software for any purpose. It is *
15 //* provided "as is" without express or implied warranty. *
16 //**************************************************************************
18 /** @file AliHLTITSClusterFinderComponent.cxx
19 @author Gaute Øvrebekk <st05886@alf.uib.no>
21 @brief Component to run offline clusterfinders
28 #include "AliHLTITSClusterFinderComponent.h"
30 #include "AliCDBEntry.h"
31 #include "AliCDBManager.h"
32 #include "AliCDBStorage.h"
33 #include "AliITSgeomTGeo.h"
34 #include "AliITSRecPoint.h"
35 #include "AliHLTITSSpacePointData.h"
36 #include "AliHLTITSClusterDataFormat.h"
37 #include <AliHLTDAQ.h>
38 #include "AliGeomManager.h"
39 #include "AliITSRecoParam.h"
40 #include "AliITSReconstructor.h"
41 #include "AliHLTITSClusterFinderSPD.h"
42 #include "AliHLTITSClusterFinderSSD.h"
49 #include "TObjString.h"
52 /** ROOT macro for the implementation of ROOT specific class methods */
53 ClassImp(AliHLTITSClusterFinderComponent);
55 AliHLTITSClusterFinderComponent::AliHLTITSClusterFinderComponent(int mode)
58 fInputDataType(kAliHLTVoidDataType),
59 fOutputDataType(kAliHLTVoidDataType),
74 fBenchmark(GetComponentID())
76 // see header file for class documentation
78 // refer to README to build package
80 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
83 case kClusterFinderSPD:
84 fInputDataType = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSPD;
85 fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSPD;
87 case kClusterFinderSDD:
88 fInputDataType = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSDD;
89 fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSDD;
91 case kClusterFinderSSD:
92 fInputDataType = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSSD;
93 fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSSD;
95 case kClusterFinderDigits:
96 fInputDataType = kAliHLTDataTypeAliTreeD|kAliHLTDataOriginITS;
97 fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITS;
100 HLTFatal("unknown cluster finder");
104 AliHLTITSClusterFinderComponent::~AliHLTITSClusterFinderComponent()
106 // see header file for class documentation
114 // Public functions to implement AliHLTComponent's interface.
115 // These functions are required for the registration process
117 const char* AliHLTITSClusterFinderComponent::GetComponentID()
119 // see header file for class documentation
121 case kClusterFinderSPD:
122 return "ITSClusterFinderSPD";
124 case kClusterFinderSDD:
125 return "ITSClusterFinderSDD";
127 case kClusterFinderSSD:
128 return "ITSClusterFinderSSD";
130 case kClusterFinderDigits:
131 return "ITSClusterFinderDigits";
137 void AliHLTITSClusterFinderComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
139 // see header file for class documentation
141 list.push_back( fInputDataType );
144 AliHLTComponentDataType AliHLTITSClusterFinderComponent::GetOutputDataType()
146 // see header file for class documentation
147 return fOutputDataType;
150 void AliHLTITSClusterFinderComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier ) {
151 // see header file for class documentation
153 inputMultiplier = 100;
156 AliHLTComponent* AliHLTITSClusterFinderComponent::Spawn() {
157 // see header file for class documentation
158 return new AliHLTITSClusterFinderComponent(fModeSwitch);
161 Int_t AliHLTITSClusterFinderComponent::DoInit( int argc, const char** argv ) {
162 // see header file for class documentation
164 fBenchmark.SetTimer(0,"total");
165 fBenchmark.SetTimer(1,"reco");
175 Int_t runNo = GetRunNo();
176 AliCDBStorage* store = AliCDBManager::Instance()->GetDefaultStorage();
183 if(store->GetLatestVersion("ITS/Calib/SPDNoisy", runNo)<0){
184 HLTError("SPDNoisy is not found in SPD/Calib");
187 if(store->GetLatestVersion("ITS/Calib/SPDDead", runNo)<0){
188 HLTError("SPDDead is not found in SPD/Calib");
191 if(store->GetLatestVersion("TRIGGER/SPD/PITConditions", runNo)<0){
192 HLTError("PITConditions is not found in TRIGGER/SPD");
197 if(store->GetLatestVersion("ITS/Calib/CalibSDD", runNo)<0){
198 HLTError("CalibSDD is not found in ITS/Calib");
201 if(store->GetLatestVersion("ITS/Calib/RespSDD", runNo)<0){
202 HLTError("RespSDD is not found in ITS/Calib");
205 if(store->GetLatestVersion("ITS/Calib/DriftSpeedSDD", runNo)<0){
206 HLTError("DriftSpeedSDD is not found in ITS/Calib");
209 if(store->GetLatestVersion("ITS/Calib/DDLMapSDD", runNo)<0){
210 HLTError("DDLMapSDD is not found in ITS/Calib");
213 if(store->GetLatestVersion("ITS/Calib/MapsTimeSDD", runNo)<0){
214 HLTError("MapsTimeSDD is not found in ITS/Calib");
219 if(store->GetLatestVersion("ITS/Calib/NoiseSSD", runNo)<0){
220 HLTError("NoiseSSD is not found in ITS/Calib");
223 if(store->GetLatestVersion("ITS/Calib/GainSSD", runNo)<0){
224 HLTError("GainSSD is not found in ITS/Calib");
227 if(store->GetLatestVersion("ITS/Calib/BadChannelsSSD", runNo)<0){
228 HLTError("BadChannelsSSD is not found in ITS/Calib");
232 //General reconstruction
233 if(store->GetLatestVersion("GRP/CTP/Scalers", runNo)<0){
234 HLTError("Scalers is not found in GRP/CTP/");
237 if(!cdbOK){return NULL;}
239 if(fModeSwitch==kClusterFinderSPD) {
240 HLTDebug("using ClusterFinder for SPD");
241 //fNModules=AliITSgeomTGeo::GetNDetectors(1)*AliITSgeomTGeo::GetNLadders(1) + AliITSgeomTGeo::GetNDetectors(2)*AliITSgeomTGeo::GetNLadders(2);
242 fId=AliHLTDAQ::DdlIDOffset("ITSSPD");
243 fNddl=AliHLTDAQ::NumberOfDdls("ITSSPD");
245 else if(fModeSwitch==kClusterFinderSDD) {
246 HLTDebug("using ClusterFinder for SDD");
247 //fNModules=AliITSgeomTGeo::GetNDetectors(3)*AliITSgeomTGeo::GetNLadders(3) + AliITSgeomTGeo::GetNDetectors(4)*AliITSgeomTGeo::GetNLadders(4);
248 fId=AliHLTDAQ::DdlIDOffset("ITSSDD");
249 fNddl=AliHLTDAQ::NumberOfDdls("ITSSDD");
251 else if(fModeSwitch==kClusterFinderSSD) {
252 HLTDebug("using ClusterFinder for SSD");
253 //fNModules=AliITSgeomTGeo::GetNDetectors(5)*AliITSgeomTGeo::GetNLadders(5) + AliITSgeomTGeo::GetNDetectors(6)*AliITSgeomTGeo::GetNLadders(6);
254 fId=AliHLTDAQ::DdlIDOffset("ITSSSD");
255 fNddl=AliHLTDAQ::NumberOfDdls("ITSSSD");
257 else if(fModeSwitch==kClusterFinderDigits) {
261 HLTFatal("No mode set for clusterfindercomponent");
264 //Removed the warning for loading default RecoParam in HLT
265 AliITSRecoParam *par = AliITSRecoParam::GetLowFluxParam();
266 AliITSReconstructor *rec = new AliITSReconstructor();
267 rec->SetRecoParam(par);
269 AliCDBManager* man = AliCDBManager::Instance();
270 if (!man->IsDefaultStorageSet()){
271 HLTError("Default CDB storage has not been set !");
275 if(AliGeomManager::GetGeometry()==NULL){
276 AliGeomManager::LoadGeometry();
279 //fgeomInit = new AliITSInitGeometry(kvSPD02,2);
280 fgeomInit = new AliITSInitGeometry(kvPPRasymmFMD,2);
281 //fgeomInit->InitAliITSgeom(fgeom);
282 fgeom = fgeomInit->CreateAliITSgeom();
284 fNModules = fgeom->GetIndexMax();
286 fClusters = new TClonesArray*[fNModules];
287 for (Int_t iModule = 0; iModule < fNModules; iModule++) {
288 fClusters[iModule] = NULL;
292 fDettype = new AliITSDetTypeRec();
293 fDettype->SetITSgeom(fgeom);
294 fDettype->SetDefaults();
295 fDettype->SetDefaultClusterFindersV2(kTRUE);
300 fRawReader = new AliRawReaderMemory();
301 fSPD = new AliHLTITSClusterFinderSPD( fDettype );
302 fSSD = new AliHLTITSClusterFinderSSD( fDettype, fRawReader );
304 TString arguments = "";
305 for ( int i = 0; i < argc; i++ ) {
306 if ( !arguments.IsNull() ) arguments += " ";
307 arguments += argv[i];
313 return Configure( arguments.Data() );
316 Int_t AliHLTITSClusterFinderComponent::DoDeinit() {
317 // see header file for class documentation
337 for (Int_t iModule = 0; iModule < fNModules; iModule++) {
338 if(fClusters[iModule] != NULL){
339 fClusters[iModule]->Delete();
340 delete fClusters[iModule];
342 fClusters[iModule] = NULL;
345 fUseOfflineFinder = 0;
350 int AliHLTITSClusterFinderComponent::DoEvent
352 const AliHLTComponentEventData& evtData,
353 const AliHLTComponentBlockData* /*blocks*/,
354 AliHLTComponentTriggerData& /*trigData*/,
355 AliHLTUInt8_t* outputPtr,
356 AliHLTUInt32_t& size,
357 vector<AliHLTComponentBlockData>& outputBlocks )
359 // see header file for class documentation
361 AliHLTUInt32_t maxBufferSize = size;
362 size = 0; // output size
364 if (!IsDataEvent()) return 0;
366 if ( evtData.fBlockCnt<=0 )
368 HLTDebug("no blocks in event" );
372 fBenchmark.StartNewEvent();
374 for( const AliHLTComponentBlockData *i= GetFirstInputBlock(fInputDataType); i!=NULL; i=GetNextInputBlock() ){
375 fBenchmark.AddInput(i->fSize);
380 if(fModeSwitch==kClusterFinderDigits) {
382 for ( const TObject *iter = GetFirstInputObject(fInputDataType); iter != NULL; iter = GetNextInputObject() ) {
383 tD = dynamic_cast<TTree*>(const_cast<TObject*>( iter ) );
385 HLTFatal("No Digit Tree found");
388 // 2010-04-17 very crude workaround: TTree objects are difficult to send
389 // The actual case: Running ITS and TPC reconstruction fails at the second event
390 // to read the ITS digits from the TreeD
392 // Reason: reading fails in TBranch::GetBasket, there a new basket is opened from
393 // a TFile object. The function TBranch::GetFile returns the file object from
394 // an internal fDirectory (TDirectory) object. This file is at the second event
395 // set to the TPC.Digits.root. The internal mismatch creates a seg fault
397 // Investigation: TBranch::Streamer uses a crude assignment after creating the
399 // fDirectory = gDirectory;
400 // gDirectory is obviously not set correctly. Setting the directory to a TFile
401 // object for the ITS digits helps to fix the internal mess. Tried also to set
402 // the Directory for the TreeD to NULL (This has only effect if ones sets it
403 // to something not NULL first, and then to NULL). But then no content, i.e.
404 // ITS clusters could be retrieved.
406 // Conclusion: TTree objects are hardly to be sent via TMessage, there are direct
407 // links to the file required anyhow.
408 TFile* dummy=new TFile("ITS.Digits.root");
409 tD->SetDirectory(dummy);
412 fDettype->SetTreeAddressD(tD);
413 fDettype->MakeBranch(tR,"R");
414 fDettype->SetTreeAddressR(tR);
417 fDettype->DigitsToRecPoints(tD,tR,0,opt,kTRUE);
419 TClonesArray * fRecPoints = NULL;
420 tR->SetBranchAddress("ITSRecPoints",&fRecPoints);
421 for(Int_t treeEntry=0;treeEntry<tR->GetEntries();treeEntry++){
422 tR->GetEntry(treeEntry);
423 for(Int_t tCloneEntry=0;tCloneEntry<fRecPoints->GetEntries();tCloneEntry++){
424 AliITSRecPoint *recpoint=(AliITSRecPoint*)fRecPoints->At(tCloneEntry);
425 fclusters.push_back(*recpoint);
435 UInt_t nClusters=fclusters.size();
437 UInt_t bufferSize = nClusters * sizeof(AliHLTITSSpacePointData) + sizeof(AliHLTITSClusterData);
438 if( size + bufferSize > maxBufferSize ){
439 HLTWarning( "Output buffer size exceed (buffer size %d, current size %d)", maxBufferSize, size+bufferSize);
445 RecPointToSpacePoint(outputPtr,size);
447 AliHLTComponentBlockData bd;
450 bd.fSize = bufferSize;
451 bd.fSpecification = 0x00000000;
452 bd.fDataType = GetOutputDataType();
453 outputBlocks.push_back( bd );
455 fBenchmark.AddOutput(bd.fSize);
462 // -- Loop over blocks
463 for( const AliHLTComponentBlockData* iter = GetFirstInputBlock(fInputDataType); iter != NULL; iter = GetNextInputBlock() ) {
465 // -- Debug output of datatype --
466 HLTDebug("Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s",
467 evtData.fEventID, evtData.fEventID,
468 DataType2Text(iter->fDataType).c_str(),
469 DataType2Text(fInputDataType).c_str());
471 // -- Check for the correct data type
472 if ( iter->fDataType != (fInputDataType) )
475 // -- Get equipment ID out of specification
476 AliHLTUInt32_t spec = iter->fSpecification;
479 for ( Int_t ii = 0; ii < fNddl ; ii++ ) { //number of ddl's
480 if ( spec & 0x00000001 ) {
487 // -- Set equipment ID to the raw reader
490 HLTWarning("The fRawReader pointer is NULL");
494 if(!fRawReader->AddBuffer((UChar_t*) iter->fPtr, iter->fSize, id)){
495 HLTWarning("Could not add buffer");
500 if(fModeSwitch==kClusterFinderSPD && !fUseOfflineFinder){ fSPD->RawdataToClusters( fRawReader, fclusters ); }
501 else if(fModeSwitch==kClusterFinderSSD && !fUseOfflineFinder){ fSSD->RawdataToClusters( fclusters ); }
503 if(fModeSwitch==kClusterFinderSPD && fUseOfflineFinder) {fDettype->DigitsToRecPoints(fRawReader,fClusters,"SPD");}
504 if(fModeSwitch==kClusterFinderSSD && fUseOfflineFinder) {fDettype->DigitsToRecPoints(fRawReader,fClusters,"SSD");}
505 if(fModeSwitch==kClusterFinderSDD) {fDettype->DigitsToRecPoints(fRawReader,fClusters,"SDD");}
506 for(int i=0;i<fNModules;i++){
507 if(fClusters[i] != NULL){
508 for(int j=0;j<fClusters[i]->GetEntriesFast();j++){
509 AliITSRecPoint *recpoint = (AliITSRecPoint*) (fClusters[i]->At(j));
510 fclusters.push_back(*recpoint);
512 fClusters[i]->Delete();
520 fRawReader->ClearBuffers();
522 UInt_t nClusters=fclusters.size();
524 UInt_t bufferSize = nClusters * sizeof(AliHLTITSSpacePointData) + sizeof(AliHLTITSClusterData);
525 if( size + bufferSize > maxBufferSize ){
526 HLTWarning( "Output buffer size exceed (buffer size %d, current size %d)", maxBufferSize, size+bufferSize);
532 RecPointToSpacePoint(outputPtr,size);
534 AliHLTComponentBlockData bd;
537 bd.fSize = bufferSize;
538 bd.fSpecification = iter->fSpecification;
539 bd.fDataType = GetOutputDataType();
540 outputBlocks.push_back( bd );
542 fBenchmark.AddOutput(bd.fSize);
550 HLTInfo(fBenchmark.GetStatistics());
554 int AliHLTITSClusterFinderComponent::Configure(const char* arguments)
559 if (!arguments) return iResult;
561 TString allArgs=arguments;
564 TObjArray* pTokens=allArgs.Tokenize(" ");
567 for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
568 argument=((TObjString*)pTokens->At(i))->GetString();
569 if (argument.IsNull()) continue;
570 if (argument.CompareTo("-use-offline-finder")==0) {
571 fUseOfflineFinder = 1;
572 HLTInfo("Off-line ClusterFinder will be used");
576 else if (argument.CompareTo("")==0) {
582 HLTError("unknown argument %s", argument.Data());
593 int AliHLTITSClusterFinderComponent::Reconfigure(const char* cdbEntry, const char* chainId)
595 // see header file for class documentation
598 const char* path="HLT/ConfigITS/ClusterFinderComponent";
599 const char* defaultNotify="";
602 defaultNotify=" (default)";
605 HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
606 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path);
608 TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
610 HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
611 iResult=Configure(pString->GetString().Data());
613 HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
616 HLTError("can not fetch object \"%s\" from CDB", path);
622 void AliHLTITSClusterFinderComponent::GetOCDBObjectDescription( TMap* const targetMap)
624 // Get a list of OCDB object description.
625 if (!targetMap) return;
627 targetMap->Add(new TObjString("ITS/Calib/SPDNoisy"),new TObjString("Calibration object for SPD" ));
628 targetMap->Add(new TObjString("ITS/Calib/SPDDead"),new TObjString("Calibration object for SPD" ));
629 targetMap->Add(new TObjString("TRIGGER/SPD/PITConditions"),new TObjString("Calibration object for SPD" ));
631 targetMap->Add(new TObjString("ITS/Calib/CalibSDD"),new TObjString("Calibration object for SDD" ));
632 targetMap->Add(new TObjString("ITS/Calib/RespSDD"),new TObjString("Calibration object for SDD" ));
633 targetMap->Add(new TObjString("ITS/Calib/DriftSpeedSDD"),new TObjString("Calibration object for SDD" ));
634 targetMap->Add(new TObjString("ITS/Calib/DDLMapSDD"),new TObjString("Calibration object for SDD" ));
635 targetMap->Add(new TObjString("ITS/Calib/MapsTimeSDD"),new TObjString("Calibration object for SDD" ));
637 targetMap->Add(new TObjString("ITS/Calib/NoiseSSD"),new TObjString("Calibration object for SSD" ));
638 targetMap->Add(new TObjString("ITS/Calib/GainSSD"),new TObjString("Calibration object for SSD" ));
639 targetMap->Add(new TObjString("ITS/Calib/BadChannelsSSD"),new TObjString("Calibration object for SSD" ));
640 //General reconstruction
641 targetMap->Add(new TObjString("GRP/CTP/Scalers"),new TObjString("General reconstruction object" ));
645 void AliHLTITSClusterFinderComponent::RecPointToSpacePoint(AliHLTUInt8_t* outputPtr,AliHLTUInt32_t& size){
646 AliHLTITSClusterData *outputClusters = reinterpret_cast<AliHLTITSClusterData*>(outputPtr + size);
647 outputClusters->fSpacePointCnt=fclusters.size();
649 for(UInt_t i=0;i<fclusters.size();i++){
650 AliITSRecPoint *recpoint = (AliITSRecPoint*) &(fclusters[i]);
651 outputClusters->fSpacePoints[clustIdx].fY=recpoint->GetY();
652 outputClusters->fSpacePoints[clustIdx].fZ=recpoint->GetZ();
653 outputClusters->fSpacePoints[clustIdx].fSigmaY2=recpoint->GetSigmaY2();
654 outputClusters->fSpacePoints[clustIdx].fSigmaZ2=recpoint->GetSigmaZ2();
655 outputClusters->fSpacePoints[clustIdx].fSigmaYZ=recpoint->GetSigmaYZ();
656 outputClusters->fSpacePoints[clustIdx].fQ=recpoint->GetQ();
657 outputClusters->fSpacePoints[clustIdx].fNy=recpoint->GetNy();
658 outputClusters->fSpacePoints[clustIdx].fNz=recpoint->GetNz();
659 outputClusters->fSpacePoints[clustIdx].fLayer=recpoint->GetLayer();
660 outputClusters->fSpacePoints[clustIdx].fIndex=recpoint->GetDetectorIndex() | recpoint->GetPindex() | recpoint->GetNindex();
661 outputClusters->fSpacePoints[clustIdx].fTracks[0]=recpoint->GetLabel(0);
662 outputClusters->fSpacePoints[clustIdx].fTracks[1]=recpoint->GetLabel(1);
663 outputClusters->fSpacePoints[clustIdx].fTracks[2]=recpoint->GetLabel(2);