]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HLT/MUON/OnlineAnalysis/AliHLTMUONFullTracker.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Protection against the zero magnetic field and cleanup in the corresponding component...
[u/mrichter/AliRoot.git] / HLT / MUON / OnlineAnalysis / AliHLTMUONFullTracker.cxx
CommitLineData
52c6d8aa 1/**************************************************************************
1ca600e9 2 * This file is property of and copyright by the ALICE HLT Project *
52c6d8aa 3 * All rights reserved. *
4 * *
5 * Primary Authors: *
6 * Indranil Das <indra.das@saha.ac.in> *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
1ca600e9 13 * about the suitability of this software for any purpose. It is *
52c6d8aa 14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17//////////////////////////////////////////////////////////////////////
18///Author : Indranil Das, SINP, INDIA
1ca600e9 19///
52c6d8aa 20///Email : indra.das@saha.ac.in || indra.ehep@gmail.com
21///
22/// This class is created to peform the a full track reconstroction
1ca600e9 23/// for online HLT for Dimuon Detector. It is based on the method of
52c6d8aa 24/// straight line tracking for slat detectors and cellular automation
25/// mehtods for finding the tracks in the quadrant detectos. The track
26/// segments of the front and back side of the spectromrter is either
1ca600e9 27/// compared using KalmanChi2 method (a slightly modified version of
52c6d8aa 28/// AliMUONTrackReconstructorK) or alternated method of simple
29/// extrapolation of tracks through dipole magnet. Finally the Track is
1ca600e9 30/// passed through Absorber to take the MCS effect and Branson Effect
52c6d8aa 31/// into account for correction of track parameters.
32/////////////////////////////////////////////////////////////////////////
33
d24a4636 34//ROOT Includes
eb672fd1 35#include "TVector3.h"
eb672fd1 36#include "TGeoGlobalMagField.h"
37
d24a4636 38//STEER Includes
eb672fd1 39#include "AliMagF.h"
40#include "AliCDBManager.h"
fd5b812e 41#include "AliCDBStorage.h"
eb672fd1 42#include "AliGeomManager.h"
eb672fd1 43
d24a4636 44//MUON Includes
eb672fd1 45#include "AliMUONTrackExtrap.h"
46#include "AliMUONTrackParam.h"
47#include "AliMUONTrackExtrap.h"
48#include "AliMUONConstants.h"
49#include "AliMUONGeometryTransformer.h"
50#include "AliMUONTrackParam.h"
51
d24a4636 52//MUON mappinf includes
fd5b812e 53#include "AliMpDEIterator.h"
eb672fd1 54
d24a4636 55//HLT Includes
fd5b812e 56#include "AliHLTLogging.h"
eb672fd1 57
d24a4636 58//HLT/MUON Includes
eb672fd1 59#include "AliHLTMUONConstants.h"
60#include "AliHLTMUONDataTypes.h"
61#include "AliHLTMUONUtils.h"
fd5b812e 62
d24a4636 63//HLT/MUON/OnlineAnalysis includes
fd5b812e 64#include "AliHLTMUONFullTracker.h"
eb672fd1 65
66using namespace std;
67
fd5b812e 68class AliRunInfo;
69class AliLog;
70class AliCDBEntry;
71class AliMpCDB;
72class AliMpSegmentation;
73class AliMpDDLStore;
74class AliGRPObject;
fd5b812e 75class TString;
76class TMap;
77
78//#define PRINT_FULL 1
79
ed5809f7 80
52c6d8aa 81#ifdef PRINT_FULL
82#define PRINT_POINTS 1
83#define PRINT_BACK 1
84#define PRINT_FRONT 1
85#define PRINT_KALMAN 1
86#define PRINT_SELECT 1
fd5b812e 87#define PRINT_OUTPUT 1
88#define PRINT_TRACKSEG 1
52c6d8aa 89///#define PRINT_DETAIL_KALMAN 1
90#endif
91
92class AliHLTMUONFullTracker;
93
eb672fd1 94const Float_t AliHLTMUONFullTracker::fgkTrackDetCoordinate[3] = {
52c6d8aa 95 155.179+20.0, 166.234+20.0,
96 (AliMUONConstants::DefaultChamberZ(4)+ AliMUONConstants::DefaultChamberZ(5))/2.0
1ca600e9 97};
52c6d8aa 98
99const Double_t AliHLTMUONFullTracker::fgkAbsoedge[4] = {92.0, 317.0,443.0,499.0};
100const Double_t AliHLTMUONFullTracker::fgkRadLen[3] = {13.875635, 11.273801, 1.765947};
101const Double_t AliHLTMUONFullTracker::fgkRho[3] = {1.750000, 2.350000, 7.880000};
102const Double_t AliHLTMUONFullTracker::fgkAtomicZ[3] = {6.000000, 11.098486,25.780000};
103const Double_t AliHLTMUONFullTracker::fgkAtomicA[3] = {12.010000, 22.334156,55.299670 };
104
105
fd5b812e 106const Int_t AliHLTMUONFullTracker::fgkMaxNofCellsPerCh = 1500;
107const Int_t AliHLTMUONFullTracker::fgkMaxNofPointsPerCh = 600;
52c6d8aa 108const Int_t AliHLTMUONFullTracker::fgkMaxNofCh = 11 ; /// 10tracking + 1 trigger
109const Int_t AliHLTMUONFullTracker::fgkMaxNofTracks = 200;
110const Int_t AliHLTMUONFullTracker::fgkMaxNofConnectedTracks = 20;
ed5809f7 111const Int_t AliHLTMUONFullTracker::fgkMaxNofTriggers = 20;
52c6d8aa 112
113
114
115AliHLTMUONFullTracker::AliHLTMUONFullTracker() :
116 AliHLTLogging(),
117 fChamberGeometryTransformer(0x0),
118 fChPoint(0x0),
119 fChPoint11(0x0),
120 fBackTrackSeg(0x0),
121 fFrontTrackSeg(0x0),
122 fExtrapSt3X(0x0),
123 fExtrapSt3Y(0x0),
124 fInclinationBack(0x0),
125 fNofConnectedfrontTrackSeg(0x0),
126 fBackToFront(0x0),
52c6d8aa 127 fNofPoints(0x0),
128 fTrackParam(0x0),
fd5b812e 129 fHalfTrack(0x0),
52c6d8aa 130 fTotNofPoints(0),
131 fTotTrackSeg(0),
52c6d8aa 132 fNofbackTrackSeg(0),
133 fNoffrontTrackSeg(0),
fd5b812e 134 fNofConnected(0),
135 fNofTracks(0),
aa205fc6 136 fDetElemList(),
ed5809f7 137 fFastTracking(kFALSE),
0fddb3cc 138 fNofInputs(0),
139 fNofTriggerInputs(0),
140 fNofTrackerInputs(0),
141 fIsMagfield(kFALSE)
52c6d8aa 142{
143
144 /// Constructor of the class
145
146 /// fChPoint = (AliHLTMUONRecHitStruct***)(malloc(10 * sizeof(AliHLTMUONRecHitStruct)));
147 /// for( Int_t ich=0;ich<10;ich++)
148 /// fChPoint[ich] = (AliHLTMUONRecHitStruct**)(malloc(300 * sizeof(AliHLTMUONRecHitStruct)));
149
150 try{
151 fChPoint = new AliHLTMUONRecHitStruct**[fgkMaxNofCh-1];
152 }catch (const std::bad_alloc&){
153 HLTError("Dynamic memory allocation failed in constructor : fChPoint");
154 throw;
155 }
156
157 for( Int_t ich=0;ich<(fgkMaxNofCh-1);ich++)
158 try{
159 fChPoint[ich] = new AliHLTMUONRecHitStruct*[fgkMaxNofPointsPerCh];
160 }catch (const std::bad_alloc&){
161 HLTError("Dynamic memory allocation failed in constructor : fChPoint[%d]",ich);
162 throw;
163 }
164
165 try{
166 fChPoint11 = new AliHLTMUONTriggerRecordStruct*[fgkMaxNofPointsPerCh];
167 }catch (const std::bad_alloc&){
168 HLTError("Dynamic memory allocation failed in constructor : fChPoint11");
169 throw;
170 }
171
172 try{
173 fBackTrackSeg = new TrackSeg[fgkMaxNofTracks];
174 }catch (const std::bad_alloc&){
175 HLTError("Dynamic memory allocation failed in constructor : fBackTrackSeg");
176 throw;
177 }
178
179 try{
180 fFrontTrackSeg = new TrackSeg[fgkMaxNofTracks];
181 }catch (const std::bad_alloc&){
182 HLTError("Dynamic memory allocation failed in constructor : fFrontTrackSeg");
183 throw;
184 }
185
186 try{
187 fExtrapSt3X = new float[fgkMaxNofTracks];
188 }catch (const std::bad_alloc&){
189 HLTError("Dynamic memory allocation failed in constructor : fExtrapSt3X");
190 throw;
191 }
192
193 try{
194 fExtrapSt3Y = new float[fgkMaxNofTracks];
195 }catch (const std::bad_alloc&){
196 HLTError("Dynamic memory allocation failed in constructor : fExtrapSt3Y");
197 throw;
198 }
199
200 try{
201 fInclinationBack = new float[fgkMaxNofTracks];
202 }catch (const std::bad_alloc&){
203 HLTError("Dynamic memory allocation failed in constructor : fInclinationBack");
204 throw;
205 }
206
207 try{
208 fNofConnectedfrontTrackSeg = new int[fgkMaxNofTracks];
209 }catch (const std::bad_alloc&){
210 HLTError("Dynamic memory allocation failed in constructor : fNofConnectedfrontTrackSeg");
211 throw;
212 }
213
214 try{
215 fBackToFront = new int*[fgkMaxNofTracks];
216 }catch (const std::bad_alloc&){
217 HLTError("Dynamic memory allocation failed in constructor : fBackToFront");
218 throw;
219 }
220
221 for( Int_t itr=0;itr<fgkMaxNofTracks;itr++)
222 try{
223 fBackToFront[itr] = new int[fgkMaxNofConnectedTracks];
224 }catch (const std::bad_alloc&){
225 HLTError("Dynamic memory allocation failed in constructor : fBackToFront[%d]",itr);
226 throw;
227 }
228
229
52c6d8aa 230 try{
231 fNofPoints = new int[fgkMaxNofCh];
232 }catch (const std::bad_alloc&){
233 HLTError("Dynamic memory allocation failed in constructor : fNofPoints");
234 throw;
235 }
fd5b812e 236
52c6d8aa 237 try{
238 fTrackParam = new AliMUONTrackParam[fgkMaxNofTracks];
239 }catch (const std::bad_alloc&){
240 HLTError("Dynamic memory allocation failed in constructor : fTrackParam");
241 throw;
242 }
243
fd5b812e 244 try{
245 fHalfTrack = new HalfTrack[fgkMaxNofTracks];
246 }catch (const std::bad_alloc&){
247 HLTError("Dynamic memory allocation failed in constructor : fHalfTrack");
248 throw;
249 }
250
52c6d8aa 251 Clear();
252
253}
254
255///__________________________________________________________________________
256
257AliHLTMUONFullTracker::~AliHLTMUONFullTracker()
258{
259 /// Destructor of the class
260
fd5b812e 261 ///delete fChamberGeometryTransformer;
52c6d8aa 262 ///free(fChPoint);
263 delete []fChPoint;
264 delete []fChPoint11;
265 delete []fBackTrackSeg;
266 delete []fFrontTrackSeg;
267 delete []fExtrapSt3X;
268 delete []fExtrapSt3Y;
269 delete []fInclinationBack;
270 delete []fNofConnectedfrontTrackSeg;
271 delete []fBackToFront;
52c6d8aa 272 delete []fNofPoints;
273 delete []fTrackParam;
fd5b812e 274 delete []fHalfTrack;
275
276 fChamberGeometryTransformer->Delete();
52c6d8aa 277
fd5b812e 278 fDetElemList.clear();
52c6d8aa 279
280}
281
282///__________________________________________________________________________
283
284Bool_t AliHLTMUONFullTracker::Clear(){
fd5b812e 285
52c6d8aa 286 /// Clear method to be called after each event
fd5b812e 287
288#ifdef PRINT_TRACKSEG
289 for(int iFrontTrackSeg=0;iFrontTrackSeg<fNoffrontTrackSeg;iFrontTrackSeg++)
290 for(int ipoint=0;ipoint<4;ipoint++)
291 if(fFrontTrackSeg[iFrontTrackSeg].fIndex[ipoint]!=-1)
ed5809f7 292 HLTImportant("FrontTrackSeg : %d\t%f\t%f\t%f\n",iFrontTrackSeg,
fd5b812e 293 fChPoint[ipoint][fFrontTrackSeg[iFrontTrackSeg].fIndex[ipoint]]->fX,
294 fChPoint[ipoint][fFrontTrackSeg[iFrontTrackSeg].fIndex[ipoint]]->fY,
295 fChPoint[ipoint][fFrontTrackSeg[iFrontTrackSeg].fIndex[ipoint]]->fZ);
ed5809f7 296 HLTImportant("\n\n");
fd5b812e 297 for(int iBackTrackSeg=0;iBackTrackSeg<fNofbackTrackSeg;iBackTrackSeg++)
298 for(int ipoint=0;ipoint<4;ipoint++)
299 if(fBackTrackSeg[iBackTrackSeg].fIndex[ipoint]!=-1)
ed5809f7 300 HLTImportant("BackTrackSeg : %d\t%f\t%f\t%f, nofFront : %d\n",iBackTrackSeg,
fd5b812e 301 fChPoint[ipoint+6][fBackTrackSeg[iBackTrackSeg].fIndex[ipoint]]->fX,
302 fChPoint[ipoint+6][fBackTrackSeg[iBackTrackSeg].fIndex[ipoint]]->fY,
303 fChPoint[ipoint+6][fBackTrackSeg[iBackTrackSeg].fIndex[ipoint]]->fZ,fNofConnectedfrontTrackSeg[iBackTrackSeg]);
304#endif
305
ed5809f7 306
307 fNofInputs = 0;
0fddb3cc 308 fNofTriggerInputs = 0;
309 fNofTrackerInputs = 0;
ed5809f7 310
52c6d8aa 311 for( Int_t ich=0;ich<fgkMaxNofCh;ich++)
312 fNofPoints[ich] = 0;
313
314 fNofbackTrackSeg = 0;
315 fNoffrontTrackSeg = 0;
316 fNofConnected = 0;
fd5b812e 317 fNofTracks = 0;
318
52c6d8aa 319 return true;
320}
321
322///__________________________________________________________________________
323
324void AliHLTMUONFullTracker::Print()
325{
326 /// Print anything mainly for debugging the codes, will be removed later
327 HLTInfo("\nPrint This--->\n");
328}
329
330///__________________________________________________________________________
331
332Bool_t AliHLTMUONFullTracker::Init()
333{
334 /// Initilation to be called once, later can be used to set/load the CDB path/entries
fd5b812e 335 HLTInfo("path : %s, run number : %d",
336 (AliCDBManager::Instance())->GetDefaultStorage()->GetURI().Data(),(AliCDBManager::Instance())->GetRun());
52c6d8aa 337 if (AliGeomManager::GetGeometry() == NULL){
338 AliGeomManager::LoadGeometry();
339 AliGeomManager::ApplyAlignObjsFromCDB("GRP MUON");
52c6d8aa 340 }
fd5b812e 341
342 Double_t b[3], x[3];
343 x[0] = 0.0 ; x[1] = 0.0 ; x[2] = -950.0;
52c6d8aa 344
fd5b812e 345 TGeoGlobalMagField::Instance()->Field(x,b);
346 //The following condition is based on the fact that at the middle of the dipole the field cannot be zero
347 if(TMath::AreEqualAbs(b[0],0.0,1.0e-5) and TMath::AreEqualAbs(b[1],0.0,1.0e-5) and TMath::AreEqualAbs(b[2],0.0,1.0e-5)){
d24a4636 348 HLTWarning("Magnetic field is not set by GRP");
0fddb3cc 349 if(not TGeoGlobalMagField::Instance()->IsLocked()){
350 TGeoGlobalMagField::Instance()->SetField(new AliMagF("Maps","Maps", -1., -1, AliMagF::k5kG));
351 TGeoGlobalMagField::Instance()->Field(x,b);
352 fIsMagfield = kTRUE;
353 }else{
354 HLTWarning("Magnetic field is not set and cannot be set since it is locked");
355 }
356 }else{
357 fIsMagfield = kTRUE;
52c6d8aa 358 }
359
fd5b812e 360
ed5809f7 361 HLTImportant("At (X,Y,Z) : (%6.2lf,%6.2lf,%6.2lf) Field (Bx,By,Bz) is (%6.2lf,%6.2lf,%6.2lf)",
fd5b812e 362 x[0],x[1],x[2],b[0],b[1],b[2]);
52c6d8aa 363
fd5b812e 364 AliMUONTrackExtrap::SetField();
365 /// see header file for class documentation
d24a4636 366 fChamberGeometryTransformer = new AliMUONGeometryTransformer();
367 if(! fChamberGeometryTransformer->LoadGeometryData()){
368 HLTError("Failed to Load Geomerty Data ");
369 }
a038ee7d 370
d24a4636 371 fDetElemList.clear();
372 for(int ich=0;ich<AliMUONConstants::NCh();ich++){
373 AliMpDEIterator it;
374 for ( it.First(ich); ! it.IsDone(); it.Next() )
375 {
376 Int_t detElemId = it.CurrentDEId();
377 fDetElemList[detElemId] = detElemId;
378 }
379 }//chamber loop
ed5809f7 380
d24a4636 381 return true;
52c6d8aa 382}
383
384///__________________________________________________________________________
385
a038ee7d 386Bool_t AliHLTMUONFullTracker::SetInput(AliHLTInt32_t /*ddl*/, const AliHLTMUONRecHitStruct *data, AliHLTInt32_t size)
52c6d8aa 387{
388 /// Set the input for rechit points
fd5b812e 389 if(int(size)>=fgkMaxNofPointsPerCh/2)
390 size = 0;
391
52c6d8aa 392 AliHLTUInt16_t detElemID;
393 AliHLTUInt8_t chamber;
394
fd5b812e 395#ifdef PRINT_POINTS
ed5809f7 396 HLTImportant("Received from DDL : %d, nofHits : %d, data : %p",ddl,size,data);
fd5b812e 397#endif
52c6d8aa 398 for( Int_t ipoint=0;ipoint<int(size);ipoint++){
a038ee7d 399 if(!data){
400 HLTError("Null Data pointer from HitRec");
401 Clear();
402 return false;
403 }
404
52c6d8aa 405 AliHLTMUONUtils::UnpackRecHitFlags(data->fFlags,chamber,detElemID);
a038ee7d 406
ed5809f7 407 if((not fDetElemList[detElemID]) or (chamber>=AliMUONConstants::NTrackingCh())){
a038ee7d 408 HLTDebug("Invalid tracking detelem : %d or chamber : %d",detElemID,chamber);
409 continue;
410 }
411
52c6d8aa 412#ifdef PRINT_POINTS
ed5809f7 413 HLTImportant("ch : %02d, detelem : %04d, (X,Y,Z) : (%8.3f,%8.3f,%8.3f)",
414 chamber,detElemID,data->fX,data->fY,data->fZ);
415
52c6d8aa 416#endif
a038ee7d 417 fChPoint[detElemID/100-1][fNofPoints[detElemID/100-1]++] = (AliHLTMUONRecHitStruct *)data;
52c6d8aa 418 data++;
419 }
ed5809f7 420 fNofInputs++;
0fddb3cc 421 fNofTrackerInputs++;
52c6d8aa 422 return true;
423}
424
425///__________________________________________________________________________
426
427Bool_t AliHLTMUONFullTracker::SetInput(AliHLTInt32_t /*ddl*/, const AliHLTMUONTriggerRecordStruct *data, AliHLTInt32_t size)
428{
429 /// Set the input for trigrecs
430
ed5809f7 431// if(int(size)>=fgkMaxNofPointsPerCh/2)
432// size = 0;
433
434 if(int(size)>fgkMaxNofTriggers){
435 HLTWarning("More triggers (%d) found than max limit : %d (not possible physics events)",int(size),fgkMaxNofTriggers);
436 Clear();
437 return false;
438 }
fd5b812e 439
52c6d8aa 440 AliHLTUInt16_t detElemID;
441 AliHLTUInt8_t chamber;
ed5809f7 442
52c6d8aa 443 for( Int_t ipoint=0;ipoint<int(size);ipoint++){
a038ee7d 444 if(!data){
445 HLTError("Null Data pointer from TrigRec");
446 Clear();
447 return false;
448 }
52c6d8aa 449 fChPoint11[fNofPoints[10]++] = (AliHLTMUONTriggerRecordStruct *)data;
450 for( Int_t ich=0;ich<4;ich++){
451 AliHLTMUONUtils::UnpackRecHitFlags((data->fHit[ich]).fFlags,chamber,detElemID);
a038ee7d 452 if((not fDetElemList[detElemID]) or (chamber<AliMUONConstants::NTrackingCh()) or (chamber > AliMUONConstants::NCh()) ){
453 HLTDebug("Invalid trigger detelem : %d or chamber : %d",detElemID,chamber);
454 continue;
455 }
52c6d8aa 456#ifdef PRINT_POINTS
ed5809f7 457 HLTImportant("size : %d, itrig : %04d, ch : %02d, detelem : %04d, (X,Y,Z) : (%8.3f,%8.3f,%8.3f)\n",
52c6d8aa 458 size,ipoint,chamber,detElemID,(data->fHit[ich]).fX,(data->fHit[ich]).fY,(data->fHit[ich]).fZ);
459#endif
460 }///ich loop
d24a4636 461 data++;
52c6d8aa 462 }///ipoint
ed5809f7 463 fNofInputs++;
0fddb3cc 464 fNofTriggerInputs++;
d24a4636 465 return true;
52c6d8aa 466}
467
d24a4636 468
52c6d8aa 469///__________________________________________________________________________
d24a4636 470
ed5809f7 471Bool_t AliHLTMUONFullTracker::CheckInput(AliHLTEventID_t iEvent)
52c6d8aa 472{
ed5809f7 473 /// Cross Check all the inputs before the starting of the tracking
474
fd5b812e 475 bool resultOk = true;
ed5809f7 476
0fddb3cc 477 //if more than expected inputs or no input from one detector, do not do anything
478 if((fNofInputs > 22) or (fNofTrackerInputs > 20) or (fNofTriggerInputs > 2)
479 or (fNofTrackerInputs == 0) or (fNofTriggerInputs == 0)){
ed5809f7 480
481 resultOk = false;
482 return resultOk;
483
484 }else{ // make double sure that no space point pointer is null
485
486 //HLTImportant("Found fNofInputs : %d less than 22",fNofInputs);
487
488 //tracker chamber test
489 for(int ich=0;ich<fgkMaxNofCh-1;ich++){
490 for(int ipt=0;ipt<fNofPoints[ich];ipt++){
491
492 if((not fChPoint[ich][ipt]) or (TMath::AreEqualAbs(fChPoint[ich][ipt]->fX,0.0,1.0e-5) and
493 TMath::AreEqualAbs(fChPoint[ich][ipt]->fY,0.0,1.0e-5) and
494 TMath::AreEqualAbs(fChPoint[ich][ipt]->fZ,0.0,1.0e-5))){
495 resultOk = false;
0fddb3cc 496 HLTError("iEvent : 0x%x, fNofTrackerInputs : %d, Nof tracker point for chamber %d, is not equal to nof valid tracker pointer",
497 iEvent,fNofTrackerInputs,ich);
ed5809f7 498 return resultOk;
499 }
500 }// tracker ch loop
501 }// // if resultOk not already false
502
503 //trigger chamber test
504 if(fNofPoints[10] == 0){
505 resultOk = false;
506 return resultOk;
507 }
508
509 for(int ipt=0;ipt<fNofPoints[10];ipt++){
510 if(not fChPoint11[ipt]){
511 resultOk = false;
0fddb3cc 512 HLTError("iEvent : 0x%x, fNofTriggerInputs : %d, Nof trigger points, is not equal to nof valid tracker pointer",
513 iEvent,fNofTriggerInputs);
ed5809f7 514 return resultOk;
515
516 }
517 }// for loop over points
fd5b812e 518
ed5809f7 519 }// if less than expected blocks found
520
521 return resultOk;
522
523}
524
525
526///__________________________________________________________________________
527
528Bool_t AliHLTMUONFullTracker::Run( AliHLTEventID_t iEvent,AliHLTMUONTrackStruct *data, AliHLTUInt32_t& size)
529{
530 /// Main Run call of the class
531
532 bool resultOk = true;
533
534 HLTDebug("Processing iEvent : 0x%X, nof triggers : %d, fNofInputs : %d",iEvent,fNofPoints[10],fNofInputs);
535
536 // for(int ich=0;ich<fgkMaxNofCh;ich++)
537// HLTDebug("\tNof hits in ich [%d] : %d\t",ich,fNofPoints[ich]);
538
539 resultOk = CheckInput(iEvent);
540
0fddb3cc 541 if((not fIsMagfield) and resultOk)
542 resultOk = false;
543
ed5809f7 544 if(resultOk){
545 resultOk = SlatTrackSeg();
546 if(not resultOk){
547 HLTDebug("Error happened in tracking through slat chambers, this event will be skipped");
548 }
fd5b812e 549 }
550 HLTDebug("Finishing SlatTrackSeg");
551
a038ee7d 552 if(resultOk){
553 resultOk = PrelimMomCalc();
554 if(not resultOk){
555 HLTDebug("Error happened in calculating preliminary momentum, this event will be skipped");
556 }
fd5b812e 557 }
558 HLTDebug("Finishing PrelimMomCalc");
fd5b812e 559
560 if(resultOk){
561 resultOk = QuadTrackSeg();
562 if(not resultOk){
563 HLTDebug("Error happened in tracking through quadrant chambers, this event will be skipped");
564 }
565 }
566 HLTDebug("Finishing QuadTrackSeg");
567
568
569 if(resultOk){
aa205fc6 570 if(fFastTracking)
571 resultOk = SelectFront();
572 else
573 resultOk = KalmanChi2Test();
574
fd5b812e 575 if(not resultOk){
576 HLTDebug("Error happened in tracking through in Kalman Chi2 checking, this event will be skipped");
577 }
578 }
579 HLTDebug("Finishing KalmanChi2Test");
52c6d8aa 580
fd5b812e 581 if(resultOk){
aa205fc6 582 resultOk = ExtrapolateToOrigin();
fd5b812e 583 if(not resultOk){
584 HLTDebug("Error happened in tracking extrapolation, this event will be skipped");
585 }
586 }
587 HLTDebug("Finishing ExtrapolateToOrigin");
588
589 if(resultOk){
590 resultOk = FillOutData(data,size);
591 if(not resultOk){
592 HLTDebug("Error happened in filling the output data, this event will be skipped");
593 }
594 }
595 HLTDebug("Finishing FillOutData");
596
597 if(!resultOk)
598 size = 0;
599
ed5809f7 600 HLTDebug("iEvent: 0x%X, has tracks : %d, triggers : %d, nof slat tracks : %d, quad tracks : %d, connected : %d\n",
fd5b812e 601 iEvent,size,fNofPoints[10],fNofbackTrackSeg,fNoffrontTrackSeg,fNofConnected);
52c6d8aa 602 Clear();
fd5b812e 603
604 return resultOk;
605
52c6d8aa 606}
607
608///__________________________________________________________________________
609
eb672fd1 610void AliHLTMUONFullTracker::Sub(const AliHLTMUONRecHitStruct *v1, const AliHLTMUONRecHitStruct *v2, AliHLTMUONRecHitStruct *v3) const
52c6d8aa 611{
612 /// Subtraction of position co-odinate of two space points
613
614 v3->fX = v1->fX - v2->fX;
615 v3->fY = v1->fY - v2->fY;
616 v3->fZ = v1->fZ - v2->fZ;
617};
618
619///__________________________________________________________________________
620
fd5b812e 621Double_t AliHLTMUONFullTracker::Angle(const AliHLTMUONRecHitStruct *v1, const AliHLTMUONRecHitStruct *v2)
52c6d8aa 622{
623 ///Angle of a straight line formed using v1 and v2
624
625 Float_t ptot2 = ((v1->fX * v1->fX) + (v1->fY * v1->fY) + (v1->fZ * v1->fZ))*
626 ((v2->fX * v2->fX) + (v2->fY * v2->fY) + (v2->fZ * v2->fZ));
627 if(ptot2 <= 0) {
628 return 0.0;
629 } else {
630 Float_t arg = ((v1->fX * v2->fX) + (v1->fY * v2->fY) + (v1->fZ * v2->fZ))/sqrt(ptot2);
631 if(arg > 1.0) arg = 1.0;
632 if(arg < -1.0) arg = -1.0;
633 return TMath::ACos(arg);
634 ///return acos(arg);
d24a4636 635 }
52c6d8aa 636
637}
638
639///__________________________________________________________________________
640
d24a4636 641Bool_t AliHLTMUONFullTracker::FillOutData(AliHLTMUONTrackStruct *track, AliHLTUInt32_t& size)
52c6d8aa 642{
643 ///Fill the output data pointers
644
645 size = (AliHLTUInt32_t(fNofbackTrackSeg)<size) ? AliHLTUInt32_t(fNofbackTrackSeg) : size;
d24a4636 646
1ca600e9 647 Bool_t hitset[16];
52c6d8aa 648 for( Int_t ibackTrackSeg=0;ibackTrackSeg<int(size);ibackTrackSeg++){
fd5b812e 649
650 if(fNofConnectedfrontTrackSeg[ibackTrackSeg]>0){
651
52c6d8aa 652
fd5b812e 653 if(not TMath::Finite(fTrackParam[ibackTrackSeg].Px())
654 || not TMath::Finite(fTrackParam[ibackTrackSeg].Py())
655 || not TMath::Finite(fTrackParam[ibackTrackSeg].Pz())) continue;
656
657#ifdef PRINT_OUTPUT
ed5809f7 658 HLTImportant("\nsize : %d, itrack : %04d, sign : %2d, Pt : %8.3f, (Px,Py,Pz) : (%8.3f,%8.3f,%8.3f)\n",
fd5b812e 659 size,ibackTrackSeg,Int_t(TMath::Sign(1.,fTrackParam[ibackTrackSeg].GetInverseBendingMomentum())),
1ca600e9 660 TMath::Sqrt(fTrackParam[ibackTrackSeg].Px()*fTrackParam[ibackTrackSeg].Px() +
fd5b812e 661 fTrackParam[ibackTrackSeg].Py()*fTrackParam[ibackTrackSeg].Py()),
662 fTrackParam[ibackTrackSeg].Px(),fTrackParam[ibackTrackSeg].Py(),
663 fTrackParam[ibackTrackSeg].Pz());
664#endif
665
666
d24a4636 667 // Bits 8 and 9 must be kept zero to prevent the track ID from conflicting with other tracker components.
668 track->fId = (ibackTrackSeg << 10) | (fBackTrackSeg[ibackTrackSeg].fTrigRec & 0xFF);
fd5b812e 669 track->fTrigRec = fBackTrackSeg[ibackTrackSeg].fTrigRec;
670 track->fPx = fTrackParam[ibackTrackSeg].Px();
671 track->fPy = fTrackParam[ibackTrackSeg].Py();
1ca600e9 672 track->fPz = fTrackParam[ibackTrackSeg].Pz();
673
fd5b812e 674 track->fChi2 = 0;
52c6d8aa 675
d24a4636 676 track->fInverseBendingMomentum = fTrackParam[ibackTrackSeg].GetInverseBendingMomentum();
677 track->fThetaY = TMath::Tan(fTrackParam[ibackTrackSeg].GetBendingSlope());
678 track->fThetaX = TMath::Tan(fTrackParam[ibackTrackSeg].GetNonBendingSlope());
fd5b812e 679
d24a4636 680 track->fZ = fTrackParam[ibackTrackSeg].GetZ();
681 track->fY = fTrackParam[ibackTrackSeg].GetBendingCoor();
682 track->fX = fTrackParam[ibackTrackSeg].GetNonBendingCoor();
fd5b812e 683
d24a4636 684
685 for( Int_t ipoint=15;ipoint>=0;ipoint--){
fd5b812e 686 track->fHit[ipoint] = AliHLTMUONConstants::NilRecHitStruct();
687 hitset[ipoint] = false;
1ca600e9 688 if(ipoint >= 6 and ipoint <= 9 and fBackTrackSeg[ibackTrackSeg].fIndex[ipoint-6]!=-1 ){
d24a4636 689 track->fHit[ipoint] = *(fChPoint[ipoint][fBackTrackSeg[ibackTrackSeg].fIndex[ipoint-6]]);
690 hitset[ipoint] = true;
691 }else if(ipoint <= 3 and fFrontTrackSeg[fBackToFront[ibackTrackSeg][0]].fIndex[ipoint]!=-1 ){
692 track->fHit[ipoint] = *(fChPoint[ipoint][fFrontTrackSeg[fBackToFront[ibackTrackSeg][0]].fIndex[ipoint]]);
693 hitset[ipoint] = true;
fd5b812e 694 }
52c6d8aa 695 }
fd5b812e 696 AliHLTMUONParticleSign sign = AliHLTMUONParticleSign(Int_t(TMath::Sign(1.,fTrackParam[ibackTrackSeg].GetInverseBendingMomentum())));
d24a4636 697 track->fFlags = AliHLTMUONUtils::PackTrackFlags(sign,hitset);
fd5b812e 698
699 track++;
700 fNofTracks++;
52c6d8aa 701
fd5b812e 702 }else{
703
704
1ca600e9 705 if(not TMath::Finite(fHalfTrack[ibackTrackSeg].fPx)
706 || not TMath::Finite(fHalfTrack[ibackTrackSeg].fPy)
707 || not TMath::Finite(fHalfTrack[ibackTrackSeg].fPz)) continue;
fd5b812e 708
709#ifdef PRINT_OUTPUT
ed5809f7 710 HLTImportant("\nsize : %d, itrack : %04d, sign : %2d, Pt : %8.3f, (Px,Py,Pz) : (%8.3f,%8.3f,%8.3f)\n",
1ca600e9 711 size,ibackTrackSeg,Int_t(TMath::Sign(1.,fTrackParam[ibackTrackSeg].GetInverseBendingMomentum())),
712 TMath::Sqrt(fHalfTrack[ibackTrackSeg].fPx*fHalfTrack[ibackTrackSeg].fPx +
fd5b812e 713 fHalfTrack[ibackTrackSeg].fPy*fHalfTrack[ibackTrackSeg].fPy),
714 fHalfTrack[ibackTrackSeg].fPx,fHalfTrack[ibackTrackSeg].fPy,
715 fHalfTrack[ibackTrackSeg].fPz);
716#endif
717
718
d24a4636 719 // Bits 8 and 9 must be kept zero to prevent the track ID from conflicting with other tracker components.
720 track->fId = (ibackTrackSeg << 10) | (fBackTrackSeg[ibackTrackSeg].fTrigRec & 0xFF);
fd5b812e 721 track->fTrigRec = fBackTrackSeg[ibackTrackSeg].fTrigRec;
722 track->fPx = fHalfTrack[ibackTrackSeg].fPx;
723 track->fPy = fHalfTrack[ibackTrackSeg].fPy;
1ca600e9 724 track->fPz = fHalfTrack[ibackTrackSeg].fPz;
fd5b812e 725
726 track->fChi2 = 0;
727
fd5b812e 728
d24a4636 729
730 for( Int_t ipoint=15;ipoint>=0;ipoint--){
fd5b812e 731 track->fHit[ipoint] = AliHLTMUONConstants::NilRecHitStruct();
732 hitset[ipoint] = false;
d24a4636 733
734 if(ipoint<=9 and ipoint>=6 and fBackTrackSeg[ibackTrackSeg].fIndex[ipoint-6]!=-1 ){
fd5b812e 735
d24a4636 736 track->fHit[ipoint] = *(fChPoint[ipoint][fBackTrackSeg[ibackTrackSeg].fIndex[ipoint-6]]);
fd5b812e 737 hitset[ipoint] = true;
738 }
739 }
740 AliHLTMUONParticleSign sign = AliHLTMUONParticleSign(fHalfTrack[ibackTrackSeg].fCharge);
741 track->fFlags = AliHLTMUONUtils::PackMansoTrackFlags(sign,hitset);
742
743 track++;
744 fNofTracks++;
745
746 }//if nof connected is more than zero or not
747 }//back track seg for loop
748
749 size = fNofTracks;
750
52c6d8aa 751 return true;
752}
753
754///__________________________________________________________________________
755
d24a4636 756
52c6d8aa 757Bool_t AliHLTMUONFullTracker::SlatTrackSeg()
758{
759
760 ///Find the Slat Track Segments
fd5b812e 761 if(fNofPoints[6]==0 and fNofPoints[7]==0){
762 HLTDebug("No Hits found in Stn4");
763 return false;
764 }else if(fNofPoints[8]==0 and fNofPoints[9]==0){
765 HLTDebug("No Hits found in Stn5");
766 return false;
767 }
52c6d8aa 768
769 Float_t trigX1,trigY1,trigZ1=AliMUONConstants::DefaultChamberZ(10);
770 Float_t trigX2,trigY2,trigZ2=AliMUONConstants::DefaultChamberZ(12);
771 Float_t extrapCh9X,extrapCh9Y,extrapCh9Z=AliMUONConstants::DefaultChamberZ(9);
772 Float_t extrapCh8X,extrapCh8Y,extrapCh8Z=AliMUONConstants::DefaultChamberZ(8);
773 Float_t extrapCh7X,extrapCh7Y,extrapCh7Z=AliMUONConstants::DefaultChamberZ(7);
774 Float_t extrapCh6X,extrapCh6Y,extrapCh6Z=AliMUONConstants::DefaultChamberZ(6);
775
52c6d8aa 776 Double_t distChFront,distChBack;
777 Int_t nofFrontChPoints,nofBackChPoints;
778 Int_t frontIndex[fgkMaxNofTracks], backIndex[fgkMaxNofTracks];
779
780 AliHLTMUONRecHitStruct p2,p3,pSeg1,pSeg2,pSeg3;
781 Double_t anglediff,anglediff1,anglediff2;
782 Double_t minAngle = 2.0;
783
eb672fd1 784 Bool_t st5TrackletFound = false;
785 Bool_t ch9PointFound = false;
786 Bool_t ch8PointFound = false;
787 Bool_t st4TrackletFound = false;
788 Bool_t ch7PointFound = false;
789 Bool_t ch6PointFound = false;
52c6d8aa 790
791 Int_t index1,index2,index3,index4;
eb672fd1 792 IntPair cells[2][fgkMaxNofTracks]; ///cell array for 5 stn for given trigger
52c6d8aa 793
d24a4636 794 Float_t maxXDeflectionExtrap = 10.0 + 4.0; ///simulation result 10.0
795 Float_t extrapRatio = 0.2; ///simulation result 0.2
796 Float_t circularWindow = 20.0 + 5.0 + 25.0; ///simulatiuon result 20
797 Float_t minAngleWindow = 2.0 + 1.0 + 2.0; ///simulation result 2.0
aa205fc6 798
799 if(fFastTracking){
800 maxXDeflectionExtrap = 10.0; ///simulation result 10.0
801 extrapRatio = 0.2; ///simulation result 0.2
802 circularWindow = 20.0 ; ///simulatiuon result 20
803 minAngleWindow = 2.0; ///simulation result 2.0
804 }
52c6d8aa 805
d24a4636 806 Float_t tx=0.0,ty=0.0;
fd5b812e 807
d24a4636 808 AliHLTUInt16_t detElemID,prevDetElemID=0xFFFF;
809 AliHLTUInt8_t chamber;
810 Int_t minTrgCh,maxTrgCh;
fd5b812e 811
812#ifdef PRINT_BACK
ed5809f7 813 HLTImportant("\nAliHLTMUONFullTracker::SlatTrackSeg()--Begin\n\n");
fd5b812e 814#endif
52c6d8aa 815
fd5b812e 816
d24a4636 817 for( Int_t itrig=0;itrig<fNofPoints[10];itrig++){
52c6d8aa 818
d24a4636 819 st5TrackletFound = false;
820 ch9PointFound = false;
821 ch8PointFound = false;
52c6d8aa 822
d24a4636 823 st4TrackletFound = false;
824 ch7PointFound = false;
825 ch6PointFound = false;
52c6d8aa 826
d24a4636 827 minTrgCh = -1;
828 maxTrgCh = -1;
52c6d8aa 829
d24a4636 830 fNofCells[0] = 0;
831 fNofCells[1] = 0;
52c6d8aa 832
d24a4636 833 for( Int_t ihit=0;ihit<4;ihit++){
52c6d8aa 834
d24a4636 835 AliHLTMUONUtils::UnpackRecHitFlags((fChPoint11[itrig]->fHit[ihit]).fFlags,chamber,detElemID);
52c6d8aa 836
d24a4636 837 if(ihit==0 and detElemID!=0)
838 minTrgCh = ihit;
839 if(ihit==1 and detElemID!=0 and prevDetElemID==0)
840 minTrgCh = ihit;
841 if(ihit==2 and detElemID!=0)
842 maxTrgCh = ihit;
843 if(ihit==3 and detElemID!=0 and prevDetElemID==0)
844 maxTrgCh = ihit;
52c6d8aa 845
d24a4636 846 prevDetElemID = detElemID;
847 }
fd5b812e 848
d24a4636 849 if(minTrgCh == -1 or maxTrgCh == -1){
850 HLTDebug("Trigger hits not found in both trigger station minTrgCh : %d, maxTrgCh : %d, not harmful to HLT chain",minTrgCh,maxTrgCh);
851 continue;
852 }
52c6d8aa 853
aa205fc6 854 if( not fFastTracking){
855 AliHLTMUONUtils::UnpackRecHitFlags((fChPoint11[itrig]->fHit[minTrgCh]).fFlags,chamber,detElemID);
856 if(not fDetElemList[detElemID]){
857 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
858 continue;
859 }
860 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,trigZ1);
861 AliHLTMUONUtils::UnpackRecHitFlags((fChPoint11[itrig]->fHit[maxTrgCh]).fFlags,chamber,detElemID);
862 if(not fDetElemList[detElemID]){
863 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
864 continue;
865 }
866 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,trigZ2);
d24a4636 867 }
aa205fc6 868
52c6d8aa 869
d24a4636 870 trigX1 = (fChPoint11[itrig]->fHit[minTrgCh]).fX;
871 trigY1 = (fChPoint11[itrig]->fHit[minTrgCh]).fY;
872 trigZ1 = (fChPoint11[itrig]->fHit[minTrgCh]).fZ;
52c6d8aa 873
d24a4636 874 trigX2 = (fChPoint11[itrig]->fHit[maxTrgCh]).fX;
875 trigY2 = (fChPoint11[itrig]->fHit[maxTrgCh]).fY;
876 trigZ2 = (fChPoint11[itrig]->fHit[maxTrgCh]).fZ;
52c6d8aa 877
52c6d8aa 878#ifdef PRINT_BACK
ed5809f7 879 HLTImportant("itrig : %d trig 1 : (%f,%f,%f) \n",itrig,trigX1,trigY1,trigZ1);
880 HLTImportant("itrig : %d trig 2 : (%f,%f,%f) \n",itrig,trigX2,trigY2,trigZ2);
52c6d8aa 881#endif
882
d24a4636 883 /////////////////////////////////////////////////// Stn 5///////////////////////////////////////////////////////////////
52c6d8aa 884
52c6d8aa 885
d24a4636 886 // #ifdef PRINT_BACK
ed5809f7 887 // HLTImportant("\textrap9 : (%f,%f,%f)\n",extrapCh9X,extrapCh9Y,extrapCh9Z);
888 // HLTImportant("\textrap8 : (%f,%f,%f)\n",extrapCh8X,extrapCh8Y,extrapCh8Z);
d24a4636 889 // #endif
52c6d8aa 890
d24a4636 891 nofFrontChPoints = 0; nofBackChPoints = 0;
aa205fc6 892
893 extrapCh9X = trigX1 * extrapCh9Z/trigZ1 ;
894 extrapCh9Y = trigY1 + (trigY2-trigY1) * (extrapCh9Z-trigZ1)/(trigZ2 - trigZ1) ;
d24a4636 895 for( Int_t ipointch9=0;ipointch9<fNofPoints[9];ipointch9++){
fd5b812e 896
aa205fc6 897 if(not fFastTracking){
898 AliHLTMUONUtils::UnpackRecHitFlags(fChPoint[9][ipointch9]->fFlags,chamber,detElemID);
899 if(not fDetElemList[detElemID]){
900 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
901 continue;
902 }
903 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,extrapCh9Z);
52c6d8aa 904
aa205fc6 905 extrapCh9X = trigX1 * extrapCh9Z/trigZ1 ;
906 extrapCh9Y = trigY1 + (trigY2-trigY1) * (extrapCh9Z-trigZ1)/(trigZ2 - trigZ1) ;
907 }
fd5b812e 908
d24a4636 909 if(nofBackChPoints < (fgkMaxNofTracks-1) &&
910 TMath::Abs(extrapCh9X-fChPoint[9][ipointch9]->fX) < maxXDeflectionExtrap &&
911 TMath::Abs(extrapCh9Y-fChPoint[9][ipointch9]->fY)/
912 ((fChPoint[9][ipointch9]->fX * fChPoint[9][ipointch9]->fX) +
913 (fChPoint[9][ipointch9]->fY * fChPoint[9][ipointch9]->fY)) <= extrapRatio ){
52c6d8aa 914
d24a4636 915 distChBack = sqrt((extrapCh9X-fChPoint[9][ipointch9]->fX)*(extrapCh9X-fChPoint[9][ipointch9]->fX)
916 + (extrapCh9Y-fChPoint[9][ipointch9]->fY)*(extrapCh9Y-fChPoint[9][ipointch9]->fY));
917 if(distChBack>circularWindow) continue;
52c6d8aa 918
919#ifdef PRINT_BACK
ed5809f7 920 HLTImportant("\t\tpoints selected in Ch9 : (%f,%f,%f)\n",
d24a4636 921 distChBack,fChPoint[9][ipointch9]->fX,fChPoint[9][ipointch9]->fY,fChPoint[9][ipointch9]->fZ);
52c6d8aa 922#endif
923
d24a4636 924 backIndex[nofBackChPoints++] = ipointch9;
925 }///if point found
926 }/// ch10 loop
52c6d8aa 927
52c6d8aa 928
aa205fc6 929 extrapCh8X = trigX1 * extrapCh8Z/trigZ1 ;
930 extrapCh8Y = trigY1 + (trigY2-trigY1) * (extrapCh8Z-trigZ1)/(trigZ2 - trigZ1) ;
d24a4636 931 for( Int_t ipointch8=0;ipointch8<fNofPoints[8];ipointch8++){
aa205fc6 932
933 if(not fFastTracking){
934 AliHLTMUONUtils::UnpackRecHitFlags(fChPoint[8][ipointch8]->fFlags,chamber,detElemID);
935 if(not fDetElemList[detElemID]){
936 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
937 continue;
938 }
939 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,extrapCh8Z);
52c6d8aa 940
aa205fc6 941 extrapCh8X = trigX1 * extrapCh8Z/trigZ1 ;
942 extrapCh8Y = trigY1 + (trigY2-trigY1) * (extrapCh8Z-trigZ1)/(trigZ2 - trigZ1) ;
d24a4636 943 }
aa205fc6 944
d24a4636 945 if( nofFrontChPoints < (fgkMaxNofTracks-1) &&
946 TMath::Abs(extrapCh8X-fChPoint[8][ipointch8]->fX) < maxXDeflectionExtrap &&
947 TMath::Abs(extrapCh8Y-fChPoint[8][ipointch8]->fY)/
948 ((fChPoint[8][ipointch8]->fX * fChPoint[8][ipointch8]->fX) +
949 (fChPoint[8][ipointch8]->fY * fChPoint[8][ipointch8]->fY)) <= extrapRatio ){
52c6d8aa 950
d24a4636 951 distChFront = sqrt((extrapCh8X-fChPoint[8][ipointch8]->fX)*(extrapCh8X-fChPoint[8][ipointch8]->fX)
952 + (extrapCh8Y-fChPoint[8][ipointch8]->fY)*(extrapCh8Y-fChPoint[8][ipointch8]->fY));
52c6d8aa 953
d24a4636 954 if(distChFront>circularWindow) continue;
52c6d8aa 955
956#ifdef PRINT_BACK
ed5809f7 957 HLTImportant("\t\tpoints selected in Ch8 : (%f,%f,%f)\n",
d24a4636 958 fChPoint[8][ipointch8]->fX,fChPoint[8][ipointch8]->fY,fChPoint[8][ipointch8]->fZ,distChFront);
52c6d8aa 959#endif
960
d24a4636 961 frontIndex[nofFrontChPoints++] = ipointch8;
962 }///if point found
963 }/// ch9 loop
964
965 if(nofBackChPoints==0 and nofFrontChPoints==0) continue;
966
967 minAngle = minAngleWindow;
968 p3.fX = trigX1 ; p3.fY = trigY1 ; p3.fZ = trigZ1 ;
969 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
970 Sub(&p3,fChPoint[9][backIndex[ibackpoint]],&pSeg2);
971 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
972 Sub(fChPoint[9][backIndex[ibackpoint]],fChPoint[8][frontIndex[ifrontpoint]],&pSeg1);
973 anglediff = TMath::RadToDeg()* Angle(&pSeg1,&pSeg2);
974 // #ifdef PRINT_BACK
ed5809f7 975 // HLTImportant("\t\ttracklet-check-St5 : anglediff : %lf, minAngle : %lf\n",anglediff,minAngle);
d24a4636 976 // #endif
977 if(anglediff<minAngle && fNofCells[1]<(fgkMaxNofTracks-1)){
978 st5TrackletFound = true;
979 cells[1][fNofCells[1]].fFirst = frontIndex[ifrontpoint];
980 cells[1][fNofCells[1]].fSecond = backIndex[ibackpoint];
981 fNofCells[1]++ ;
52c6d8aa 982#ifdef PRINT_BACK
ed5809f7 983 HLTImportant("\t\ttracklet-St5 : anglediff : %lf\n",anglediff);
984 HLTImportant("\t\t\tCh9 : (%f,%f,%f)\n",fChPoint[9][backIndex[ibackpoint]]->fX,
d24a4636 985 fChPoint[9][backIndex[ibackpoint]]->fY,fChPoint[9][backIndex[ibackpoint]]->fZ);
ed5809f7 986 HLTImportant("\t\t\tCh8 : (%f,%f,%f)\n",fChPoint[8][frontIndex[ifrontpoint]]->fX,
d24a4636 987 fChPoint[8][frontIndex[ifrontpoint]]->fY,fChPoint[8][frontIndex[ifrontpoint]]->fZ);
52c6d8aa 988#endif
d24a4636 989 }///anglediff condition
990 }///front
991 }///back
52c6d8aa 992
993
994
995
d24a4636 996 /// If tracklet not found, search for the single space point in Ch9 or in Ch8
997 if(!st5TrackletFound){
52c6d8aa 998
d24a4636 999 minAngle = minAngleWindow;
1000 p3.fX = trigX2 ; p3.fY = trigY2 ; p3.fZ = trigZ2 ;
1001 p2.fX = trigX1 ; p2.fY = trigY1 ; p2.fZ = trigZ1 ;
1002 Sub(&p3,&p2,&pSeg2);
52c6d8aa 1003
d24a4636 1004 ///Search in Ch9
1005 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1006 Sub(&p2,fChPoint[9][backIndex[ibackpoint]],&pSeg1);
1007 anglediff = TMath::RadToDeg()*Angle(&pSeg1,&pSeg2);
1008 if(anglediff<minAngle && fNofCells[1]<(fgkMaxNofTracks-1)){
1009 ch9PointFound = true;
1010 cells[1][fNofCells[1]].fFirst = -1;
1011 cells[1][fNofCells[1]].fSecond = backIndex[ibackpoint];
1012 fNofCells[1]++ ;
52c6d8aa 1013#ifdef PRINT_BACK
ed5809f7 1014 HLTImportant("\t\tno st tracklet and single point-Ch9 : anglediff : %lf\n",anglediff);
52c6d8aa 1015#endif
d24a4636 1016 }
1017 }///back
1018
1019 ///Search in Ch8
1020 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
1021 Sub(&p2,fChPoint[8][frontIndex[ifrontpoint]],&pSeg1);
1022 anglediff = TMath::RadToDeg()*Angle(&pSeg1,&pSeg2);
1023 if(anglediff<minAngle && fNofCells[1]<(fgkMaxNofTracks-1)){
1024 ch8PointFound = true;
1025 cells[1][fNofCells[1]].fFirst = frontIndex[ifrontpoint];
1026 cells[1][fNofCells[1]].fSecond = -1;
1027 fNofCells[1]++ ;
52c6d8aa 1028#ifdef PRINT_BACK
ed5809f7 1029 HLTImportant("\t\tno st tracklet and single point-Ch8 : anglediff : %lf\n",anglediff);
52c6d8aa 1030#endif
d24a4636 1031 }
1032 }///front
1033 }///if no tracklets found condition
52c6d8aa 1034
1035#ifdef PRINT_BACK
ed5809f7 1036 HLTImportant("\tnofTracks found after stn 5 : %d\n",fNofCells[1]);
52c6d8aa 1037#endif
1038
d24a4636 1039 if(!st5TrackletFound && !ch9PointFound && !ch8PointFound) continue;
52c6d8aa 1040
d24a4636 1041 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
52c6d8aa 1042
1043
d24a4636 1044 /////////////////////////////////////////////////// Stn 4///////////////////////////////////////////////////////////////
52c6d8aa 1045
d24a4636 1046 // #ifdef PRINT_BACK
ed5809f7 1047 // HLTImportant("\textrap7 : (%f,%f,%f)\n",extrapCh7X,extrapCh7Y,extrapCh7Z);
1048 // HLTImportant("\textrap6 : (%f,%f,%f)\n",extrapCh6X,extrapCh6Y,extrapCh6Z);
d24a4636 1049 // #endif
52c6d8aa 1050
d24a4636 1051 nofFrontChPoints = 0; nofBackChPoints = 0;
aa205fc6 1052
1053 extrapCh7X = trigX1 * extrapCh7Z/trigZ1 ;
1054 extrapCh7Y = trigY1 + (trigY2-trigY1) * (extrapCh7Z-trigZ1)/(trigZ2 - trigZ1) ;
d24a4636 1055 for( Int_t ipointch7=0;ipointch7<fNofPoints[7];ipointch7++){
aa205fc6 1056
1057 if(not fFastTracking){
1058 AliHLTMUONUtils::UnpackRecHitFlags(fChPoint[7][ipointch7]->fFlags,chamber,detElemID);
1059 if(not fDetElemList[detElemID]){
1060 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
1061 continue;
1062 }
1063 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,extrapCh7Z);
fd5b812e 1064
aa205fc6 1065 extrapCh7X = trigX1 * extrapCh7Z/trigZ1 ;
1066 extrapCh7Y = trigY1 + (trigY2-trigY1) * (extrapCh7Z-trigZ1)/(trigZ2 - trigZ1) ;
1067 }
fd5b812e 1068
d24a4636 1069 if( nofBackChPoints < (fgkMaxNofTracks-1) &&
1070 TMath::Abs(extrapCh7X-fChPoint[7][ipointch7]->fX) < maxXDeflectionExtrap &&
1071 TMath::Abs(extrapCh7Y-fChPoint[7][ipointch7]->fY)/
1072 ((fChPoint[7][ipointch7]->fX * fChPoint[7][ipointch7]->fX) +
1073 (fChPoint[7][ipointch7]->fY * fChPoint[7][ipointch7]->fY)) <= extrapRatio ){
52c6d8aa 1074
d24a4636 1075 distChBack = sqrt((extrapCh7X-fChPoint[7][ipointch7]->fX)*(extrapCh7X-fChPoint[7][ipointch7]->fX)
1076 + (extrapCh7Y-fChPoint[7][ipointch7]->fY)*(extrapCh7Y-fChPoint[7][ipointch7]->fY));
52c6d8aa 1077
d24a4636 1078 if(distChBack>circularWindow) continue;
52c6d8aa 1079#ifdef PRINT_BACK
ed5809f7 1080 HLTImportant("\t\tpoints selected in Ch7 : (%f,%f,%f)\n",
d24a4636 1081 fChPoint[7][ipointch7]->fX,fChPoint[7][ipointch7]->fY,fChPoint[7][ipointch7]->fZ);
52c6d8aa 1082#endif
1083
d24a4636 1084 backIndex[nofBackChPoints++] = ipointch7;
1085 }///if point found
1086 }///ch8 loop
52c6d8aa 1087
aa205fc6 1088 extrapCh6X = trigX1 * extrapCh6Z/trigZ1 ;
1089 extrapCh6Y = trigY1 + (trigY2-trigY1) * (extrapCh6Z-trigZ1)/(trigZ2 - trigZ1) ;
d24a4636 1090 for( Int_t ipointch6=0;ipointch6<fNofPoints[6];ipointch6++){
fd5b812e 1091
aa205fc6 1092 if(not fFastTracking){
1093 AliHLTMUONUtils::UnpackRecHitFlags(fChPoint[6][ipointch6]->fFlags,chamber,detElemID);
1094 if(not fDetElemList[detElemID]){
1095 HLTDebug("Invalid detection element : %d, not harmful to HLT chain",detElemID);
1096 continue;
1097 }
1098 fChamberGeometryTransformer->Local2Global(detElemID,0.0,0.0,0.0,tx,ty,extrapCh6Z);
52c6d8aa 1099
52c6d8aa 1100
aa205fc6 1101 extrapCh6X = trigX1 * extrapCh6Z/trigZ1 ;
1102 extrapCh6Y = trigY1 + (trigY2-trigY1) * (extrapCh6Z-trigZ1)/(trigZ2 - trigZ1) ;
1103 }
fd5b812e 1104
d24a4636 1105 if(nofFrontChPoints < (fgkMaxNofTracks-1) &&
1106 TMath::Abs(extrapCh6X-fChPoint[6][ipointch6]->fX) < maxXDeflectionExtrap &&
1107 TMath::Abs(extrapCh6Y-fChPoint[6][ipointch6]->fY)/
1108 ((fChPoint[6][ipointch6]->fX * fChPoint[6][ipointch6]->fX) +
1109 (fChPoint[6][ipointch6]->fY * fChPoint[6][ipointch6]->fY)) <= extrapRatio ){
52c6d8aa 1110
d24a4636 1111 distChFront = sqrt((extrapCh6X-fChPoint[6][ipointch6]->fX)*(extrapCh6X-fChPoint[6][ipointch6]->fX)
1112 + (extrapCh6Y-fChPoint[6][ipointch6]->fY)*(extrapCh6Y-fChPoint[6][ipointch6]->fY));
1113 if(distChFront>circularWindow) continue;
52c6d8aa 1114
1115#ifdef PRINT_BACK
ed5809f7 1116 HLTImportant("\t\tpoints selected in Ch6 : (%f,%f,%f)\n",
d24a4636 1117 fChPoint[6][ipointch6]->fX,fChPoint[6][ipointch6]->fY,fChPoint[6][ipointch6]->fZ);
52c6d8aa 1118#endif
d24a4636 1119 frontIndex[nofFrontChPoints++] = ipointch6;
1120 }///if point found
1121 }/// ch7 loop
1122
1123 if(nofBackChPoints==0 and nofFrontChPoints==0) continue;
1124
1125 minAngle = minAngleWindow;
1126 p3.fX = trigX1 ; p3.fY = trigY1 ; p3.fZ = trigZ1 ;
1127 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1128 Sub(&p3,fChPoint[7][backIndex[ibackpoint]],&pSeg2);
1129 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
1130 Sub(fChPoint[7][backIndex[ibackpoint]],fChPoint[6][frontIndex[ifrontpoint]],&pSeg1);
1131 anglediff = TMath::RadToDeg() * Angle(&pSeg1,&pSeg2);
1132 if(anglediff<minAngle && fNofCells[0]<(fgkMaxNofTracks-1)){
1133 st4TrackletFound = true;
1134 cells[0][fNofCells[0]].fFirst = frontIndex[ifrontpoint];
1135 cells[0][fNofCells[0]].fSecond = backIndex[ibackpoint];
1136 fNofCells[0]++ ;
52c6d8aa 1137#ifdef PRINT_BACK
ed5809f7 1138 HLTImportant("\t\ttracklet-St4 : anglediff : %lf\n",anglediff);
1139 HLTImportant("\t\t\tCh7 : (%f,%f,%f)\n",fChPoint[7][backIndex[ibackpoint]]->fX,
d24a4636 1140 fChPoint[7][backIndex[ibackpoint]]->fY,fChPoint[7][backIndex[ibackpoint]]->fZ);
ed5809f7 1141 HLTImportant("\t\t\tCh6 : (%f,%f,%f)\n",fChPoint[6][frontIndex[ifrontpoint]]->fX,
d24a4636 1142 fChPoint[6][frontIndex[ifrontpoint]]->fY,fChPoint[6][frontIndex[ifrontpoint]]->fZ);
52c6d8aa 1143#endif
d24a4636 1144 }///anglediff condn
1145 }///front
1146 }///back
52c6d8aa 1147
1148
1149
1150
d24a4636 1151 /// If tracklet not found search for the single space point in Ch7 or in Ch6
1152 if(!st4TrackletFound){
52c6d8aa 1153
d24a4636 1154 minAngle = minAngleWindow;
1155 p3.fX = trigX2 ; p3.fY = trigY2 ; p3.fZ = trigZ2 ;
1156 p2.fX = trigX1 ; p2.fY = trigY1 ; p2.fZ = trigZ1 ;
1157 Sub(&p3,&p2,&pSeg2);
52c6d8aa 1158
d24a4636 1159 ///Search in Ch7
1160 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1161 Sub(&p2,fChPoint[7][backIndex[ibackpoint]],&pSeg1);
1162 anglediff = TMath::RadToDeg()*Angle(&pSeg1,&pSeg2);
1163 if(anglediff<minAngle && fNofCells[0]<(fgkMaxNofTracks-1)){
1164 ch7PointFound = true;
1165 cells[0][fNofCells[0]].fFirst = -1;
1166 cells[0][fNofCells[0]].fSecond = backIndex[ibackpoint];
1167 fNofCells[0]++ ;
52c6d8aa 1168#ifdef PRINT_BACK
ed5809f7 1169 HLTImportant("\t\tno st tracklet and single point-Ch7 : anglediff : %lf\n",anglediff);
52c6d8aa 1170#endif
d24a4636 1171 }
1172 }///back
1173
1174 ///Search in Ch6
1175 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
1176 Sub(&p2,fChPoint[6][frontIndex[ifrontpoint]],&pSeg1);
1177 anglediff = TMath::RadToDeg()*Angle(&pSeg1,&pSeg2);
1178 if(anglediff<minAngle && fNofCells[0]<(fgkMaxNofTracks-1)){
1179 ch6PointFound = true;
1180 cells[0][fNofCells[0]].fFirst = frontIndex[ifrontpoint];
1181 cells[0][fNofCells[0]].fSecond = -1;
1182 fNofCells[0]++ ;
52c6d8aa 1183#ifdef PRINT_BACK
ed5809f7 1184 HLTImportant("\t\tno st tracklet and single point-Ch6 : anglediff : %lf\n",anglediff);
52c6d8aa 1185#endif
d24a4636 1186 }
1187 }///front
1188 }///if no tracklets found condition
52c6d8aa 1189
1190#ifdef PRINT_BACK
ed5809f7 1191 HLTImportant("\tnofTracks found after stn 4 : %d\n",fNofCells[0]);
52c6d8aa 1192#endif
1193
d24a4636 1194 if(!st4TrackletFound && !ch7PointFound && !ch6PointFound) continue;
52c6d8aa 1195
d24a4636 1196 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
1197 ;
52c6d8aa 1198
d24a4636 1199 ////////////////////////////////////////////// Analyse and fill trackseg array////////////////////////////////////////
1200 ;
52c6d8aa 1201#ifdef PRINT_BACK
ed5809f7 1202 HLTImportant("\tfNofbackTrackSeg : %d, st5TrackletFound : %d, st4TrackletFound : %d\n",fNofbackTrackSeg,st5TrackletFound,st4TrackletFound);
52c6d8aa 1203#endif
1204
d24a4636 1205 if(st5TrackletFound && st4TrackletFound){
1206
1207 minAngle = minAngleWindow;
1208
1209 for( Int_t itrackletfront=0;itrackletfront<fNofCells[0];itrackletfront++){
1210 index1 = cells[0][itrackletfront].fFirst ;
1211 index2 = cells[0][itrackletfront].fSecond ;
1212 Sub(fChPoint[7][index2],fChPoint[6][index1],&pSeg1);
1213 for( Int_t itrackletback=0;itrackletback<fNofCells[1];itrackletback++){
1214 index3 = cells[1][itrackletback].fFirst ;
1215 index4 = cells[1][itrackletback].fSecond ;
1216 Sub(fChPoint[8][index3],fChPoint[7][index2],&pSeg2);
1217 Sub(fChPoint[9][index4],fChPoint[8][index3],&pSeg3);
1218 anglediff = Angle(&pSeg1,&pSeg2) + Angle(&pSeg2,&pSeg3);
1219 if(anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1220 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1221 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;
1222 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1223 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1224 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1225 minAngle = anglediff;
fd5b812e 1226#ifdef PRINT_BACK
ed5809f7 1227 HLTImportant("\t\ttracklet-St4 and St5 : anglediff : %lf\n",anglediff);
1228 HLTImportant("\t\t\tCh9 : (%f,%f,%f)\n",fChPoint[9][index4]->fX,
d24a4636 1229 fChPoint[9][index4]->fY,fChPoint[9][index4]->fZ);
ed5809f7 1230 HLTImportant("\t\t\tCh8 : (%f,%f,%f)\n",fChPoint[8][index3]->fX,
d24a4636 1231 fChPoint[8][index3]->fY,fChPoint[8][index3]->fZ);
ed5809f7 1232 HLTImportant("\t\t\tCh7 : (%f,%f,%f)\n",fChPoint[7][index2]->fX,
d24a4636 1233 fChPoint[7][index2]->fY,fChPoint[7][index2]->fZ);
ed5809f7 1234 HLTImportant("\t\t\tCh6 : (%f,%f,%f)\n",fChPoint[6][index1]->fX,
d24a4636 1235 fChPoint[6][index1]->fY,fChPoint[6][index1]->fZ);
fd5b812e 1236#endif
d24a4636 1237 }///if minangle
1238 }///for of front ch
1239 }///for loop of back ch
1240
1241 if(minAngle<minAngleWindow)
1242 fNofbackTrackSeg++;
1243
1244 }else if(st5TrackletFound && (ch7PointFound || ch6PointFound)){
1245
1246
1247 nofFrontChPoints = 0; nofBackChPoints = 0;
1248 for( Int_t ifrontpoint=0;ifrontpoint<fNofCells[0];ifrontpoint++){
1249 if(cells[0][ifrontpoint].fFirst==-1 && nofBackChPoints<(fgkMaxNofTracks-1))
1250 backIndex[nofBackChPoints++] = cells[0][ifrontpoint].fSecond;
1251 else if(cells[0][ifrontpoint].fSecond==-1 && nofFrontChPoints<(fgkMaxNofTracks-1))
1252 frontIndex[nofFrontChPoints++] = cells[0][ifrontpoint].fFirst;
1253 }
1254
1255 minAngle = minAngleWindow;
1256 if(nofFrontChPoints>0 && nofBackChPoints>0){
1257 for( Int_t itrackletback=0;itrackletback<fNofCells[1];itrackletback++){
1258 index3 = cells[1][itrackletback].fFirst ;
1259 index4 = cells[1][itrackletback].fSecond ;
1260 Sub(fChPoint[9][index4],fChPoint[8][index3],&pSeg3);
1261 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1262 Sub(fChPoint[8][index3],fChPoint[7][backIndex[ibackpoint]],&pSeg2);
1263 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
1264 Sub(fChPoint[7][backIndex[ibackpoint]],fChPoint[6][frontIndex[ifrontpoint]],&pSeg1);
1265 anglediff = TMath::RadToDeg()*(Angle(&pSeg1,&pSeg2) + Angle(&pSeg2,&pSeg3))/2.0;
1266
1267 if(anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1268 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = frontIndex[ifrontpoint];
1269 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = backIndex[ibackpoint] ;
1270 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1271 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1272 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1273 minAngle = anglediff;
1274 continue;
1275 }
1276
1277 Sub(fChPoint[8][index3],fChPoint[6][frontIndex[ifrontpoint]],&pSeg1);
1278 anglediff1 = TMath::RadToDeg() * Angle(&pSeg1,&pSeg3);
1279 anglediff2 = TMath::RadToDeg() * Angle(&pSeg2,&pSeg3);
1280
1281 if( anglediff1 < anglediff2 && anglediff1<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1282 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = frontIndex[ifrontpoint];
1283 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = -1;
1284 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1285 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1286 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1287 minAngle = anglediff1;
1288 continue;
1289 }
1290
1291 if( anglediff2 < anglediff1 && anglediff2<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1292 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = -1;
1293 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = backIndex[ibackpoint] ;
1294 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1295 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1296 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1297 minAngle = anglediff2;
1298 }
52c6d8aa 1299
d24a4636 1300 }///loop of ifrontpoint
1301 }///loop of ibackpoint
1302 }/// for loop of St5 cells
1303 }else if(nofFrontChPoints>0){
52c6d8aa 1304
d24a4636 1305 for( Int_t itrackletback=0;itrackletback<fNofCells[1];itrackletback++){
1306 index3 = cells[1][itrackletback].fFirst ;
1307 index4 = cells[1][itrackletback].fSecond ;
1308 Sub(fChPoint[9][index4],fChPoint[8][index3],&pSeg3);
52c6d8aa 1309
d24a4636 1310 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
1311 Sub(fChPoint[8][index3],fChPoint[6][frontIndex[ifrontpoint]],&pSeg2);
52c6d8aa 1312
d24a4636 1313 anglediff = TMath::RadToDeg() * Angle(&pSeg2,&pSeg3);
1314 if( anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1315 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = frontIndex[ifrontpoint];
1316 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = -1;
1317 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1318 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1319 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1320 minAngle = anglediff;
1321 }///if anglediff
1322 }///backch loop
1323 }///tracklet loop
1324
1325 }else{ /// if(nofBackChPoints>0){
1326 for( Int_t itrackletback=0;itrackletback<fNofCells[1];itrackletback++){
1327 index3 = cells[1][itrackletback].fFirst ;
1328 index4 = cells[1][itrackletback].fSecond ;
1329 Sub(fChPoint[9][index4],fChPoint[8][index3],&pSeg3);
52c6d8aa 1330
d24a4636 1331 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1332 Sub(fChPoint[8][index3],fChPoint[7][backIndex[ibackpoint]],&pSeg2);
52c6d8aa 1333
d24a4636 1334 anglediff = TMath::RadToDeg() * Angle(&pSeg2,&pSeg3);
52c6d8aa 1335
d24a4636 1336 if( anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1337 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = -1;
1338 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = backIndex[ibackpoint] ;
1339 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = index3;
1340 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = index4;
1341 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1342 minAngle = anglediff;
1343 }///if anglediff
1344 }///backch loop
1345 }///tracklet loop
1346
1347 }///condn for if(nofFrontChPoints>0)
1348
1349 if(minAngle<minAngleWindow)
1350 fNofbackTrackSeg++;
1351
1352 }else if((ch9PointFound || ch8PointFound) && st4TrackletFound){
1353
1354 nofFrontChPoints = 0; nofBackChPoints = 0;
1355 for( Int_t ibackpoint=0;ibackpoint<fNofCells[1];ibackpoint++){
1356 if(cells[1][ibackpoint].fFirst==-1 && nofBackChPoints<(fgkMaxNofTracks-1))
1357 backIndex[nofBackChPoints++] = cells[1][ibackpoint].fSecond;
1358 else if(nofFrontChPoints<(fgkMaxNofTracks-1))
1359 frontIndex[nofFrontChPoints++] = cells[1][ibackpoint].fFirst;
1360 }
1361
1362 minAngle = minAngleWindow;
1363 if(nofFrontChPoints>0 && nofBackChPoints>0){
1364
1365 for( Int_t itrackletfront=0;itrackletfront<fNofCells[0];itrackletfront++){
1366 index1 = cells[0][itrackletfront].fFirst ;
1367 index2 = cells[0][itrackletfront].fSecond ;
52c6d8aa 1368
d24a4636 1369 Sub(fChPoint[7][index2],fChPoint[6][index1],&pSeg1);
52c6d8aa 1370
d24a4636 1371 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
52c6d8aa 1372
d24a4636 1373 Sub(fChPoint[8][frontIndex[ifrontpoint]],fChPoint[7][index2],&pSeg2);
52c6d8aa 1374
d24a4636 1375 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
52c6d8aa 1376
d24a4636 1377 Sub(fChPoint[9][backIndex[ibackpoint]],fChPoint[8][frontIndex[ifrontpoint]],&pSeg3);
52c6d8aa 1378
d24a4636 1379 anglediff = TMath::RadToDeg()*(Angle(&pSeg1,&pSeg2) + Angle(&pSeg2,&pSeg3))/2.0;
1380 if(anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1381 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1382 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;
1383 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = frontIndex[ifrontpoint];
1384 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = backIndex[ibackpoint] ;
1385 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1386 minAngle = anglediff;
1387 continue;
1388 }
52c6d8aa 1389
d24a4636 1390 Sub(fChPoint[9][backIndex[ibackpoint]],fChPoint[7][index2],&pSeg3);
52c6d8aa 1391
d24a4636 1392 anglediff1 = TMath::RadToDeg() * Angle(&pSeg1,&pSeg2);
1393 anglediff2 = TMath::RadToDeg() * Angle(&pSeg1,&pSeg3);
52c6d8aa 1394
d24a4636 1395 if( anglediff1 < anglediff2 && anglediff1<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1396 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1397 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;
1398 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = frontIndex[ifrontpoint];
1399 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = -1;
1400 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1401 minAngle = anglediff1;
1402 continue;
1403 }
52c6d8aa 1404
d24a4636 1405 if( anglediff2 < anglediff1 && anglediff2<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1406 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1407 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;
1408 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = -1;
1409 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = backIndex[ibackpoint] ;
1410 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1411 minAngle = anglediff2;
1412 }
52c6d8aa 1413
d24a4636 1414 }///loop of ifrontpoint
1415 }///loop of ibackpoint
1416 }/// for loop of St5 cells
1417 }else if(nofFrontChPoints>0){
1418
1419 for( Int_t itrackletfront=0;itrackletfront<fNofCells[0];itrackletfront++){
1420 index1 = cells[0][itrackletfront].fFirst ;
1421 index2 = cells[0][itrackletfront].fSecond ;
52c6d8aa 1422
d24a4636 1423 Sub(fChPoint[7][index2],fChPoint[6][index1],&pSeg1);
52c6d8aa 1424
d24a4636 1425 for( Int_t ifrontpoint=0;ifrontpoint<nofFrontChPoints;ifrontpoint++){
52c6d8aa 1426
d24a4636 1427 Sub(fChPoint[8][frontIndex[ifrontpoint]],fChPoint[7][index2],&pSeg2);
52c6d8aa 1428
d24a4636 1429 anglediff = TMath::RadToDeg() * Angle(&pSeg1,&pSeg2);
1430 if( anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1431 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1432 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;;
1433 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = frontIndex[ifrontpoint];
1434 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = -1;
1435 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1436 minAngle = anglediff;
1437 }///if anglediff
1438 }///point loop
1439 }///tracklet loop
1440
1441 }else{ /// if(nofBackChPoints>0){
1442
1443 for( Int_t itrackletfront=0;itrackletfront<fNofCells[0];itrackletfront++){
1444 index1 = cells[0][itrackletfront].fFirst ;
1445 index2 = cells[0][itrackletfront].fSecond ;
1446
1447 Sub(fChPoint[7][index2],fChPoint[6][index1],&pSeg1);
1448
1449 for( Int_t ibackpoint=0;ibackpoint<nofBackChPoints;ibackpoint++){
1450 Sub(fChPoint[9][backIndex[ibackpoint]],fChPoint[6][index1],&pSeg3);
1451 anglediff = TMath::RadToDeg()* Angle(&pSeg1,&pSeg3);
1452 if(anglediff<minAngle && fNofbackTrackSeg<(fgkMaxNofTracks-1)){
1453 fBackTrackSeg[fNofbackTrackSeg].fIndex[0] = index1;
1454 fBackTrackSeg[fNofbackTrackSeg].fIndex[1] = index2;
1455 fBackTrackSeg[fNofbackTrackSeg].fIndex[2] = -1;
1456 fBackTrackSeg[fNofbackTrackSeg].fIndex[3] = backIndex[ibackpoint] ;
1457 fBackTrackSeg[fNofbackTrackSeg].fTrigRec = fChPoint11[itrig]->fId;
1458 minAngle = anglediff;
1459 }
1460 }///backch loop
1461 }///tracklet loop
1462
1463 }///condn for if(nofFrontChPoints>0)
1464
1465 if(minAngle<minAngleWindow)
1466 fNofbackTrackSeg++;
1467
1468 }else if((ch9PointFound || ch8PointFound) && (ch7PointFound || ch6PointFound)){
1469 ///To Do : To be analysed for two points out of four slat chambers
1470 }
1471
1472 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
fd5b812e 1473#ifdef PRINT_BACK
ed5809f7 1474 HLTImportant("\n");
fd5b812e 1475#endif
52c6d8aa 1476
d24a4636 1477 }///trigger loop
52c6d8aa 1478
d24a4636 1479 Float_t meanX1,meanX2,meanY1,meanY2,meanZ1,meanZ2;
fd5b812e 1480
d24a4636 1481 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
52c6d8aa 1482
52c6d8aa 1483#ifdef PRINT_BACK
ed5809f7 1484 HLTImportant("Index : (%d,%d,%d,%d) \n",fBackTrackSeg[ibacktrackseg].fIndex[0],
d24a4636 1485 fBackTrackSeg[ibacktrackseg].fIndex[1],fBackTrackSeg[ibacktrackseg].fIndex[2],
1486 fBackTrackSeg[ibacktrackseg].fIndex[3]);
52c6d8aa 1487#endif
1488
d24a4636 1489 if(fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1 && fBackTrackSeg[ibacktrackseg].fIndex[1]!=-1 ){
1490 meanX1 = (fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fX
1491 + fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fX)/2.0 ;
1492 meanY1 = (fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fY
1493 + fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fY)/2.0 ;
1494 meanZ1 = (fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fZ
1495 + fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fZ)/2.0 ;
1496 }else if(fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1 && fBackTrackSeg[ibacktrackseg].fIndex[1]==-1 ){
1497 meanX1 = fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fX ;
1498 meanY1 = fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fY ;
1499 meanZ1 = fChPoint[6][fBackTrackSeg[ibacktrackseg].fIndex[0]]->fZ ;
1500 }else{
1501 meanX1 = fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fX ;
1502 meanY1 = fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fY ;
1503 meanZ1 = fChPoint[7][fBackTrackSeg[ibacktrackseg].fIndex[1]]->fZ ;
1504 }
1505
1506 if(fBackTrackSeg[ibacktrackseg].fIndex[2]!=-1 && fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1 ){
1507 meanX2 = (fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fX
1508 + fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fX)/2.0 ;
1509 meanY2 = (fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fY
1510 + fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fY)/2.0 ;
1511 meanZ2 = (fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fZ
1512 + fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fZ)/2.0 ;
1513 }else if(fBackTrackSeg[ibacktrackseg].fIndex[2]!=-1 && fBackTrackSeg[ibacktrackseg].fIndex[3]==-1 ){
1514 meanX2 = fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fX ;
1515 meanY2 = fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fY ;
1516 meanZ2 = fChPoint[8][fBackTrackSeg[ibacktrackseg].fIndex[2]]->fZ ;
1517 }else{
1518 meanX2 = fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fX ;
1519 meanY2 = fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fY ;
1520 meanZ2 = fChPoint[9][fBackTrackSeg[ibacktrackseg].fIndex[3]]->fZ ;
1521 }
1522 fExtrapSt3X[ibacktrackseg] = meanX1 + (fgkTrackDetCoordinate[2]-meanZ1)*(meanX2-meanX1)/(meanZ2-meanZ1);
1523 fExtrapSt3Y[ibacktrackseg] = meanY1 + (fgkTrackDetCoordinate[2]-meanZ1)*(meanY2-meanY1)/(meanZ2-meanZ1);
1524 fInclinationBack[ibacktrackseg] = (meanX2-meanX1)/(meanZ2-meanZ1) ;
1525 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
1526 }///backtrigseg loop
fd5b812e 1527
1528#ifdef PRINT_BACK
ed5809f7 1529 HLTImportant("AliHLTMUONFullTracker::SlatTrackSeg()--End\n");
1530 HLTImportant("\n\n");
fd5b812e 1531#endif
1532
d24a4636 1533 return true;
fd5b812e 1534}
1535///__________________________________________________________________________
1536
1537Bool_t AliHLTMUONFullTracker::PrelimMomCalc()
1538{
1539 /// momentum calculation for half tracks
1540
1541 Cluster clus1,clus2;
1542 Float_t xf,yf,thetaDev,zf = 0.5*(AliMUONConstants::DefaultChamberZ(4) + AliMUONConstants::DefaultChamberZ(5));
1543 Float_t p,px,py,pz;
52c6d8aa 1544
1545 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
1546
1547
fd5b812e 1548 Int_t maxIndex = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?3:2;
1549 Int_t maxCh = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?9:8;
1550
1551 Int_t minIndex = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?0:1;
1552 Int_t minCh = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?6:7;
52c6d8aa 1553
fd5b812e 1554 clus2.fX = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fX ;
1555 clus2.fY = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fY ;
1556 clus2.fZ = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fZ ;
52c6d8aa 1557
fd5b812e 1558 clus1.fX = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fX ;
1559 clus1.fY = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fY ;
1560 clus1.fZ = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fZ ;
52c6d8aa 1561
fd5b812e 1562 thetaDev= (1/zf)*(clus1.fY*clus2.fZ - clus2.fY*clus1.fZ)/(clus2.fZ - clus1.fZ);
1563 xf = clus1.fX*zf/clus1.fZ;
1564 yf = clus2.fY - ((clus2.fY - clus1.fY)*(clus2.fZ-zf))/(clus2.fZ - clus1.fZ);
1565 p = 3.0*0.3/thetaDev;
1566 px = p*xf/zf;
1567 py = p*yf/zf;;
1568 pz = sqrt((p*p-(px*px + py*py)));
1569 if(zf<0)pz = -pz;
1570
1571 fHalfTrack[ibacktrackseg].fCharge = (p>0)?-1:+1;
1572 fHalfTrack[ibacktrackseg].fPx = p*xf/zf;
1573 fHalfTrack[ibacktrackseg].fPy = p*yf/zf;
1574 fHalfTrack[ibacktrackseg].fPz = pz;
1575
1576 }
52c6d8aa 1577
1578 return true;
1579}
1580
52c6d8aa 1581///__________________________________________________________________________
1582
1583Bool_t AliHLTMUONFullTracker::QuadTrackSeg()
1584{
1585 ///Find the Track Segments in the Quadrant chambers
fd5b812e 1586 if(fNofPoints[0]==0 and fNofPoints[1]==0){
1587 HLTDebug("No Hits found in Stn4");
1588 return false;
1589 }else if(fNofPoints[2]==0 and fNofPoints[3]==0){
1590 HLTDebug("No Hits found in Stn5");
1591 return false;
1592 }else if(fNofbackTrackSeg==0){
1593 HLTDebug("No Hits found in Stn5 and Stn4 so no tracking is done in quadrants");
1594 return false;
1595 }
1596
52c6d8aa 1597
1598 Float_t meanX1,meanX2,meanY1,meanY2,meanZ1,meanZ2;
1599 Float_t expectSt3X,expectSt3Y,inclinationFront;
1600
1601 AliHLTMUONRecHitStruct pSeg1,pSeg2,pSeg3;
1602 Double_t anglediff;///,anglediff1,anglediff2;
d24a4636 1603 Double_t minAngle = -1.0;
1604 Int_t indexMinAngleFront = -1;
1605 Int_t indexMinAngleBack = -1;
1606 Int_t backIndex = -1;
1607
1608 Int_t ch3CellPoint[fgkMaxNofCellsPerCh],ch2CellPoint[fgkMaxNofCellsPerCh],nofSt2Cells=0;
1609 Int_t ch1CellPoint[fgkMaxNofCellsPerCh],ch0CellPoint[fgkMaxNofCellsPerCh],nofSt1Cells=0;
1610 Bool_t isConnected[fgkMaxNofTracks];
1611
1612 Float_t distDiffX = 4.0; ///simulation result 4.0
1613 Float_t distDiffY = 10.0 ; ///simulation result 4.0
1614 ///float closeToY0 = 10.0; ///simulation result 10.0
1615 Float_t minAngleWindow = 2.0 ; ///simulation result 2.0
1616 Float_t diffDistStX = 25.0; ///simulation result 25.0
1617 Float_t diffDistStY = 75.0; ///simulation result 25.0
1618 Float_t st3WindowX = 40.0 ; ///simulation result 40.0
1619 Float_t st3WindowY = 10.0; ///simulation result 10.0
1620
aa205fc6 1621 if(fFastTracking){
1622 distDiffX = 4.0; ///simulation result 4.0
1623 distDiffY = 4.0 ; ///simulation result 4.0
1624 ///float closeToY0 = 10.0; ///simulation result 10.0
1625 minAngleWindow = 2.0 ; ///simulation result 2.0
1626 diffDistStX = 25.0; ///simulation result 25.0
1627 diffDistStY = 25.0; ///simulation result 25.0
1628 st3WindowX = 40.0 ; ///simulation result 40.0
1629 st3WindowY = 10.0; ///simulation result 10.0
1630 }
1631
d24a4636 1632 /// Float_t inclinationWindow = 0.04; ///inclination window
1633 /// Float_t st3WindowXOp2 = 40.0 ; ///simulation result 40.0
1634 /// Float_t st3WindowYOp2 = 10.0; ///simulation result 10.0
52c6d8aa 1635
1636
1637#ifdef PRINT_FRONT
ed5809f7 1638 HLTImportant("\nAliHLTMUONFullTracker::QuadTrackSeg()--Begin\n\n");
52c6d8aa 1639#endif
52c6d8aa 1640
d24a4636 1641 for( Int_t ibackpoint=0;ibackpoint<fNofPoints[3];ibackpoint++){
1642 for( Int_t ifrontpoint=0;ifrontpoint<fNofPoints[2];ifrontpoint++){
52c6d8aa 1643
d24a4636 1644 if(TMath::Abs(fChPoint[3][ibackpoint]->fX - fChPoint[2][ifrontpoint]->fX)<distDiffX
1645 && TMath::Abs(fChPoint[3][ibackpoint]->fY - fChPoint[2][ifrontpoint]->fY)<distDiffY){
fd5b812e 1646
52c6d8aa 1647
d24a4636 1648 /// if((TMath::Abs(fChPoint[3][ibackpoint]->fY) > closeToY0 &&
1649 /// TMath::Abs(fChPoint[3][ibackpoint]->fY) < TMath::Abs(fChPoint[2][ifrontpoint]->fY)) ||
1650 /// nofSt2Cells >= (fgkMaxNofCellsPerCh-1)) continue;
52c6d8aa 1651
d24a4636 1652 if(nofSt2Cells >= (fgkMaxNofCellsPerCh-1)) continue;
52c6d8aa 1653
d24a4636 1654 // #ifdef PRINT_FRONT
ed5809f7 1655 // HLTImportant("\t\t\tCh3 : %d, (%f,%f,%f)\n",
d24a4636 1656 // nofSt2Cells,fChPoint[3][ibackpoint]->fX,fChPoint[3][ibackpoint]->fY,fChPoint[3][ibackpoint]->fZ);
ed5809f7 1657 // HLTImportant("\t\t\tCh2 :(%f,%f,%f)\n\n",
d24a4636 1658 // fChPoint[2][ifrontpoint]->fX,fChPoint[2][ifrontpoint]->fY,fChPoint[2][ifrontpoint]->fZ);
1659 // #endif
52c6d8aa 1660
d24a4636 1661 ch3CellPoint[nofSt2Cells] = ibackpoint;
1662 ch2CellPoint[nofSt2Cells] = ifrontpoint;
1663 nofSt2Cells++;
52c6d8aa 1664
d24a4636 1665 }///if point found
1666 }///frontch
1667 }///backch
fd5b812e 1668
d24a4636 1669 if(nofSt2Cells==0){
1670 HLTDebug("No Hits found in Stn2");
1671 return false;
1672 }
fd5b812e 1673
d24a4636 1674 for( Int_t ibackpoint=0;ibackpoint<fNofPoints[1];ibackpoint++){
1675 for( Int_t ifrontpoint=0;ifrontpoint<fNofPoints[0];ifrontpoint++){
fd5b812e 1676
d24a4636 1677 if(TMath::Abs(fChPoint[1][ibackpoint]->fX - fChPoint[0][ifrontpoint]->fX)< distDiffX
1678 && TMath::Abs(fChPoint[1][ibackpoint]->fY - fChPoint[0][ifrontpoint]->fY)<distDiffY){
fd5b812e 1679
52c6d8aa 1680
d24a4636 1681 /// if((TMath::Abs(fChPoint[1][ibackpoint]->fY) > closeToY0 &&
1682 /// TMath::Abs(fChPoint[1][ibackpoint]->fY) < TMath::Abs(fChPoint[0][ifrontpoint]->fY)) ||
1683 /// nofSt1Cells >= (fgkMaxNofCellsPerCh-1)) continue;
52c6d8aa 1684
d24a4636 1685 if(nofSt1Cells >= (fgkMaxNofCellsPerCh-1)) continue;
52c6d8aa 1686
1687
d24a4636 1688 // #ifdef PRINT_FRONT
ed5809f7 1689 // HLTImportant("\t\t\tCh1 : %d, (%f,%f,%f)\n",
d24a4636 1690 // nofSt1Cells,fChPoint[1][ibackpoint]->fX,fChPoint[1][ibackpoint]->fY,fChPoint[1][ibackpoint]->fZ);
ed5809f7 1691 // HLTImportant("\t\t\tCh0 :(%f,%f,%f)\n\n",
d24a4636 1692 // fChPoint[0][ifrontpoint]->fX,fChPoint[0][ifrontpoint]->fY,fChPoint[0][ifrontpoint]->fZ);
1693 // #endif
1694 ch1CellPoint[nofSt1Cells] = ibackpoint;
1695 ch0CellPoint[nofSt1Cells] = ifrontpoint;
1696 nofSt1Cells++;
1697 }///if point found
1698 }///frontch
1699 }///backch
1700 if(nofSt1Cells==0){
1701 HLTDebug("No Hits found in Stn1");
1702 return false;
1703 }
52c6d8aa 1704
1705#ifdef PRINT_FRONT
ed5809f7 1706 HLTImportant("\tnofSt1Cells : %d, nofSt2Cells : %d\n",nofSt1Cells,nofSt2Cells);
52c6d8aa 1707#endif
fd5b812e 1708
d24a4636 1709 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++)
1710 isConnected[ibacktrackseg] = false;
1711
1712
1713 ///First Check for Tracklets in two front stations
1714 for( Int_t itrackletfront=0;itrackletfront<nofSt1Cells;itrackletfront++){
1715
1716 Sub(fChPoint[1][ch1CellPoint[itrackletfront]],fChPoint[0][ch0CellPoint[itrackletfront]],&pSeg1);
1717
1718 minAngle = minAngleWindow;
1719 indexMinAngleBack = -1;
1720 indexMinAngleFront = -1;
1721
1722 meanX1 = (fChPoint[0][ch0CellPoint[itrackletfront]]->fX
1723 + fChPoint[1][ch1CellPoint[itrackletfront]]->fX)/2.0 ;
1724 meanY1 = (fChPoint[0][ch0CellPoint[itrackletfront]]->fY
1725 + fChPoint[1][ch1CellPoint[itrackletfront]]->fY)/2.0 ;
1726 meanZ1 = (fChPoint[0][ch0CellPoint[itrackletfront]]->fZ
1727 + fChPoint[1][ch1CellPoint[itrackletfront]]->fZ)/2.0 ;
1728
1729 for( Int_t itrackletback=0;itrackletback<nofSt2Cells;itrackletback++){
1730 // #ifdef PRINT_FRONT
1731 // cout<<"\tBefore "
1732 // <<TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fX - fChPoint[1][ch1CellPoint[itrackletfront]]->fX)
1733 // <<"\t"
1734 // <<TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fY - fChPoint[1][ch1CellPoint[itrackletfront]]->fY)
1735 // <<endl;
1736 // #endif
1737 if(TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fX -
1738 fChPoint[1][ch1CellPoint[itrackletfront]]->fX) > diffDistStX ||
1739 TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fY -
1740 fChPoint[1][ch1CellPoint[itrackletfront]]->fY) > diffDistStY ) continue;
1741
1742 meanX2 = (fChPoint[2][ch2CellPoint[itrackletback]]->fX
1743 + fChPoint[3][ch3CellPoint[itrackletback]]->fX)/2.0 ;
1744 meanY2 = (fChPoint[2][ch2CellPoint[itrackletback]]->fY
1745 + fChPoint[3][ch3CellPoint[itrackletback]]->fY)/2.0 ;
1746 meanZ2 = (fChPoint[2][ch2CellPoint[itrackletback]]->fZ
1747 + fChPoint[3][ch3CellPoint[itrackletback]]->fZ)/2.0 ;
1748
1749 expectSt3X = meanX2 + (fgkTrackDetCoordinate[2]-meanZ2)*(meanX2-meanX1)/(meanZ2-meanZ1);
1750 expectSt3Y = meanY2 + (fgkTrackDetCoordinate[2]-meanZ2)*(meanY2-meanY1)/(meanZ2-meanZ1);
1751 inclinationFront = (meanX2-meanX1)/(meanZ2-meanZ1) ;
1752
1753 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
52c6d8aa 1754
d24a4636 1755 if( /// TMath::Abs(inclinationBack[ibacktrackseg]-inclinationFront)<0.04 &&
fd5b812e 1756 TMath::Abs((expectSt3X-fExtrapSt3X[ibacktrackseg])) < st3WindowX
1757 && TMath::Abs((expectSt3Y-fExtrapSt3Y[ibacktrackseg])) < st3WindowY){
52c6d8aa 1758
d24a4636 1759 Sub(fChPoint[2][ch2CellPoint[itrackletback]],fChPoint[1][ch1CellPoint[itrackletfront]],&pSeg2);
1760 Sub(fChPoint[3][ch3CellPoint[itrackletback]],fChPoint[2][ch2CellPoint[itrackletback]],&pSeg3);
52c6d8aa 1761
d24a4636 1762 anglediff = TMath::RadToDeg()* (Angle(&pSeg1,&pSeg2) + Angle(&pSeg2,&pSeg3));
1763
1764 // #ifdef PRINT_FRONT
ed5809f7 1765 // HLTImportant("\t\t\tanglediff : %lf\n",anglediff);
d24a4636 1766 // #endif
1767 if(anglediff<minAngle){
1768 minAngle = anglediff;
1769 indexMinAngleBack = itrackletback;
1770 indexMinAngleFront = itrackletfront;
1771 backIndex = ibacktrackseg;
1772 isConnected[ibacktrackseg] = true;
1773 }
1774 }///matching tracklet
1775 }///for loop of backtrackseg
52c6d8aa 1776
d24a4636 1777 }///for of back ch
52c6d8aa 1778
d24a4636 1779 if(minAngle < minAngleWindow && indexMinAngleFront!=-1
1780 && indexMinAngleBack!=-1 && fNoffrontTrackSeg<(fgkMaxNofTracks-1)){
52c6d8aa 1781
d24a4636 1782 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[0] = ch0CellPoint[indexMinAngleFront];
1783 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[1] = ch1CellPoint[indexMinAngleFront];
1784 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[2] = ch2CellPoint[indexMinAngleBack];
1785 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[3] = ch3CellPoint[indexMinAngleBack];
52c6d8aa 1786
d24a4636 1787 fBackToFront[backIndex][fNofConnectedfrontTrackSeg[backIndex]++] = fNoffrontTrackSeg;
1788 fNoffrontTrackSeg++;
1789 }///condition to find valid tracklet
52c6d8aa 1790
d24a4636 1791 }///for loop of front ch
52c6d8aa 1792
d24a4636 1793 Int_t nofNCfBackTrackSeg = 0;
1794 Int_t ncfBackTrackSeg[fgkMaxNofTracks];
52c6d8aa 1795
d24a4636 1796 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++)
1797 if(!isConnected[ibacktrackseg])
1798 ncfBackTrackSeg[nofNCfBackTrackSeg++] = ibacktrackseg;
1799 else
1800 fNofConnected++;
52c6d8aa 1801
1802
1803#ifdef PRINT_FRONT
ed5809f7 1804 HLTImportant("\tfNofConnected : %d, nofNCfBackTrackSeg : %d\n",fNofConnected,nofNCfBackTrackSeg);
1805 HLTImportant("\tfNofPoints[3] : %d, fNofPoints[2] : %d\n",fNofPoints[3],fNofPoints[2]);
d24a4636 1806 if(nofNCfBackTrackSeg==0){
ed5809f7 1807 HLTImportant("All fBackTrackSegs are connected with fFrontTrackSegs, no need to search further\n");
1808 HLTImportant("AliHLTMUONFullTracker::QuadTrackSeg()--End\n\n");
d24a4636 1809 }
52c6d8aa 1810#endif
1811
d24a4636 1812 if(nofNCfBackTrackSeg==0) return true;
fd5b812e 1813
1814
d24a4636 1815 ///Next Check for tracklet in Stn1 and space point in Stn2
1816 Bool_t isbackpoint=false,isfrontpoint=false;
1817 for( Int_t itrackletfront=0;itrackletfront<nofSt1Cells;itrackletfront++){
1818 Sub(fChPoint[1][ch1CellPoint[itrackletfront]],fChPoint[0][ch0CellPoint[itrackletfront]],&pSeg1);
1819 minAngle = minAngleWindow;
1820 indexMinAngleBack = -1;
1821 indexMinAngleFront = -1;
fd5b812e 1822
d24a4636 1823 for( Int_t ibackpoint=0;ibackpoint<fNofPoints[3];ibackpoint++){
1824 if(/// hasCh3Cells[ibackpoint] == true &&
fd5b812e 1825 TMath::Abs(fChPoint[3][ibackpoint]->fX -
1826 fChPoint[1][ch1CellPoint[itrackletfront]]->fX) > diffDistStX ||
1827 TMath::Abs(fChPoint[3][ibackpoint]->fY -
1828 fChPoint[1][ch1CellPoint[itrackletfront]]->fY) > diffDistStY ) continue;
1829
d24a4636 1830 expectSt3X = fChPoint[3][ibackpoint]->fX + (fgkTrackDetCoordinate[2] - fChPoint[3][ibackpoint]->fZ)*
1831 (fChPoint[3][ibackpoint]->fX - fChPoint[1][ch1CellPoint[itrackletfront]]->fX)/
1832 (fChPoint[3][ibackpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ);
1833 expectSt3Y = fChPoint[3][ibackpoint]->fY + (fgkTrackDetCoordinate[2] - fChPoint[3][ibackpoint]->fZ)*
1834 (fChPoint[3][ibackpoint]->fY - fChPoint[1][ch1CellPoint[itrackletfront]]->fY)/
1835 (fChPoint[3][ibackpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ);
1836 inclinationFront = (fChPoint[3][ibackpoint]->fX - fChPoint[1][ch1CellPoint[itrackletfront]]->fX)/
1837 (fChPoint[3][ibackpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ) ;
fd5b812e 1838
d24a4636 1839 for( Int_t ibacktrackseg=0;ibacktrackseg<nofNCfBackTrackSeg;ibacktrackseg++){
52c6d8aa 1840
d24a4636 1841 if(/// TMath::Abs(inclinationBack[ncfBackTrackSeg[ibacktrackseg]]-inclinationFront)< inclinationWindow &&
fd5b812e 1842 TMath::Abs((expectSt3X-fExtrapSt3X[ncfBackTrackSeg[ibacktrackseg]])) < st3WindowX
1843 && TMath::Abs((expectSt3Y-fExtrapSt3Y[ncfBackTrackSeg[ibacktrackseg]])) < st3WindowY){
52c6d8aa 1844
d24a4636 1845 Sub(fChPoint[3][ibackpoint],fChPoint[1][ch1CellPoint[itrackletfront]],&pSeg2);
52c6d8aa 1846
d24a4636 1847 anglediff = TMath::RadToDeg()* Angle(&pSeg1,&pSeg2) ;
1848 // #ifdef PRINT_FRONT
ed5809f7 1849 // HLTImportant("\t\t annglediff(Ch4) : %lf\n",anglediff);
d24a4636 1850 // #endif
1851 if(anglediff<minAngle){
1852 minAngle = anglediff;
1853 indexMinAngleBack = ibackpoint;
1854 indexMinAngleFront = itrackletfront;
1855 backIndex = ncfBackTrackSeg[ibacktrackseg];
1856 isConnected[ncfBackTrackSeg[ibacktrackseg]] = true;
1857 isbackpoint = true;
1858 isfrontpoint = false;
1859 }///angle min condn
1860 }///matching tracklet
1861 }///loop on not found back trackseg
1862 }///backpoint loop
1863
1864 for( Int_t ifrontpoint=0;ifrontpoint<fNofPoints[2];ifrontpoint++){
1865 if(/// hasCh2Cells[ifrontpoint] == true &&
fd5b812e 1866 TMath::Abs(fChPoint[2][ifrontpoint]->fX -
1867 fChPoint[1][ch1CellPoint[itrackletfront]]->fX) > diffDistStX ||
1868 TMath::Abs(fChPoint[2][ifrontpoint]->fY -
1869 fChPoint[1][ch1CellPoint[itrackletfront]]->fY) > diffDistStY ) continue;
1870
d24a4636 1871 expectSt3X = fChPoint[2][ifrontpoint]->fX + (fgkTrackDetCoordinate[2] - fChPoint[2][ifrontpoint]->fZ)*
1872 (fChPoint[2][ifrontpoint]->fX - fChPoint[1][ch1CellPoint[itrackletfront]]->fX)/
1873 (fChPoint[2][ifrontpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ);
1874 expectSt3Y = fChPoint[2][ifrontpoint]->fY + (fgkTrackDetCoordinate[2] - fChPoint[2][ifrontpoint]->fZ)*
1875 (fChPoint[2][ifrontpoint]->fY - fChPoint[1][ch1CellPoint[itrackletfront]]->fY)/
1876 (fChPoint[2][ifrontpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ);
1877 inclinationFront = (fChPoint[2][ifrontpoint]->fX - fChPoint[1][ch1CellPoint[itrackletfront]]->fX)/
1878 (fChPoint[2][ifrontpoint]->fZ - fChPoint[1][ch1CellPoint[itrackletfront]]->fZ) ;
1879 // #ifdef PRINT_FRONT
ed5809f7 1880 // HLTImportant("\t\texpectSt3X : %f, expectSt3Y : %f, inclinationFront : %f\n",expectSt3X,expectSt3Y,inclinationFront);
d24a4636 1881 // #endif
1882
1883 for( Int_t ibacktrackseg=0;ibacktrackseg<nofNCfBackTrackSeg;ibacktrackseg++){
52c6d8aa 1884
d24a4636 1885 if( /// TMath::Abs(inclinationBack[ncfBackTrackSeg[ibacktrackseg]]-inclinationFront)< inclinationWindow &&
fd5b812e 1886 TMath::Abs((expectSt3X-fExtrapSt3X[ncfBackTrackSeg[ibacktrackseg]])) < st3WindowX
1887 && TMath::Abs((expectSt3Y-fExtrapSt3Y[ncfBackTrackSeg[ibacktrackseg]])) < st3WindowY){
52c6d8aa 1888
d24a4636 1889 Sub(fChPoint[2][ifrontpoint],fChPoint[1][ch1CellPoint[itrackletfront]],&pSeg2);
1890
1891 anglediff = TMath::RadToDeg()* Angle(&pSeg1,&pSeg2) ;
1892 // #ifdef PRINT_FRONT
ed5809f7 1893 // HLTImportant("\t\t annglediff(Ch3) : %lf\n",anglediff);
d24a4636 1894 // #endif
1895 if(anglediff<minAngle){
1896 minAngle = anglediff;
1897 indexMinAngleBack = ifrontpoint;
1898 indexMinAngleFront = itrackletfront;
1899 backIndex = ncfBackTrackSeg[ibacktrackseg];
1900 isConnected[ncfBackTrackSeg[ibacktrackseg]] = true;
1901 isbackpoint = false;
1902 isfrontpoint = true;
1903
1904 }///angle min condn
1905 }///matching tracklet
1906 }///loop on not found back trackseg
1907 }///backpoint loop
1908
1909 if(minAngle < minAngleWindow && indexMinAngleFront!=-1 &&
1910 indexMinAngleBack!=-1 && fNoffrontTrackSeg<(fgkMaxNofTracks-1)){
1911
1912 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[0] = ch0CellPoint[indexMinAngleFront];
1913 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[1] = ch1CellPoint[indexMinAngleFront];
1914 if(isfrontpoint && !isbackpoint){
1915 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[2] = indexMinAngleBack;
1916 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[3] = -1;
1917 }else{
1918 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[2] = -1;
1919 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[3] = indexMinAngleBack;
1920 }
1921 fBackToFront[backIndex][fNofConnectedfrontTrackSeg[backIndex]++] = fNoffrontTrackSeg;
1922 fNoffrontTrackSeg++;
1923 }///condition to find valid tracklet
1924
1925 }///front ch
1926
1927
1928 Int_t nofSNCfBackTrackSeg = 0;
1929 Int_t sncfBackTrackSeg[fgkMaxNofTracks];
1930
1931 for( Int_t ibacktrackseg=0;ibacktrackseg<nofNCfBackTrackSeg;ibacktrackseg++)
1932 if(!isConnected[ncfBackTrackSeg[ibacktrackseg]])
1933 sncfBackTrackSeg[nofSNCfBackTrackSeg++] = ncfBackTrackSeg[ibacktrackseg];
1934 else
1935 fNofConnected++;
52c6d8aa 1936
1937#ifdef PRINT_FRONT
ed5809f7 1938 HLTImportant("\tfNofConnected : %d, nofSNCfBackTrackSeg : %d\n",fNofConnected,nofSNCfBackTrackSeg);
d24a4636 1939 if(nofSNCfBackTrackSeg==0){
ed5809f7 1940 HLTImportant("All fBackTrackSegs are connected with fFrontTrackSegs, no need to search further\n");
1941 HLTImportant("AliHLTMUONFullTracker::QuadTrackSeg()--End\n\n");
d24a4636 1942 }
52c6d8aa 1943#endif
1944
d24a4636 1945 if(nofSNCfBackTrackSeg==0) return true;
fd5b812e 1946
d24a4636 1947 ///Last Check for tracklet in Stn2 and space point in Stn1
1948 for( Int_t itrackletback=0;itrackletback<nofSt2Cells;itrackletback++){
1949 Sub(fChPoint[3][ch3CellPoint[itrackletback]],fChPoint[2][ch2CellPoint[itrackletback]],&pSeg1);
1950 minAngle = minAngleWindow ;
1951 indexMinAngleBack = -1;
1952 indexMinAngleFront = -1;
fd5b812e 1953
d24a4636 1954 for( Int_t ibackpoint=0;ibackpoint<fNofPoints[1];ibackpoint++){
1955 if(/// hasCh1Cells[ibackpoint] == true &&
fd5b812e 1956 TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fX -
1957 fChPoint[1][ibackpoint]->fX) > diffDistStX ||
1958 TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fY -
1959 fChPoint[1][ibackpoint]->fY) > diffDistStY) continue;
1960
d24a4636 1961 expectSt3X = fChPoint[2][ch2CellPoint[itrackletback]]->fX + (fgkTrackDetCoordinate[2] - fChPoint[2][ch2CellPoint[itrackletback]]->fZ)*
1962 (fChPoint[2][ch2CellPoint[itrackletback]]->fX - fChPoint[1][ibackpoint]->fX)/
1963 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[1][ibackpoint]->fZ);
1964 expectSt3Y = fChPoint[2][ch2CellPoint[itrackletback]]->fY + (fgkTrackDetCoordinate[2] - fChPoint[2][ch2CellPoint[itrackletback]]->fZ)*
1965 (fChPoint[2][ch2CellPoint[itrackletback]]->fY - fChPoint[1][ibackpoint]->fY)/
1966 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[1][ibackpoint]->fZ);
1967 inclinationFront = (fChPoint[2][ch2CellPoint[itrackletback]]->fX - fChPoint[1][ibackpoint]->fX)/
1968 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[1][ibackpoint]->fZ) ;
fd5b812e 1969
d24a4636 1970 for( Int_t ibacktrackseg=0;ibacktrackseg<nofSNCfBackTrackSeg;ibacktrackseg++){
52c6d8aa 1971
d24a4636 1972 if(/// TMath::Abs(inclinationBack[sncfBackTrackSeg[ibacktrackseg]]-inclinationFront)<inclinationWindow &&
fd5b812e 1973 TMath::Abs((expectSt3X-fExtrapSt3X[sncfBackTrackSeg[ibacktrackseg]])) < st3WindowX
1974 && TMath::Abs((expectSt3Y-fExtrapSt3Y[sncfBackTrackSeg[ibacktrackseg]])) < st3WindowY ){
52c6d8aa 1975
d24a4636 1976 Sub(fChPoint[2][ch2CellPoint[itrackletback]],fChPoint[1][ibackpoint],&pSeg2);
52c6d8aa 1977
d24a4636 1978 anglediff = TMath::RadToDeg()* Angle(&pSeg1,&pSeg2) ;
1979 if(anglediff<minAngle){
1980 minAngle = anglediff;
1981 indexMinAngleBack = itrackletback;
1982 indexMinAngleFront = ibackpoint;
1983 backIndex = sncfBackTrackSeg[ibacktrackseg];
1984 isConnected[sncfBackTrackSeg[ibacktrackseg]] = true;
1985 isbackpoint = true;
1986 isfrontpoint = false;
1987 }///angle min condn
1988 }///matching tracklet
1989 }///loop on not found back trackseg
1990 }///backpoint loop
1991
1992 for( Int_t ifrontpoint=0;ifrontpoint<fNofPoints[0];ifrontpoint++){
1993 if(/// hasCh0Cells[ifrontpoint] == true &&
fd5b812e 1994 TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fX -
1995 fChPoint[0][ifrontpoint]->fX) > diffDistStX ||
1996 TMath::Abs(fChPoint[2][ch2CellPoint[itrackletback]]->fY -
1997 fChPoint[0][ifrontpoint]->fY) > diffDistStY ) continue;
1998
d24a4636 1999 expectSt3X = fChPoint[2][ch2CellPoint[itrackletback]]->fX + (fgkTrackDetCoordinate[2] - fChPoint[2][ch2CellPoint[itrackletback]]->fZ)*
2000 (fChPoint[2][ch2CellPoint[itrackletback]]->fX - fChPoint[0][ifrontpoint]->fX)/
2001 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[0][ifrontpoint]->fZ);
2002 expectSt3Y = fChPoint[2][ch2CellPoint[itrackletback]]->fY + (fgkTrackDetCoordinate[2] - fChPoint[2][ch2CellPoint[itrackletback]]->fZ)*
2003 (fChPoint[2][ch2CellPoint[itrackletback]]->fY - fChPoint[0][ifrontpoint]->fY)/
2004 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[0][ifrontpoint]->fZ);
2005 inclinationFront = (fChPoint[2][ch2CellPoint[itrackletback]]->fX - fChPoint[0][ifrontpoint]->fX)/
2006 (fChPoint[2][ch2CellPoint[itrackletback]]->fZ - fChPoint[0][ifrontpoint]->fZ) ;
fd5b812e 2007
d24a4636 2008 for( Int_t ibacktrackseg=0;ibacktrackseg<nofSNCfBackTrackSeg;ibacktrackseg++){
fd5b812e 2009
d24a4636 2010 if(/// TMath::Abs(inclinationBack[sncfBackTrackSeg[ibacktrackseg]]-inclinationFront)<inclinationWindow &&
fd5b812e 2011 TMath::Abs((expectSt3X-fExtrapSt3X[sncfBackTrackSeg[ibacktrackseg]])) < st3WindowX
2012 && TMath::Abs((expectSt3Y-fExtrapSt3Y[sncfBackTrackSeg[ibacktrackseg]])) < st3WindowY ){
52c6d8aa 2013
d24a4636 2014 Sub(fChPoint[2][ch2CellPoint[itrackletback]],fChPoint[0][ifrontpoint],&pSeg2);
52c6d8aa 2015
d24a4636 2016 anglediff = TMath::RadToDeg() * Angle(&pSeg1,&pSeg2) ;
2017 if(anglediff<minAngle){
2018 minAngle = anglediff;
2019 indexMinAngleBack = itrackletback;
2020 indexMinAngleFront = ifrontpoint;
2021 backIndex = sncfBackTrackSeg[ibacktrackseg];
2022 isConnected[sncfBackTrackSeg[ibacktrackseg]] = true;
2023 isbackpoint = false;
2024 isfrontpoint = true;
52c6d8aa 2025
d24a4636 2026 }///angle min condn
2027 }///matching tracklet
2028 }///loop on not found back trackseg
2029 }///backpoint loop
2030
2031 if(minAngle < minAngleWindow && indexMinAngleFront!=-1 &&
2032 indexMinAngleBack!=-1 && fNoffrontTrackSeg<(fgkMaxNofTracks-1)){
2033
2034 if(isfrontpoint){
2035 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[0] = indexMinAngleFront;
2036 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[1] = -1;
2037 }else{
2038 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[0] = -1;
2039 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[1] = indexMinAngleFront;
2040 }
52c6d8aa 2041
d24a4636 2042 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[2] = ch2CellPoint[indexMinAngleBack];
2043 fFrontTrackSeg[fNoffrontTrackSeg].fIndex[3] = ch3CellPoint[indexMinAngleBack];
52c6d8aa 2044
d24a4636 2045 fBackToFront[backIndex][fNofConnectedfrontTrackSeg[backIndex]++] = fNoffrontTrackSeg;
2046 fNoffrontTrackSeg++;
2047 }///condition to find valid tracklet
52c6d8aa 2048
d24a4636 2049 }///front ch
52c6d8aa 2050
d24a4636 2051 for( Int_t ibacktrackseg=0;ibacktrackseg<nofSNCfBackTrackSeg;ibacktrackseg++)
2052 if(isConnected[sncfBackTrackSeg[ibacktrackseg]])
2053 fNofConnected++;
52c6d8aa 2054
2055#ifdef PRINT_FRONT
ed5809f7 2056 HLTImportant("\tfNofConnected : %d\n",fNofConnected);
2057 HLTImportant("Three spacepoints are found in fFrontTrackSegs\n");
2058 HLTImportant("AliHLTMUONFullTracker::QuadTrackSeg()--End\n\n");
52c6d8aa 2059#endif
2060
2061
d24a4636 2062 return true;
52c6d8aa 2063}
2064
2065///__________________________________________________________________________
2066
2067Double_t AliHLTMUONFullTracker::KalmanFilter(AliMUONTrackParam &trackParamAtCluster, Cluster *cluster)
2068{
2069 //// Compute new track parameters and their covariances including new cluster using kalman filter
2070 //// return the additional track chi2
2071
2072#ifdef PRINT_DETAIL_KALMAN
ed5809f7 2073 HLTImportant("AliHLTMUONFullTracker::KalmanFilter()--Begin\n\n");
52c6d8aa 2074#endif
2075
2076
2077 /// Get actual track parameters (p)
d24a4636 2078 TMatrixD param(trackParamAtCluster.GetParameters());
52c6d8aa 2079#ifdef PRINT_DETAIL_KALMAN
d24a4636 2080 Info("\tKalmanFilter","param.Print() [p]");
2081 param.Print();
ed5809f7 2082 HLTImportant("GetZ : %lf\n",trackParamAtCluster.GetZ());
52c6d8aa 2083#endif
2084
2085
2086
d24a4636 2087 /// Get new cluster parameters (m)
2088 ///AliMUONVCluster *cluster = trackParamAtCluster.GetClusterPtr();
2089 TMatrixD clusterParam(5,1);
2090 clusterParam.Zero();
2091 /// clusterParam(0,0) = cluster->GetX();
2092 /// clusterParam(2,0) = cluster->GetY();
2093 clusterParam(0,0) = cluster->fX;
2094 clusterParam(2,0) = cluster->fY;
52c6d8aa 2095#ifdef PRINT_DETAIL_KALMAN
d24a4636 2096 Info("\tKalmanFilter","clusterParam.Print() [m]");
2097 clusterParam.Print();
52c6d8aa 2098#endif
2099
2100
2101
d24a4636 2102 /// Compute the actual parameter weight (W)
2103 TMatrixD paramWeight(trackParamAtCluster.GetCovariances());
52c6d8aa 2104#ifdef PRINT_DETAIL_KALMAN
d24a4636 2105 Info("\tKalmanFilter","Covariance : [C]");
2106 paramWeight.Print();
52c6d8aa 2107#endif
2108
d24a4636 2109 if (paramWeight.Determinant() != 0) {
2110 paramWeight.Invert();
2111 } else {
2112 Warning("KalmanFilter"," Determinant = 0");
a038ee7d 2113 return 1.0e10;
d24a4636 2114 }
52c6d8aa 2115
2116#ifdef PRINT_DETAIL_KALMAN
d24a4636 2117 Info("\tKalmanFilter","Weight Matrix inverse of Covariance [W = C^-1]");
2118 paramWeight.Print();
52c6d8aa 2119#endif
2120
2121
d24a4636 2122 /// Compute the new cluster weight (U)
2123 TMatrixD clusterWeight(5,5);
2124 clusterWeight.Zero();
2125 clusterWeight(0,0) = 1. / cluster->fErrX2;
2126 clusterWeight(2,2) = 1. / cluster->fErrY2;
52c6d8aa 2127#ifdef PRINT_DETAIL_KALMAN
d24a4636 2128 Info("\tKalmanFilter","clusterWeight.Print() [U]");
ed5809f7 2129 HLTImportant("\tErrX2 : %lf, ErrY2 : %lf\n",cluster->fErrX2,cluster->fErrY2);
d24a4636 2130 clusterWeight.Print();
52c6d8aa 2131#endif
2132
2133
2134
2135
d24a4636 2136 /// Compute the new parameters covariance matrix ( (W+U)^-1 )
2137 TMatrixD newParamCov(paramWeight,TMatrixD::kPlus,clusterWeight);
52c6d8aa 2138#ifdef PRINT_DETAIL_KALMAN
d24a4636 2139 Info("\tKalmanFilter","newParamCov.Print() [(W+U)]");
2140 newParamCov.Print();
52c6d8aa 2141#endif
d24a4636 2142 if (newParamCov.Determinant() != 0) {
2143 newParamCov.Invert();
2144 } else {
2145 Warning("RunKalmanFilter"," Determinant = 0");
a038ee7d 2146 return 1.0e10;
d24a4636 2147 }
52c6d8aa 2148#ifdef PRINT_DETAIL_KALMAN
d24a4636 2149 Info("\tKalmanFilter","newParamCov.Print() [(W+U)^-1] (new covariances[W] for trackParamAtCluster)");
2150 newParamCov.Print();
52c6d8aa 2151#endif
2152
d24a4636 2153 /// Save the new parameters covariance matrix
2154 trackParamAtCluster.SetCovariances(newParamCov);
52c6d8aa 2155
d24a4636 2156 /// Compute the new parameters (p' = ((W+U)^-1)U(m-p) + p)
2157 TMatrixD tmp(clusterParam,TMatrixD::kMinus,param);
2158 TMatrixD tmp2(clusterWeight,TMatrixD::kMult,tmp); /// U(m-p)
2159 TMatrixD newParam(newParamCov,TMatrixD::kMult,tmp2); /// ((W+U)^-1)U(m-p)
2160 newParam += param; /// ((W+U)^-1)U(m-p) + p
52c6d8aa 2161#ifdef PRINT_DETAIL_KALMAN
d24a4636 2162 Info("\tKalmanFilter","newParam.Print() [p' = ((W+U)^-1)U(m-p) + p] (new parameters[p] for trackParamAtCluster)");
2163 newParam.Print();
52c6d8aa 2164#endif
2165
d24a4636 2166 /// Save the new parameters
2167 trackParamAtCluster.SetParameters(newParam);
ed5809f7 2168 /// HLTImportant(Form("Pt : %lf\n",TMath::Sqrt(trackParamAtCluster.Px()*trackParamAtCluster.Px() +
d24a4636 2169 /// trackParamAtCluster.Py()*trackParamAtCluster.Py())));
52c6d8aa 2170
2171
d24a4636 2172 /// Compute the additional chi2 (= ((p'-p)^-1)W(p'-p) + ((p'-m)^-1)U(p'-m))
2173 tmp = newParam; /// p'
2174 tmp -= param; /// (p'-p)
2175 TMatrixD tmp3(paramWeight,TMatrixD::kMult,tmp); /// W(p'-p)
2176 TMatrixD addChi2Track(tmp,TMatrixD::kTransposeMult,tmp3); /// ((p'-p)^-1)W(p'-p)
2177 tmp = newParam; /// p'
2178 tmp -= clusterParam; /// (p'-m)
2179 TMatrixD tmp4(clusterWeight,TMatrixD::kMult,tmp); /// U(p'-m)
2180 addChi2Track += TMatrixD(tmp,TMatrixD::kTransposeMult,tmp4); /// ((p'-p)^-1)W(p'-p) + ((p'-m)^-1)U(p'-m)
52c6d8aa 2181
2182#ifdef PRINT_DETAIL_KALMAN
d24a4636 2183 Info("\tKalmanFilter","addChi2Track.Print() [additional chi2 = ((p'-p)^-1)W(p'-p) + ((p'-m)^-1)U(p'-m)))]");
2184 addChi2Track.Print();
ed5809f7 2185 HLTImportant("AliHLTMUONFullTracker::KalmanFilter()--End\n\n");
52c6d8aa 2186#endif
2187
2188
d24a4636 2189 return addChi2Track(0,0);
52c6d8aa 2190}
2191
2192///__________________________________________________________________________
2193
2194Double_t AliHLTMUONFullTracker::TryOneCluster(const AliMUONTrackParam &trackParam, Cluster* cluster,
fd5b812e 2195 AliMUONTrackParam &trackParamAtCluster, Bool_t updatePropagator)
52c6d8aa 2196{
2197 //// Test the compatibility between the track and the cluster (using trackParam's covariance matrix):
2198 //// return the corresponding Chi2
2199 //// return trackParamAtCluster
2200
2201 /// extrapolate track parameters and covariances at the z position of the tested cluster
2202 trackParamAtCluster = trackParam;
2203 AliMUONTrackExtrap::ExtrapToZCov(&trackParamAtCluster, cluster->fZ, updatePropagator);
2204
2205 /// set pointer to cluster into trackParamAtCluster
d24a4636 2206 ///trackParamAtCluster.SetClusterPtr(cluster);
52c6d8aa 2207
d24a4636 2208 /// Set differences between trackParam and cluster in the bending and non bending directions
2209 Double_t dX = cluster->fX - trackParamAtCluster.GetNonBendingCoor();
2210 Double_t dY = cluster->fY - trackParamAtCluster.GetBendingCoor();
52c6d8aa 2211
d24a4636 2212 /// Calculate errors and covariances
2213 const TMatrixD& kParamCov = trackParamAtCluster.GetCovariances();
2214 Double_t sigmaX2 = kParamCov(0,0) + cluster->fErrX2;
2215 Double_t sigmaY2 = kParamCov(2,2) + cluster->fErrY2;
52c6d8aa 2216
d24a4636 2217 /// Compute chi2
2218 return dX * dX / sigmaX2 + dY * dY / sigmaY2;
52c6d8aa 2219
2220}
2221
2222///__________________________________________________________________________
2223
fd5b812e 2224Bool_t AliHLTMUONFullTracker::TryOneClusterFast(const AliMUONTrackParam &trackParam, const Cluster* cluster)
52c6d8aa 2225{
2226 //// Test the compatibility between the track and the cluster
2227 //// given the track resolution + the maximum-distance-to-track value
2228 //// and assuming linear propagation of the track:
2229 //// return kTRUE if they are compatibles
2230
eb672fd1 2231 Float_t sigmaCutForTracking = 6.0;
2232 Float_t maxNonBendingDistanceToTrack = 1.0;
2233 Float_t maxBendingDistanceToTrack = 1.0;
52c6d8aa 2234
2235 Double_t dZ = cluster->fZ - trackParam.GetZ();
2236 Double_t dX = cluster->fX - (trackParam.GetNonBendingCoor() + trackParam.GetNonBendingSlope() * dZ);
2237 Double_t dY = cluster->fY - (trackParam.GetBendingCoor() + trackParam.GetBendingSlope() * dZ);
2238 const TMatrixD& kParamCov = trackParam.GetCovariances();
2239 Double_t errX2 = kParamCov(0,0) + dZ * dZ * kParamCov(1,1) + 2. * dZ * kParamCov(0,1);
2240 Double_t errY2 = kParamCov(2,2) + dZ * dZ * kParamCov(3,3) + 2. * dZ * kParamCov(2,3);
2241
eb672fd1 2242 Double_t dXmax = sigmaCutForTracking * TMath::Sqrt(errX2) +
2243 maxNonBendingDistanceToTrack;
2244 Double_t dYmax = sigmaCutForTracking * TMath::Sqrt(errY2) +
2245 maxBendingDistanceToTrack;
52c6d8aa 2246
2247 if (TMath::Abs(dX) > dXmax || TMath::Abs(dY) > dYmax) return kFALSE;
2248
2249 return kTRUE;
2250
2251}
2252///__________________________________________________________________________
2253
2254void AliHLTMUONFullTracker::PropagateTracks(Double_t charge, Float_t& px, Float_t& py, Float_t& pz,
fd5b812e 2255 Float_t& xr, Float_t& yr, Float_t& zr, Float_t zprop)
52c6d8aa 2256{
2257 ///
2258 /// propagate in magnetic field between hits of indices i1 and i2
2259 ///
2260
2261 Double_t vect[7], vout[7];
2262 Double_t step = -5.0;
2263 Double_t zMax = zprop+.5;
2264
2265 vect[0] = xr;
2266 vect[1] = yr;
2267 vect[2] = zr;
2268 vect[6] = sqrt((px)*(px) + (py)*(py) + (pz)*(pz));
2269 vect[3] = px/vect[6];
2270 vect[4] = py/vect[6];
2271 vect[5] = pz/vect[6];
2272 /// cout<<"vec[2] : "<<vect[2]<<", zMax : "<<zMax<<endl;
2273
d24a4636 2274 Int_t nSteps = 0;
2275 while ((vect[2] < zMax) && (nSteps < 10000)) {
2276 nSteps++;
2277 /// OneStepRungekutta(charge, step, vect, vout);
2278 OneStepHelix3(charge,step,vect,vout);
2279 ///SetPoint(fCount,vout[0],vout[1],vout[2]);
2280 ///fCount++;
ed5809f7 2281 /// HLTImportant("(x,y,z) : (%f,%f,%f)\n",vout[0],vout[1],vout[2]);
d24a4636 2282 for (Int_t i = 0; i < 7; i++) {
2283 vect[i] = vout[i];
52c6d8aa 2284 }
d24a4636 2285 }
52c6d8aa 2286
d24a4636 2287 xr = vect[0];
2288 yr = vect[1];
2289 zr = vect[2];
52c6d8aa 2290
d24a4636 2291 px = vect[3]*vect[6];
2292 py = vect[4]*vect[6];
2293 pz = vect[5]*vect[6];
2294 return;
52c6d8aa 2295}
2296
2297///__________________________________________________________________________
2298
eb672fd1 2299void AliHLTMUONFullTracker::OneStepHelix3(Double_t field, Double_t step, const Double_t *vect, Double_t *vout) const
52c6d8aa 2300{
2301 //// <pre>
2302 //// ******************************************************************
2303 //// * *
2304 //// * Tracking routine in a constant field oriented *
2305 //// * along axis 3 *
2306 //// * Tracking is performed with a conventional *
2307 //// * helix step method *
2308 //// * *
2309 //// * ==>Called by : USER, GUSWIM *
2310 //// * Authors R.Brun, M.Hansroul ********* *
2311 //// * Rewritten V.Perevoztchikov
2312 //// * *
2313 //// ******************************************************************
2314 //// </pre>
2315
2316 Double_t hxp[3];
2317 Double_t h4, hp, rho, tet;
2318 Double_t sint, sintt, tsint, cos1t;
2319 Double_t f1, f2, f3, f4, f5, f6;
2320
2321 const Int_t kix = 0;
2322 const Int_t kiy = 1;
2323 const Int_t kiz = 2;
2324 const Int_t kipx = 3;
2325 const Int_t kipy = 4;
2326 const Int_t kipz = 5;
2327 const Int_t kipp = 6;
2328
2329 const Double_t kec = 2.9979251e-4;
2330
2331 ///
d24a4636 2332 /// ------------------------------------------------------------------
2333 ///
2334 /// units are kgauss,centimeters,gev/c
2335 ///
2336 vout[kipp] = vect[kipp];
2337 h4 = field * kec;
2338
2339 hxp[0] = - vect[kipy];
2340 hxp[1] = + vect[kipx];
52c6d8aa 2341
d24a4636 2342 hp = vect[kipz];
2343
2344 rho = -h4/vect[kipp];
2345 tet = rho * step;
2346 if (TMath::Abs(tet) > 0.15) {
2347 sint = TMath::Sin(tet);
2348 sintt = (sint/tet);
2349 tsint = (tet-sint)/tet;
2350 cos1t = 2.* TMath::Sin(0.5*tet) * TMath::Sin(0.5*tet)/tet;
2351 } else {
2352 tsint = tet*tet/36.;
2353 sintt = (1. - tsint);
2354 sint = tet*sintt;
2355 cos1t = 0.5*tet;
2356 }
52c6d8aa 2357
d24a4636 2358 f1 = step * sintt;
2359 f2 = step * cos1t;
2360 f3 = step * tsint * hp;
2361 f4 = -tet*cos1t;
2362 f5 = sint;
2363 f6 = tet * cos1t * hp;
52c6d8aa 2364
d24a4636 2365 vout[kix] = vect[kix] + f1*vect[kipx] + f2*hxp[0];
2366 vout[kiy] = vect[kiy] + f1*vect[kipy] + f2*hxp[1];
2367 vout[kiz] = vect[kiz] + f1*vect[kipz] + f3;
52c6d8aa 2368
d24a4636 2369 vout[kipx] = vect[kipx] + f4*vect[kipx] + f5*hxp[0];
2370 vout[kipy] = vect[kipy] + f4*vect[kipy] + f5*hxp[1];
2371 vout[kipz] = vect[kipz] + f4*vect[kipz] + f6;
52c6d8aa 2372
d24a4636 2373 return;
52c6d8aa 2374}
2375
2376///__________________________________________________________________________
2377
2378///______________________________________________________________________________
2379
2380void AliHLTMUONFullTracker::OneStepRungekutta(Double_t charge, Double_t step,
fd5b812e 2381 const Double_t* vect, Double_t* vout)
52c6d8aa 2382{
2383 //// ******************************************************************
2384 //// * *
2385 //// * Runge-Kutta method for tracking a particle through a magnetic *
2386 //// * field. Uses Nystroem algorithm (See Handbook Nat. Bur. of *
2387 //// * Standards, procedure 25.5.20) *
2388 //// * *
2389 //// * Input parameters *
2390 //// * CHARGE Particle charge *
2391 //// * STEP Step size *
2392 //// * VECT Initial co-ords,direction cosines,momentum *
2393 //// * Output parameters *
2394 //// * VOUT Output co-ords,direction cosines,momentum *
2395 //// * User routine called *
2396 //// * CALL GUFLD(X,F) *
2397 //// * *
2398 //// * ==>Called by : <USER>, GUSWIM *
2399 //// * Authors R.Brun, M.Hansroul ********* *
2400 //// * V.Perevoztchikov (CUT STEP implementation) *
2401 //// * *
2402 //// * *
2403 //// ******************************************************************
2404
2405 Double_t h2, h4, f[4];
2406 Double_t xyzt[3], a, b, c, ph,ph2;
2407 Double_t secxs[4],secys[4],seczs[4],hxp[3];
2408 Double_t g1, g2, g3, g4, g5, g6, ang2, dxt, dyt, dzt;
2409 Double_t est, at, bt, ct, cba;
2410 Double_t f1, f2, f3, f4, rho, tet, hnorm, hp, rho1, sint, cost;
2411
2412 Double_t x;
2413 Double_t y;
2414 Double_t z;
2415
2416 Double_t xt;
2417 Double_t yt;
2418 Double_t zt;
2419
2420 Double_t maxit = 1992;
2421 Double_t maxcut = 11;
2422
2423 const Double_t kdlt = 1e-4;
2424 const Double_t kdlt32 = kdlt/32.;
2425 const Double_t kthird = 1./3.;
2426 const Double_t khalf = 0.5;
2427 const Double_t kec = 2.9979251e-4;
2428
2429 const Double_t kpisqua = 9.86960440109;
2430 const Int_t kix = 0;
2431 const Int_t kiy = 1;
2432 const Int_t kiz = 2;
2433 const Int_t kipx = 3;
2434 const Int_t kipy = 4;
2435 const Int_t kipz = 5;
2436
2437 /// *.
d24a4636 2438 /// *. ------------------------------------------------------------------
2439 /// *.
2440 /// * this constant is for units cm,gev/c and kgauss
2441 /// *
2442 Int_t iter = 0;
2443 Int_t ncut = 0;
2444 for(Int_t j = 0; j < 7; j++)
2445 vout[j] = vect[j];
2446
2447 Double_t pinv = kec * charge / vect[6];
2448 Double_t tl = 0.;
2449 Double_t h = step;
2450 Double_t rest;
2451
2452
2453 do {
2454 rest = step - tl;
2455 if (TMath::Abs(h) > TMath::Abs(rest)) h = rest;
2456 ///cmodif: call gufld(vout,f) changed into:
2457 TGeoGlobalMagField::Instance()->Field(vout,f);
2458
2459 /// *
2460 /// * start of integration
2461 /// *
2462 x = vout[0];
2463 y = vout[1];
2464 z = vout[2];
2465 a = vout[3];
2466 b = vout[4];
2467 c = vout[5];
2468
2469 h2 = khalf * h;
2470 h4 = khalf * h2;
2471 ph = pinv * h;
2472 ph2 = khalf * ph;
2473 secxs[0] = (b * f[2] - c * f[1]) * ph2;
2474 secys[0] = (c * f[0] - a * f[2]) * ph2;
2475 seczs[0] = (a * f[1] - b * f[0]) * ph2;
2476 ang2 = (secxs[0]*secxs[0] + secys[0]*secys[0] + seczs[0]*seczs[0]);
2477 if (ang2 > kpisqua) break;
2478
2479 dxt = h2 * a + h4 * secxs[0];
2480 dyt = h2 * b + h4 * secys[0];
2481 dzt = h2 * c + h4 * seczs[0];
2482 xt = x + dxt;
2483 yt = y + dyt;
2484 zt = z + dzt;
2485 /// *
2486 /// * second intermediate point
52c6d8aa 2487 /// *
52c6d8aa 2488
d24a4636 2489 est = TMath::Abs(dxt) + TMath::Abs(dyt) + TMath::Abs(dzt);
2490 if (est > h) {
2491 if (ncut++ > maxcut) break;
2492 h *= khalf;
2493 continue;
2494 }
2495
2496 xyzt[0] = xt;
2497 xyzt[1] = yt;
2498 xyzt[2] = zt;
2499
2500 ///cmodif: call gufld(xyzt,f) changed into:
2501 TGeoGlobalMagField::Instance()->Field(xyzt,f);
2502
2503 at = a + secxs[0];
2504 bt = b + secys[0];
2505 ct = c + seczs[0];
2506
2507 secxs[1] = (bt * f[2] - ct * f[1]) * ph2;
2508 secys[1] = (ct * f[0] - at * f[2]) * ph2;
2509 seczs[1] = (at * f[1] - bt * f[0]) * ph2;
2510 at = a + secxs[1];
2511 bt = b + secys[1];
2512 ct = c + seczs[1];
2513 secxs[2] = (bt * f[2] - ct * f[1]) * ph2;
2514 secys[2] = (ct * f[0] - at * f[2]) * ph2;
2515 seczs[2] = (at * f[1] - bt * f[0]) * ph2;
2516 dxt = h * (a + secxs[2]);
2517 dyt = h * (b + secys[2]);
2518 dzt = h * (c + seczs[2]);
2519 xt = x + dxt;
2520 yt = y + dyt;
2521 zt = z + dzt;
2522 at = a + 2.*secxs[2];
2523 bt = b + 2.*secys[2];
2524 ct = c + 2.*seczs[2];
2525
2526 est = TMath::Abs(dxt)+TMath::Abs(dyt)+TMath::Abs(dzt);
2527 if (est > 2.*TMath::Abs(h)) {
2528 if (ncut++ > maxcut) break;
2529 h *= khalf;
2530 continue;
2531 }
2532
2533 xyzt[0] = xt;
2534 xyzt[1] = yt;
2535 xyzt[2] = zt;
2536
2537 ///cmodif: call gufld(xyzt,f) changed into:
2538 TGeoGlobalMagField::Instance()->Field(xyzt,f);
2539
2540 z = z + (c + (seczs[0] + seczs[1] + seczs[2]) * kthird) * h;
2541 y = y + (b + (secys[0] + secys[1] + secys[2]) * kthird) * h;
2542 x = x + (a + (secxs[0] + secxs[1] + secxs[2]) * kthird) * h;
2543
2544 secxs[3] = (bt*f[2] - ct*f[1])* ph2;
2545 secys[3] = (ct*f[0] - at*f[2])* ph2;
2546 seczs[3] = (at*f[1] - bt*f[0])* ph2;
2547 a = a+(secxs[0]+secxs[3]+2. * (secxs[1]+secxs[2])) * kthird;
2548 b = b+(secys[0]+secys[3]+2. * (secys[1]+secys[2])) * kthird;
2549 c = c+(seczs[0]+seczs[3]+2. * (seczs[1]+seczs[2])) * kthird;
52c6d8aa 2550
d24a4636 2551 est = TMath::Abs(secxs[0]+secxs[3] - (secxs[1]+secxs[2]))
2552 + TMath::Abs(secys[0]+secys[3] - (secys[1]+secys[2]))
2553 + TMath::Abs(seczs[0]+seczs[3] - (seczs[1]+seczs[2]));
2554
2555 if (est > kdlt && TMath::Abs(h) > 1.e-4) {
2556 if (ncut++ > maxcut) break;
2557 h *= khalf;
2558 continue;
2559 }
2560
2561 ncut = 0;
2562 /// * if too many iterations, go to helix
2563 if (iter++ > maxit) break;
2564
2565 tl += h;
2566 if (est < kdlt32)
2567 h *= 2.;
2568 cba = 1./ TMath::Sqrt(a*a + b*b + c*c);
2569 vout[0] = x;
2570 vout[1] = y;
2571 vout[2] = z;
2572 vout[3] = cba*a;
2573 vout[4] = cba*b;
2574 vout[5] = cba*c;
2575 rest = step - tl;
2576 if (step < 0.) rest = -rest;
2577 if (rest < 1.e-5*TMath::Abs(step)) return;
2578
2579 } while(1);
2580
2581 /// angle too big, use helix
2582
2583 f1 = f[0];
2584 f2 = f[1];
2585 f3 = f[2];
2586 f4 = TMath::Sqrt(f1*f1+f2*f2+f3*f3);
2587 rho = -f4*pinv;
2588 tet = rho * step;
2589
2590 hnorm = 1./f4;
2591 f1 = f1*hnorm;
2592 f2 = f2*hnorm;
2593 f3 = f3*hnorm;
2594
2595 hxp[0] = f2*vect[kipz] - f3*vect[kipy];
2596 hxp[1] = f3*vect[kipx] - f1*vect[kipz];
2597 hxp[2] = f1*vect[kipy] - f2*vect[kipx];
2598
2599 hp = f1*vect[kipx] + f2*vect[kipy] + f3*vect[kipz];
2600
2601 rho1 = 1./rho;
2602 sint = TMath::Sin(tet);
2603 cost = 2.*TMath::Sin(khalf*tet)*TMath::Sin(khalf*tet);
2604
2605 g1 = sint*rho1;
2606 g2 = cost*rho1;
2607 g3 = (tet-sint) * hp*rho1;
2608 g4 = -cost;
2609 g5 = sint;
2610 g6 = cost * hp;
2611
2612 vout[kix] = vect[kix] + g1*vect[kipx] + g2*hxp[0] + g3*f1;
2613 vout[kiy] = vect[kiy] + g1*vect[kipy] + g2*hxp[1] + g3*f2;
2614 vout[kiz] = vect[kiz] + g1*vect[kipz] + g2*hxp[2] + g3*f3;
2615
2616 vout[kipx] = vect[kipx] + g4*vect[kipx] + g5*hxp[0] + g6*f1;
2617 vout[kipy] = vect[kipy] + g4*vect[kipy] + g5*hxp[1] + g6*f2;
2618 vout[kipz] = vect[kipz] + g4*vect[kipz] + g5*hxp[2] + g6*f3;
2619
2620 return;
52c6d8aa 2621}
2622
2623///______________________________________________________________________________
2624
2625
2626Bool_t AliHLTMUONFullTracker::SelectFront()
2627{
fd5b812e 2628 /// Select the front trackseg connected with back trackseg
eb672fd1 2629
52c6d8aa 2630 Cluster clus1,clus2;
2631 Int_t minIndex=0,maxIndex=0;
2632 Int_t minCh=0,maxCh=0;
2633 Int_t ifronttrackseg=0;
2634
2635 Float_t xf,yf,zf,thetaDev;
2636 Float_t p,spx,spy,spz,px,py,pz,sx,sy,sz,x,y,z;
2637 Double_t charge;
2638 Float_t dist,mindist;
2639 Int_t frontsegIndex ;
2640
2641 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
2642
2643 if(fNofConnectedfrontTrackSeg[ibacktrackseg]<=0) continue;
2644
2645 /// if(fBackTrackSeg[ibacktrackseg].fIndex[2]==-1 || fBackTrackSeg[ibacktrackseg].fIndex[3]==-1) continue;
2646
d24a4636 2647 ifronttrackseg = fBackToFront[ibacktrackseg][fNofConnectedfrontTrackSeg[ibacktrackseg]-1];
52c6d8aa 2648
2649#ifdef PRINT_SELECT
ed5809f7 2650 HLTImportant("AliHLTMUONFullTracker::SelectFront()--Begin\n\n");
2651 HLTImportant("\tbacktrack : %d is connected with : %d front tracks\n",
d24a4636 2652 ibacktrackseg,fNofConnectedfrontTrackSeg[ibacktrackseg]);
52c6d8aa 2653#endif
2654
d24a4636 2655 maxIndex = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?3:2;
2656 maxCh = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?9:8;
52c6d8aa 2657
d24a4636 2658 minIndex = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?0:1;
2659 minCh = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?6:7;
52c6d8aa 2660
2661
2662
d24a4636 2663 clus2.fX = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fX ;
2664 clus2.fY = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fY ;
2665 clus2.fZ = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fZ ;
52c6d8aa 2666
d24a4636 2667 clus1.fX = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fX ;
2668 clus1.fY = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fY ;
2669 clus1.fZ = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fZ ;
52c6d8aa 2670
d24a4636 2671 zf = 0.5*(AliMUONConstants::DefaultChamberZ(4) + AliMUONConstants::DefaultChamberZ(5));
2672 thetaDev= (1/zf)*(clus1.fY*clus2.fZ - clus2.fY*clus1.fZ)/(clus2.fZ - clus1.fZ);
2673 xf = clus1.fX*zf/clus1.fZ;
2674 yf = clus2.fY - ((clus2.fY - clus1.fY)*(clus2.fZ-zf))/(clus2.fZ - clus1.fZ);
2675 p = 3.0*0.3/thetaDev;
2676 charge = (p>0)?-1:+1;
52c6d8aa 2677
d24a4636 2678 spx = p*xf/zf;
2679 spy = p*yf/zf;
2680 spz = sqrt((p*p-(spx*spx + spy*spy)));
52c6d8aa 2681
d24a4636 2682 if(zf<0)spz = -spz;
52c6d8aa 2683
d24a4636 2684 sx = clus1.fX; sy = clus1.fY; sz = clus1.fZ;
2685 mindist = 200000.0;
2686 frontsegIndex = -1;
2687 for( Int_t iconnected=0;iconnected<fNofConnectedfrontTrackSeg[ibacktrackseg];iconnected++){
52c6d8aa 2688
d24a4636 2689 ifronttrackseg = fBackToFront[ibacktrackseg][iconnected];
52c6d8aa 2690
d24a4636 2691 minIndex = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
2692 minCh = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
2693 clus1.fX = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fX ;
2694 clus1.fY = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fY ;
2695 clus1.fZ = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fZ ;
52c6d8aa 2696
d24a4636 2697 x = sx; y = sy; z = sz;
2698 px = spx; py = spy; pz = spz;
2699 PropagateTracks(charge,px,py,pz,x,y,z,clus1.fZ);
52c6d8aa 2700
d24a4636 2701 dist = sqrt((clus1.fX-x)*(clus1.fX-x) +
2702 (clus1.fY-y)*(clus1.fY-y));
2703 if(dist<mindist){
2704 mindist = dist;
2705 frontsegIndex = ifronttrackseg;
2706 }
2707 }///for loop on all connected front track segs
52c6d8aa 2708
d24a4636 2709 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
2710 ///have to check later
2711 if(frontsegIndex==-1) continue;
52c6d8aa 2712
d24a4636 2713 fBackToFront[ibacktrackseg][fNofConnectedfrontTrackSeg[ibacktrackseg]++] = frontsegIndex;
fd5b812e 2714
d24a4636 2715 /// fTrackParam[ibacktrackseg] = trackParam;
52c6d8aa 2716
52c6d8aa 2717#ifdef PRINT_SELECT
ed5809f7 2718 HLTImportant("\tbacktrack : %d is connected with : %d front tracks\n",
d24a4636 2719 ibacktrackseg,fNofConnectedfrontTrackSeg[ibacktrackseg]);
ed5809f7 2720 HLTImportant("AliHLTMUONFullTracker::SelectFront()--End\n\n");
52c6d8aa 2721
2722#endif
2723
2724 }///backtrackSeg loop
2725
2726
d24a4636 2727 return true;
52c6d8aa 2728}
2729
2730///__________________________________________________________________________
2731
2732Bool_t AliHLTMUONFullTracker::KalmanChi2Test()
2733{
2734
2735 /// Kalman Chi2 test for trach segments selection
2736
2737 Cluster clus1,clus2;
2738 Int_t minIndex=0,maxIndex=0;
2739 Int_t minCh=0,maxCh=0;
2740 Int_t ifronttrackseg=0;
2741 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
2742
2743 if(fNofConnectedfrontTrackSeg[ibacktrackseg]<=0) continue;
2744
2745 /// if(fBackTrackSeg[ibacktrackseg].fIndex[2]==-1 || fBackTrackSeg[ibacktrackseg].fIndex[3]==-1) continue;
2746
d24a4636 2747 ifronttrackseg = fBackToFront[ibacktrackseg][fNofConnectedfrontTrackSeg[ibacktrackseg]-1];
52c6d8aa 2748
2749#ifdef PRINT_KALMAN
ed5809f7 2750 HLTImportant("AliHLTMUONFullTracker::KalmanChi2Test()--Begin\n\n");
2751 HLTImportant("\tbacktrack : %d is connected with : %d front tracks, front track index %d\n",
d24a4636 2752 ibacktrackseg,fNofConnectedfrontTrackSeg[ibacktrackseg],ifronttrackseg);
52c6d8aa 2753#endif
d24a4636 2754 maxIndex = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?3:2;
2755 maxCh = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?9:8;
52c6d8aa 2756
d24a4636 2757 minIndex = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?0:1;
2758 minCh = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?6:7;
52c6d8aa 2759
2760
d24a4636 2761 clus2.fX = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fX ;
2762 clus2.fY = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fY ;
2763 clus2.fZ = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fZ ;
2764 clus2.fErrX2 = 0.020736;
2765 clus2.fErrY2 = 0.000100;
52c6d8aa 2766
d24a4636 2767 clus1.fX = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fX ;
2768 clus1.fY = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fY ;
2769 clus1.fZ = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fZ ;
2770 clus1.fErrX2 = 0.020736;
2771 clus1.fErrY2 = 0.000100;
2772
2773
2774 AliMUONTrackParam trackParam;
2775 Double_t dZ = clus2.fZ - clus1.fZ;
2776 trackParam.SetNonBendingCoor(clus2.fX);
2777 trackParam.SetBendingCoor(clus2.fY);
2778 trackParam.SetZ(clus2.fZ);
2779 trackParam.SetNonBendingSlope((clus2.fX - clus1.fX) / dZ);
2780 trackParam.SetBendingSlope((clus2.fY - clus1.fY) / dZ);
2781 Double_t bendingImpact = clus1.fY - clus1.fZ * trackParam.GetBendingSlope();
2782 Double_t inverseBendingMomentum = 1.
2783 / AliMUONTrackExtrap::GetBendingMomentumFromImpactParam(bendingImpact);
2784 trackParam.SetInverseBendingMomentum(inverseBendingMomentum);
52c6d8aa 2785
2786#ifdef PRINT_KALMAN
ed5809f7 2787 HLTImportant("\t\tCh%d : (%f,%f,%f)\n",maxCh,clus2.fX,clus2.fY,clus2.fZ);
2788 HLTImportant("\t\tCh%d : (%f,%f,%f)\n",minCh,clus1.fX,clus1.fY,clus1.fZ);
2789 HLTImportant("\t\tFor minCh : %d, GetBeMom : %lf\n",minCh,trackParam.GetInverseBendingMomentum());
52c6d8aa 2790#endif
2791
d24a4636 2792 /// trackParam->SetClusterPtr(clus[8]);
52c6d8aa 2793
d24a4636 2794 /// => Reset track parameter covariances at last clus (as if the other cluss did not exist)
2795 TMatrixD lastParamCov(5,5);
2796 lastParamCov.Zero();
2797 lastParamCov(0,0) = clus1.fErrX2;
2798 lastParamCov(1,1) = 100. * ( clus2.fErrX2 + clus1.fErrX2 ) / dZ / dZ;
2799 lastParamCov(2,2) = clus1.fErrY2;
2800 lastParamCov(3,3) = 100. * ( clus2.fErrY2 + clus1.fErrY2 ) / dZ / dZ;
2801 lastParamCov(4,4) = ((10.0*10.0 + (clus2.fZ * clus2.fZ * clus1.fErrY2 +
2802 clus1.fZ * clus1.fZ * clus2.fErrY2)
2803 / dZ / dZ) /bendingImpact / bendingImpact +
2804 0.1 * 0.1) * inverseBendingMomentum * inverseBendingMomentum;
2805 trackParam.SetCovariances(lastParamCov);
52c6d8aa 2806
d24a4636 2807 AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(minCh),1.);
2808 /// AliMUONTrackExtrap::ExtrapToZCov(trackParam, AliMUONConstants::DefaultChamberZ(minCh),kTRUE);
2809 ///AliMUONTrackExtrap::ExtrapToZCov(&trackParam, clus1.fZ,kTRUE);
52c6d8aa 2810
d24a4636 2811 LinearExtrapToZ(&trackParam, clus1.fZ);
52c6d8aa 2812
d24a4636 2813 trackParam.SetTrackChi2(0.);
52c6d8aa 2814
d24a4636 2815 Double_t chi2 = 0.0;
52c6d8aa 2816
d24a4636 2817 chi2 = KalmanFilter(trackParam,&clus1);
a038ee7d 2818
2819 if( chi2 > 1.0e9 /* is order to check TMath::AreEqualAbs(chi2,,1.0e-5)*/ ) {
2820 HLTWarning("Kalman Chi2 calculation cannot be completed...skipping slat track %d (c.p = 1)",ibacktrackseg);
2821 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
2822 continue;
2823 }
52c6d8aa 2824
2825#ifdef PRINT_KALMAN
ed5809f7 2826 HLTImportant("\t\tFor minCh : %d, Chi2 = %lf, GetBeMom : %lf\n",minCh,chi2,trackParam.GetInverseBendingMomentum());
d24a4636 2827 /// TMatrixD para2(trackParam->GetParameters());
2828 /// para2.Print();
52c6d8aa 2829#endif
2830
d24a4636 2831 AliMUONTrackParam extrapTrackParamAtCluster2,minChi2Param;
2832 Double_t chi2OfCluster;
2833 Bool_t tryonefast;
2834 Double_t minChi2=1000000.0;
2835 Int_t frontsegIndex = -1;
2836 extrapTrackParamAtCluster2 = trackParam ;
2837 minChi2Param = trackParam ;
2838
2839 for( Int_t iconnected=0;iconnected<fNofConnectedfrontTrackSeg[ibacktrackseg];iconnected++){
2840
2841 ifronttrackseg = fBackToFront[ibacktrackseg][iconnected];
2842 AliMUONTrackParam extrapTrackParam;
2843 trackParam = extrapTrackParamAtCluster2;
2844
2845 minIndex = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
2846 minCh = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
2847 clus1.fX = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fX ;
2848 clus1.fY = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fY ;
2849 clus1.fZ = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fZ ;
2850 clus1.fErrX2 = 0.020736;
2851 clus1.fErrY2 = 0.000100;
fd5b812e 2852
2853
d24a4636 2854 AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(minCh),1.);
2855 trackParam.ResetPropagator();
2856 AliMUONTrackExtrap::ExtrapToZCov(&trackParam,clus1.fZ, kTRUE);
52c6d8aa 2857
d24a4636 2858 tryonefast = TryOneClusterFast(trackParam, &clus1);
2859 ///if (!tryonefast) continue;
52c6d8aa 2860
d24a4636 2861 /// try to add the current cluster accuratly
2862 chi2OfCluster = TryOneCluster(trackParam, &clus1 , extrapTrackParam,kTRUE);
fd5b812e 2863
d24a4636 2864 extrapTrackParam.SetExtrapParameters(extrapTrackParam.GetParameters());
2865 extrapTrackParam.SetExtrapCovariances(extrapTrackParam.GetCovariances());
fd5b812e 2866
d24a4636 2867 chi2 = KalmanFilter(extrapTrackParam,&clus1);
a038ee7d 2868 if( chi2 > 1.0e9 /* is order to check TMath::AreEqualAbs(chi2,,1.0e-5)*/ ) {
2869 HLTWarning("Kalman Chi2 calculation cannot be completed...skipping slat track %d (c.p = 2)",ibacktrackseg);
2870 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
2871 continue;
2872 }
d24a4636 2873 if(chi2<minChi2){
2874 minChi2 = chi2;
2875 minChi2Param = extrapTrackParam;
2876 frontsegIndex = ifronttrackseg;
2877 }
2878 }///for loop on all connected front track segs
52c6d8aa 2879
d24a4636 2880 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
2881 ///have to check later
2882 if(frontsegIndex==-1) continue;
52c6d8aa 2883
d24a4636 2884 fBackToFront[ibacktrackseg][fNofConnectedfrontTrackSeg[ibacktrackseg]++] = frontsegIndex;
2885 ifronttrackseg = frontsegIndex;
52c6d8aa 2886
d24a4636 2887 trackParam = minChi2Param ;
52c6d8aa 2888
2889#ifdef PRINT_KALMAN
ed5809f7 2890 HLTImportant("\t\t\tCh%d : (%f,%f,%f)\n",minCh,clus1.fX,clus1.fY,clus1.fZ);
2891 HLTImportant("\t\t\tFor minCh : %d, Chi2 = %lf, GetBeMom : %lf\t",minCh,chi2,trackParam.GetInverseBendingMomentum());
2892 HLTImportant(Form("Pt : %lf\n",TMath::Sqrt(trackParam.Px()*trackParam.Px() +
d24a4636 2893 trackParam.Py()*trackParam.Py())));
52c6d8aa 2894#endif
2895
2896
d24a4636 2897 minIndex = (fFrontTrackSeg[ifronttrackseg].fIndex[0]!=-1)?0:1;
2898 minCh = (fFrontTrackSeg[ifronttrackseg].fIndex[0]!=-1)?0:1;
2899 clus1.fX = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fX ;
2900 clus1.fY = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fY ;
2901 clus1.fZ = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fZ ;
2902 clus1.fErrX2 = 0.020736;
2903 clus1.fErrY2 = 0.000100;
2904
2905 AliMUONTrackExtrap::AddMCSEffect(&trackParam,AliMUONConstants::ChamberThicknessInX0(minCh),1.);
2906 trackParam.ResetPropagator();
2907 ///AliMUONTrackExtrap::ExtrapToZCov(&trackParam,clus1.fZ, kTRUE);
2908 LinearExtrapToZ(&trackParam, clus1.fZ);
2909
2910 tryonefast = TryOneClusterFast(trackParam, &clus1);
2911 ///if (!tryonefast) continue;
52c6d8aa 2912
d24a4636 2913 chi2OfCluster = TryOneCluster(trackParam, &clus1 , extrapTrackParamAtCluster2,kTRUE);
52c6d8aa 2914
d24a4636 2915 extrapTrackParamAtCluster2.SetExtrapParameters(extrapTrackParamAtCluster2.GetParameters());
2916 extrapTrackParamAtCluster2.SetExtrapCovariances(extrapTrackParamAtCluster2.GetCovariances());
52c6d8aa 2917
d24a4636 2918 chi2 = KalmanFilter(extrapTrackParamAtCluster2,&clus1);
a038ee7d 2919 if(chi2 > 1.0e9 /* is order to check TMath::AreEqualAbs(chi2,,1.0e-5)*/ ) {
2920 HLTWarning("Kalman Chi2 calculation cannot be completed...skipping slat track %d (c.p = 3)",ibacktrackseg);
2921 fNofConnectedfrontTrackSeg[ibacktrackseg] = 0;
2922 continue;
2923 }
2924
d24a4636 2925 trackParam = extrapTrackParamAtCluster2;
52c6d8aa 2926
2927#ifdef PRINT_KALMAN
ed5809f7 2928 HLTImportant("\t\tCh%d : (%f,%f,%f)\n",minCh,clus1.fX,clus1.fY,clus1.fZ);
2929 HLTImportant("\t\tFor minCh : %d, Chi2 = %lf, GetBeMom : %lf\t",minCh,chi2,trackParam.GetInverseBendingMomentum());
2930 HLTImportant(Form("Pt : %lf\n",TMath::Sqrt(extrapTrackParamAtCluster2.Px()*extrapTrackParamAtCluster2.Px() +
d24a4636 2931 extrapTrackParamAtCluster2.Py()*extrapTrackParamAtCluster2.Py())));
2932 trackParam.Print();
ed5809f7 2933 HLTImportant("AliHLTMUONFullTracker::KalmanChi2Test()--End\n\n");
52c6d8aa 2934#endif
d24a4636 2935 ///AliMUONTrackExtrap::ExtrapToVertex(&trackParam, 0., 0., 0., 0., 0.);
2936 //trackParam.SetInverseBendingMomentum(trackParam.GetCharge());
2937 fTrackParam[ibacktrackseg] = trackParam;
52c6d8aa 2938
2939
2940
2941 }///trig loop
2942
2943
2944
d24a4636 2945 return true;
52c6d8aa 2946}
2947
2948///__________________________________________________________________________
2949
2950
2951Double_t AliHLTMUONFullTracker::BetheBloch(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ)
2952{
2953 //// Returns the mean total momentum energy loss of muon with total momentum='pTotal'
2954 //// in the absorber layer of lenght='pathLength', density='rho', A='atomicA' and Z='atomicZ'
2955 Double_t muMass = 0.105658369; /// GeV
d24a4636 2956 Double_t eMass = 0.510998918e-3; /// GeV
2957 Double_t k = 0.307075e-3; /// GeV.g^-1.cm^2
2958 Double_t i = 9.5e-9; /// mean exitation energy per atomic Z (GeV)
2959 Double_t p2=pTotal*pTotal;
2960 Double_t beta2=p2/(p2 + muMass*muMass);
52c6d8aa 2961
d24a4636 2962 Double_t w = k * rho * pathLength * atomicZ / atomicA / beta2;
52c6d8aa 2963
d24a4636 2964 if (beta2/(1-beta2)>3.5*3.5)
2965 return w * (log(2.*eMass*3.5/(i*atomicZ)) + 0.5*log(beta2/(1-beta2)) - beta2);
52c6d8aa 2966
d24a4636 2967 return w * (log(2.*eMass*beta2/(1-beta2)/(i*atomicZ)) - beta2);
52c6d8aa 2968}
2969
2970
2971///__________________________________________________________________________
2972
2973Double_t AliHLTMUONFullTracker::EnergyLossFluctuation2(Double_t pTotal, Double_t pathLength, Double_t rho, Double_t atomicA, Double_t atomicZ)
2974{
2975 //// Returns the total momentum energy loss fluctuation of muon with total momentum='pTotal'
2976 //// in the absorber layer of lenght='pathLength', density='rho', A='atomicA' and Z='atomicZ'
2977 Double_t muMass = 0.105658369; /// GeV
d24a4636 2978 ///Double_t eMass = 0.510998918e-3; /// GeV
2979 Double_t k = 0.307075e-3; /// GeV.g^-1.cm^2
2980 Double_t p2=pTotal*pTotal;
2981 Double_t beta2=p2/(p2 + muMass*muMass);
52c6d8aa 2982
d24a4636 2983 Double_t fwhm = 2. * k * rho * pathLength * atomicZ / atomicA / beta2; /// FWHM of the energy loss Landau distribution
2984 Double_t sigma2 = fwhm * fwhm / (8.*log(2.)); /// gaussian: fwmh = 2 * srqt(2*ln(2)) * sigma (i.e. fwmh = 2.35 * sigma)
52c6d8aa 2985
d24a4636 2986 ///sigma2 = k * rho * pathLength * atomicZ / atomicA * eMass; /// sigma2 of the energy loss gaussian distribution
52c6d8aa 2987
d24a4636 2988 return sigma2;
52c6d8aa 2989}
2990
2991
2992///__________________________________________________________________________
2993
fd5b812e 2994void AliHLTMUONFullTracker::LinearExtrapToZ(AliMUONTrackParam* trackParam, Double_t zEnd)
52c6d8aa 2995{
2996 //// Track parameters (and their covariances if any) linearly extrapolated to the plane at "zEnd".
2997 //// On return, results from the extrapolation are updated in trackParam.
2998
fd5b812e 2999 if ( TMath::AreEqualAbs(trackParam->GetZ(),zEnd,1.0e-5) ) return; /// nothing to be done if same z
52c6d8aa 3000
d24a4636 3001 Double_t dZ = zEnd - trackParam->GetZ();
3002 trackParam->SetNonBendingCoor(trackParam->GetNonBendingCoor() + trackParam->GetNonBendingSlope() * dZ);
3003 trackParam->SetBendingCoor(trackParam->GetBendingCoor() + trackParam->GetBendingSlope() * dZ);
3004 trackParam->SetZ(zEnd);
52c6d8aa 3005}
3006
3007void AliHLTMUONFullTracker::CorrectELossEffectInAbsorber(AliMUONTrackParam* param, Double_t eLoss)
3008{
fd5b812e 3009 /// Energy loss correction
52c6d8aa 3010 Double_t nonBendingSlope = param->GetNonBendingSlope();
3011 Double_t bendingSlope = param->GetBendingSlope();
3012 param->SetInverseBendingMomentum(param->GetCharge() / (param->P() + eLoss) *
3013 TMath::Sqrt(1.0 + nonBendingSlope*nonBendingSlope + bendingSlope*bendingSlope) /
3014 TMath::Sqrt(1.0 + bendingSlope*bendingSlope));
3015}
3016
3017///__________________________________________________________________________
3018
3019void AliHLTMUONFullTracker::Cov2CovP(const TMatrixD &param, TMatrixD &cov)
3020{
3021 //// change coordinate system: (X, SlopeX, Y, SlopeY, q/Pyz) -> (X, SlopeX, Y, SlopeY, q*PTot)
3022 //// parameters (param) are given in the (X, SlopeX, Y, SlopeY, q/Pyz) coordinate system
3023
3024 /// charge * total momentum
3025 Double_t qPTot = TMath::Sqrt(1. + param(1,0)*param(1,0) + param(3,0)*param(3,0)) /
3026 TMath::Sqrt(1. + param(3,0)*param(3,0)) / param(4,0);
3027
3028 /// Jacobian of the opposite transformation
d24a4636 3029 TMatrixD jacob(5,5);
3030 jacob.UnitMatrix();
3031 jacob(4,1) = qPTot * param(1,0) / (1. + param(1,0)*param(1,0) + param(3,0)*param(3,0));
3032 jacob(4,3) = - qPTot * param(1,0) * param(1,0) * param(3,0) /
3033 (1. + param(3,0)*param(3,0)) / (1. + param(1,0)*param(1,0) + param(3,0)*param(3,0));
3034 jacob(4,4) = - qPTot / param(4,0);
3035
3036 /// compute covariances in new coordinate system
3037 TMatrixD tmp(5,5);
3038 tmp.MultT(cov,jacob);
3039 cov.Mult(jacob,tmp);
52c6d8aa 3040
3041}
3042
3043///__________________________________________________________________________
3044
3045void AliHLTMUONFullTracker::CovP2Cov(const TMatrixD &param, TMatrixD &covP)
3046{
3047 //// change coordinate system: (X, SlopeX, Y, SlopeY, q*PTot) -> (X, SlopeX, Y, SlopeY, q/Pyz)
3048 //// parameters (param) are given in the (X, SlopeX, Y, SlopeY, q/Pyz) coordinate system
3049
3050 /// charge * total momentum
3051 Double_t qPTot = TMath::Sqrt(1. + param(1,0)*param(1,0) + param(3,0)*param(3,0)) /
3052 TMath::Sqrt(1. + param(3,0)*param(3,0)) / param(4,0);
3053
3054 /// Jacobian of the transformation
d24a4636 3055 TMatrixD jacob(5,5);
3056 jacob.UnitMatrix();
3057 jacob(4,1) = param(4,0) * param(1,0) / (1. + param(1,0)*param(1,0) + param(3,0)*param(3,0));
3058 jacob(4,3) = - param(4,0) * param(1,0) * param(1,0) * param(3,0) /
3059 (1. + param(3,0)*param(3,0)) / (1. + param(1,0)*param(1,0) + param(3,0)*param(3,0));
3060 jacob(4,4) = - param(4,0) / qPTot;
3061
3062 /// compute covariances in new coordinate system
3063 TMatrixD tmp(5,5);
3064 tmp.MultT(covP,jacob);
3065 covP.Mult(jacob,tmp);
52c6d8aa 3066
3067}
3068
3069///__________________________________________________________________________
3070
3071
3072void AliHLTMUONFullTracker::CorrectMCSEffectInAbsorber(AliMUONTrackParam* param,
fd5b812e 3073 Double_t xVtx, Double_t yVtx, Double_t zVtx,
3074 Double_t absZBeg, Double_t f1, Double_t f2)
52c6d8aa 3075{
3076 //// Correct parameters and corresponding covariances using Branson correction
3077 //// - input param are parameters and covariances at the end of absorber
3078 //// - output param are parameters and covariances at vertex
3079 //// Absorber correction parameters are supposed to be calculated at the current track z-position
3080
3081 /// Position of the Branson plane (spectro. (z<0))
3082 Double_t zB = (f1>0.) ? absZBeg - f2/f1 : 0.;
3083
3084 LinearExtrapToZ(param,zB);
3085
3086 /// compute track parameters at vertex
d24a4636 3087 TMatrixD newParam(5,1);
3088 newParam.Zero();
3089 newParam(0,0) = xVtx;
3090 newParam(1,0) = (param->GetNonBendingCoor() - xVtx) / (zB - zVtx);
3091 newParam(2,0) = yVtx;
3092 newParam(3,0) = (param->GetBendingCoor() - yVtx) / (zB - zVtx);
3093 newParam(4,0) = param->GetCharge() / param->P() *
3094 TMath::Sqrt(1.0 + newParam(1,0)*newParam(1,0) + newParam(3,0)*newParam(3,0)) /
3095 TMath::Sqrt(1.0 + newParam(3,0)*newParam(3,0));
3096
3097 TMatrixD paramCovP(5,5);
3098
3099 /// Get covariances in (X, SlopeX, Y, SlopeY, q*PTot) coordinate system
3100 paramCovP.Use(param->GetCovariances());
3101
3102 Cov2CovP(param->GetParameters(),paramCovP);
3103
3104 /// Get the covariance matrix in the (XVtx, X, YVtx, Y, q*PTot) coordinate system
3105 /// TMatrixD paramCovVtx(5,5);
3106 Double_t element44 = paramCovP(4,4);
3107 paramCovP.Zero();
3108 paramCovP(4,4) = element44;
3109
3110 CovP2Cov(newParam,paramCovP);
3111
3112 /// Set parameters and covariances at vertex
3113 param->SetParameters(newParam);
3114 param->SetZ(zVtx);
3115 param->SetCovariances(paramCovP);
52c6d8aa 3116}
3117
3118///__________________________________________________________________________
3119
aa205fc6 3120Bool_t AliHLTMUONFullTracker::ExtrapolateToOrigin()
52c6d8aa 3121{
3122 /// Extrapolation to origin through absorber
3123
3124 Int_t minIndex=0,maxIndex=0;
3125 Int_t minCh=0,maxCh=0;
3126 Int_t ifronttrackseg = -1;
3127 AliMUONTrackParam trackP;
fd5b812e 3128 Double_t bSlope, nbSlope;
52c6d8aa 3129 AliHLTMUONRecHitStruct p1,p2,pSeg1,pSeg2;
eb672fd1 3130 Double_t pyz = -1.0;
52c6d8aa 3131 TVector3 v1,v2,v3,v4;
3132 Double_t eLoss1,eLoss2,eLoss3;
3133 Double_t b;
3134 Double_t zE,zB,dzE,dzB;
eb672fd1 3135 Double_t f0,f1,f2;
3136 Double_t f0Sum,f1Sum,f2Sum;
52c6d8aa 3137 Double_t fXVertex=0.0,fYVertex=0.0,fZVertex=0.0;
fd5b812e 3138 Double_t thetaDev;
52c6d8aa 3139
fd5b812e 3140 for( Int_t ibacktrackseg=0;ibacktrackseg<fNofbackTrackSeg;ibacktrackseg++){
3141
52c6d8aa 3142 if(fNofConnectedfrontTrackSeg[ibacktrackseg]<=0) continue;
3143
3144 maxIndex = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?3:2;
3145 maxCh = (fBackTrackSeg[ibacktrackseg].fIndex[3]!=-1)?9:8;
3146
3147 minIndex = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?0:1;
3148 minCh = (fBackTrackSeg[ibacktrackseg].fIndex[0]!=-1)?6:7;
3149
3150 p2.fX = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fX ;
3151 p2.fY = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fY ;
3152 p2.fZ = fChPoint[maxCh][fBackTrackSeg[ibacktrackseg].fIndex[maxIndex]]->fZ ;
3153
3154 p1.fX = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fX ;
3155 p1.fY = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fY ;
3156 p1.fZ = fChPoint[minCh][fBackTrackSeg[ibacktrackseg].fIndex[minIndex]]->fZ ;
3157
fd5b812e 3158 thetaDev= (p1.fY*p2.fZ - p2.fY*p1.fZ)/(p2.fZ - p1.fZ);
3159
52c6d8aa 3160 Sub(&p2,&p1,&pSeg1);
3161
3162 ifronttrackseg = fBackToFront[ibacktrackseg][0];
3163
3164 maxIndex = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
3165 maxCh = (fFrontTrackSeg[ifronttrackseg].fIndex[3]!=-1)?3:2;
3166
3167 minIndex = (fFrontTrackSeg[ifronttrackseg].fIndex[0]!=-1)?0:1;
3168 minCh = (fFrontTrackSeg[ifronttrackseg].fIndex[0]!=-1)?0:1;
3169
3170 p2.fX = fChPoint[maxCh][fFrontTrackSeg[ifronttrackseg].fIndex[maxIndex]]->fX ;
3171 p2.fY = fChPoint[maxCh][fFrontTrackSeg[ifronttrackseg].fIndex[maxIndex]]->fY ;
3172 p2.fZ = fChPoint[maxCh][fFrontTrackSeg[ifronttrackseg].fIndex[maxIndex]]->fZ ;
3173
3174 p1.fX = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fX ;
3175 p1.fY = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fY ;
3176 p1.fZ = fChPoint[minCh][fFrontTrackSeg[ifronttrackseg].fIndex[minIndex]]->fZ ;
3177
3178 Sub(&p2,&p1,&pSeg2);
fd5b812e 3179
3180 if(thetaDev > 0)
3181 pyz = (3.0*0.3/sin(Angle(&pSeg1,&pSeg2)));/// * sqrt(x3*x3 + y3*y3)/z3 ;
3182 else
3183 pyz = -(3.0*0.3/sin(Angle(&pSeg1,&pSeg2)));/// * sqrt(x3*x3 + y3*y3)/z3 ;
52c6d8aa 3184
d24a4636 3185 nbSlope = (p2.fX - p1.fX)/(p2.fZ - p1.fZ);
3186 bSlope = (p2.fY - p1.fY)/(p2.fZ - p1.fZ);
52c6d8aa 3187
d24a4636 3188 trackP.SetZ(p1.fZ);
3189 trackP.SetNonBendingCoor(p1.fX);
3190 trackP.SetNonBendingSlope(nbSlope);
3191 trackP.SetBendingCoor(p1.fY);
3192 trackP.SetBendingSlope(bSlope);
3193 trackP.SetInverseBendingMomentum(1.0/pyz) ;
52c6d8aa 3194
3195
aa205fc6 3196
3197
3198 if(not fFastTracking){
3199 trackP = fTrackParam[ibacktrackseg] ;
3200 }
52c6d8aa 3201
aa205fc6 3202 LinearExtrapToZ(&trackP,fgkAbsoedge[3]);
3203 v4.SetXYZ(trackP.GetNonBendingCoor(),trackP.GetBendingCoor(),trackP.GetZ());
3204 LinearExtrapToZ(&trackP,fgkAbsoedge[2]);
3205 v3.SetXYZ(trackP.GetNonBendingCoor(),trackP.GetBendingCoor(),trackP.GetZ());
3206 LinearExtrapToZ(&trackP,fgkAbsoedge[1]);
3207 v2.SetXYZ(trackP.GetNonBendingCoor(),trackP.GetBendingCoor(),trackP.GetZ());
3208 LinearExtrapToZ(&trackP,fgkAbsoedge[0]);
3209 v1.SetXYZ(trackP.GetNonBendingCoor(),trackP.GetBendingCoor(),trackP.GetZ());
52c6d8aa 3210
aa205fc6 3211 eLoss1 = BetheBloch(trackP.P(), (v4-v3).Mag(), fgkRho[2], fgkAtomicA[2], fgkAtomicZ[2]);
3212 eLoss2 = BetheBloch(trackP.P(), (v3-v2).Mag(), fgkRho[1], fgkAtomicA[1], fgkAtomicZ[1]);
3213 eLoss3 = BetheBloch(trackP.P(), (v2-v1).Mag(), fgkRho[0], fgkAtomicA[0], fgkAtomicZ[0]);
3214
3215 /// sigmaELoss1 = EnergyLossFluctuation2(trackP.P(), (v4-v3).Mag(), rho[2], atomicA[2], atomicZ[2]);
3216 /// sigmaELoss2 = EnergyLossFluctuation2(trackP.P(), (v3-v2).Mag(), rho[1], atomicA[1], atomicZ[1]);
3217 /// sigmaELoss3 = EnergyLossFluctuation2(trackP.P(), (v2-v1).Mag(), rho[0], atomicA[0], atomicZ[0]);
52c6d8aa 3218
aa205fc6 3219 /// eDiff = totELoss-(eLoss1+eLoss2+eLoss3);
3220 /// sigmaELossDiff = sigmaTotELoss ;///- (sigmaELoss1+sigmaELoss2+sigmaELoss3);
52c6d8aa 3221
aa205fc6 3222 /// CorrectELossEffectInAbsorber(&trackP, 0.5*(eLoss1+eLoss2+eLoss3), 0.5*(sigmaELoss1+sigmaELoss2+sigmaELoss3));
52c6d8aa 3223
3224
aa205fc6 3225 ///CorrectELossEffectInAbsorber(&trackP, totELoss,sigmaTotELoss);
52c6d8aa 3226
aa205fc6 3227 CorrectELossEffectInAbsorber(&trackP, 0.7*(eLoss1+eLoss2+eLoss3));
52c6d8aa 3228
aa205fc6 3229 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
3230 f0Sum = 0.0; f1Sum = 0.0; f2Sum = 0.0;
52c6d8aa 3231
aa205fc6 3232 b = (v4.Z()-v1.Z())/((v4-v1).Mag());
52c6d8aa 3233
aa205fc6 3234 zB = v1.Z();
3235 zE = b*((v2-v1).Mag()) + zB;
3236 dzB = zB - v1.Z();
3237 dzE = zE - v1.Z();
52c6d8aa 3238
aa205fc6 3239 f0 = ((v2-v1).Mag())/fgkRadLen[0];
3240 f1 = (dzE*dzE - dzB*dzB) / b / b / fgkRadLen[0] /2.;
3241 f2 = (dzE*dzE*dzE - dzB*dzB*dzB) / b / b / b / fgkRadLen[0] / 3.;
52c6d8aa 3242
aa205fc6 3243 f0Sum += f0;
3244 f1Sum += f1;
3245 f2Sum += f2;
52c6d8aa 3246
aa205fc6 3247 zB = zE;
3248 zE = b*((v3-v2).Mag()) + zB;
3249 dzB = zB - v1.Z();
3250 dzE = zE - v1.Z();
52c6d8aa 3251
aa205fc6 3252 f0 = ((v3-v2).Mag())/fgkRadLen[1];
3253 f1 = (dzE*dzE - dzB*dzB) / b / b / fgkRadLen[1] /2.;
3254 f2 = (dzE*dzE*dzE - dzB*dzB*dzB) / b / b / b / fgkRadLen[1] / 3.;
52c6d8aa 3255
aa205fc6 3256 f0Sum += f0;
3257 f1Sum += f1;
3258 f2Sum += f2;
52c6d8aa 3259
aa205fc6 3260 zB = zE;
3261 zE = b*((v4-v3).Mag()) + zB;
3262 dzB = zB - v1.Z();
3263 dzE = zE - v1.Z();
52c6d8aa 3264
aa205fc6 3265 f0 = ((v4-v3).Mag())/fgkRadLen[2];
3266 f1 = (dzE*dzE - dzB*dzB) / b / b / fgkRadLen[2] /2.;
3267 f2 = (dzE*dzE*dzE - dzB*dzB*dzB) / b / b / b / fgkRadLen[2] / 3.;
52c6d8aa 3268
aa205fc6 3269 f0Sum += f0;
3270 f1Sum += f1;
3271 f2Sum += f2;
52c6d8aa 3272
3273
aa205fc6 3274 ///AddMCSEffectInAbsorber(&trackP,(v4-v1).Mag(),f0Sum,f1Sum,f2Sum);
3275
3276 CorrectMCSEffectInAbsorber(&trackP,fXVertex,fYVertex, fZVertex,AliMUONConstants::AbsZBeg(),f1Sum,f2Sum);
3277 CorrectELossEffectInAbsorber(&trackP, 0.5*(eLoss1+eLoss2+eLoss3));
52c6d8aa 3278
d24a4636 3279
3280 ///AliMUONTrackExtrap::ExtrapToVertex(&trackP, 0., 0., 0., 0., 0.);
3281 trackP.SetZ(p1.fZ);
3282 trackP.SetNonBendingCoor(p1.fX);
3283 trackP.SetBendingCoor(p1.fY);
3284 LinearExtrapToZ(&trackP, 0.0);
52c6d8aa 3285
d24a4636 3286 fTrackParam[ibacktrackseg] = trackP;
52c6d8aa 3287
3288 }///backtrackseg for loop
3289
d24a4636 3290 return true;
52c6d8aa 3291}
3292