194ae0beb3f3ea224c04c82096c75e7568d9c52b
[u/mrichter/AliRoot.git] / HLT / src / AliHLTReconstructor.cxx
1 // $Id$
2
3 ///////////////////////////////////////////////////////////////////////////////
4 //                                                                           //
5 // class for HLT reconstruction                                              //
6 // <Cvetan.Cheshkov@cern.ch>                                                 //
7 // <loizides@ikf.uni-frankfurt.de>                                           //
8 ///////////////////////////////////////////////////////////////////////////////
9
10 // very ugly but it has to work fast
11 #ifdef use_reconstruction
12
13 #include <Riostream.h>
14 #include <TSystem.h>
15 #include <TArrayF.h>
16 #include <TObjString.h>
17
18 #include <AliRunLoader.h>
19 #include <AliHeader.h>
20 #include <AliGenEventHeader.h>
21 #include <AliESD.h>
22 #include <AliESDHLTtrack.h>
23
24 #include "AliHLTStandardIncludes.h"
25 #include "AliHLTLogging.h"
26 //#include "AliLevel3.h"
27 //#include "AliHLTEvaluate.h"
28 #include "AliHLTReconstructor.h"
29 #include "AliHLTTransform.h"
30 #include "AliHLTHough.h"
31 #include "AliHLTFileHandler.h"
32 #include "AliHLTTrack.h"
33 #include "AliHLTHoughTrack.h"
34 #include "AliHLTTrackArray.h"
35
36 #include "AliLog.h"
37 #include "AliRun.h"
38 #include "AliITS.h"
39 #include "AliHLTITStracker.h"
40 #include "AliHLTTPCtracker.h"
41 #include "MUON/src/AliRoot/AliHLTMUONTracker.h"
42 #include "MUON/src/AliRoot/AliHLTMUONHitReconstructor.h"
43 #include "AliRawReader.h"
44 #include "AliHLTSystem.h"
45
46 #if __GNUC__== 3
47 using namespace std;
48 #endif
49
50 const char* kHLTDefaultLibs[]= {
51   "libAliHLTUtil.so", 
52   "libAliHLTTPC.so", 
53   //  "libAliHLTSample.so",
54   "libAliHLTPHOS.so",
55   NULL
56 };
57
58 ClassImp(AliHLTReconstructor)
59
60 AliHLTReconstructor::AliHLTReconstructor()
61   : 
62   AliReconstructor(),
63   fDoHough(0),
64   fDoTracker(1),
65   fDoBench(0),
66   fDoCleanUp(0),
67   fpSystem(NULL)
68
69   //constructor
70 #ifndef use_logging
71   AliHLTLog::fgLevel=AliHLTLog::kWarning;
72 #endif
73 }
74
75 AliHLTReconstructor::AliHLTReconstructor(Bool_t doTracker, Bool_t doHough)
76   : 
77   AliReconstructor(),
78   fDoHough(doHough),
79   fDoTracker(doTracker),
80   fDoBench(0),
81   fDoCleanUp(0),
82   fpSystem(new AliHLTSystem)
83 {
84   //constructor
85 #ifndef use_logging
86   AliHLTLog::fgLevel=AliHLTLog::kWarning;
87 #endif
88 }
89
90 AliHLTReconstructor::AliHLTReconstructor(const AliHLTReconstructor&)
91   :
92   AliReconstructor(),
93   fDoHough(0),
94   fDoTracker(0),
95   fDoBench(0),
96   fDoCleanUp(0),
97   fpSystem(NULL)
98 {
99   // not a valid copy constructor
100 }
101
102 AliHLTReconstructor& AliHLTReconstructor::operator=(const AliHLTReconstructor&)
103 {
104   // not a valid assignment operator
105   fDoHough=0;
106   fDoTracker=0;
107   fDoBench=0;
108   fDoCleanUp=0;
109   fpSystem=NULL;
110   return *this;
111 }
112
113 AliHLTReconstructor::~AliHLTReconstructor()
114
115   //destructor
116   if(fDoCleanUp){
117     char name[256];
118     gSystem->Exec("rm -rf hlt");
119     sprintf(name, "rm -f confmap_*.root confmap_*.dat");
120     gSystem->Exec(name);
121     gSystem->Exec("rm -rf hough");
122     sprintf(name, "rm -f hough_*.root hough_*.dat");
123     gSystem->Exec(name);
124   }
125   if (fpSystem) {
126     delete fpSystem;
127   }
128   fpSystem=NULL;
129 }
130
131 void AliHLTReconstructor::Init(AliRunLoader* runLoader)
132 {
133   // init the reconstructor
134   if(!runLoader) {
135     AliError("Missing RunLoader! 0x0");
136     return;
137   }
138
139   if (!fpSystem) fpSystem=new AliHLTSystem;
140   if (!fpSystem) {
141     AliError("can not create AliHLTSystem object");
142     return;
143   }
144   if (fpSystem->CheckStatus(AliHLTSystem::kError)) {
145     AliError("HLT system in error state");
146     return;
147   }
148
149   TString libs("");
150   TString option = GetOption();
151   TObjArray* pTokens=option.Tokenize(" ");
152   if (pTokens) {
153     int iEntries=pTokens->GetEntries();
154     for (int i=0; i<iEntries; i++) {
155       TString token=(((TObjString*)pTokens->At(i))->GetString());
156       if (token.Contains("loglevel=")) {
157         TString param=token.ReplaceAll("loglevel=", "");
158         if (param.IsDigit()) {
159           fpSystem->SetGlobalLoggingLevel((AliHLTComponentLogSeverity)param.Atoi());
160         } else if (param.BeginsWith("0x") &&
161                    param.Replace(0,2,"",0).IsHex()) {
162           int severity=0;
163           sscanf(param.Data(),"%x", &severity);
164           fpSystem->SetGlobalLoggingLevel((AliHLTComponentLogSeverity)severity);
165         } else {
166           AliWarning("wrong parameter for option \'loglevel=\', (hex) number expected");
167         }
168       } else if (token.Contains("alilog=off")) {
169         fpSystem->SwitchAliLog(0);
170       } else if (token.BeginsWith("lib") && token.EndsWith(".so")) {
171         libs+=token;
172         libs+=" ";
173       } else {
174         AliWarning(Form("unknown option: %s", token.Data()));
175       }
176     }
177     delete pTokens;
178   }
179   
180   Bool_t bForceLibLoad=0;
181   if (bForceLibLoad=(libs.IsNull())) {
182     const char** deflib=kHLTDefaultLibs;
183     while (*deflib) {
184       libs+=*deflib++;
185       libs+=" ";
186     }
187   }
188   if ((bForceLibLoad || !fpSystem->CheckStatus(AliHLTSystem::kLibrariesLoaded)) &&
189       (fpSystem->LoadComponentLibraries(libs.Data())<0)) {
190     AliError("error while loading HLT libraries");
191     return;
192   }
193   if (!fpSystem->CheckStatus(AliHLTSystem::kReady) &&
194       (fpSystem->Configure(runLoader))<0) {
195     AliError("error during HLT system configuration");
196     return;
197   }
198 }
199
200 void AliHLTReconstructor::Reconstruct(AliRunLoader* runLoader) const
201 {
202   // reconstruction of simulated data
203   Reconstruct(runLoader, NULL);
204 }
205
206 void AliHLTReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawReader) const 
207 {
208   // reconstruction of real data if rawReader!=NULL
209   if(!runLoader) {
210     AliError("Missing RunLoader! 0x0");
211     return;
212   }
213
214   Int_t nEvents = runLoader->GetNumberOfEvents();
215   int iResult=0;
216
217   if (fpSystem) {
218     if (fpSystem->CheckStatus(AliHLTSystem::kError)) {
219       AliError("HLT system in error state");
220       return;
221     }
222     if ((iResult=fpSystem->Reconstruct(nEvents, runLoader, rawReader))>=0) {
223     }
224   }
225 }
226
227 void AliHLTReconstructor::ReconstructWithConformalMapping(AliRunLoader* runLoader,Int_t iEvent) const
228 {
229   // reconstruct with conformal mapper
230   /*
231   AliLevel3 *fHLT = new AliLevel3(runLoader);
232   fHLT->Init("./", AliLevel3::kRunLoader, 1);
233
234   Int_t phiSegments = 50;
235   Int_t etaSegments = 100;
236   Int_t trackletlength = 3;
237   Int_t tracklength = 10;
238   Int_t rowscopetracklet = 2;
239   Int_t rowscopetrack = 10;
240   Double_t minPtFit = 0;
241   Double_t maxangle = 0.1745;
242   Double_t goodDist = 5;
243   Double_t maxphi = 0.1;
244   Double_t maxeta = 0.1;
245   Double_t hitChi2Cut = 20;
246   Double_t goodHitChi2 = 5;
247   Double_t trackChi2Cut = 10;
248   Double_t xyerror = -1;
249   Double_t zerror =  -1;
250   
251   fHLT->SetClusterFinderParam(xyerror,zerror,kTRUE);
252   fHLT->SetTrackerParam(phiSegments, etaSegments, 
253                         trackletlength, tracklength,
254                         rowscopetracklet, rowscopetrack,
255                         minPtFit, maxangle, goodDist, hitChi2Cut,
256                         goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kTRUE);
257   fHLT->SetTrackerParam(phiSegments, etaSegments, 
258                         trackletlength, tracklength,
259                         rowscopetracklet, rowscopetrack,
260                         minPtFit, maxangle, goodDist, hitChi2Cut,
261                         goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kFALSE);
262   fHLT->SetMergerParameters(2,3,0.003,0.1,0.05);
263   fHLT->DoMc();
264   fHLT->DoNonVertexTracking(); // 2 tracking passes, last without vertex contraint.
265   fHLT->WriteFiles("./hlt/");  
266   fHLT->ProcessEvent(0, 35, iEvent);
267   if(fDoBench){
268     char filename[256];
269     sprintf(filename, "confmap_%d",iEvent);
270     fHLT->DoBench(filename);
271   }
272
273   delete fHLT;
274   */
275 }
276
277 void AliHLTReconstructor::ReconstructWithHoughTransform(AliRunLoader* runLoader,Int_t iEvent) const
278 {
279   //reconstruct with hough
280   //not used anymore, Hough tracking is moved out of the local
281   //reconstruction chain
282   Float_t ptmin = 0.1*AliHLTTransform::GetSolenoidField();
283
284   Float_t zvertex = 0;
285   TArrayF mcVertex(3); 
286   AliHeader * header = runLoader->GetHeader();
287   if (header) {
288     AliGenEventHeader * genHeader = header->GenEventHeader();
289     if (genHeader) genHeader->PrimaryVertex(mcVertex);
290   }
291   zvertex = mcVertex[2];
292
293   AliInfo(Form("Hough Transform will run with ptmin=%f and zvertex=%f", ptmin, zvertex));
294
295   AliHLTHough *hough = new AliHLTHough();
296     
297   hough->SetThreshold(4);
298   hough->CalcTransformerParams(ptmin);
299   hough->SetPeakThreshold(70,-1);
300   hough->SetRunLoader(runLoader);
301   hough->Init("./", kFALSE, 100, kFALSE,4,0,0,zvertex);
302   hough->SetAddHistograms();
303
304   for(Int_t slice=0; slice<=35; slice++)
305     {
306       //cout<<"Processing slice "<<slice<<endl;
307       hough->ReadData(slice,iEvent);
308       hough->Transform();
309       hough->AddAllHistogramsRows();
310       hough->FindTrackCandidatesRow();
311       //hough->WriteTracks(slice,"./hough");
312       hough->AddTracks();
313     }
314   hough->WriteTracks("./hough");
315   
316   if(fDoBench){
317     char filename[256];
318     sprintf(filename, "hough_%d",iEvent);
319     hough->DoBench(filename);
320   }
321   delete hough;
322 }
323
324 void AliHLTReconstructor::FillESD(AliRunLoader* runLoader, 
325                                   AliESD* esd) const
326 {
327   //fill the esd file with found tracks
328   if(!runLoader) {
329     AliError("Missing RunLoader! 0x0");
330     return;
331   }
332   Int_t iEvent = runLoader->GetEventNumber();
333   if (fpSystem) {
334     if (fpSystem->CheckStatus(AliHLTSystem::kError)) {
335       AliError("HLT system in error state");
336       return;
337     }
338     if (!fpSystem->CheckStatus(AliHLTSystem::kReady)) {
339       AliError("HLT system in wrong state");
340       return;
341     }
342     fpSystem->FillESD(iEvent, runLoader, esd);
343   }
344   /*
345   if(fDoTracker) FillESDforConformalMapping(esd,iEvent);
346   if(fDoHough) FillESDforHoughTransform(esd,iEvent);
347   */
348 }
349
350 void AliHLTReconstructor::FillESDforConformalMapping(AliESD* esd,Int_t iEvent) const
351 {
352   //fill esd with tracks from conformal mapping
353   /*
354   Int_t slicerange[2]={0,35};
355   Int_t good = (int)(0.4*AliHLTTransform::GetNRows());
356   Int_t nclusters = (int)(0.4*AliHLTTransform::GetNRows());
357   Int_t nminpointsontracks = (int)(0.3*AliHLTTransform::GetNRows());
358   Float_t ptmin = 0.;
359   Float_t ptmax = 0.;
360   Float_t maxfalseratio = 0.1;
361   
362   AliHLTEvaluate *fHLTEval = new AliHLTEvaluate("./hlt",nclusters,good,ptmin,ptmax,slicerange);
363   fHLTEval->SetMaxFalseClusters(maxfalseratio);
364   fHLTEval->LoadData(iEvent,kTRUE);
365   fHLTEval->AssignPIDs();
366   fHLTEval->AssignIDs();
367   AliHLTTrackArray *fTracks = fHLTEval->GetTracks();
368   if(!fTracks){
369     delete fHLTEval;
370     return;
371   }
372   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
373     {
374       AliHLTTrack *tpt = (AliHLTTrack *)fTracks->GetCheckedTrack(i);
375       if(!tpt) continue; 
376       if(tpt->GetNumberOfPoints() < nminpointsontracks) continue;
377       
378       AliESDHLTtrack *esdtrack = new AliESDHLTtrack() ; 
379
380       esdtrack->SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
381       esdtrack->SetNHits(tpt->GetNHits());
382       esdtrack->SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ());
383       esdtrack->SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ());
384       esdtrack->SetPt(tpt->GetPt());
385       esdtrack->SetPsi(tpt->GetPsi());
386       esdtrack->SetTgl(tpt->GetTgl());
387       esdtrack->SetCharge(tpt->GetCharge());
388       esdtrack->SetMCid(tpt->GetMCid());
389       esdtrack->SetSector(tpt->GetSector());
390       esdtrack->SetPID(tpt->GetPID());
391       esdtrack->ComesFromMainVertex(tpt->ComesFromMainVertex());
392
393       esd->AddHLTConfMapTrack(esdtrack);
394       delete esdtrack;
395     }
396   delete fHLTEval;
397   */
398 }
399
400 void AliHLTReconstructor::FillESDforHoughTransform(AliESD* esd,Int_t iEvent) const
401 {
402   //fill esd with tracks from hough
403   char filename[256];
404   sprintf(filename,"./hough/tracks_%d.raw",iEvent);
405   
406   AliHLTFileHandler *tfile = new AliHLTFileHandler();
407   if(!tfile->SetBinaryInput(filename)){
408     AliError(Form("Missing file %s", filename));
409     return;
410   }
411   
412   AliHLTTrackArray *fTracks = new AliHLTTrackArray("AliHLTHoughTrack");
413   tfile->Binary2TrackArray(fTracks);
414   tfile->CloseBinaryInput();
415   delete tfile;
416   if(!fTracks) return; 
417   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
418     {
419       AliHLTHoughTrack *tpt = (AliHLTHoughTrack *)fTracks->GetCheckedTrack(i);
420       if(!tpt) continue; 
421       
422       AliESDHLTtrack *esdtrack = new AliESDHLTtrack() ; 
423
424       esdtrack->SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
425       esdtrack->SetNHits(tpt->GetNHits());
426       esdtrack->SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ());
427       esdtrack->SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ());
428       esdtrack->SetPt(tpt->GetPt());
429       esdtrack->SetPsi(tpt->GetPsi());
430       esdtrack->SetTgl(tpt->GetTgl());
431       esdtrack->SetCharge(tpt->GetCharge());
432       esdtrack->SetMCid(tpt->GetMCid());
433       esdtrack->SetWeight(tpt->GetWeight());
434       esdtrack->SetSector(tpt->GetSector());
435       esdtrack->SetBinXY(tpt->GetBinX(),tpt->GetBinY(),tpt->GetSizeX(),tpt->GetSizeY());
436       esdtrack->SetPID(tpt->GetPID());
437       esdtrack->ComesFromMainVertex(tpt->ComesFromMainVertex());
438
439       esd->AddHLTHoughTrack(esdtrack);
440       delete esdtrack;
441     }
442
443   delete fTracks;
444 }
445
446 AliTracker* AliHLTReconstructor::CreateTracker(AliRunLoader* runLoader) const
447 {
448   //Create HLT trackers for TPC and ITS
449
450   TString opt = GetOption();
451   if(!opt.CompareTo("TPC")) {
452     // Create Hough tracker for TPC
453     return new AliHLTTPCtracker(runLoader);
454   }
455   if(!opt.CompareTo("ITS")) {
456     // Create ITS tracker
457     return new AliHLTITStracker(0);
458   }
459   if(!opt.CompareTo("MUON")) {
460     return new AliHLTMUONTracker(runLoader);
461   }
462
463   return NULL;
464 }
465
466 void AliHLTReconstructor::FillDHLTRecPoint(AliRawReader* rawReader, Int_t nofEvent, Int_t dcCut = 0) const
467 {
468   // Hit recontruction for dimuon-HLT
469   AliHLTMUONHitReconstructor kdHLTRec(rawReader);
470
471   Int_t iEvent = 0 ;
472   kdHLTRec.SetDCCut(dcCut);
473   TString lookupTablePath = getenv("ALICE_ROOT");
474   lookupTablePath += "/HLT/MUON/src/AliRoot/Lut";
475   kdHLTRec.Init(lookupTablePath.Data(),lookupTablePath.Data());
476
477    while(rawReader->NextEvent() && iEvent < nofEvent){
478      AliInfo(Form("Event : %d",iEvent));
479      kdHLTRec.WriteDHLTRecHits(iEvent);  
480      iEvent++;
481    }
482
483 }
484
485 #endif