HLT/TPCLib/AliHLTTPCGlobalMerger.cxx

   1 // @(#) $Id$
   2 // Original: AliHLTGlobalMerger.cxx,v 1.13 2005/06/14 10:55:21 cvetan
   3
   4 /**************************************************************************
   5  * This file is property of and copyright by the ALICE HLT Project        *
   6  * ALICE Experiment at CERN, All rights reserved.                         *
   7  *                                                                        *
   8  * Primary Authors: Uli Frankenfeld, maintained by                          *
   9  *                  Matthias Richter <Matthias.Richter@ift.uib.no>        *
  10  *                  for The ALICE HLT Project.                            *
  11  *                                                                        *
  12  * Permission to use, copy, modify and distribute this software and its   *
  13  * documentation strictly for non-commercial purposes is hereby granted   *
  14  * without fee, provided that the above copyright notice appears in all   *
  15  * copies and that both the copyright notice and this permission notice   *
  16  * appear in the supporting documentation. The authors make no claims     *
  17  * about the suitability of this software for any purpose. It is          *
  18  * provided "as is" without express or implied warranty.                  *
  19  **************************************************************************/
  20
  21 /** @file   AliHLTTPCGlobalMerger.cxx
  22     @author Uli Frankenfeld, maintained by Matthias Richter
  23     @date
  24     @brief  The HLT TPC slice merger
  25 */
  26
  27 #include "AliHLTTPCLogging.h"
  28 #include "AliHLTTPCGlobalMerger.h"
  29 #include "AliHLTTPCTrack.h"
  30 #include "AliHLTTPCTransform.h"
  31 #include "AliHLTTPCTrackArray.h"
  32
  33 #if __GNUC__ >= 3
  34 using namespace std;
  35 #endif
  36
  37 ClassImp(AliHLTTPCGlobalMerger)
  38
  39 AliHLTTPCGlobalMerger::AliHLTTPCGlobalMerger()
  40   :
  41   fNSlices(0),
  42   fFirst(0),
  43   fLast(0)
  44 {
  45   //Default constructor. Use Setup to specify and setup the necessary parameters and arrays.
  46   Is2Global(kTRUE);
  47   SetParameter(0,0,0,0,0);
  48 }
  49
  50
  51 AliHLTTPCGlobalMerger::~AliHLTTPCGlobalMerger()
  52 {
  53   //Destructor
  54 }
  55
  56 void AliHLTTPCGlobalMerger::Setup(Int_t first,Int_t last)
  57 {
  58   //Used to setup the arrays and everything
  59
  60   fNSlices = last-first+1;
  61   fFirst = first;
  62   fLast = last;
  63   InitMerger(last-first+1);
  64 }
  65
  66 void AliHLTTPCGlobalMerger::InitSlice(Int_t slice)
  67 {
  68   //
  69   // Select Sector The following FillTracks call will
  70   // fill this Sector
  71   //
  72   fSlice = slice;
  73   fCurrentTracks = fSlice - fFirst;
  74 }
  75
  76 Double_t AliHLTTPCGlobalMerger::CheckTracks(AliHLTTPCTrack *innertrack,AliHLTTPCTrack *outertrack,Int_t slice)
  77 {
  78   //Compare the tracks by propagating the outermost track to the last and first point plane
  79   //of the innermost track. This plane is defined by the padrow plane where these points
  80   //are.
  81
  82   if(innertrack->GetCharge()!=outertrack->GetCharge()) return -1;
  83
  84   Float_t angle = 0;//perpendicular to the padrowplane (in local system)
  85   AliHLTTPCTransform::Local2GlobalAngle(&angle,slice);
  86   Double_t dx[2],dy[2],dz[2];
  87   Double_t diff =-1;
  88   AliHLTTPCTrack *tracks[2];
  89   tracks[0] = innertrack;
  90   tracks[1] = outertrack;
  91   SortGlobalTracks(tracks,2);
  92   innertrack = tracks[0];
  93   outertrack = tracks[1];
  94
  95   Float_t point[3];
  96
  97   point[0]=innertrack->GetLastPointX();
  98   point[1]=innertrack->GetLastPointY();
  99   point[2]=innertrack->GetLastPointZ();
 100   AliHLTTPCTransform::Global2LocHLT(point,slice);
 101
 102   outertrack->CalculateReferencePoint(angle,point[0]);//local x = global distance to padrowplane
 103   if(!outertrack->IsPoint()) return diff;
 104   dx[0] = fabs(outertrack->GetPointX()-innertrack->GetLastPointX());
 105   dy[0] = fabs(outertrack->GetPointY()-innertrack->GetLastPointY());
 106   dz[0] = fabs(outertrack->GetPointZ()-innertrack->GetLastPointZ());
 107
 108   point[0]=innertrack->GetFirstPointX();
 109   point[1]=innertrack->GetFirstPointY();
 110   point[2]=innertrack->GetFirstPointZ();
 111   AliHLTTPCTransform::Global2LocHLT(point,slice);
 112
 113   outertrack->CalculateReferencePoint(angle,point[0]);//local x = global distance to padrowplane
 114   if(!outertrack->IsPoint()) return diff;
 115   dx[1] = fabs(outertrack->GetPointX()-innertrack->GetFirstPointX());
 116   dy[1] = fabs(outertrack->GetPointY()-innertrack->GetFirstPointY());
 117   dz[1] = fabs(outertrack->GetPointZ()-innertrack->GetFirstPointZ());
 118
 119   diff=0;//This was a tough bug to find....
 120   for(Int_t i=0; i<2; i++)
 121     diff += sqrt(dx[i]*dx[i] + dy[i]*dy[i] + dz[i]*dz[i]);
 122   return diff;
 123 }
 124
 125 void AliHLTTPCGlobalMerger::SlowMerge(const Char_t *path)
 126 {
 127   //Tuning of parameters. This matches _all_ tracks between two neighbouring
 128   //slices, and merges the ones which are closest in space. The difference
 129   //is written to a ntuppel, which can be used as a input for the SetParameters
 130   //when using normal Merge function.
 131
 132
 133   void* ntuple=GetNtuple();
 134   AliHLTTPCTrack *track[2];
 135   AliHLTTPCTrackArray *tout = GetOutTracks();
 136   if(fNSlices<2)
 137     {
 138       LOG(AliHLTTPCLog::kWarning,"AliHLTTPCGlobalMerger::SlowMerge","Slice Number")
 139         <<"Need more than one Slice!"<<ENDLOG;
 140       return;
 141     }
 142
 143   for(Int_t i=0; i<fNSlices; i++)
 144     {
 145       //if(fNSlices!=18 && i+1 == fNSlices) continue;
 146       Int_t slice = fFirst + i;
 147       AliHLTTPCTrackArray *ttt0=GetInTracks(i);
 148       Int_t slice2 = i+1;
 149       //if(slice2==fNSlices) slice2 =0;
 150
 151       //Make sure slices are on the same side of the TPC
 152       if(slice2 == 18) slice2=0;
 153       else if(slice2 == 36) slice2=18;
 154       AliHLTTPCTrackArray *ttt1=GetInTracks(slice2);
 155       //10 degrees -> the border of the slices in local coordinates
 156       Float_t angle = AliHLTTPCTransform::Pi()/18;
 157       AliHLTTPCTransform::Local2GlobalAngle(&angle,slice);
 158
 159       //In the two following cases, the angle is 2*pi, so set it back to 0 in order for
 160       //the calculation of crossing points to be correct.
 161       if(slice==17 || slice==35)
 162         angle=0;
 163       if(i==0)
 164         ttt0->QSort();
 165       ttt1->QSort();
 166       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
 167         {
 168           AliHLTTPCTrack *track0=ttt0->GetCheckedTrack(s0);
 169           if(!track0) continue;
 170           track0->CalculateHelix();
 171           track0->CalculateEdgePoint(angle);
 172           //      if(track0->IsPoint()) AddTrack(tout,track0);
 173         }
 174       for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
 175         {
 176           AliHLTTPCTrack *track1=ttt1->GetCheckedTrack(s1);
 177           if(!track1) continue;
 178           track1->CalculateHelix();
 179           track1->CalculateEdgePoint(angle);
 180           //      if(track1->IsPoint())  AddTrack(tout,track1);
 181         }
 182       Bool_t merge = kTRUE;
 183       while(merge)
 184         {
 185           Int_t min0=-1,min1=-1;
 186           Double_t min=5;
 187           Int_t n0=ttt0->GetNTracks(),n1=ttt1->GetNTracks();
 188           for(Int_t s0=0;s0<n0;s0++)
 189             {
 190               AliHLTTPCTrack *track0=ttt0->GetCheckedTrack(s0);
 191               if(!track0) continue;
 192               if(!track0->IsPoint()) continue;
 193               for(Int_t s1=0;s1<n1;s1++)
 194                 {
 195                   AliHLTTPCTrack *track1=ttt1->GetCheckedTrack(s1);
 196                   if(!track1) continue;
 197                   if(!track1->IsPoint()) continue;
 198
 199                   //Double_t diff = TrackDiff(track0,track1,angle);
 200                   Double_t diff = CheckTracks(track0,track1,slice);
 201                   //PrintDiff(track0,track1);
 202                   if(diff>=0&&diff<min)
 203                     {
 204                       min=diff;
 205                       min0=s0;
 206                       min1=s1;
 207                     }
 208                 }
 209             }
 210           if(min0>=0&&min1>=0)
 211             {
 212               AliHLTTPCTrack *track0=ttt0->GetTrack(min0);
 213               AliHLTTPCTrack *track1=ttt1->GetTrack(min1);
 214               track[0] = track0;
 215               track[1] = track1;
 216               SortGlobalTracks(track,2);
 217               track1->CalculateEdgePoint((angle+AliHLTTPCTransform::Pi()/9));
 218               if(track1->IsPoint())//Check if the track will cross the boundary of yet another slice.
 219                 MultiMerge(ttt1,track,2);
 220               else
 221                 MultiMerge(tout,track,2);
 222               track0->CalculateReferencePoint(angle);
 223               track1->CalculateReferencePoint(angle);
 224               //PrintDiff(track0,track1);
 225               FillNtuple(ntuple,track0,track1);
 226               ttt0->Remove(min0);
 227               ttt1->Remove(min1);
 228
 229             }
 230           else merge = kFALSE;
 231         }
 232       ttt0->Compress();
 233       ttt1->Compress();
 234       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCGlobalMerger::SlowMerge","Result")
 235         <<AliHLTTPCLog::kDec<<"Merged Tracks: "<<tout->GetNTracks()<<" at:"
 236         <<angle<<ENDLOG;
 237     }
 238   const int namelen=1024;
 239   Char_t fname[namelen];
 240   snprintf(fname,namelen,"%s/merge_parameters.root",path);
 241   WriteNtuple(fname,ntuple);
 242 }
 243
 244 void AliHLTTPCGlobalMerger::Merge()
 245 {
 246   //Normal merging procedure. Matches tracks which are within limits
 247   //set by SetParameters. Parameters can be tuned by SlowMerge.
 248
 249   AliHLTTPCTrack *track[2];
 250   AliHLTTPCTrackArray *tout = GetOutTracks();
 251   if(fNSlices<2)
 252     {
 253       LOG(AliHLTTPCLog::kWarning,"AliHLTTPCGlobalMerger::Merge","Slice Number")
 254         <<"Need more than one Slice!"<<ENDLOG;
 255       return;
 256     }
 257   for(Int_t i=0; i<fNSlices; i++)
 258     {
 259       //if(fNSlices!=18 && i+1 == fNSlices) continue;
 260       Int_t slice = fFirst + i;
 261       AliHLTTPCTrackArray *ttt0=GetInTracks(i);
 262       Int_t slice2 = i+1;
 263       //if(slice2==fNSlices) slice2 =0;
 264
 265       //Make sure slices are on the same side of the TPC
 266       if(slice2 == 18) slice2=0;
 267       else if(slice2 == 36) slice2=18;
 268       AliHLTTPCTrackArray *ttt1=GetInTracks(slice2);
 269       //10 degrees -> the border of the slices in local coordinates
 270       Float_t angle = AliHLTTPCTransform::Pi()/18;
 271       AliHLTTPCTransform::Local2GlobalAngle(&angle,slice);
 272
 273       //In the two following cases, the angle is 2*pi, so set it back to 0 in order for
 274       //the calculation of crossing points to be correct.
 275       if(slice==17 || slice==35)
 276         angle=0;
 277       if(i==0)
 278         ttt0->QSort();
 279       ttt1->QSort();
 280       Bool_t *ismatched0  = new Bool_t[ttt0->GetNTracks()];
 281       Bool_t *ismatched1  = new Bool_t[ttt1->GetNTracks()];
 282       Int_t n0=0,n1=0;
 283       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
 284         {
 285           ismatched0[s0]=kFALSE;
 286           AliHLTTPCTrack *track0=ttt0->GetCheckedTrack(s0);
 287           if(!track0) continue;
 288           //track0->CalculateHelix();    //This it done in TrackArray
 289           track0->CalculateEdgePoint(angle);
 290           if(track0->IsPoint())
 291             {
 292               n0++;
 293               track0->CalculateReferencePoint(angle);
 294             }
 295         }
 296       for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
 297         {
 298           ismatched1[s1]=kFALSE;
 299           AliHLTTPCTrack *track1=ttt1->GetCheckedTrack(s1);
 300           if(!track1) continue;
 301           //track1->CalculateHelix();   //This is done in TrackArray
 302           track1->CalculateEdgePoint(angle);
 303           if(track1->IsPoint())
 304             {
 305               n1++;
 306               track1->CalculateReferencePoint(angle);
 307             }
 308         }
 309       for(Int_t s0=0;s0<ttt0->GetNTracks();s0++)
 310         {
 311           if(ismatched0[s0]) continue;
 312           AliHLTTPCTrack *track0=ttt0->GetCheckedTrack(s0);
 313           if(!track0) continue;
 314           if(!track0->IsPoint()) continue;
 315           for(Int_t s1=0;s1<ttt1->GetNTracks();s1++)
 316             {
 317               if(ismatched1[s1]) continue;
 318               AliHLTTPCTrack *track1=ttt1->GetCheckedTrack(s1);
 319               if(!track1) continue;
 320               if(!track1->IsPoint()) continue;
 321               if(IsRTrack(track0,track1))
 322                 {
 323                   track[0] = track0;
 324                   track[1] = track1;
 325                   SortGlobalTracks(track,2);
 326                   Double_t r0 = pow(track[0]->GetLastPointX(),2)+
 327                     pow(track[0]->GetLastPointY(),2);
 328                   Double_t r1 = pow(track[1]->GetFirstPointX(),2)+
 329                     pow(track[1]->GetFirstPointY(),2);
 330                   if(r0<r1)
 331                     {
 332                       MultiMerge(tout,track,2);
 333                       ismatched0[s0]=kTRUE;
 334                       ismatched1[s1]=kTRUE;
 335                       ttt0->Remove(s0);
 336                       ttt1->Remove(s1);
 337                       break;
 338                       /*
 339                         The track is merged, so we will _not_ look for more matches.
 340                         Because there could easily be more matches, if a track is being
 341                         split within the sector....
 342                       */
 343                     }
 344                 }
 345             }
 346         }
 347       LOG(AliHLTTPCLog::kInformational,"AliHLTTPCGlobalMerger::Merge","Result")
 348         <<AliHLTTPCLog::kDec<<"slice0: "<<n0<<" slice1: "<<n1
 349         <<" Merged Tracks: "<<tout->GetNTracks()<<ENDLOG;
 350       delete [] ismatched0;
 351       delete [] ismatched1;
 352     }
 353 }