]>
Commit | Line | Data |
---|---|---|
31f19f09 | 1 | // $Id$ |
2 | ||
4aa41877 | 3 | #include "AliHLTStandardIncludes.h" |
4 | #include "AliHLTRootTypes.h" | |
3e87ef69 | 5 | #include <sys/time.h> |
31f19f09 | 6 | #include <TNtuple.h> |
7 | #include <TTimer.h> | |
4aa41877 | 8 | #include "AliHLTTrack.h" |
9 | #include "AliHLTKalmanTrack.h" | |
10 | #include "AliHLTBenchmark.h" | |
11 | #include "AliHLTMemHandler.h" | |
12 | #include "AliHLTFileHandler.h" | |
13 | #include "AliHLTDataHandler.h" | |
14 | #include "AliHLTTransform.h" | |
15 | #include "AliHLTSpacePointData.h" | |
16 | #include "AliHLTDigitData.h" | |
17 | #include "AliHLTLogging.h" | |
18 | #include "AliHLTTrackArray.h" | |
19 | #include "AliHLTTrackSegmentData.h" | |
20 | #include "AliHLTInterMerger.h" | |
21 | #include "AliHLTTrackMerger.h" | |
22 | #include "AliHLTKalman.h" | |
de3c3890 | 23 | |
3e87ef69 | 24 | /* |
4aa41877 | 25 | AliHLTKalman |
3e87ef69 | 26 | */ |
27 | ||
4aa41877 | 28 | ClassImp(AliHLTKalman) |
3e87ef69 | 29 | |
4aa41877 | 30 | AliHLTKalman::AliHLTKalman(Char_t *datapath, Int_t *slice, Int_t min_clusters = 0){ |
3e87ef69 | 31 | // Constructor |
3e87ef69 | 32 | if (slice) |
33 | { | |
34 | fMinSlice = slice[0]; | |
35 | fMaxSlice = slice[1]; | |
36 | } | |
37 | else | |
38 | { | |
39 | fMinSlice = 0; | |
40 | fMaxSlice = 35; | |
41 | } | |
42 | ||
43 | sprintf(fPath,"%s",datapath); | |
44 | fMinPointsOnTrack = 0; //min_clusters; | |
45 | fTracks = 0; | |
46 | fKalmanTracks = 0; | |
47 | ||
48 | // NB! fNrow under only for single-patch, must include other possibilities | |
49 | // later on. ?? Maybe better also to put it in an Init-function | |
50 | fRow[0][0] = 0; | |
4aa41877 | 51 | fRow[0][1] = AliHLTTransform::GetLastRow(-1); |
0a86fbb7 | 52 | fWriteOut = kTRUE; |
3e87ef69 | 53 | fBenchmark = 0; |
b2a02bce | 54 | |
3e87ef69 | 55 | } |
56 | ||
4aa41877 | 57 | AliHLTKalman::~AliHLTKalman() |
3e87ef69 | 58 | { |
59 | // Destructor | |
60 | if (fBenchmark) delete fBenchmark; | |
61 | } | |
62 | ||
4aa41877 | 63 | void AliHLTKalman::Init() |
3e87ef69 | 64 | { |
4aa41877 | 65 | fBenchmark = new AliHLTBenchmark(); |
3e87ef69 | 66 | } |
67 | ||
4aa41877 | 68 | void AliHLTKalman::LoadTracks(Int_t event, Bool_t sp) |
3e87ef69 | 69 | { |
0a86fbb7 | 70 | // Load space points and tracks from conformal tracker |
3e87ef69 | 71 | // Must also be possible to take seeds (and clusters) from Hough-transform?? |
72 | ||
3e87ef69 | 73 | Double_t initTime,cpuTime; |
74 | initTime = GetCpuTime(); | |
75 | fBenchmark->Start("Load tracks"); | |
76 | ||
0a86fbb7 | 77 | // Load space points |
3e87ef69 | 78 | Char_t fname[1024]; |
4aa41877 | 79 | AliHLTFileHandler *clusterfile[36][6]; |
3e87ef69 | 80 | for(Int_t s=fMinSlice; s<=fMaxSlice; s++) |
81 | { | |
4aa41877 | 82 | for(Int_t p=0; p<AliHLTTransform::GetNPatches(); p++) |
3e87ef69 | 83 | { |
84 | Int_t patch; | |
85 | if(sp==kTRUE) | |
86 | patch=-1; | |
87 | else | |
88 | patch=p; | |
89 | ||
90 | fClusters[s][p] = 0; | |
4aa41877 | 91 | clusterfile[s][p] = new AliHLTFileHandler(); |
3e87ef69 | 92 | if(event<0) |
93 | sprintf(fname,"%s/points_%d_%d.raw",fPath,s,patch); | |
94 | else | |
95 | sprintf(fname,"%s/points_%d_%d_%d.raw",fPath,event,s,patch); | |
96 | if(!clusterfile[s][p]->SetBinaryInput(fname)) | |
97 | { | |
4aa41877 | 98 | LOG(AliHLTLog::kError,"AliHLTEvaluation::Setup","File Open") |
3e87ef69 | 99 | <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; |
100 | delete clusterfile[s][p]; | |
101 | clusterfile[s][p] = 0; | |
102 | continue; | |
103 | } | |
4aa41877 | 104 | fClusters[s][p] = (AliHLTSpacePointData*)clusterfile[s][p]->Allocate(); |
3e87ef69 | 105 | clusterfile[s][p]->Binary2Memory(fNcl[s][p],fClusters[s][p]); |
106 | clusterfile[s][p]->CloseBinaryInput(); | |
107 | if(sp==kTRUE) | |
108 | break; | |
109 | delete clusterfile[s][p]; | |
110 | } | |
111 | } | |
112 | ||
113 | // Load tracks | |
de3c3890 | 114 | //sprintf(fname,"%s/kalmantracks_%d.raw",fPath,event); |
3e87ef69 | 115 | sprintf(fname,"%s/tracks_%d.raw",fPath,event); |
4aa41877 | 116 | AliHLTFileHandler *tfile = new AliHLTFileHandler(); |
3e87ef69 | 117 | if(!tfile->SetBinaryInput(fname)){ |
4aa41877 | 118 | LOG(AliHLTLog::kError,"AliHLTEvaluation::Setup","File Open") |
3e87ef69 | 119 | <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; |
120 | return; | |
121 | } | |
122 | if(fTracks) | |
123 | delete fTracks; | |
4aa41877 | 124 | fTracks = new AliHLTTrackArray(); |
3e87ef69 | 125 | tfile->Binary2TrackArray(fTracks); |
126 | tfile->CloseBinaryInput(); | |
127 | delete tfile; | |
3e87ef69 | 128 | fBenchmark->Stop("Load tracks"); |
129 | cpuTime = GetCpuTime() - initTime; | |
4aa41877 | 130 | LOG(AliHLTLog::kInformational,"AliHLTKalman::LoadTracks()","Timing") |
3e87ef69 | 131 | <<"Loading tracks in "<<cpuTime*1000<<" ms"<<ENDLOG; |
132 | ||
133 | } | |
134 | ||
4aa41877 | 135 | void AliHLTKalman::ProcessTracks() |
3e87ef69 | 136 | { |
0a86fbb7 | 137 | // Run the Kalman filter algorithm on the loaded tracks. |
138 | // If the track is OK, the loaded track is saved in file kalmantracks_0.raw | |
139 | // The kalman filter variables (that is the state vector, covariance matrix | |
140 | // and chi2) is written to root-file kalmantracks.root. | |
b2a02bce | 141 | // Should be extended to correct the track parameters of the loaded track |
142 | // or save the kalmantrack instead of the loaded track.?? | |
143 | UInt_t fEvent = 0; | |
3e87ef69 | 144 | |
145 | Double_t initTime,cpuTime; | |
146 | initTime = GetCpuTime(); | |
147 | ||
148 | fBenchmark->Start("Process tracks"); | |
149 | ||
3e87ef69 | 150 | fTracks->QSort(); |
151 | ||
4aa41877 | 152 | fKalmanTracks = new AliHLTTrackArray(); |
3e87ef69 | 153 | |
0bd0c1ef | 154 | // Make a ntuple to store state vector, covariance matrix and chisquare |
0a86fbb7 | 155 | // Will eventually not need a TTree?? |
0bd0c1ef | 156 | TNtuple *kalmanTree = new TNtuple("kalmanTree","kalmantracks","x0:x1:x2:x3:x4:c0:c1:c2:c3:c4:c5:c6:c7:c8:c9:c10:c11:c12:c13:c14:chisq"); |
3e87ef69 | 157 | Float_t meas[21]; |
b2a02bce | 158 | |
3e87ef69 | 159 | // Go through the tracks from conformal or hough tracker |
160 | for (Int_t iTrack = 0; iTrack < fTracks->GetNTracks(); iTrack++) | |
161 | { | |
b2a02bce | 162 | /*Double_t initTime,cpuTime; |
163 | initTime = GetCpuTime(); | |
164 | fBenchmark->Start("Process tracks");*/ | |
0a86fbb7 | 165 | |
4aa41877 | 166 | AliHLTTrack *track = (AliHLTTrack*)fTracks->GetCheckedTrack(iTrack); |
3e87ef69 | 167 | if (!track) continue; |
168 | if (track->GetNumberOfPoints() < fMinPointsOnTrack) continue; | |
169 | ||
4aa41877 | 170 | AliHLTKalmanTrack *kalmantrack = new AliHLTKalmanTrack(); |
3e87ef69 | 171 | |
b2a02bce | 172 | Bool_t save = kTRUE; |
0bd0c1ef | 173 | |
de3c3890 | 174 | /*if (InitKalmanTrack(kalmantrack, track) == 0) |
b2a02bce | 175 | { |
0bd0c1ef | 176 | save = kFALSE; |
177 | continue; | |
de3c3890 | 178 | }*/ |
0bd0c1ef | 179 | |
de3c3890 | 180 | if (MakeKalmanSeed(kalmantrack,track) ==0) |
181 | { | |
182 | save = kFALSE; | |
183 | continue; | |
184 | } | |
3e87ef69 | 185 | |
186 | if (Propagate(kalmantrack, track) == 0) | |
187 | { | |
188 | save = kFALSE; | |
189 | } | |
de3c3890 | 190 | |
191 | if (save) {// cout << track->GetPt() << endl; | |
3e87ef69 | 192 | Float_t x[5]; |
193 | kalmantrack->GetStateVector(x); | |
194 | Float_t c[15]; | |
195 | kalmantrack->GetCovariance(c); | |
196 | Float_t chisq = kalmantrack->GetChisq(); | |
197 | meas[0] = x[0]; | |
198 | meas[1] = x[1]; | |
199 | meas[2] = x[2]; | |
200 | meas[3] = x[3]; | |
201 | meas[4] = x[4]; | |
202 | meas[5] = c[0]; | |
203 | meas[6] = c[1]; | |
204 | meas[7] = c[2]; | |
205 | meas[8] = c[3]; | |
206 | meas[9] = c[4]; | |
207 | meas[10] = c[5]; | |
208 | meas[11] = c[6]; | |
209 | meas[12] = c[7]; | |
210 | meas[13] = c[8]; | |
211 | meas[14] = c[9]; | |
212 | meas[15] = c[10]; | |
213 | meas[16] = c[11]; | |
214 | meas[17] = c[12]; | |
215 | meas[18] = c[13]; | |
216 | meas[19] = c[14]; | |
217 | meas[20] = chisq; | |
0a86fbb7 | 218 | |
3e87ef69 | 219 | // Add the track to the trackarray |
4aa41877 | 220 | AliHLTTrack *outtrack = (AliHLTTrack*)fKalmanTracks->NextTrack(); |
0a86fbb7 | 221 | outtrack->Set(track); |
de3c3890 | 222 | // SET THE PARAMETERS ACCORDING TO KALMAN FILTER |
223 | outtrack->SetTgl(x[3]); | |
224 | // The factor 2 in the expression for Pt is not included in the similar offline expression. However | |
225 | // it should be like this if I use a factor 1/2 in the calculation of par4?? | |
4aa41877 | 226 | //outtrack->SetPt(1/(2*TMath::Abs(1e-9*TMath::Abs(x[4])/x[4] + x[4])/(0.0029980*AliHLTTransform::GetBField()))); |
227 | //outtrack->SetPt(1/(TMath::Abs(x[4])/0.0029980*AliHLTTransform::GetBField())); | |
de3c3890 | 228 | outtrack->SetPsi(x[2]); |
229 | //outtrack->Set(track); | |
0a86fbb7 | 230 | |
3e87ef69 | 231 | // Fill the ntuple with the state vector, covariance matrix and |
0bd0c1ef | 232 | // chisquare |
233 | kalmanTree->Fill(meas); | |
3e87ef69 | 234 | } |
235 | ||
236 | delete track; | |
237 | delete kalmantrack; | |
0a86fbb7 | 238 | |
b2a02bce | 239 | /*fBenchmark->Stop("Process tracks"); |
240 | cpuTime = GetCpuTime() - initTime; | |
4aa41877 | 241 | LOG(AliHLTLog::kInformational,"AliHLTKalman::ProcessTracks()","Timing") |
b2a02bce | 242 | <<"Processed track "<<iTrack<<" in "<<cpuTime*1000<<" ms"<<ENDLOG;*/ |
3e87ef69 | 243 | } |
244 | ||
245 | fBenchmark->Stop("Process tracks"); | |
246 | cpuTime = GetCpuTime() - initTime; | |
4aa41877 | 247 | LOG(AliHLTLog::kInformational,"AliHLTKalman::ProcessTracks()","Timing") |
3e87ef69 | 248 | <<"Process tracks in "<<cpuTime*1000<<" ms"<<ENDLOG; |
0a86fbb7 | 249 | |
3e87ef69 | 250 | if (fWriteOut) |
251 | { | |
252 | Char_t tname[80]; | |
253 | sprintf(tname,"%s/kalmantracks_%d.raw",fWriteOutPath,fEvent); | |
4aa41877 | 254 | AliHLTMemHandler *mem = new AliHLTMemHandler(); |
3e87ef69 | 255 | mem->SetBinaryOutput(tname); |
256 | mem->TrackArray2Binary(fKalmanTracks); | |
257 | mem->CloseBinaryOutput(); | |
258 | delete mem; | |
259 | } | |
b2a02bce | 260 | |
0a86fbb7 | 261 | // This will be removed?? |
262 | TFile *out = new TFile("kalmantracks.root","recreate"); | |
3e87ef69 | 263 | kalmanTree->Write(); |
264 | out->Close(); | |
265 | ||
266 | delete kalmanTree; | |
0a86fbb7 | 267 | |
3e87ef69 | 268 | } |
269 | ||
4aa41877 | 270 | Int_t AliHLTKalman::MakeKalmanSeed(AliHLTKalmanTrack *kalmantrack, AliHLTTrack *track) |
de3c3890 | 271 | { |
272 | Int_t num_of_clusters = track->GetNumberOfPoints(); | |
273 | ||
274 | UInt_t *hitnum = track->GetHitNumbers(); | |
275 | UInt_t id; | |
276 | ||
277 | id = hitnum[0]; | |
278 | Int_t slice0 = (id>>25) & 0x7f; | |
279 | Int_t patch0 = (id>>22) & 0x7; | |
280 | UInt_t pos0 = id&0x3fffff; | |
4aa41877 | 281 | AliHLTSpacePointData *points0 = fClusters[slice0][patch0]; |
de3c3890 | 282 | |
283 | id = hitnum[Int_t(num_of_clusters/2)]; | |
284 | Int_t slice1 = (id>>25) & 0x7f; | |
285 | Int_t patch1 = (id>>22) & 0x7; | |
286 | UInt_t pos1 = id&0x3fffff; | |
4aa41877 | 287 | AliHLTSpacePointData *points1 = fClusters[slice1][patch1]; |
de3c3890 | 288 | |
289 | id = hitnum[num_of_clusters-1]; | |
290 | Int_t slice2 = (id>>25) & 0x7f; | |
291 | Int_t patch2 = (id>>22) & 0x7; | |
292 | UInt_t pos2 = id&0x3fffff; | |
4aa41877 | 293 | AliHLTSpacePointData *points2 = fClusters[slice2][patch2]; |
de3c3890 | 294 | |
295 | return kalmantrack->MakeSeed(track, points0, pos0, slice0, points1, pos1, slice1, points2, pos2, slice2); | |
296 | } | |
297 | ||
4aa41877 | 298 | Int_t AliHLTKalman::InitKalmanTrack(AliHLTKalmanTrack *kalmantrack, AliHLTTrack *track) |
b2a02bce | 299 | { |
300 | UInt_t *hitnum = track->GetHitNumbers(); | |
301 | UInt_t id; | |
302 | ||
303 | id = hitnum[0]; | |
304 | Int_t slice = (id>>25) & 0x7f; | |
305 | Int_t patch = (id>>22) & 0x7; | |
306 | UInt_t pos = id&0x3fffff; | |
4aa41877 | 307 | AliHLTSpacePointData *points = fClusters[slice][patch]; |
b2a02bce | 308 | |
0bd0c1ef | 309 | return kalmantrack->Init(track, points, pos, slice); |
3e87ef69 | 310 | } |
311 | ||
4aa41877 | 312 | Int_t AliHLTKalman::Propagate(AliHLTKalmanTrack *kalmantrack, AliHLTTrack *track) |
3e87ef69 | 313 | { |
0a86fbb7 | 314 | // This function propagtes the kalmantrack to the next cluster of the loaded |
315 | // track | |
3e87ef69 | 316 | Int_t num_of_clusters = track->GetNumberOfPoints(); |
3e87ef69 | 317 | UInt_t *hitnum = track->GetHitNumbers(); |
318 | UInt_t id; | |
0a86fbb7 | 319 | UInt_t badpoint = 0; |
3e87ef69 | 320 | |
321 | for (Int_t icl = 1; icl < num_of_clusters; icl++) | |
322 | { | |
0a86fbb7 | 323 | |
3e87ef69 | 324 | id = hitnum[icl]; |
325 | Int_t slice = (id>>25) & 0x7f; | |
b2a02bce | 326 | Int_t patch = (id>>22) & 0x7; |
3e87ef69 | 327 | UInt_t pos = id&0x3fffff; |
b2a02bce | 328 | |
4aa41877 | 329 | AliHLTSpacePointData *points = fClusters[slice][patch]; |
3e87ef69 | 330 | if (!points) continue; |
0bd0c1ef | 331 | if (kalmantrack->Propagate(points,pos,slice) == 0) |
0a86fbb7 | 332 | { |
333 | badpoint++; | |
334 | continue; | |
335 | } | |
3e87ef69 | 336 | } |
b2a02bce | 337 | |
338 | // If too many clusters are missing, the track is no good | |
339 | if (badpoint >= UInt_t(num_of_clusters*0.8)) return 0; | |
340 | return 1; | |
3e87ef69 | 341 | } |
342 | ||
4aa41877 | 343 | Double_t AliHLTKalman::GetCpuTime() |
3e87ef69 | 344 | { |
345 | //Return the Cputime in seconds. | |
346 | struct timeval tv; | |
347 | gettimeofday( &tv, NULL ); | |
348 | return tv.tv_sec+(((Double_t)tv.tv_usec)/1000000.); | |
349 | //return (Double_t)(clock()) / CLOCKS_PER_SEC; | |
350 | } |