Fit mathematics improved, obsollete GBTracker cleaned up
[u/mrichter/AliRoot.git] / HLT / TPCLib / tracking-ca / AliHLTTPCCATrackletConstructor.cxx
CommitLineData
00d07bcd 1// @(#) $Id: AliHLTTPCCATrackletConstructor.cxx 27042 2008-07-02 12:06:02Z richterm $
ce565086 2// **************************************************************************
00d07bcd 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: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
7// Ivan Kisel <kisel@kip.uni-heidelberg.de> *
8// for The ALICE HLT Project. *
9// *
10// Permission to use, copy, modify and distribute this software and its *
11// documentation strictly for non-commercial purposes is hereby granted *
12// without fee, provided that the above copyright notice appears in all *
13// copies and that both the copyright notice and this permission notice *
14// appear in the supporting documentation. The authors make no claims *
15// about the suitability of this software for any purpose. It is *
16// provided "as is" without express or implied warranty. *
ce565086 17// *
00d07bcd 18//***************************************************************************
19
ce565086 20
00d07bcd 21#include "AliHLTTPCCATracker.h"
22#include "AliHLTTPCCATrackParam.h"
4687b8fc 23#include "AliHLTTPCCATrackParam.h"
00d07bcd 24#include "AliHLTTPCCAGrid.h"
25#include "AliHLTTPCCAHitArea.h"
26#include "AliHLTTPCCAMath.h"
27#include "AliHLTTPCCADef.h"
ce565086 28#include "AliHLTTPCCATracklet.h"
00d07bcd 29#include "AliHLTTPCCATrackletConstructor.h"
693d2443 30//#include "AliHLTTPCCAPerformance.h"
31//#include "TH1D.h"
32
33//#define DRAW
34
35#ifdef DRAW
36#include "AliHLTTPCCADisplay.h"
37#endif
00d07bcd 38
00d07bcd 39
40GPUd() void AliHLTTPCCATrackletConstructor::Step0
4687b8fc 41( Int_t nBlocks, Int_t /*nThreads*/, Int_t iBlock, Int_t iThread,
42 AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, AliHLTTPCCATrackParam &/*tParam*/ )
00d07bcd 43{
44 // reconstruction of tracklets, step 0
45
4687b8fc 46 r.fIsMemThread = ( iThread<NMemThreads() );
47 if( iThread==0 ){
48 Int_t nTracks = *tracker.NTracklets();
49 Int_t nTrPerBlock = nTracks/nBlocks+1;
00d07bcd 50 s.fNRows = tracker.Param().NRows();
51 s.fItr0 = nTrPerBlock*iBlock;
52 s.fItr1 = s.fItr0 + nTrPerBlock;
53 if( s.fItr1> nTracks ) s.fItr1 = nTracks;
4687b8fc 54 s.fUsedHits = tracker.HitWeights();
00d07bcd 55 s.fMinStartRow = 158;
56 s.fMaxStartRow = 0;
57 }
58 if( iThread<32 ){
59 s.fMinStartRow32[iThread] = 158;
60 s.fMaxStartRow32[iThread] = 0;
61 }
62}
63
64
65GPUd() void AliHLTTPCCATrackletConstructor::Step1
4687b8fc 66( Int_t /*nBlocks*/, Int_t /*nThreads*/, Int_t /*iBlock*/, Int_t iThread,
67 AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, AliHLTTPCCATrackParam &tParam )
00d07bcd 68{
4687b8fc 69 // reconstruction of tracklets, step 1
00d07bcd 70
4687b8fc 71 r.fItr= s.fItr0 + ( iThread - NMemThreads() );
00d07bcd 72 r.fGo = (!r.fIsMemThread) && ( r.fItr<s.fItr1 );
73 r.fSave = r.fGo;
74 r.fNHits=0;
75
76 if( !r.fGo ) return;
77
78 r.fStage = 0;
00d07bcd 79
ce565086 80 AliHLTTPCCATracklet &tracklet = tracker.Tracklets()[r.fItr];
81
00d07bcd 82 UInt_t kThread = iThread %32;//& 00000020;
ce565086 83 if( SAVE() ) for( Int_t i=0; i<160; i++ ) tracklet.SetRowHit(i,-1);
4687b8fc 84
85 Int_t id = tracker.TrackletStartHits()[r.fItr];
ce565086 86 r.fStartRow = AliHLTTPCCATracker::ID2IRow(id);
87 r.fEndRow = r.fStartRow;
88 r.fFirstRow = r.fStartRow;
89 r.fLastRow = r.fFirstRow;
00d07bcd 90 r.fCurrIH = AliHLTTPCCATracker::ID2IHit(id);
4687b8fc 91
ce565086 92 CAMath::AtomicMin( &s.fMinStartRow32[kThread], r.fStartRow);
93 CAMath::AtomicMax( &s.fMaxStartRow32[kThread], r.fStartRow);
693d2443 94 tParam.SetSinPhi(0);
95 tParam.SetDzDs(0);
15d2e9cf 96 tParam.SetQPt(0);
97 tParam.SetSignCosPhi(1);
693d2443 98 tParam.SetChi2(0);
99 tParam.SetNDF(-3);
100 tParam.SetCov( 0,1);
101 tParam.SetCov( 1,0);
102 tParam.SetCov( 2,1);
103 tParam.SetCov( 3,0);
104 tParam.SetCov( 4,0);
105 tParam.SetCov( 5,1);
106 tParam.SetCov( 6,0);
107 tParam.SetCov( 7,0);
108 tParam.SetCov( 8,0);
109 tParam.SetCov( 9,1);
110 tParam.SetCov(10,0);
111 tParam.SetCov(11,0);
112 tParam.SetCov(12,0);
113 tParam.SetCov(13,0);
15d2e9cf 114 tParam.SetCov(14,10.);
325a2bc4 115
00d07bcd 116}
117
118GPUd() void AliHLTTPCCATrackletConstructor::Step2
4687b8fc 119( Int_t /*nBlocks*/, Int_t nThreads, Int_t /*iBlock*/, Int_t iThread,
120 AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &/*r*/, AliHLTTPCCATracker &/*tracker*/, AliHLTTPCCATrackParam &/*tParam*/ )
00d07bcd 121{
122 // reconstruction of tracklets, step 2
4687b8fc 123
00d07bcd 124 if( iThread==0 ){
ce565086 125 //CAMath::AtomicMinGPU(&s.fMinRow, s.fMinRow32[iThread]);
4687b8fc 126 Int_t minStartRow = 158;
127 Int_t maxStartRow = 0;
128 Int_t n = (nThreads>32 ) ?32 :nThreads;
129 for( Int_t i=0; i<n; i++ ){
00d07bcd 130 if( s.fMinStartRow32[i]<minStartRow ) minStartRow = s.fMinStartRow32[i];
131 if( s.fMaxStartRow32[i]>maxStartRow ) maxStartRow = s.fMaxStartRow32[i];
132 }
133 s.fMinStartRow = minStartRow;
134 s.fMaxStartRow = maxStartRow;
4687b8fc 135 }
00d07bcd 136}
137
138GPUd() void AliHLTTPCCATrackletConstructor::ReadData
139( Int_t iThread, AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, Int_t iRow )
4687b8fc 140{
00d07bcd 141 // reconstruction of tracklets, read data step
4687b8fc 142
00d07bcd 143 if( r.fIsMemThread ){
693d2443 144 const AliHLTTPCCARow &row = tracker.Row(iRow);
4687b8fc 145 Bool_t jr = !r.fCurrentData;
00d07bcd 146 Int_t n = row.FullSize();
693d2443 147 const uint4* gMem = tracker.RowData() + row.FullOffset();
00d07bcd 148 uint4 *sMem = s.fData[jr];
4687b8fc 149 for( Int_t i=iThread; i<n; i+=NMemThreads() ) sMem[i] = gMem[i];
150 }
00d07bcd 151}
152
153
154GPUd() void AliHLTTPCCATrackletConstructor::StoreTracklet
4687b8fc 155( Int_t /*nBlocks*/, Int_t /*nThreads*/, Int_t /*iBlock*/, Int_t /*iThread*/,
156 AliHLTTPCCASharedMemory &/*s*/, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, AliHLTTPCCATrackParam &tParam )
00d07bcd 157{
158 // reconstruction of tracklets, tracklet store step
4687b8fc 159
00d07bcd 160 if( !r.fSave ) return;
161
693d2443 162 //AliHLTTPCCAPerformance::Instance().HNHitsPerTrackCand()->Fill(r.fNHits);
163
00d07bcd 164 do{
693d2443 165 {
166 //std::cout<<"tracklet to store: "<<r.fItr<<", nhits = "<<r.fNHits<<std::endl;
167 }
168
169 if( r.fNHits<5 ){
00d07bcd 170 r.fNHits = 0;
171 break;
172 }
ce565086 173
15d2e9cf 174 if(0){
175 if( 1./.5 < CAMath::Abs(tParam.QPt()) ){ //SG!!!
ce565086 176 r.fNHits = 0;
177 break;
178 }
179 }
00d07bcd 180
181 {
182 Bool_t ok=1;
693d2443 183 const Float_t *c = tParam.Cov();
4687b8fc 184 for( Int_t i=0; i<15; i++ ) ok = ok && CAMath::Finite(c[i]);
185 for( Int_t i=0; i<5; i++ ) ok = ok && CAMath::Finite(tParam.Par()[i]);
00d07bcd 186 ok = ok && (tParam.X()>50);
187
188 if( c[0]<=0 || c[2]<=0 || c[5]<=0 || c[9]<=0 || c[14]<=0 ) ok = 0;
189
190 if(!ok){
191 r.fNHits = 0;
192 break;
193 }
194 }
195 }while(0);
196
197 if( !SAVE() ) return;
198
ce565086 199 AliHLTTPCCATracklet &tracklet = tracker.Tracklets()[r.fItr];
200
201 tracklet.SetNHits(r.fNHits);
00d07bcd 202
203 if( r.fNHits>0 ){
693d2443 204#ifdef DRAW
205 if(0){
206 std::cout<<"store tracklet "<<r.fItr<<", nhits = "<<r.fNHits<<std::endl;
207 if( AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kBlue,1.) ){
208 AliHLTTPCCADisplay::Instance().Ask();
209 }
210 }
211#endif
693d2443 212 if( CAMath::Abs(tParam.Par()[4])<1.e-8 ) tParam.SetPar( 4, 1.e-8);
ce565086 213 tracklet.SetFirstRow( r.fFirstRow );
214 tracklet.SetLastRow( r.fLastRow );
215 tracklet.SetParam( tParam );
4687b8fc 216 Int_t w = (r.fNHits<<16)+r.fItr;
217 for( Int_t iRow=0; iRow<160; iRow++ ){
ce565086 218 Int_t ih = tracklet.RowHit( iRow );
00d07bcd 219 if( ih>=0 ){
693d2443 220 Int_t ihTot = tracker.Row(iRow).FirstHit() + ih;
ce565086 221 CAMath::AtomicMax( tracker.HitWeights() + ihTot, w );
00d07bcd 222 }
223 }
224 }
225}
226
227GPUd() void AliHLTTPCCATrackletConstructor::UpdateTracklet
4687b8fc 228( Int_t /*nBlocks*/, Int_t /*nThreads*/, Int_t /*iBlock*/, Int_t /*iThread*/,
229 AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, AliHLTTPCCATrackParam &tParam,Int_t iRow )
230{
00d07bcd 231 // reconstruction of tracklets, tracklets update step
232
693d2443 233 //std::cout<<"Update tracklet: "<<r.fItr<<" "<<r.fGo<<" "<<r.fStage<<" "<<iRow<<std::endl;
15d2e9cf 234 Bool_t drawSearch = 0;//r.fItr==2;
235 Bool_t drawFit = 0;//r.fItr==2;
236 Bool_t drawFitted = drawFit ;//|| 1;//r.fItr==16;
693d2443 237
00d07bcd 238 if( !r.fGo ) return;
693d2443 239
240 const Int_t kMaxRowGap = 4;
00d07bcd 241
ce565086 242 AliHLTTPCCATracklet &tracklet = tracker.Tracklets()[r.fItr];
243
693d2443 244 const AliHLTTPCCARow &row = tracker.Row(iRow);
00d07bcd 245
4687b8fc 246 Float_t y0 = row.Grid().YMin();
247 Float_t stepY = row.HstepY();
248 Float_t z0 = row.Grid().ZMin();
249 Float_t stepZ = row.HstepZ();
250 Float_t stepYi = row.HstepYi();
251 Float_t stepZi = row.HstepZi();
00d07bcd 252
253 if( r.fStage == 0 ){ // fitting part
254 do{
255
ce565086 256 if( iRow<r.fStartRow || r.fCurrIH<0 ) break;
257
258 if( (iRow-r.fStartRow)%2!=0 ) break; // SG!!! - jump over the row
259
00d07bcd 260 uint4 *tmpint4 = s.fData[r.fCurrentData];
261 ushort2 hh = reinterpret_cast<ushort2*>(tmpint4)[r.fCurrIH];
262
263 Int_t oldIH = r.fCurrIH;
264 r.fCurrIH = reinterpret_cast<Short_t*>(tmpint4)[row.FullLinkOffset() + r.fCurrIH];
265
4687b8fc 266 Float_t x = row.X();
267 Float_t y = y0 + hh.x*stepY;
268 Float_t z = z0 + hh.y*stepZ;
693d2443 269 if( drawFit ) std::cout<<" fit tracklet: new hit "<<oldIH<<", xyz="<<x<<" "<<y<<" "<<z<<std::endl;
270
ce565086 271 if( iRow==r.fStartRow ){
693d2443 272 tParam.SetX( x );
273 tParam.SetY( y );
274 tParam.SetZ( z );
275 r.fLastY = y;
276 r.fLastZ = z;
277 //#ifdef DRAW
278 if( drawFit ) std::cout<<" fit tracklet "<<r.fItr<<", row "<<iRow<<" first row"<<std::endl;
279 //#endif
00d07bcd 280 }else{
693d2443 281
282 Float_t err2Y, err2Z;
283 Float_t dx = x - tParam.X();
284 Float_t dy = y - r.fLastY;//tParam.Y();
285 Float_t dz = z - r.fLastZ;//tParam.Z();
286 r.fLastY = y;
287 r.fLastZ = z;
288
289 Float_t ri = 1./CAMath::Sqrt(dx*dx+dy*dy);
ce565086 290 if( iRow==r.fStartRow+1 ){
693d2443 291 tParam.SetSinPhi( dy*ri );
15d2e9cf 292 tParam.SetSignCosPhi( dx );
693d2443 293 tParam.SetDzDs( dz*ri );
15d2e9cf 294 std::cout<<"Init. errors... "<<r.fItr<<std::endl;
693d2443 295 tracker.GetErrors2( iRow, tParam, err2Y, err2Z );
15d2e9cf 296 std::cout<<"Init. errors = "<<err2Y<<" "<<err2Z<<std::endl;
693d2443 297 tParam.SetCov( 0, err2Y );
298 tParam.SetCov( 2, err2Z );
299 }
300 if( drawFit ){
301 //#ifdef DRAW
302 std::cout<<" fit tracklet "<<r.fItr<<", row "<<iRow<<" transporting.."<<std::endl;
303 std::cout<<" params before transport="<<std::endl;
304 tParam.Print();
305 //#endif
306 }
307 Float_t sinPhi, cosPhi;
308 if( r.fNHits>=10 && CAMath::Abs(tParam.SinPhi())<.99 ){
309 sinPhi = tParam.SinPhi();
310 cosPhi = CAMath::Sqrt(1-sinPhi*sinPhi);
311 }else{
312 sinPhi = dy*ri;
313 cosPhi = dx*ri;
314 }
315 //#ifdef DRAW
316 if( drawFit ) std::cout<<"sinPhi0 = "<<sinPhi<<", cosPhi0 = "<<cosPhi<<std::endl;
317 //#endif
15d2e9cf 318 if( !tParam.TransportToX( x, sinPhi, cosPhi, tracker.Param().Bz(),-1 ) ){
693d2443 319 //#ifdef DRAW
320 if( drawFit ) std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<": can not transport!!"<<std::endl;
321//#endif
ce565086 322 if( SAVE() ) tracklet.SetRowHit( iRow, -1 );
00d07bcd 323 break;
324 }
15d2e9cf 325 //std::cout<<"mark1 "<<r.fItr<<std::endl;
326 //tParam.Print();
693d2443 327 tracker.GetErrors2( iRow, tParam.GetZ(), sinPhi, cosPhi, tParam.GetDzDs(), err2Y, err2Z );
15d2e9cf 328 //std::cout<<"mark2"<<std::endl;
693d2443 329
330 if( drawFit ){
331 //#ifdef DRAW
332 std::cout<<" params after transport="<<std::endl;
333 tParam.Print();
334 std::cout<<"fit tracklet before filter: "<<r.fItr<<", row "<<iRow<<" errs="<<err2Y<<" "<<err2Z<<std::endl;
335 //#endif
336#ifdef DRAW
337 AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kBlue, 2.,1);
338 AliHLTTPCCADisplay::Instance().Ask();
339#endif
340 }
15d2e9cf 341 if( !tParam.Filter( y, z, err2Y, err2Z, .99 ) ) {
693d2443 342 //#ifdef DRAW
343 if( drawFit ) std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<": can not filter!!"<<std::endl;
344 //#endif
ce565086 345 if( SAVE() ) tracklet.SetRowHit( iRow, -1 );
00d07bcd 346 break;
347 }
348 }
ce565086 349 if( SAVE() ) tracklet.SetRowHit( iRow, oldIH );
693d2443 350 if( drawFit ){
351 //#ifdef DRAW
352 std::cout<<"fit tracklet after filter "<<r.fItr<<", row "<<iRow<<std::endl;
353 tParam.Print();
354 //#endif
355#ifdef DRAW
356 AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kGreen,2.);
357 AliHLTTPCCADisplay::Instance().Ask();
358#endif
359 }
00d07bcd 360 r.fNHits++;
361 r.fLastRow = iRow;
ce565086 362 r.fEndRow = iRow;
00d07bcd 363 break;
364 } while(0);
693d2443 365
366 if( r.fCurrIH<0 ){
367 //#ifdef DRAW
368 if( drawFitted ) std::cout<<"fitted tracklet "<<r.fItr<<", nhits="<<r.fNHits<<std::endl;
369 //#endif
370 r.fStage = 1;
371 //AliHLTTPCCAPerformance::Instance().HNHitsPerSeed()->Fill(r.fNHits);
ce565086 372 if( r.fNHits<3 ){ r.fNHits=0; r.fGo = 0;}//SG!!!
693d2443 373 if( CAMath::Abs(tParam.SinPhi())>.999 ){
374 //#ifdef DRAW
375 if( drawFitted ) std::cout<<" fitted tracklet error: sinPhi="<<tParam.SinPhi()<<std::endl;
376 //#endif
377 r.fNHits=0; r.fGo = 0;
378 }else{
15d2e9cf 379 //tParam.SetCosPhi( CAMath::Sqrt(1-tParam.SinPhi()*tParam.SinPhi()) );
693d2443 380 }
381 if( drawFitted ){
382 //#ifdef DRAW
383 std::cout<<"fitted tracklet "<<r.fItr<<" miss="<<r.fNMissed<<" go="<<r.fGo<<std::endl;
384 tParam.Print();
385#ifdef DRAW
386 AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kBlue);
387 AliHLTTPCCADisplay::Instance().Ask();
388#endif
389 }
390 }
00d07bcd 391 }
392 else // forward/backward searching part
693d2443 393 {
00d07bcd 394 do{
693d2443 395 if( drawSearch ){
396 //#ifdef DRAW
397 std::cout<<"search tracklet "<<r.fItr<<" row "<<iRow<<" miss="<<r.fNMissed<<" go="<<r.fGo<<" stage="<<r.fStage<<std::endl;
398 //#endif
399 }
ce565086 400
401 if( r.fStage == 2 && ( (iRow>=r.fEndRow) ||
402 (iRow>=r.fStartRow && (iRow-r.fStartRow)%2==0 )
403 ) ) break;
693d2443 404 if( r.fNMissed>kMaxRowGap ){
00d07bcd 405 break;
406 }
693d2443 407
00d07bcd 408 r.fNMissed++;
693d2443 409
4687b8fc 410 Float_t x = row.X();
411 Float_t err2Y, err2Z;
693d2443 412 if( drawSearch ){
413 //#ifdef DRAW
414 std::cout<<"tracklet "<<r.fItr<<" before transport to row "<<iRow<<" : "<<std::endl;
415 tParam.Print();
416 //#endif
417 }
15d2e9cf 418 if( !tParam.TransportToX( x, tParam.SinPhi(), tParam.GetCosPhi(), tracker.Param().Bz(), .99 ) ){
693d2443 419 //#ifdef DRAW
420 if( drawSearch ) std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<": can not transport!!"<<std::endl;
421 //#endif
422 break;
423 }
424 if( drawSearch ){
425 //#ifdef DRAW
426 std::cout<<"tracklet "<<r.fItr<<" after transport to row "<<iRow<<" : "<<std::endl;
427 tParam.Print();
428#ifdef DRAW
429 AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kBlue,2.,1);
430 AliHLTTPCCADisplay::Instance().Ask();
431#endif
432 }
00d07bcd 433 uint4 *tmpint4 = s.fData[r.fCurrentData];
434
435 ushort2 *hits = reinterpret_cast<ushort2*>(tmpint4);
4687b8fc 436
437 Float_t fY = tParam.GetY();
438 Float_t fZ = tParam.GetZ();
00d07bcd 439 Int_t best = -1;
440
441 { // search for the closest hit
442
443 Int_t ds;
444 Int_t fY0 = (Int_t) ((fY - y0)*stepYi);
445 Int_t fZ0 = (Int_t) ((fZ - z0)*stepZi);
446 Int_t ds0 = ( ((int)1)<<30);
447 ds = ds0;
448
449 UInt_t fIndYmin;
450 UInt_t fHitYfst=1, fHitYlst=0, fHitYfst1=1, fHitYlst1=0;
451
452 fIndYmin = row.Grid().GetBin( (float)(fY-1.), (float)(fZ-1.) );
693d2443 453 if( drawSearch ){
454#ifdef DRAW
455 std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<": grid N="<<row.Grid().N() <<std::endl;
456 std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<": minbin="<<fIndYmin <<std::endl;
457#endif
458 }
4687b8fc 459 {
460 Int_t nY = row.Grid().Ny();
693d2443 461
4687b8fc 462 UShort_t *sGridP = (reinterpret_cast<UShort_t*>(tmpint4)) + row.FullGridOffset();
463 fHitYfst = sGridP[fIndYmin];
464 fHitYlst = sGridP[fIndYmin+2];
465 fHitYfst1 = sGridP[fIndYmin+nY];
466 fHitYlst1 = sGridP[fIndYmin+nY+2];
693d2443 467 if( drawSearch ){
468#ifdef DRAW
469 std::cout<<" Grid, row "<<iRow<<": nHits="<<row.NHits()<<", grid n="<<row.Grid().N()<<", c[n]="<<sGridP[row.Grid().N()]<<std::endl;
470 std::cout<<"hit steps = "<<stepY<<" "<<stepZ<<std::endl;
471 std::cout<<" Grid bins:"<<std::endl;
472 for( unsigned int i=0; i<row.Grid().N(); i++ ){
473 std::cout<<" bin "<<i<<": ";
474 for( int j=sGridP[i]; j<sGridP[i+1]; j++){
475 ushort2 hh = hits[j];
476 Float_t y = y0 + hh.x*stepY;
477 Float_t z = z0 + hh.y*stepZ;
478 std::cout<<"["<<j<<"|"<<y<<","<<z<<"] ";
479 }
480 std::cout<<std::endl;
481 }
482#endif
483 }
484 if( sGridP[row.Grid().N()]!=row.NHits() ){
485#ifdef DRAW
486 std::cout<<" grid, row "<<iRow<<": nHits="<<row.NHits()<<", grid n="<<row.Grid().N()<<", c[n]="<<sGridP[row.Grid().N()]<<std::endl;
487 //exit(0);
488#endif
489 }
490 }
491 if( drawSearch ){
492 //#ifdef DRAW
493 std::cout<<" tracklet "<<r.fItr<<", row "<<iRow<<", yz= "<<fY<<","<<fZ<<": search hits="<<fHitYfst<<" "<<fHitYlst<<" / "<< fHitYfst1<<" "<<fHitYlst1<<std::endl;
494 std::cout<<" hit search :"<<std::endl;
495 //#endif
4687b8fc 496 }
00d07bcd 497 for( UInt_t fIh = fHitYfst; fIh<fHitYlst; fIh++ ){
498 ushort2 hh = hits[fIh];
499 Int_t ddy = (Int_t)(hh.x) - fY0;
500 Int_t ddz = (Int_t)(hh.y) - fZ0;
693d2443 501 Int_t dds = CAMath::Abs(ddy) + CAMath::Abs(ddz);
502 if( drawSearch ){
503 //#ifdef DRAW
504 std::cout<<fIh<<": hityz= "<<hh.x<<" "<<hh.y<<"("<<hh.x*stepY<<" "<<hh.y*stepZ<<"), trackyz="<<fY0<<" "<<fZ0<<"("<<fY0*stepY<<" "<<fZ0*stepZ<<"), dy,dz,ds= "<<ddy<<" "<<ddz<<" "<<dds<<"("<<ddy*stepY<<" "<<ddz*stepZ<<std::endl;
505 //#endif
506 }
00d07bcd 507 if( dds<ds ){
508 ds = dds;
509 best = fIh;
510 }
511 }
512
513 for( UInt_t fIh = fHitYfst1; fIh<fHitYlst1; fIh++ ){
514 ushort2 hh = hits[fIh];
515 Int_t ddy = (Int_t)(hh.x) - fY0;
516 Int_t ddz = (Int_t)(hh.y) - fZ0;
693d2443 517 Int_t dds = CAMath::Abs(ddy) + CAMath::Abs(ddz);
518 if( drawSearch ){
519 //#ifdef DRAW
520 std::cout<<fIh<<": hityz= "<<hh.x<<" "<<hh.y<<"("<<hh.x*stepY<<" "<<hh.y*stepZ<<"), trackyz="<<fY0<<" "<<fZ0<<"("<<fY0*stepY<<" "<<fZ0*stepZ<<"), dy,dz,ds= "<<ddy<<" "<<ddz<<" "<<dds<<"("<<ddy*stepY<<" "<<ddz*stepZ<<std::endl;
521 //#endif
522 }
00d07bcd 523 if( dds<ds ){
524 ds = dds;
525 best = fIh;
526 }
527 }
528 }// end of search for the closest hit
529
530 if( best<0 ) break;
693d2443 531 if( drawSearch ){
532 //#ifdef DRAW
533 std::cout<<"hit search "<<r.fItr<<", row "<<iRow<<" hit "<<best<<" found"<<std::endl;
534#ifdef DRAW
535 AliHLTTPCCADisplay::Instance().DrawSliceHit( iRow, best, kRed, 1. );
536 AliHLTTPCCADisplay::Instance().Ask();
537 AliHLTTPCCADisplay::Instance().DrawSliceHit( iRow, best, kWhite, 1 );
538 AliHLTTPCCADisplay::Instance().DrawSliceHit( iRow, best );
539#endif
540 }
541
00d07bcd 542 ushort2 hh = hits[best];
15d2e9cf 543
544 //std::cout<<"mark 3, "<<r.fItr<<std::endl;
545 //tParam.Print();
00d07bcd 546 tracker.GetErrors2( iRow, *((AliHLTTPCCATrackParam*)&tParam), err2Y, err2Z );
15d2e9cf 547 //std::cout<<"mark 4"<<std::endl;
548
4687b8fc 549 Float_t y = y0 + hh.x*stepY;
550 Float_t z = z0 + hh.y*stepZ;
551 Float_t dy = y - fY;
552 Float_t dz = z - fZ;
00d07bcd 553
ce565086 554 const Float_t kFactor = tracker.Param().HitPickUpFactor()*tracker.Param().HitPickUpFactor()*3.5*3.5;
00d07bcd 555 Float_t sy2 = kFactor*( tParam.GetErr2Y() + err2Y );
556 Float_t sz2 = kFactor*( tParam.GetErr2Z() + err2Z );
ce565086 557 if( sy2 > 2. ) sy2 = 2.;
693d2443 558 if( sz2 > 2. ) sz2 = 2.;
ce565086 559
693d2443 560 if( drawSearch ){
561 //#ifdef DRAW
562 std::cout<<"dy,sy= "<< dy <<" "<<CAMath::Sqrt(sy2)<<", dz,sz= "<< dz<<" "<<CAMath::Sqrt(sz2)<<std::endl;
563 std::cout<<"dy,dz= "<< dy <<" "<<dz<<", sy,sz= "<<CAMath::Sqrt(sy2)<<" "<<CAMath::Sqrt(sz2)<<", sy,sz= "<<CAMath::Sqrt(kFactor*( tParam.GetErr2Y() + err2Y ))<<" "<<CAMath::Sqrt(kFactor*( tParam.GetErr2Z() + err2Z ))<<std::endl;
564 //#endif
565 }
ce565086 566 if( CAMath::FMulRZ(dy,dy)>sy2 || CAMath::FMulRZ(dz,dz)>sz2 ){
693d2443 567 if( drawSearch ){
568 //#ifdef DRAW
00d07bcd 569
693d2443 570 std::cout<<"found hit is out of the chi2 window\n "<<std::endl;
571 //#endif
572 }
573 break;
574 }
575#ifdef DRAW
576 //if( SAVE() ) hitstore[ iRow ] = best;
577 //std::cout<<"hit search before filter: "<<r.fItr<<", row "<<iRow<<std::endl;
578 //AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kBlue);
579 //AliHLTTPCCADisplay::Instance().Ask();
580#endif
15d2e9cf 581 if( !tParam.Filter( y, z, err2Y, err2Z, .99 ) ){
693d2443 582 if( drawSearch ){
583 //#ifdef DRAW
584 std::cout<<"tracklet "<<r.fItr<<" at row "<<iRow<<" : can not filter!!!! "<<std::endl;
585 //#endif
586 }
587 break;
588 }
ce565086 589 if( SAVE() ) tracklet.SetRowHit( iRow, best );
693d2443 590 if( drawSearch ){
591 //#ifdef DRAW
592 std::cout<<"tracklet "<<r.fItr<<" after filter at row "<<iRow<<" : "<<std::endl;
593 tParam.Print();
594#ifdef DRAW
595 AliHLTTPCCADisplay::Instance().DrawTracklet(tParam, hitstore, kRed);
596 AliHLTTPCCADisplay::Instance().Ask();
597#endif
598 }
599 r.fNHits++;
00d07bcd 600 r.fNMissed=0;
693d2443 601 if( r.fStage ==1 ) r.fLastRow = iRow;
602 else r.fFirstRow = iRow;
00d07bcd 603 }while(0);
604 }
605}
606
607
608
609GPUd() void AliHLTTPCCATrackletConstructor::Thread
610( Int_t nBlocks, Int_t nThreads, Int_t iBlock, Int_t iThread, Int_t iSync,
4687b8fc 611 AliHLTTPCCASharedMemory &s, AliHLTTPCCAThreadMemory &r, AliHLTTPCCATracker &tracker, AliHLTTPCCATrackParam &tParam )
00d07bcd 612{
613
614 // reconstruction of tracklets
615 if( iSync==0 )
616 {
4687b8fc 617 Step0( nBlocks, nThreads, iBlock, iThread, s, r, tracker, tParam );
00d07bcd 618 }
619 else if( iSync==1 )
620 {
4687b8fc 621 Step1( nBlocks, nThreads, iBlock, iThread, s, r, tracker, tParam );
00d07bcd 622 }
623 else if( iSync==2 )
624 {
4687b8fc 625 Step2( nBlocks, nThreads, iBlock, iThread, s, r, tracker, tParam );
00d07bcd 626 }
627
628 else if( iSync==3 )
629
630 {
631 r.fCurrentData = 1;
632 ReadData( iThread, s, r, tracker, s.fMinStartRow );
633 r.fCurrentData = 0;
634 r.fNMissed = 0;
635 }
4687b8fc 636 else if( iSync==3+159 +1 )
00d07bcd 637 {
638 r.fCurrentData = 1;
639 Int_t nextRow = s.fMaxStartRow-1;
640 if( nextRow<0 ) nextRow = 0;
641 ReadData( iThread, s, r, tracker, nextRow );
642 r.fCurrentData = 0;
643 r.fNMissed = 0;
644 r.fStage = 2;
693d2443 645 if( r.fGo ){
ce565086 646 const AliHLTTPCCARow &row = tracker.Row(r.fEndRow);
693d2443 647 Float_t x = row.X();
15d2e9cf 648 if( !tParam.TransportToX( x, tracker.Param().Bz(),.999 ) ) r.fGo = 0;
693d2443 649 }
00d07bcd 650 }
651
4687b8fc 652 else if( iSync<=3+159+1+159 )
00d07bcd 653 {
4687b8fc 654 Int_t iRow, nextRow;
655 if( iSync<=3+159 ){
656 iRow = iSync - 4;
00d07bcd 657 if( iRow < s.fMinStartRow ) return;
658 nextRow = iRow+1;
659 if( nextRow>158 ) nextRow = 158;
660 }else{
4687b8fc 661 iRow = 158 - (iSync - 4-159-1);
00d07bcd 662 if( iRow >= s.fMaxStartRow ) return;
663 nextRow = iRow-1;
664 if( nextRow<0 ) nextRow = 0;
665 }
666
4687b8fc 667 if( r.fIsMemThread ){
668 ReadData( iThread, s, r, tracker, nextRow );
669 }else{
670 UpdateTracklet( nBlocks, nThreads, iBlock, iThread,
671 s, r, tracker, tParam, iRow );
672 }
673 r.fCurrentData = !r.fCurrentData;
00d07bcd 674 }
675
4687b8fc 676 else if( iSync== 4+159*2 +1+1 ) //
00d07bcd 677 {
4687b8fc 678 StoreTracklet( nBlocks, nThreads, iBlock, iThread,
00d07bcd 679 s, r, tracker, tParam );
680 }
681}
682