Added two pass tracking for confmapper.
[u/mrichter/AliRoot.git] / HLT / src / AliHLTReconstructor.cxx
1 // $Id$
2
3 ///////////////////////////////////////////////////////////////////////////////
4 //                                                                           //
5 // class for HLT reconstruction                                              //
6 // <Cvetan.Cheshkov@cern.ch>                                                 //
7 ///////////////////////////////////////////////////////////////////////////////
8
9 // very ugly but it has to work fast
10 #ifdef use_reconstruction
11
12 #include <Riostream.h>
13 #include <TSystem.h>
14 #include <TArrayF.h>
15
16 #include "AliL3StandardIncludes.h"
17 #include "AliL3Logging.h"
18 #include "AliLevel3.h"
19 #include "AliL3Evaluate.h"
20 #include "AliHLTReconstructor.h"
21 #include "AliL3Transform.h"
22 #include "AliL3Hough.h"
23 #include "AliL3FileHandler.h"
24 #include "AliL3Track.h"
25 #include "AliL3HoughTrack.h"
26 #include "AliL3TrackArray.h"
27 #include "AliRunLoader.h"
28 #include "AliHeader.h"
29 #include "AliGenEventHeader.h"
30 #include "AliESD.h"
31 #include "AliESDHLTtrack.h"
32
33 #if __GNUC__== 3
34 using namespace std;
35 #endif
36
37
38 ClassImp(AliHLTReconstructor)
39
40 void AliHLTReconstructor::Reconstruct(AliRunLoader* runLoader) const
41 {
42   if(!runLoader) {
43     LOG(AliL3Log::kFatal,"AliHLTReconstructor::Reconstruct","RunLoader")
44       <<" Missing RunLoader! 0x0"<<ENDLOG;
45     return;
46   }
47   gSystem->Exec("rm -rf hlt");
48   gSystem->MakeDirectory("hlt");
49   gSystem->Exec("rm -rf hough");
50   gSystem->MakeDirectory("hough");
51
52   Bool_t isinit=AliL3Transform::Init(runLoader);
53   if(!isinit){
54     LOG(AliL3Log::kError,"AliHLTReconstructor::Reconstruct","Transformer")
55      << "Could not create transform settings, please check log for error messages!" << ENDLOG;
56     return;
57   }
58
59   Int_t nEvents = runLoader->GetNumberOfEvents();
60
61   for(Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
62     runLoader->GetEvent(iEvent);
63
64     if(fDoTracker) ReconstructWithConformalMapping(runLoader,iEvent);
65     if(fDoHough) ReconstructWithHoughTransform(runLoader,iEvent);
66   }
67 }
68
69 void AliHLTReconstructor::ReconstructWithConformalMapping(AliRunLoader* runLoader,Int_t iEvent) const
70 {
71   AliLevel3 *fHLT = new AliLevel3(runLoader);
72   fHLT->Init("./", AliLevel3::kRunLoader, 1);
73
74   Int_t phiSegments = 50;
75   Int_t etaSegments = 100;
76   Int_t trackletlength = 3;
77   Int_t tracklength = 10;
78   Int_t rowscopetracklet = 2;
79   Int_t rowscopetrack = 10;
80   Double_t minPtFit = 0;
81   Double_t maxangle = 0.1745;
82   Double_t goodDist = 5;
83   Double_t maxphi = 0.1;
84   Double_t maxeta = 0.1;
85   Double_t hitChi2Cut = 20;
86   Double_t goodHitChi2 = 5;
87   Double_t trackChi2Cut = 10;
88   Double_t xyerror = -1;
89   Double_t zerror =  -1;
90   
91   fHLT->SetClusterFinderParam(xyerror,zerror,kTRUE);
92   fHLT->SetTrackerParam(phiSegments, etaSegments, 
93                         trackletlength, tracklength,
94                         rowscopetracklet, rowscopetrack,
95                         minPtFit, maxangle, goodDist, hitChi2Cut,
96                         goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kTRUE);
97   fHLT->SetTrackerParam(phiSegments, etaSegments, 
98                         trackletlength, tracklength,
99                         rowscopetracklet, rowscopetrack,
100                         minPtFit, maxangle, goodDist, hitChi2Cut,
101                         goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kFALSE);
102   fHLT->SetMergerParameters(2,3,0.003,0.1,0.05);
103   fHLT->DoMc();
104   fHLT->DoNonVertexTracking(); /*2 tracking passes, last without vertex contraint.*/
105   fHLT->WriteFiles("./hlt/");  
106   fHLT->ProcessEvent(0, 35, iEvent);
107
108   char filename[256];
109   sprintf(filename, "confmap_%d",iEvent);
110   fHLT->DoBench(filename);
111
112   delete fHLT;
113 }
114
115 void AliHLTReconstructor::ReconstructWithHoughTransform(AliRunLoader* runLoader,Int_t iEvent) const
116 {
117   Float_t ptmin = 0.1*AliL3Transform::GetSolenoidField();
118
119   Float_t zvertex = 0;
120   TArrayF mcVertex(3); 
121   runLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
122   zvertex = mcVertex[2];
123
124   LOG(AliL3Log::kInformational,"AliHLTReconstructor::Reconstruct","HoughTransform")
125     <<" Hough Transform will run with ptmin="<<ptmin<<" and zvertex="<<zvertex<<ENDLOG;
126
127   AliL3Hough *hough = new AliL3Hough();
128     
129   hough->SetThreshold(4);
130   hough->SetTransformerParams(76,140,ptmin,-1);
131   hough->SetPeakThreshold(70,-1);
132   hough->SetRunLoader(runLoader);
133   hough->Init("./", kFALSE, 100, kFALSE,4,0,0,zvertex);
134   hough->SetAddHistograms();
135
136   for(Int_t slice=0; slice<=35; slice++)
137     {
138       //cout<<"Processing slice "<<slice<<endl;
139       hough->ReadData(slice,iEvent);
140       hough->Transform();
141       hough->AddAllHistogramsRows();
142       hough->FindTrackCandidatesRow();
143       //hough->WriteTracks(slice,"./hough");
144       hough->AddTracks();
145     }
146   hough->WriteTracks("./hough");
147   
148   char filename[256];
149   sprintf(filename, "hough_%d",iEvent);
150   hough->DoBench(filename);
151
152   delete hough;
153 }
154
155 void AliHLTReconstructor::FillESD(AliRunLoader* runLoader, 
156                                   AliESD* esd) const
157 {
158   Int_t iEvent = runLoader->GetEventNumber();
159
160   if(fDoTracker) FillESDforConformalMapping(esd,iEvent);
161   if(fDoHough) FillESDforHoughTransform(esd,iEvent);
162
163 #if 0
164   char name[256];
165   gSystem->Exec("rm -rf hlt");
166   sprintf(name, "rm -f confmap_%d.root confmap_%d.dat",iEvent,iEvent);
167   gSystem->Exec(name);
168   gSystem->Exec("rm -rf hough");
169   sprintf(name, "rm -f hough_%d.root hough_%d.dat",iEvent,iEvent);
170   gSystem->Exec(name);
171 #endif
172 }
173
174 void AliHLTReconstructor::FillESDforConformalMapping(AliESD* esd,Int_t iEvent) const
175 {
176   //Assign MC labels for found tracks
177   Int_t slicerange[2]={0,35};
178   Int_t good = (int)(0.4*AliL3Transform::GetNRows());
179   Int_t nclusters = (int)(0.4*AliL3Transform::GetNRows());
180   Float_t ptmin = 0.;
181   Float_t ptmax = 0.;
182   Float_t maxfalseratio = 0.1;
183   
184   AliL3Evaluate *fHLTEval = new AliL3Evaluate("./hlt",nclusters,good,ptmin,ptmax,slicerange);
185   fHLTEval->SetMaxFalseClusters(maxfalseratio);
186
187   fHLTEval->LoadData(iEvent,kTRUE);
188   fHLTEval->AssignPIDs();
189   fHLTEval->AssignIDs();
190   AliL3TrackArray *fTracks = fHLTEval->GetTracks();
191   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
192     {
193       AliL3Track *tpt = (AliL3Track *)fTracks->GetCheckedTrack(i);
194       if(!tpt) continue; 
195       
196       AliESDHLTtrack *esdtrack = new AliESDHLTtrack() ; 
197
198       esdtrack->SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
199       esdtrack->SetNHits(tpt->GetNHits());
200       esdtrack->SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ());
201       esdtrack->SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ());
202       esdtrack->SetPt(tpt->GetPt());
203       esdtrack->SetPsi(tpt->GetPsi());
204       esdtrack->SetTgl(tpt->GetTgl());
205       esdtrack->SetCharge(tpt->GetCharge());
206       esdtrack->SetMCid(tpt->GetMCid());
207       esdtrack->SetSector(tpt->GetSector());
208       esdtrack->SetPID(tpt->GetPID());
209       esdtrack->ComesFromMainVertex(tpt->ComesFromMainVertex());
210
211       esd->AddHLTConfMapTrack(esdtrack);
212       delete esdtrack;
213     }
214   delete fHLTEval;
215 }
216
217 void AliHLTReconstructor::FillESDforHoughTransform(AliESD* esd,Int_t iEvent) const
218 {
219   char filename[256];
220   sprintf(filename,"./hough/tracks_%d.raw",iEvent);
221   
222   AliL3FileHandler *tfile = new AliL3FileHandler();
223   if(!tfile->SetBinaryInput(filename)){
224     Error("FillESD","Inputfile ",filename," does not exist");
225     return;
226   }
227   
228   AliL3TrackArray *fTracks = new AliL3TrackArray("AliL3HoughTrack");
229   tfile->Binary2TrackArray(fTracks);
230   tfile->CloseBinaryInput();
231   delete tfile;
232   
233   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
234     {
235       AliL3HoughTrack *tpt = (AliL3HoughTrack *)fTracks->GetCheckedTrack(i);
236       if(!tpt) continue; 
237       
238       AliESDHLTtrack *esdtrack = new AliESDHLTtrack() ; 
239
240       esdtrack->SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow());
241       esdtrack->SetNHits(tpt->GetNHits());
242       esdtrack->SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ());
243       esdtrack->SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ());
244       esdtrack->SetPt(tpt->GetPt());
245       esdtrack->SetPsi(tpt->GetPsi());
246       esdtrack->SetTgl(tpt->GetTgl());
247       esdtrack->SetCharge(tpt->GetCharge());
248       esdtrack->SetMCid(tpt->GetMCid());
249       esdtrack->SetWeight(tpt->GetWeight());
250       esdtrack->SetSector(tpt->GetSector());
251       esdtrack->SetBinXY(tpt->GetBinX(),tpt->GetBinY(),tpt->GetSizeX(),tpt->GetSizeY());
252       esdtrack->SetPID(tpt->GetPID());
253       esdtrack->ComesFromMainVertex(tpt->ComesFromMainVertex());
254
255       esd->AddHLTHoughTrack(esdtrack);
256       delete esdtrack;
257     }
258
259   delete fTracks;
260 }
261
262 #endif