PHOS/AliPHOSTrackSegmentMakerv1.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15 /* $Id$ */
  16
  17 /* History of cvs commits:
  18  *
  19  * $Log$
  20  * Revision 1.93  2007/10/10 09:05:10  schutz
  21  * Changing name QualAss to QA
  22  *
  23  * Revision 1.92  2007/08/28 12:55:08  policheh
  24  * Loaders removed from the reconstruction code (C.Cheshkov)
  25  *
  26  * Revision 1.91  2007/08/07 14:12:03  kharlov
  27  * Quality assurance added (Yves Schutz)
  28  *
  29  * Revision 1.90  2007/07/11 13:43:30  hristov
  30  * New class AliESDEvent, backward compatibility with the old AliESD (Christian)
  31  *
  32  * Revision 1.89  2007/07/03 08:13:04  kharlov
  33  * Bug fix in CPV local coordinates
  34  *
  35  * Revision 1.88  2007/06/27 09:11:07  kharlov
  36  * Bug fix for CPV-EMC distance
  37  *
  38  * Revision 1.87  2007/05/04 14:49:29  policheh
  39  * AliPHOSRecPoint inheritance from AliCluster
  40  *
  41  * Revision 1.86  2007/04/02 15:00:16  cvetan
  42  * No more calls to gAlice in the reconstruction
  43  *
  44  * Revision 1.85  2007/03/28 19:18:15  kharlov
  45  * RecPoints recalculation in TSM removed
  46  *
  47  * Revision 1.84  2007/03/07 07:01:21  hristov
  48  * Fixing copy/paste erro. Additional protections
  49  *
  50  * Revision 1.83  2007/03/06 21:07:37  kharlov
  51  * DP: xz CPV-EMC distance filled to TS
  52  *
  53  * Revision 1.82  2007/03/06 06:54:48  kharlov
  54  * DP:Calculation of cluster properties dep. on vertex added
  55  *
  56  * Revision 1.81  2007/02/05 10:02:40  kharlov
  57  * Module numbering is corrected
  58  *
  59  * Revision 1.80  2006/08/28 10:01:56  kharlov
  60  * Effective C++ warnings fixed (Timur Pocheptsov)
  61  *
  62  * Revision 1.79  2006/04/25 12:41:15  hristov
  63  * Moving non-persistent data to AliESDfriend (Yu.Belikov)
  64  *
  65  * Revision 1.78  2005/11/18 13:04:51  hristov
  66  * Bug fix
  67  *
  68  * Revision 1.77  2005/11/17 23:34:36  hristov
  69  * Corrected logics
  70  *
  71  * Revision 1.76  2005/11/17 22:29:12  hristov
  72  * Faster version, no attempt to match tracks outside the PHOS acceptance
  73  *
  74  * Revision 1.75  2005/11/17 12:35:27  hristov
  75  * Use references instead of objects. Avoid to create objects when they are not really needed
  76  *
  77  * Revision 1.74  2005/07/08 14:01:36  hristov
  78  * Tracking in non-uniform nmagnetic field (Yu.Belikov)
  79  *
  80  * Revision 1.73  2005/05/28 14:19:05  schutz
  81  * Compilation warnings fixed by T.P.
  82  *
  83  */
  84
  85 //_________________________________________________________________________
  86 // Implementation version 1 of algorithm class to construct PHOS track segments
  87 // Track segment for PHOS is list of
  88 //        EMC RecPoint + (possibly) CPV RecPoint
  89 // To find TrackSegments we do the following:
  90 //  for each EMC RecPoints we look at
  91 //   CPV RecPoints in the radious fRcpv.
  92 //  If there is such a CPV RecPoint,
  93 //   we make "Link" it is just indexes of EMC and CPV RecPoint and distance
  94 //   between them in the PHOS plane.
  95 //  Then we sort "Links" and starting from the
  96 //   least "Link" pointing to the unassined EMC and CPV RecPoints assing them to
  97 //   new TrackSegment.
  98 // If there is no CPV RecPoint we make TrackSegment
  99 // consisting from EMC alone. There is no TrackSegments without EMC RecPoint.
 100 //// In principle this class should be called from AliPHOSReconstructor, but
 101 // one can use it as well in standalone mode.
 102 // Use  case:
 103 //  root [0] AliPHOSTrackSegmentMakerv1 * t = new AliPHOSTrackSegmentMaker("galice.root", "tracksegmentsname", "recpointsname")
 104 //  Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
 105 //               // reads gAlice from header file "galice.root", uses recpoints stored in the branch names "recpointsname" (default = "Default")
 106 //               // and saves recpoints in branch named "tracksegmentsname" (default = "recpointsname")
 107 //  root [1] t->ExecuteTask()
 108 //  root [3] t->SetTrackSegmentsBranch("max distance 5 cm")
 109 //  root [4] t->ExecuteTask("deb all time")
 110 //
 111 //*-- Author: Dmitri Peressounko (RRC Ki & SUBATECH) & Yves Schutz (SUBATECH)
 112 //
 113
 114 // --- ROOT system ---
 115 #include "TVector3.h"
 116 #include "TTree.h"
 117 #include "TBenchmark.h"
 118
 119 // --- Standard library ---
 120 #include "Riostream.h"
 121 // --- AliRoot header files ---
 122 #include "AliPHOSGeometry.h"
 123 #include "AliPHOSTrackSegmentMakerv1.h"
 124 #include "AliPHOSTrackSegment.h"
 125 #include "AliPHOSLink.h"
 126 #include "AliESDEvent.h"
 127 #include "AliESDtrack.h"
 128 #include "AliPHOSEmcRecPoint.h"
 129 #include "AliPHOSCpvRecPoint.h"
 130 #include "AliLog.h"
 131 #include "AliMagF.h"
 132 #include "AliMagF.h"
 133 #include "AliTracker.h"
 134 #include "AliGeomManager.h"
 135 #include "AliCluster.h"
 136 #include "AliKalmanTrack.h"
 137 #include "AliGlobalQADataMaker.h"
 138 #include "AliVParticle.h"
 139
 140
 141 ClassImp( AliPHOSTrackSegmentMakerv1)
 142
 143
 144 //____________________________________________________________________________
 145 AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1() :
 146   AliPHOSTrackSegmentMaker(),
 147   fDefaultInit(kTRUE),
 148   fWrite(kFALSE),
 149   fNTrackSegments(0),
 150   fRcpv(0.f),
 151   fRtpc(0.f),
 152   fVtx(0.f),
 153   fLinkUpArray(0),
 154   fEmcFirst(0),
 155   fEmcLast(0),
 156   fCpvFirst(0),
 157   fCpvLast(0),
 158   fModule(0),
 159   fTrackSegments(NULL)
 160 {
 161   // default ctor (to be used mainly by Streamer)
 162   InitParameters() ;
 163 }
 164
 165 //____________________________________________________________________________
 166 AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(AliPHOSGeometry *geom) :
 167   AliPHOSTrackSegmentMaker(geom),
 168   fDefaultInit(kFALSE),
 169   fWrite(kFALSE),
 170   fNTrackSegments(0),
 171   fRcpv(0.f),
 172   fRtpc(0.f),
 173   fVtx(0.f),
 174   fLinkUpArray(0),
 175   fEmcFirst(0),
 176   fEmcLast(0),
 177   fCpvFirst(0),
 178   fCpvLast(0),
 179   fModule(0),
 180   fTrackSegments(NULL)
 181 {
 182   // ctor
 183   InitParameters() ;
 184   Init() ;
 185   fESD = 0;
 186 }
 187
 188
 189 AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(const AliPHOSTrackSegmentMakerv1 & tsm) :
 190   AliPHOSTrackSegmentMaker(tsm),
 191   fDefaultInit(kFALSE),
 192   fWrite(kFALSE),
 193   fNTrackSegments(0),
 194   fRcpv(0.f),
 195   fRtpc(0.f),
 196   fVtx(0.f),
 197   fLinkUpArray(0),
 198   fEmcFirst(0),
 199   fEmcLast(0),
 200   fCpvFirst(0),
 201   fCpvLast(0),
 202   fModule(0),
 203   fTrackSegments(NULL)
 204 {
 205   // cpy ctor: no implementation yet
 206   // requested by the Coding Convention
 207   Fatal("cpy ctor", "not implemented") ;
 208 }
 209
 210
 211 //____________________________________________________________________________
 212  AliPHOSTrackSegmentMakerv1::~AliPHOSTrackSegmentMakerv1()
 213 {
 214   // dtor
 215   // fDefaultInit = kTRUE if TrackSegmentMaker created by default ctor (to get just the parameters)
 216   if (!fDefaultInit)
 217     delete fLinkUpArray ;
 218   if (fTrackSegments) {
 219     fTrackSegments->Delete();
 220     delete fTrackSegments;
 221   }
 222 }
 223
 224 //____________________________________________________________________________
 225 void  AliPHOSTrackSegmentMakerv1::FillOneModule()
 226 {
 227   // Finds first and last indexes between which
 228   // clusters from one PHOS module are
 229
 230   //First EMC clusters
 231   Int_t totalEmc = fEMCRecPoints->GetEntriesFast() ;
 232   for(fEmcFirst = fEmcLast; (fEmcLast < totalEmc) &&
 233         ((dynamic_cast<AliPHOSRecPoint *>(fEMCRecPoints->At(fEmcLast)))->GetPHOSMod() == fModule );
 234       fEmcLast ++)  ;
 235
 236   //Now CPV clusters
 237   Int_t totalCpv = fCPVRecPoints->GetEntriesFast() ;
 238
 239     for(fCpvFirst = fCpvLast; (fCpvLast < totalCpv) &&
 240          ((dynamic_cast<AliPHOSRecPoint *>(fCPVRecPoints->At(fCpvLast)))->GetPHOSMod() == fModule );
 241        fCpvLast ++) ;
 242
 243 }
 244
 245 //____________________________________________________________________________
 246 void  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint * emcClu,
 247                                                          AliPHOSCpvRecPoint * cpvClu,
 248                                                          Int_t &trackindex,
 249                                                          Float_t &dx, Float_t &dz) const
 250 {
 251   // Calculates the distance between the EMC RecPoint and the CPV RecPoint
 252   // If no CPV, calculates the distance between the EMC RecPoint and the track
 253   // prolongation to the PHOS module plane.
 254   // Clusters are sorted in "rows" and "columns" of width 1 cm
 255
 256 //  Float_t delta = 1 ;  // Width of the rows in sorting of RecPoints (in cm)
 257 //                       // if you change this value, change it as well in xxxRecPoint::Compare()
 258
 259   if(!cpvClu) {
 260
 261     trackindex = -1;
 262     dx=999.;
 263     dz=999.;
 264
 265     if(!emcClu) {
 266       return;
 267     }
 268
 269     // *** Start the matching
 270     Int_t nt=fESD->GetNumberOfTracks();
 271     Int_t iPHOSMod = emcClu->GetPHOSMod()  ;
 272     //Calculate actual distance to PHOS module
 273     TVector3 globaPos ;
 274     fGeom->Local2Global(iPHOSMod, 0.,0., globaPos) ;
 275     const Double_t rPHOS = globaPos.Pt() ; //Distance to PHOS module
 276     const Double_t kYmax = 72.+10. ; //Size of the module (with some reserve) in phi direction
 277     const Double_t kZmax = 64.+10. ; //Size of the module (with some reserve) in z direction
 278     const Double_t kAlpha= 20./180.*TMath::Pi() ; //TPC sector angular size
 279     Double_t minDistance = 1.e6;
 280
 281     TMatrixF gmat;
 282     TVector3 gposRecPoint; // global (in ALICE frame) position of rec. point
 283     emcClu->GetGlobalPosition(gposRecPoint,gmat);
 284     Double_t gposTrack[3] ;
 285
 286     Double_t bz = AliTracker::GetBz() ; //B-Field for approximate matching
 287     Double_t b[3];
 288     for (Int_t i=0; i<nt; i++) {
 289       AliESDtrack *esdTrack=fESD->GetTrack(i);
 290
 291 //     // Skip the tracks having "wrong" status (has to be checked/tuned)
 292 //     ULong_t status = esdTrack->GetStatus();
 293 //     if ((status & AliESDtrack::kTRDout)   == 0) continue;
 294 //     if ((status & AliESDtrack::kTRDrefit) == 1) continue;
 295
 296       AliExternalTrackParam t(*esdTrack);
 297       Int_t isec=Int_t(t.GetAlpha()/kAlpha);
 298       Int_t imod=-isec-2; // PHOS module
 299
 300       Double_t y;                       // Some tracks do not reach the PHOS
 301       if (!t.GetYAt(rPHOS,bz,y)) continue; //    because of the bending
 302
 303       Double_t z; t.GetZAt(rPHOS,bz,z);
 304       if (TMath::Abs(z) > kZmax) continue; // Some tracks miss the PHOS in Z
 305
 306       Bool_t ok=kTRUE;
 307       while (TMath::Abs(y) > kYmax) {   // Find the matching module
 308         Double_t alp=t.GetAlpha();
 309         if (y > kYmax) {
 310           if (!t.Rotate(alp+kAlpha)) {ok=kFALSE; break;}
 311           imod--;
 312         } else if (y < -kYmax) {
 313           if (!t.Rotate(alp-kAlpha)) {ok=kFALSE; break;}
 314           imod++;
 315         }
 316         if (!t.GetYAt(rPHOS,bz,y)) {ok=kFALSE; break;}
 317       }
 318       if (!ok) continue; // Track rotation failed
 319
 320
 321       if(imod!= iPHOSMod-1)
 322         continue; //not even approximate coincidence
 323
 324       //t.CorrectForMaterial(...); // Correct for the TOF material, if needed
 325       t.GetBxByBz(b) ;
 326       t.PropagateToBxByBz(rPHOS,b);        // Propagate to the matching module
 327       t.GetXYZ(gposTrack) ;
 328
 329       Double_t ddx = gposTrack[0] - gposRecPoint.X(), ddy = gposTrack[1] - gposRecPoint.Y(), ddz = gposTrack[2] - gposRecPoint.Z();
 330       Double_t d2 = ddx*ddx + ddy*ddy + ddz*ddz;
 331       if(d2 < minDistance) {
 332         dx = TMath::Sign(TMath::Sqrt(ddx*ddx + ddy*ddy),ddx) ;
 333         dz = ddz ;
 334         trackindex=i;
 335         minDistance=d2 ;
 336       }
 337     }
 338     return ;
 339   }
 340
 341
 342   TVector3 emcGlobal;
 343   fGeom->GetGlobalPHOS((AliPHOSRecPoint*)emcClu,emcGlobal);
 344
 345   // Radius from IP to current point
 346   Double_t rEMC = TMath::Abs(emcGlobal.Pt());
 347
 348   // Extrapolate the global track direction to EMC
 349   // and find the closest track
 350
 351   Int_t nTracks = fESD->GetNumberOfTracks();
 352
 353   AliESDtrack *track;
 354   Double_t xyz[] = {-1,-1,-1};
 355   Double_t pxyz[3];
 356   Double_t zEMC,xEMC;
 357   Int_t module;
 358   TVector3 vecP;
 359   TVector3 locClu;
 360
 361   Float_t minDistance = 1.e6;
 362   Float_t dr;
 363
 364   for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
 365     track = fESD->GetTrack(iTrack);
 366     if (!track->GetXYZAt(rEMC, fESD->GetMagneticField(), xyz)) continue;
 367
 368     AliDebug(1,Form("Event %d, iTrack: %d, (%.3f,%.3f,%.3f)",
 369                     fESD->GetEventNumberInFile(),iTrack,xyz[0],xyz[1],xyz[2]));
 370
 371     if (track->GetPxPyPzAt(rEMC,fESD->GetMagneticField(),pxyz)) {
 372
 373       vecP.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
 374       fGeom->ImpactOnEmc(xyz,vecP.Theta(),vecP.Phi(),module,zEMC,xEMC) ;
 375
 376       if(!module) continue;
 377       AliDebug(1,Form("\t\tTrack hit PHOS! Module: %d, (x,z)=(%.3f,%.3f)",module,xEMC,zEMC));
 378
 379       if(emcClu->GetPHOSMod() != module) continue;
 380
 381       // match track to EMC cluster
 382       emcClu->GetLocalPosition(locClu);
 383
 384       Float_t delta_x = xEMC - locClu.X();
 385       Float_t delta_z = zEMC - locClu.Z();
 386       dr = TMath::Sqrt(delta_x*delta_x + delta_z*delta_z);
 387       AliDebug(1,Form("\tMatch iTrack=%d: (dx,dz)=(%.3f,%.3f)",iTrack,delta_x,delta_z));
 388
 389       if(dr<minDistance) {
 390         trackindex = iTrack;
 391         minDistance = dr;
 392         dx = delta_x;
 393         dz = delta_z;
 394       }
 395     }
 396
 397   }
 398
 399   if(trackindex>=0)
 400     AliDebug(1,Form("\t\tBest match for (xClu,zClu,eClu)=(%.3f,%.3f,%.3f): iTrack=%d, dR=%.3f",
 401                     locClu.X(),locClu.Z(),emcClu->GetEnergy(),
 402                     trackindex,TMath::Sqrt(dx*dx+dz*dz)));
 403   return;
 404
 405   Float_t distance2Track = fRtpc ;
 406
 407   trackindex = -1 ; // closest track within fRCpv
 408
 409   TVector3 vecEmc ;   // Local position of EMC recpoint
 410   TVector3 vecPloc ;     // Momentum direction at CPV plain
 411
 412   //toofar = kTRUE ;
 413   if(emcClu->GetPHOSMod() != cpvClu->GetPHOSMod()){
 414     dx=999. ;
 415     dz=999. ;
 416     return ;
 417   }
 418
 419   emcClu->GetLocalPosition(vecEmc) ;
 420
 421   Double_t xCPV,zCPV ; //EMC-projected coordinates of CPV cluster
 422   TVector3 cpvGlobal; // Global position of the CPV recpoint
 423   fGeom->GetGlobalPHOS((AliPHOSRecPoint*)cpvClu,cpvGlobal);
 424   Double_t vtxCPV[3]={cpvGlobal.X(),cpvGlobal.Y(),cpvGlobal.Z()} ;
 425   Int_t dummyMod ;
 426
 427   if (fESD == 0x0) {
 428     //if no track information available, assume straight line from IP to emcal
 429     fGeom->ImpactOnEmc(vtxCPV,cpvGlobal.Theta(),cpvGlobal.Phi(),dummyMod,zCPV,xCPV) ;
 430     dx=xCPV - vecEmc.X() ;
 431     dz=zCPV - vecEmc.Z() ;
 432     return ;
 433   }
 434
 435   //if there is ESD try to correct distance using TPC information on particle direct in CPV
 436   if (fESD != 0x0) {
 437
 438     Double_t rCPV = cpvGlobal.Pt() ;// Radius from IP to current point
 439
 440     // Extrapolate the global track direction if any to CPV and find the closest track
 441     Int_t iClosestTrack = -1;
 442     TVector3 inPHOS ;
 443
 444     for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
 445       track = fESD->GetTrack(iTrack);
 446       if (!track->GetXYZAt(rCPV, fESD->GetMagneticField(), xyz))
 447         continue; //track coord on the cylinder of PHOS radius
 448       if ((TMath::Abs(xyz[0])+TMath::Abs(xyz[1])+TMath::Abs(xyz[2]))<=0)
 449         continue;
 450       //Check if this track hits PHOS
 451       inPHOS.SetXYZ(xyz[0],xyz[1],xyz[2]);
 452       distance2Track = inPHOS.Angle(cpvGlobal) ;
 453       // Find the closest track to the CPV recpoint
 454       if (distance2Track < minDistance) {
 455         minDistance = distance2Track;
 456         iClosestTrack = iTrack;
 457       }
 458     }
 459
 460     if (iClosestTrack != -1) {
 461       track = fESD->GetTrack(iClosestTrack);
 462       if (track->GetPxPyPzAt(rCPV, fESD->GetMagneticField(), pxyz)) { // track momentum ibid.
 463         vecP.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
 464         fGeom->ImpactOnEmc(vtxCPV,vecP.Theta(),vecP.Phi(),dummyMod,zCPV,xCPV) ;
 465       }
 466     }
 467
 468     if(minDistance < fRtpc ){
 469       trackindex = iClosestTrack ;
 470     }
 471   }
 472   if(trackindex!=-1){
 473     // If the closest global track is found, calculate EMC-CPV distance from it
 474     dx=xCPV - vecEmc.X() ;
 475     dz=zCPV - vecEmc.Z() ;
 476   }
 477   else{
 478     // If no global track was found, just take the nearest CPV point
 479     fGeom->ImpactOnEmc(vtxCPV,cpvGlobal.Theta(),cpvGlobal.Phi(),dummyMod,zCPV,xCPV) ;
 480     dx=xCPV - vecEmc.X() ;
 481     dz=zCPV - vecEmc.Z() ;
 482   }
 483   return ;
 484 }
 485 //____________________________________________________________________________
 486 void  AliPHOSTrackSegmentMakerv1::Init()
 487 {
 488   // Make all memory allocations that are not possible in default constructor
 489
 490   fLinkUpArray  = new TClonesArray("AliPHOSLink", 1000);
 491   fTrackSegments = new TClonesArray("AliPHOSTrackSegment",100);
 492   fTrackSegments->SetName("TRACKS");
 493 }
 494
 495 //____________________________________________________________________________
 496 void  AliPHOSTrackSegmentMakerv1::InitParameters()
 497 {
 498   //Initializes parameters
 499   fRcpv      = 10. ;
 500   fRtpc      = 4. ;
 501   fEmcFirst  = 0 ;
 502   fEmcLast   = 0 ;
 503   fCpvFirst  = 0 ;
 504   fCpvLast   = 0 ;
 505   fLinkUpArray = 0 ;
 506   fWrite                   = kTRUE ;
 507 }
 508
 509
 510 //____________________________________________________________________________
 511 void  AliPHOSTrackSegmentMakerv1::MakeLinks()const
 512 {
 513   // Finds distances (links) between all EMC and CPV clusters,
 514   // which are not further apart from each other than fRcpv
 515   // and sort them in accordance with this distance
 516
 517   fLinkUpArray->Clear() ;
 518
 519   AliPHOSCpvRecPoint * cpv ;
 520   AliPHOSEmcRecPoint * emcclu ;
 521
 522   Int_t iLinkUp  = 0 ;
 523
 524   Int_t iEmcRP;
 525   for(iEmcRP = fEmcFirst; iEmcRP < fEmcLast; iEmcRP++ ) {
 526     emcclu = dynamic_cast<AliPHOSEmcRecPoint *>(fEMCRecPoints->At(iEmcRP)) ;
 527
 528     //Bool_t toofar ;
 529     Int_t iCpv = 0 ;
 530     for(iCpv = fCpvFirst; iCpv < fCpvLast;iCpv++ ) {
 531
 532       cpv = dynamic_cast<AliPHOSCpvRecPoint *>(fCPVRecPoints->At(iCpv)) ;
 533       Int_t track = -1 ;
 534       Float_t dx,dz ;
 535       GetDistanceInPHOSPlane(emcclu, cpv, track,dx,dz) ;
 536       if(TMath::Sqrt(dx*dx+dz*dz) < fRcpv ){
 537         new ((*fLinkUpArray)[iLinkUp++])  AliPHOSLink(dx, dz, iEmcRP, iCpv, track) ;
 538       }
 539     }
 540   }
 541
 542   fLinkUpArray->Sort() ;  //first links with smallest distances
 543 }
 544
 545 //____________________________________________________________________________
 546 void  AliPHOSTrackSegmentMakerv1::MakePairs()
 547 {
 548   // Using the previously made list of "links", we found the smallest link - i.e.
 549   // link with the least distance between EMC and CPV and pointing to still
 550   // unassigned RecParticles. We assign these RecPoints to TrackSegment and
 551   // remove them from the list of "unassigned".
 552
 553   //Make arrays to mark clusters already chosen
 554   Int_t * emcExist = 0;
 555   if(fEmcLast > fEmcFirst)
 556     emcExist = new Int_t[fEmcLast-fEmcFirst] ;
 557
 558   Int_t index;
 559   for(index = 0; index <fEmcLast-fEmcFirst; index ++)
 560     emcExist[index] = 1 ;
 561
 562   Bool_t * cpvExist = 0;
 563   if(fCpvLast > fCpvFirst)
 564     cpvExist = new Bool_t[fCpvLast-fCpvFirst] ;
 565   for(index = 0; index <fCpvLast-fCpvFirst; index ++)
 566     cpvExist[index] = kTRUE ;
 567
 568
 569   // Finds the smallest links and makes pairs of CPV and EMC clusters with smallest distance
 570   TIter nextUp(fLinkUpArray) ;
 571
 572   AliPHOSLink * linkUp ;
 573
 574   AliPHOSCpvRecPoint * nullpointer = 0 ;
 575
 576   while ( (linkUp =  static_cast<AliPHOSLink *>(nextUp()) ) ){
 577
 578     if(emcExist[linkUp->GetEmc()-fEmcFirst] != -1){
 579
 580       if(cpvExist[linkUp->GetCpv()-fCpvFirst]){ //CPV still exist
 581          Float_t dx,dz ;
 582          linkUp->GetXZ(dx,dz) ;
 583          new ((* fTrackSegments)[fNTrackSegments])
 584            AliPHOSTrackSegment(dynamic_cast<AliPHOSEmcRecPoint *>(fEMCRecPoints->At(linkUp->GetEmc())) ,
 585                                dynamic_cast<AliPHOSCpvRecPoint *>(fCPVRecPoints->At(linkUp->GetCpv())) ,
 586                                linkUp->GetTrack(),dx,dz) ;
 587
 588        (dynamic_cast<AliPHOSTrackSegment *>(fTrackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
 589        fNTrackSegments++ ;
 590        emcExist[linkUp->GetEmc()-fEmcFirst] = -1 ; //Mark emc  that Cpv was found
 591        //mark CPV recpoint as already used
 592        cpvExist[linkUp->GetCpv()-fCpvFirst] = kFALSE ;
 593       } //if CpvUp still exist
 594     }
 595   }
 596
 597   //look through emc recPoints left without CPV
 598   if(emcExist){ //if there is emc rec point
 599     Int_t iEmcRP ;
 600     for(iEmcRP = 0; iEmcRP < fEmcLast-fEmcFirst  ; iEmcRP++ ){
 601       if(emcExist[iEmcRP] > 0 ){
 602         Int_t track = -1 ;
 603         Float_t dx,dz ;
 604         AliPHOSEmcRecPoint *emcclu = dynamic_cast<AliPHOSEmcRecPoint *>(fEMCRecPoints->At(iEmcRP+fEmcFirst));
 605         GetDistanceInPHOSPlane(emcclu, 0, track,dx,dz);
 606         if(track<0)
 607           new ((*fTrackSegments)[fNTrackSegments]) AliPHOSTrackSegment(emcclu,nullpointer) ;
 608         else
 609           new ((*fTrackSegments)[fNTrackSegments]) AliPHOSTrackSegment(emcclu,0,track,dx,dz);
 610         (dynamic_cast<AliPHOSTrackSegment *>(fTrackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
 611         fNTrackSegments++;
 612       }
 613     }
 614   }
 615   delete [] emcExist ;
 616   delete [] cpvExist ;
 617 }
 618
 619 //____________________________________________________________________________
 620 void AliPHOSTrackSegmentMakerv1::Clusters2TrackSegments(Option_t *option)
 621 {
 622   // Steering method to perform track segment construction for the current event
 623   // Returns an array with the found track-segments.
 624
 625   if(strstr(option,"tim"))
 626     gBenchmark->Start("PHOSTSMaker");
 627
 628   if(strstr(option,"print")) {
 629     Print() ;
 630     return ;
 631   }
 632
 633   //Make some initializations
 634   fNTrackSegments = 0 ;
 635   fEmcFirst = 0 ;
 636   fEmcLast  = 0 ;
 637   fCpvFirst = 0 ;
 638   fCpvLast  = 0 ;
 639
 640   fTrackSegments->Clear();
 641
 642   //   if(!ReadRecPoints(ievent))   continue; //reads RecPoints for event ievent
 643
 644   for(fModule = 1; fModule <= fGeom->GetNModules() ; fModule++ ) {
 645     FillOneModule() ;
 646     MakeLinks() ;
 647     MakePairs() ;
 648   }
 649
 650   if(strstr(option,"deb"))
 651     PrintTrackSegments(option);
 652
 653   if(strstr(option,"tim")){
 654     gBenchmark->Stop("PHOSTSMaker");
 655     Info("Exec", "took %f seconds for making TS",
 656          gBenchmark->GetCpuTime("PHOSTSMaker"));
 657   }
 658 }
 659
 660 //____________________________________________________________________________
 661 void AliPHOSTrackSegmentMakerv1::Print(const Option_t *)const
 662 {
 663   //  Print TrackSegmentMaker parameters
 664
 665   TString message("") ;
 666   if( strcmp(GetName(), "") != 0 ) {
 667     message = "\n======== AliPHOSTrackSegmentMakerv1 ========\n" ;
 668     message += "Making Track segments\n" ;
 669     message += "with parameters:\n" ;
 670     message += "     Maximal EMC - CPV distance (cm) %f\n" ;
 671     message += "============================================\n" ;
 672     Info("Print", message.Data(),fRcpv) ;
 673   }
 674   else
 675     Info("Print", "AliPHOSTrackSegmentMakerv1 not initialized ") ;
 676 }
 677
 678 //____________________________________________________________________________
 679 void AliPHOSTrackSegmentMakerv1::PrintTrackSegments(Option_t * option)
 680 {
 681   // option deb - prints # of found TrackSegments
 682   // option deb all - prints as well indexed of found RecParticles assigned to the TS
 683
 684   Info("PrintTrackSegments", "Results from TrackSegmentMaker:") ;
 685   printf("        Found %d TrackSegments\n", fTrackSegments->GetEntriesFast() );
 686
 687   if(strstr(option,"all")) {  // printing found TS
 688     printf("TrackSegment #  EMC RP#  CPV RP#\n") ;
 689     Int_t index;
 690     for (index = 0 ; index <fTrackSegments->GetEntriesFast() ; index++) {
 691       AliPHOSTrackSegment * ts = (AliPHOSTrackSegment * )fTrackSegments->At(index) ;
 692       printf("   %d           %d        %d \n", ts->GetIndexInList(), ts->GetEmcIndex(), ts->GetCpvIndex() ) ;
 693     }
 694   }
 695 }
 696