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