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