3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project *
5 //* ALICE Experiment at CERN, All rights reserved. *
8 //* for The ALICE HLT Project. *
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 //**************************************************************************
19 /** @file AliHLTTRDTrackerV1Component.cxx
20 @author Theodor Rascanu
22 @brief A TRDTrackerV1 processing component for the HLT.
29 #include "AliHLTTRDTrackerV1Component.h"
30 #include "AliHLTTRDDefinitions.h"
31 #include "AliHLTTRDTrack.h"
32 #include "AliHLTTRDUtils.h"
37 #include "AliGeomManager.h"
38 #include "AliCDBManager.h"
39 #include "AliCDBStorage.h"
40 #include "AliCDBEntry.h"
41 #include "AliESDEvent.h"
42 #include "AliESDfriend.h"
44 #include "AliTRDcalibDB.h"
45 #include "AliTRDReconstructor.h"
46 #include "AliTRDtrackerV1.h"
47 #include "AliTRDrecoParam.h"
53 ClassImp(AliHLTTRDTrackerV1Component)
55 void AliHLTTRDTrackerV1Component::AliHLTTRDESDEvent::CreateStdContent()
57 TClonesArray* tracksArray = new TClonesArray("AliESDtrack",0);
58 tracksArray->SetName(AliESDEvent::fgkESDListName[AliESDEvent::kTracks]);
59 AddObject(tracksArray);
63 void AliHLTTRDTrackerV1Component::AliHLTTRDESDEvent::Streamer(TBuffer &/*R__b*/)
65 AliFatal("class is for internal us only and not for streaming");
68 AliHLTTRDTrackerV1Component::AliHLTTRDTrackerV1Component():
70 fOutputPercentage(100), // By default we copy to the output exactly what we got as input
79 fgeometryFileName(""),
81 fOutputV1Tracks(kTRUE),
82 fHighLevelOutput(kFALSE),
83 fEmulateHLTTracks(kFALSE),
84 fImproveTracklets(kFALSE)
86 // Default constructor
90 AliHLTTRDTrackerV1Component::~AliHLTTRDTrackerV1Component()
95 const char* AliHLTTRDTrackerV1Component::GetComponentID()
97 // Return the component ID const char *
98 return "TRDTrackerV1"; // The ID of this component
101 void AliHLTTRDTrackerV1Component::GetInputDataTypes( vector<AliHLTComponent_DataType>& list)
103 // Get the list of input data
104 list.clear(); // We do not have any requirements for our input data type(s).
105 list.push_back(AliHLTTRDDefinitions::fgkClusterDataType);
108 AliHLTComponentDataType AliHLTTRDTrackerV1Component::GetOutputDataType()
110 // Get the output data type
111 return kAliHLTMultipleDataType;
114 int AliHLTTRDTrackerV1Component::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
116 // Get the output data types
118 tgtList.push_back(kAliHLTDataTypeTrack | kAliHLTDataOriginTRD);
119 tgtList.push_back(AliHLTTRDDefinitions::fgkTracksDataType);
120 return tgtList.size();
123 void AliHLTTRDTrackerV1Component::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
125 // Get the output data size
127 inputMultiplier = fOutputV1Tracks ? 2*((double)fOutputPercentage)/100.0 : 0.5*((double)fOutputPercentage)/100.0;
130 // Spawn function, return new instance of this class
131 AliHLTComponent* AliHLTTRDTrackerV1Component::Spawn()
133 return new AliHLTTRDTrackerV1Component;
137 int AliHLTTRDTrackerV1Component::DoInit( int argc, const char** argv )
139 // perform initialization. We check whether our relative output size is specified in the arguments.
142 fReconstructor = new AliTRDReconstructor();
143 HLTDebug("TRDReconstructor at 0x%x", fReconstructor);
144 fESD = new AliHLTTRDESDEvent();
145 fESD->CreateStdContent();
147 TString configuration="";
149 for (int i=0; i<argc && iResult>=0; i++) {
151 if (!configuration.IsNull()) configuration+=" ";
152 configuration+=argument;
155 if (!configuration.IsNull()) {
156 iResult=Configure(configuration.Data());
158 iResult=Reconfigure(NULL, NULL);
161 if(iResult<0) return iResult;
163 fTracker = new AliTRDtrackerV1();
164 HLTDebug("TRDTracker at 0x%x", fTracker);
165 fTracker->SetReconstructor(fReconstructor);
167 fClusterArray = new TClonesArray("AliTRDcluster"); // would be nice to allocate memory for all clusters here.
172 int AliHLTTRDTrackerV1Component::DoDeinit()
174 // Deinitialization of the component
176 fTracker->SetClustersOwner(kFALSE);
180 fClusterArray->Delete();
181 delete fClusterArray;
182 fClusterArray = NULL;
184 // We need to set clusters in Reconstructor to null to prevent from
185 // double deleting, since we delete TClonesArray by ourself.
186 fReconstructor->SetClusters(0x0);
187 delete fReconstructor;
188 fReconstructor = NULL;
192 AliTRDcalibDB::Terminate();
197 int AliHLTTRDTrackerV1Component::DoEvent( const AliHLTComponentEventData& evtData,
198 const AliHLTComponentBlockData* blocks,
199 AliHLTComponent_TriggerData& /*trigData*/,
200 AliHLTUInt8_t* outputPtr,
201 AliHLTUInt32_t& size,
202 vector<AliHLTComponent_BlockData>& outputBlocks )
206 HLTDebug("NofBlocks %i", evtData.fBlockCnt );
208 AliHLTUInt32_t totalSize = 0, offset = 0;
210 AliHLTComponentDataType expectedDataType = AliHLTTRDDefinitions::fgkClusterDataType;
211 for ( unsigned long iBlock = 0; iBlock < evtData.fBlockCnt; iBlock++ )
213 const AliHLTComponentBlockData &block = blocks[iBlock];
214 AliHLTComponentDataType inputDataType = block.fDataType;
216 if(inputDataType != expectedDataType)
218 HLTDebug( "Block # %i/%i; Event 0x%08LX (%Lu) Wrong received datatype: %s - Skipping",
219 iBlock, evtData.fBlockCnt-1,
220 evtData.fEventID, evtData.fEventID,
221 DataType2Text(inputDataType).c_str());
225 HLTDebug("We get the right data type: Block # %i/%i; Event 0x%08LX (%Lu) Received datatype: %s; Block Size: %i",
226 iBlock, evtData.fBlockCnt-1,
227 evtData.fEventID, evtData.fEventID,
228 DataType2Text(inputDataType).c_str(),
233 unsigned long constBase;
234 double inputMultiplier;
235 GetOutputDataSize(constBase,inputMultiplier);
236 if(size<(constBase+block.fSize*inputMultiplier)){
237 HLTWarning("Memory Block given might be too small: %i < %i; Event %Lu", size, constBase+block.fSize*inputMultiplier, evtData.fEventID);
242 //fESD->SetMagneticField(GetBz());
244 AliHLTTRDUtils::ReadClusters(fClusterArray, block.fPtr, block.fSize, &fNtimeBins);
245 HLTDebug("Reading number of time bins from input block. Setting number of timebins to %d", fNtimeBins);
246 AliTRDtrackerV1::SetNTimeBins(fNtimeBins);
248 HLTDebug("TClonesArray of clusters: nbEntries = %i", fClusterArray->GetEntriesFast());
249 fTracker->LoadClusters(fClusterArray);
251 fTracker->Clusters2Tracks(fESD);
253 Int_t nTracks = fESD->GetNumberOfTracks();
254 HLTInfo("Number of tracks == %d ==", nTracks);
256 TClonesArray* trdTracks;
257 trdTracks = fTracker->GetListOfTracks();
259 if(fHighLevelOutput){
260 if(fEmulateHLTTracks && trdTracks){
261 // TClonesArray* oldArr = trdTracks;
262 trdTracks = new TClonesArray(*trdTracks);
263 AliHLTTRDUtils::EmulateHLTTracks(trdTracks);
264 // if(oldArr->At(0)){
265 // HLTInfo("Old Track:");
266 // ((AliTRDtrackV1*)oldArr->At(0))->Print("a");
267 // HLTInfo("\nNew Track:");
268 // ((AliTRDtrackV1*)trdTracks->At(0))->Print("a");
273 strg.String() += fNtimeBins;
275 PushBack(trdTracks, AliHLTTRDDefinitions::fgkHiLvlTracksDataType, block.fSpecification);
277 TClonesArray temp("AliTRDtrackV1");
278 PushBack(&temp, AliHLTTRDDefinitions::fgkHiLvlTracksDataType, block.fSpecification);
280 PushBack(&strg, AliHLTTRDDefinitions::fgkHiLvlTracksDataType, block.fSpecification);
282 if(fEmulateHLTTracks && trdTracks){
288 HLTDebug("We have an output ESDEvent: 0x%x with %i tracks", fESD, nTracks);
289 AliHLTUInt32_t addedSize = AliHLTTRDUtils::AddESDToOutput(fESD, outputPtr+offset);
290 totalSize += addedSize;
293 AliHLTComponentBlockData bd;
295 //bd.fPtr = outputPtr;
297 bd.fSize = addedSize;
298 bd.fSpecification = block.fSpecification;
299 bd.fDataType = kAliHLTDataTypeTrack | kAliHLTDataOriginTRD;
300 outputBlocks.push_back( bd );
301 HLTDebug("BD ptr 0x%x, offset %i, size %i, datav1Type %s, spec 0x%x ", bd.fPtr, bd.fOffset, bd.fSize, DataType2Text(bd.fDataType).c_str(), bd.fSpecification);
304 if (fOutputV1Tracks && trdTracks){
305 HLTDebug("We have an output array: pointer to trdTracks = 0x%x, nbEntries = %i", trdTracks, trdTracks->GetEntriesFast());
307 addedSize = AliHLTTRDUtils::AddTracksToOutput(trdTracks, outputPtr+offset, fNtimeBins);
308 totalSize += addedSize;
312 //bd.fPtr = outputPtr;
314 bd.fSize = addedSize;
315 bd.fSpecification = block.fSpecification;
316 bd.fDataType = AliHLTTRDDefinitions::fgkTracksDataType;
317 outputBlocks.push_back( bd );
318 HLTDebug("BD ptr 0x%x, offset %i, size %i, dataType %s, spec 0x%x ", bd.fPtr, bd.fOffset, bd.fSize, DataType2Text(bd.fDataType).c_str(), bd.fSpecification);
323 HLTDebug("totalSize: %i", totalSize);
325 // if ( totalSize > allocSize )
327 // HLTError("Too much data; Data written over allowed buffer. Amount written: %lu, allowed amount: %lu.",
328 // totalSize, size );
332 //here we are deleting clusters (but not the TClonesArray itself)
333 fTracker->UnloadClusters();
334 AliTRDReconstructor::SetClusters(0x0);
335 fClusterArray->Delete();
340 HLTDebug("Event is done. size written to the output is %i", size);
344 int AliHLTTRDTrackerV1Component::Configure(const char* arguments){
346 if (!arguments) return iResult;
348 TString allArgs=arguments;
352 TObjArray* pTokens=allArgs.Tokenize(" ");
354 for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
355 argument=((TObjString*)pTokens->At(i))->GetString();
356 if (argument.IsNull()) continue;
358 if (argument.CompareTo("output_percentage")==0) {
359 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
360 HLTInfo("Setting output percentage to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
361 fOutputPercentage=((TObjString*)pTokens->At(i))->GetString().Atoi();
364 else if (argument.CompareTo("-solenoidBz")==0) {
365 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
366 HLTWarning("argument -solenoidBz is deprecated, magnetic field set up globally (%f)", GetBz());
369 else if (argument.CompareTo("-geometry")==0) {
370 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
371 HLTInfo("Setting geometry to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
372 fgeometryFileName=((TObjString*)pTokens->At(i))->GetString();
375 else if (argument.CompareTo("-lowflux")==0) {
377 HLTInfo("Low flux reconstruction selected");
380 else if (argument.CompareTo("-highflux")==0) {
382 HLTInfo("High flux reconstruction selected");
385 else if (argument.CompareTo("-cosmics")==0) {
387 HLTInfo("Cosmics reconstruction selected");
390 else if (argument.CompareTo("-HLTflag")==0) {
391 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
392 TString toCompareTo=((TObjString*)pTokens->At(i))->GetString();
393 if (toCompareTo.CompareTo("yes")==0){
394 HLTInfo("Setting HLTflag to: %s", toCompareTo.Data());
397 else if (toCompareTo.CompareTo("no")==0){
398 HLTInfo("Setting HLTflag to: %s", toCompareTo.Data());
402 HLTError("unknown argument for HLTflag: %s", toCompareTo.Data());
408 else if (argument.CompareTo("-outputV1Tracks")==0) {
409 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
410 TString toCompareTo=((TObjString*)pTokens->At(i))->GetString();
411 if (toCompareTo.CompareTo("yes")==0){
412 HLTInfo("Setting OutputV1Tracks to: %s", toCompareTo.Data());
413 fOutputV1Tracks=kTRUE;
415 else if (toCompareTo.CompareTo("no")==0){
416 HLTInfo("Setting OutputV1Tracks to: %s", toCompareTo.Data());
417 fOutputV1Tracks=kFALSE;
420 HLTError("unknown argument for OutputV1Tracks: %s", toCompareTo.Data());
426 else if (argument.CompareTo("-highLevelOutput")==0) {
427 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
428 TString toCompareTo=((TObjString*)pTokens->At(i))->GetString();
429 if (toCompareTo.CompareTo("yes")==0){
430 HLTWarning("Setting highLevelOutput to: %s", toCompareTo.Data());
431 fHighLevelOutput=kTRUE;
433 else if (toCompareTo.CompareTo("no")==0){
434 HLTInfo("Setting highLevelOutput to: %s", toCompareTo.Data());
435 fHighLevelOutput=kFALSE;
438 HLTError("unknown argument for highLevelOutput: %s", toCompareTo.Data());
444 else if (argument.CompareTo("-emulateHLToutput")==0) {
445 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
446 TString toCompareTo=((TObjString*)pTokens->At(i))->GetString();
447 if (toCompareTo.CompareTo("yes")==0){
448 HLTWarning("Setting emulateHLToutput to: %s", toCompareTo.Data());
449 fEmulateHLTTracks=kTRUE;
451 else if (toCompareTo.CompareTo("no")==0){
452 HLTInfo("Setting emulateHLToutput to: %s", toCompareTo.Data());
453 fEmulateHLTTracks=kFALSE;
456 HLTError("unknown argument for emulateHLToutput: %s", toCompareTo.Data());
462 else if (argument.CompareTo("-PIDmethod")==0) {
463 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
464 TString toCompareTo=((TObjString*)pTokens->At(i))->GetString();
465 if (toCompareTo.CompareTo("LH")==0){
466 HLTInfo("Setting PID method to: %s", toCompareTo.Data());
469 else if (toCompareTo.CompareTo("NN")==0){
470 HLTInfo("Setting PID method to: %s", toCompareTo.Data());
473 else if (toCompareTo.CompareTo("TM")==0){
474 HLTInfo("Setting PID method to: %s", toCompareTo.Data());
478 HLTError("unknown argument for PID method: %s", toCompareTo.Data());
486 HLTError("unknown argument: %s", argument.Data());
494 HLTError("missing parameter for argument %s", argument.Data());
503 int AliHLTTRDTrackerV1Component::SetParams()
506 if(!AliCDBManager::Instance()->IsDefaultStorageSet()){
507 HLTError("DefaultStorage is not set in CDBManager");
510 if(AliCDBManager::Instance()->GetRun()<0){
511 HLTError("Run Number is not set in CDBManager");
514 HLTInfo("CDB default storage: %s; RunNo: %i", (AliCDBManager::Instance()->GetDefaultStorage()->GetBaseFolder()).Data(), AliCDBManager::Instance()->GetRun());
516 if(!AliGeomManager::GetGeometry()){
517 if(fgeometryFileName.CompareTo("")==0 || !TFile::Open(fgeometryFileName.Data())){
518 HLTInfo("Loading standard geometry file");
519 AliGeomManager::LoadGeometry();
521 HLTWarning("Loading NON-standard geometry file");
522 AliGeomManager::LoadGeometry(fgeometryFileName.Data());
524 if(!AliGeomManager::GetGeometry()){
525 HLTError("Could not load geometry");
528 HLTInfo("Applying Alignment from CDB object");
529 AliGeomManager::ApplyAlignObjsFromCDB("TRD");
532 HLTInfo("Geometry Already Loaded!");
535 if(fReconstructor->GetRecoParam()){
536 fRecoParam = new AliTRDrecoParam(*fReconstructor->GetRecoParam());
537 HLTInfo("RecoParam already set!");
539 if(fRecoParamType == 0){
540 HLTDebug("Low flux params init.");
541 fRecoParam = AliTRDrecoParam::GetLowFluxParam();
543 if(fRecoParamType == 1){
544 HLTDebug("High flux params init.");
545 fRecoParam = AliTRDrecoParam::GetHighFluxParam();
547 if(fRecoParamType == 2){
548 HLTDebug("Cosmic Test params init.");
549 fRecoParam = AliTRDrecoParam::GetCosmicTestParam();
555 HLTError("No reco params initialized. Sniffing big trouble!");
560 case 0: fRecoParam->SetPIDNeuralNetwork(kFALSE); break;
561 case 1: fRecoParam->SetPIDNeuralNetwork(kTRUE); break;
562 case 2: fRecoParam->SetPIDNeuralNetwork(kFALSE); break;
565 fRecoParam->SetImproveTracklets(fImproveTracklets);
567 fRecoParam->SetStreamLevel(AliTRDrecoParam::kTracker, 0);
568 fReconstructor->SetRecoParam(fRecoParam);
570 TString recoOptions="sa,!cw";
573 recoOptions += ",hlt";
575 HLTDebug("Reconstructor options are: %s",recoOptions.Data());
576 fReconstructor->SetOption(recoOptions.Data());
581 int AliHLTTRDTrackerV1Component::Reconfigure(const char* cdbEntry, const char* chainId)
583 // see header file for class documentation
586 const char* path="HLT/ConfigTRD/TrackerV1Component";
587 const char* defaultNotify="";
590 defaultNotify=" (default)";
593 HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
594 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path/*,GetRunNo()*/);
596 TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
598 HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
599 iResult=Configure(pString->GetString().Data());
601 HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
604 HLTError("cannot fetch object \"%s\" from CDB", path);
612 int AliHLTTRDTrackerV1Component::ReadPreprocessorValues(const char* modules)
614 // see header file for class documentation
617 TString str(modules);
618 if(str.Contains("HLT") || str.Contains("TRD") || str.Contains("GRP")){
624 void AliHLTTRDTrackerV1Component::GetOCDBObjectDescription(TMap* const targetMap){
625 // Get a list of OCDB object description needed for the particular component
626 if (!targetMap) return;
627 targetMap->Add(new TObjString("HLT/ConfigTRD/TrackerV1Component"), new TObjString("component arguments"));
628 targetMap->Add(new TObjString("TRD/Calib/ChamberGainFactor"), new TObjString("gain factor of chambers"));
629 targetMap->Add(new TObjString("TRD/Calib/ChamberT0"), new TObjString("T0 of chambers"));
630 targetMap->Add(new TObjString("TRD/Calib/ChamberVdrift"), new TObjString("drift velocity of chambers"));
631 targetMap->Add(new TObjString("TRD/Calib/DetNoise"), new TObjString("noise of chambers"));
632 targetMap->Add(new TObjString("TRD/Calib/LocalGainFactor"), new TObjString("per pad gain factor"));
633 targetMap->Add(new TObjString("TRD/Calib/LocalT0"), new TObjString("per pad T0"));
634 targetMap->Add(new TObjString("TRD/Calib/LocalVdrift"), new TObjString("per pad drift velocity"));
635 targetMap->Add(new TObjString("TRD/Calib/PadNoise"), new TObjString("per pad noise"));
636 targetMap->Add(new TObjString("TRD/Calib/PadStatus"), new TObjString("pad status"));
637 targetMap->Add(new TObjString("TRD/Calib/PRFWidth"), new TObjString("pad response function"));
638 targetMap->Add(new TObjString("TRD/Calib/ChamberStatus"), new TObjString("status of chambers"));
639 targetMap->Add(new TObjString("TRD/Calib/PIDLQ"), new TObjString("likelyhood PID"));
640 targetMap->Add(new TObjString("TRD/Calib/PIDNN"), new TObjString("neuronal network PID"));
641 targetMap->Add(new TObjString("TRD/Calib/PIDThresholds"), new TObjString("threshold for PID"));