]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/ITS/AliHLTITSClusterFinderComponent.cxx
MC labels added for ITS tracks ?\127 ?\127?\127
[u/mrichter/AliRoot.git] / HLT / ITS / AliHLTITSClusterFinderComponent.cxx
1 // $Id$
2 //**************************************************************************
3 //* This file is property of and copyright by the ALICE HLT Project        * 
4 //* ALICE Experiment at CERN, All rights reserved.                         *
5 //*                                                                        *
6 //* Primary Authors: Gaute Øvrebekk <st05886@alf.uib.no>                   *
7 //*                  for The ALICE HLT Project.                            *
8 //*                                                                        *
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 //**************************************************************************
17
18 /** @file   AliHLTITSClusterFinderComponent.cxx
19     @author Gaute Øvrebekk <st05886@alf.uib.no>
20     @date   
21     @brief  Component to run offline clusterfinders
22 */
23
24 #if __GNUC__>= 3
25 using namespace std;
26 #endif
27
28 #include "AliHLTITSClusterFinderComponent.h" 
29
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"
43 #include "TMap.h"
44
45 #include <cstdlib>
46 #include <cerrno>
47 #include "TFile.h"
48 #include "TString.h"
49 #include "TObjString.h"
50 #include <sys/time.h>
51
52 /** ROOT macro for the implementation of ROOT specific class methods */
53 ClassImp(AliHLTITSClusterFinderComponent);
54
55 AliHLTITSClusterFinderComponent::AliHLTITSClusterFinderComponent(int mode)
56   :
57   fModeSwitch(mode),
58   fInputDataType(kAliHLTVoidDataType),
59   fOutputDataType(kAliHLTVoidDataType),
60   fUseOfflineFinder(0),
61   fNModules(0),
62   fId(0),
63   fNddl(0),
64   fClusters(NULL),
65   fRawReader(NULL),
66   fDettype(NULL),
67   fgeom(NULL),
68   fgeomInit(NULL),
69   fSPD(NULL),
70   fSSD(NULL),
71   tD(NULL),
72   tR(NULL),
73   fclusters(),
74   fBenchmark(GetComponentID())
75
76   // see header file for class documentation
77   // or
78   // refer to README to build package
79   // or
80   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
81
82   switch(fModeSwitch){
83   case kClusterFinderSPD:
84     fInputDataType  = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSPD;
85     fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSPD;
86     break;
87   case kClusterFinderSDD:        
88     fInputDataType  = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSDD;
89     fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSDD;
90     break;
91   case kClusterFinderSSD:
92     fInputDataType  = kAliHLTDataTypeDDLRaw | kAliHLTDataOriginITSSSD;
93     fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITSSSD;
94     break;
95   case kClusterFinderDigits:
96     fInputDataType  = kAliHLTDataTypeAliTreeD|kAliHLTDataOriginITS;
97     fOutputDataType = kAliHLTDataTypeClusters|kAliHLTDataOriginITS;
98     break;
99   default:
100     HLTFatal("unknown cluster finder");
101   }
102 }
103
104 AliHLTITSClusterFinderComponent::~AliHLTITSClusterFinderComponent() 
105 {
106   // see header file for class documentation
107   delete fRawReader;
108   delete fDettype;
109   delete fgeomInit;  
110   delete fSPD;
111   delete fSSD;
112 }
113
114 // Public functions to implement AliHLTComponent's interface.
115 // These functions are required for the registration process
116
117 const char* AliHLTITSClusterFinderComponent::GetComponentID()
118 {
119   // see header file for class documentation
120   switch(fModeSwitch){
121   case kClusterFinderSPD:
122     return "ITSClusterFinderSPD";
123     break;
124   case kClusterFinderSDD:        
125     return "ITSClusterFinderSDD";
126     break;
127   case kClusterFinderSSD:
128     return "ITSClusterFinderSSD";
129     break;
130   case kClusterFinderDigits:
131     return "ITSClusterFinderDigits";
132     break;
133   }
134   return "";
135 }
136
137 void AliHLTITSClusterFinderComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list) 
138 {
139   // see header file for class documentation
140   list.clear(); 
141   list.push_back( fInputDataType );
142 }
143
144 AliHLTComponentDataType AliHLTITSClusterFinderComponent::GetOutputDataType() 
145 {
146   // see header file for class documentation
147   return fOutputDataType;
148 }
149
150 void AliHLTITSClusterFinderComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier ) {
151   // see header file for class documentation
152   constBase = 0;
153   inputMultiplier = 100;
154 }
155
156 AliHLTComponent* AliHLTITSClusterFinderComponent::Spawn() {
157   // see header file for class documentation
158   return new AliHLTITSClusterFinderComponent(fModeSwitch);
159 }
160         
161 Int_t AliHLTITSClusterFinderComponent::DoInit( int argc, const char** argv ) {
162   // see header file for class documentation
163   fBenchmark.Reset();
164   fBenchmark.SetTimer(0,"total");
165   fBenchmark.SetTimer(1,"reco");
166
167   /*
168   fStatTime = 0;
169   fStatTimeAll = 0;
170   fStatTimeC = 0;
171   fStatTimeAllC = 0;
172   fStatNEv = 0;
173   */
174   
175   Int_t runNo = GetRunNo();
176   AliCDBStorage* store = AliCDBManager::Instance()->GetDefaultStorage();
177   if (!store) {
178     return NULL;
179   }
180
181   bool cdbOK = true;
182   //OCDB for SPD
183   if(store->GetLatestVersion("ITS/Calib/SPDNoisy", runNo)<0){
184     HLTError("SPDNoisy is not found in SPD/Calib");
185     cdbOK = false;
186   }
187   if(store->GetLatestVersion("ITS/Calib/SPDDead", runNo)<0){
188     HLTError("SPDDead is not found in SPD/Calib");
189     cdbOK = false;
190   }
191   if(store->GetLatestVersion("TRIGGER/SPD/PITConditions", runNo)<0){
192     HLTError("PITConditions is not found in TRIGGER/SPD");
193     cdbOK = false;
194   }
195   
196   //OCDB for SDD
197   if(store->GetLatestVersion("ITS/Calib/CalibSDD", runNo)<0){
198     HLTError("CalibSDD is not found in ITS/Calib");
199     cdbOK = false;
200   }
201   if(store->GetLatestVersion("ITS/Calib/RespSDD", runNo)<0){
202     HLTError("RespSDD is not found in ITS/Calib");
203     cdbOK = false;
204   }
205   if(store->GetLatestVersion("ITS/Calib/DriftSpeedSDD", runNo)<0){
206     HLTError("DriftSpeedSDD is not found in ITS/Calib");
207     cdbOK = false;
208   }
209   if(store->GetLatestVersion("ITS/Calib/DDLMapSDD", runNo)<0){
210     HLTError("DDLMapSDD is not found in ITS/Calib");
211     cdbOK = false;
212   }
213   if(store->GetLatestVersion("ITS/Calib/MapsTimeSDD", runNo)<0){
214     HLTError("MapsTimeSDD is not found in ITS/Calib");
215     cdbOK = false;
216   }
217
218   //OCDB for SSD
219   if(store->GetLatestVersion("ITS/Calib/NoiseSSD", runNo)<0){
220     HLTError("NoiseSSD is not found in ITS/Calib");
221     cdbOK = false;
222   }
223   if(store->GetLatestVersion("ITS/Calib/GainSSD", runNo)<0){
224     HLTError("GainSSD is not found in ITS/Calib");
225     cdbOK = false;
226   }
227   if(store->GetLatestVersion("ITS/Calib/BadChannelsSSD", runNo)<0){
228     HLTError("BadChannelsSSD is not found in ITS/Calib");
229     cdbOK = false;
230   }
231   
232   //General reconstruction
233   if(store->GetLatestVersion("GRP/CTP/Scalers", runNo)<0){
234     HLTError("Scalers is not found in GRP/CTP/");
235     cdbOK = false;
236   }
237   if(!cdbOK){return NULL;}
238
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");
244   }
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");
250   }
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");
256   }
257   else if(fModeSwitch==kClusterFinderDigits) {
258     //tR = new TTree();
259   }
260   else{
261      HLTFatal("No mode set for clusterfindercomponent");
262   }
263
264   //Removed the warning for loading default RecoParam in HLT
265   AliITSRecoParam *par = AliITSRecoParam::GetLowFluxParam();
266   AliITSReconstructor *rec = new AliITSReconstructor();
267   rec->SetRecoParam(par);
268   
269   AliCDBManager* man = AliCDBManager::Instance();
270   if (!man->IsDefaultStorageSet()){
271     HLTError("Default CDB storage has not been set !");
272     return -ENOENT;
273   }
274
275   if(AliGeomManager::GetGeometry()==NULL){
276     AliGeomManager::LoadGeometry();
277   }
278
279  //fgeomInit = new AliITSInitGeometry(kvSPD02,2);
280   fgeomInit = new AliITSInitGeometry(kvPPRasymmFMD,2);
281   //fgeomInit->InitAliITSgeom(fgeom);
282   fgeom = fgeomInit->CreateAliITSgeom();
283  
284   fNModules = fgeom->GetIndexMax();
285
286   fClusters = new TClonesArray*[fNModules]; 
287   for (Int_t iModule = 0; iModule < fNModules; iModule++) {
288     fClusters[iModule] = NULL;
289   } 
290
291   //set dettype
292   fDettype = new AliITSDetTypeRec();
293   fDettype->SetITSgeom(fgeom); 
294   fDettype->SetDefaults();
295   fDettype->SetDefaultClusterFindersV2(kTRUE); 
296
297   if ( fRawReader )
298     return -EINPROGRESS;
299
300   fRawReader = new AliRawReaderMemory();
301   fSPD = new AliHLTITSClusterFinderSPD( fDettype );
302   fSSD = new AliHLTITSClusterFinderSSD( fDettype, fRawReader );
303
304   TString arguments = "";
305   for ( int i = 0; i < argc; i++ ) {
306     if ( !arguments.IsNull() ) arguments += " ";
307     arguments += argv[i];
308   }
309
310   tD = NULL;
311   tR = NULL;
312
313   return Configure( arguments.Data() );
314 }
315
316 Int_t AliHLTITSClusterFinderComponent::DoDeinit() {
317   // see header file for class documentation
318
319   if ( fRawReader )
320     delete fRawReader;
321   fRawReader = NULL;
322
323   if ( fDettype )
324     delete fDettype;
325   fDettype = NULL;
326
327   if ( fgeomInit )
328     delete fgeomInit;
329   fgeomInit = NULL;
330   
331   delete fSPD;
332   fSPD = 0;
333
334   delete fSSD;
335   fSSD = 0;
336
337   for (Int_t iModule = 0; iModule < fNModules; iModule++) {
338     if(fClusters[iModule] != NULL){
339       fClusters[iModule]->Delete();
340       delete fClusters[iModule];
341     }
342     fClusters[iModule] = NULL;
343   } 
344   
345   fUseOfflineFinder = 0;
346
347   return 0;
348 }
349
350 int AliHLTITSClusterFinderComponent::DoEvent
351 (
352  const AliHLTComponentEventData& evtData,
353  const AliHLTComponentBlockData* /*blocks*/,
354  AliHLTComponentTriggerData& /*trigData*/,
355  AliHLTUInt8_t* outputPtr,
356  AliHLTUInt32_t& size,
357   vector<AliHLTComponentBlockData>& outputBlocks )
358 {
359  // see header file for class documentation
360   
361   AliHLTUInt32_t maxBufferSize = size;
362   size = 0; // output size
363
364   if (!IsDataEvent()) return 0;
365   
366   if ( evtData.fBlockCnt<=0 )
367     {
368       HLTDebug("no blocks in event" );
369       return 0;
370     }
371
372   fBenchmark.StartNewEvent();
373   fBenchmark.Start(0);
374   for( const AliHLTComponentBlockData *i= GetFirstInputBlock(fInputDataType); i!=NULL; i=GetNextInputBlock() ){
375     fBenchmark.AddInput(i->fSize);
376   }
377
378   Int_t ret = 0;
379
380   if(fModeSwitch==kClusterFinderDigits) {
381
382     for ( const TObject *iter = GetFirstInputObject(fInputDataType); iter != NULL; iter = GetNextInputObject() ) {  
383       tD = dynamic_cast<TTree*>(const_cast<TObject*>( iter ) );
384       if(!tD){
385         HLTFatal("No Digit Tree found");
386         return -1;
387       }
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
391       //
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
396       //
397       // Investigation: TBranch::Streamer uses a crude assignment after creating the
398       // TBranch object
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.
405       //
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);
410       tR = new TTree();
411       tR->SetDirectory(0);
412       fDettype->SetTreeAddressD(tD);
413       fDettype->MakeBranch(tR,"R");
414       fDettype->SetTreeAddressR(tR);
415       Option_t *opt="All";
416       fBenchmark.Start(1);
417       fDettype->DigitsToRecPoints(tD,tR,0,opt,kTRUE);
418       fBenchmark.Stop(1);
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);
426         }
427       }
428       
429       if(tR){
430         tR->Delete();
431       }
432
433       tD->SetDirectory(0);
434       delete dummy;
435       UInt_t nClusters=fclusters.size();
436       
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);
440         ret = -ENOSPC;      
441         break;          
442       }
443       if( nClusters>0 ){
444         fBenchmark.Start(1);
445         RecPointToSpacePoint(outputPtr,size);
446         fBenchmark.Stop(1);
447         AliHLTComponentBlockData bd;
448         FillBlockData( bd );
449         bd.fOffset = size;
450         bd.fSize = bufferSize;
451         bd.fSpecification = 0x00000000;
452         bd.fDataType = GetOutputDataType();
453         outputBlocks.push_back( bd );
454         size += bufferSize;
455         fBenchmark.AddOutput(bd.fSize);
456         fclusters.clear();      
457       }
458     }
459   }
460   else{
461       
462     // -- Loop over blocks
463     for( const AliHLTComponentBlockData* iter = GetFirstInputBlock(fInputDataType); iter != NULL; iter = GetNextInputBlock() ) {
464       
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());
470       
471       // -- Check for the correct data type
472       if ( iter->fDataType != (fInputDataType) )  
473         continue;
474       
475       // -- Get equipment ID out of specification
476       AliHLTUInt32_t spec = iter->fSpecification;
477       
478       Int_t id = fId;
479       for ( Int_t ii = 0; ii < fNddl ; ii++ ) {   //number of ddl's
480         if ( spec & 0x00000001 ) {
481           id += ii;
482           break;
483         }
484         spec = spec >> 1 ;
485       }
486       
487       // -- Set equipment ID to the raw reader
488       
489       if(!fRawReader){
490         HLTWarning("The fRawReader pointer is NULL");
491         continue;
492       }
493
494       if(!fRawReader->AddBuffer((UChar_t*) iter->fPtr, iter->fSize, id)){
495         HLTWarning("Could not add buffer");
496       }
497
498       fBenchmark.Start(1);
499
500       if(fModeSwitch==kClusterFinderSPD && !fUseOfflineFinder){ fSPD->RawdataToClusters( fRawReader, fclusters ); }
501       else if(fModeSwitch==kClusterFinderSSD && !fUseOfflineFinder){ fSSD->RawdataToClusters( fclusters ); }
502       else{
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);
511             }
512             fClusters[i]->Delete();
513             delete fClusters[i];
514           }
515           fClusters[i] = NULL;
516         }     
517       }
518       fBenchmark.Stop(1);
519       
520       fRawReader->ClearBuffers();    
521          
522       UInt_t nClusters=fclusters.size();
523       
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);
527         ret = -ENOSPC;      
528         break;          
529       }
530       if( nClusters>0 ){
531
532         RecPointToSpacePoint(outputPtr,size);
533
534         AliHLTComponentBlockData bd;
535         FillBlockData( bd );
536         bd.fOffset = size;
537         bd.fSize = bufferSize;
538         bd.fSpecification = iter->fSpecification;
539         bd.fDataType = GetOutputDataType();
540         outputBlocks.push_back( bd );
541         size += bufferSize;
542         fBenchmark.AddOutput(bd.fSize);
543         fclusters.clear();      
544       }
545       
546     } // input blocks
547   }
548
549   fBenchmark.Stop(0);
550   HLTInfo(fBenchmark.GetStatistics());
551   return ret;
552 }
553
554 int AliHLTITSClusterFinderComponent::Configure(const char* arguments)
555 {
556   
557   int iResult=0;
558   
559   if (!arguments) return iResult;
560   
561   TString allArgs=arguments;
562   TString argument;
563   
564   TObjArray* pTokens=allArgs.Tokenize(" ");
565   
566   if (pTokens) {
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");
573         continue;
574       }
575       /*
576       else if (argument.CompareTo("")==0) {
577         HLTInfo("");
578         continue;
579       }
580       */
581       else {
582         HLTError("unknown argument %s", argument.Data());
583         iResult=-EINVAL;
584         break;
585       }
586     }
587     delete pTokens;
588   }
589   
590   return iResult;
591 }
592
593 int AliHLTITSClusterFinderComponent::Reconfigure(const char* cdbEntry, const char* chainId)
594 {
595   // see header file for class documentation
596   int iResult=0;
597   
598   const char* path="HLT/ConfigITS/ClusterFinderComponent";
599   const char* defaultNotify="";
600   if (cdbEntry) {
601     path=cdbEntry;
602     defaultNotify=" (default)";
603   }
604   if (path) {
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);
607     if (pEntry) {
608       TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
609       if (pString) {
610         HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
611         iResult=Configure(pString->GetString().Data());
612       } else {
613         HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
614       }
615     } else {
616       HLTError("can not fetch object \"%s\" from CDB", path);
617     }
618   }
619   
620   return iResult;
621 }
622 void AliHLTITSClusterFinderComponent::GetOCDBObjectDescription( TMap* const targetMap)
623 {
624   // Get a list of OCDB object description.
625   if (!targetMap) return;
626   //SPD
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" ));
630   //SDD
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" ));
636   //SSD
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" ));
642 }
643
644
645 void AliHLTITSClusterFinderComponent::RecPointToSpacePoint(AliHLTUInt8_t* outputPtr,AliHLTUInt32_t& size){
646   AliHLTITSClusterData *outputClusters = reinterpret_cast<AliHLTITSClusterData*>(outputPtr + size);
647   outputClusters->fSpacePointCnt=fclusters.size();    
648   int clustIdx=0;
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);      
664     clustIdx++;
665   }
666 }