HBTAN/AliHBTClusterMap.cxx

   1 #include "AliHBTClusterMap.h"
   2
   3 //_________________________________________________
   4 ///////////////////////////////////////////////////
   5 //
   6 // class AliHBTClusterMap
   7 //
   8 // class that describes cluster occupation at TPC
   9 // Each padraw has a corresponding bit in fPadRawMap
  10 //
  11 //
  12 // more info: http://aliweb.cern.ch/people/skowron/analyzer/index.html
  13 // Piotr.Skowronski@cern.ch
  14 //
  15 ///////////////////////////////////////////////////////////////////////////
  16
  17
  18 #include "AliESDtrack.h"
  19 #include "AliTPCtrack.h"
  20 #include "AliHBTParticle.h"
  21 #include <TString.h>
  22 const Int_t AliHBTClusterMap::fNPadRows = 159;
  23
  24 AliHBTClusterMap::AliHBTClusterMap():
  25  fPadRawMap(fNPadRows)
  26 {
  27 //ctor
  28
  29 }
  30 /***********************************************************************/
  31 AliHBTClusterMap::AliHBTClusterMap(AliESDtrack* track):
  32  fPadRawMap(fNPadRows)
  33 {
  34  //cotor
  35  track->Print();//to shut up compiler warning
  36
  37 }
  38 /***********************************************************************/
  39
  40 AliHBTClusterMap::AliHBTClusterMap(AliTPCtrack* track):
  41  fPadRawMap(fNPadRows)
  42 {
  43  //cotor
  44
  45  //Does not work since indeces in the claster index array
  46  //in the TPC track does not correspond to the padraw segmatation
  47  if (AliHBTParticle::GetDebug() > 9)
  48    Info("AliHBTClusterMap",
  49       "#####################################################################");
  50  if (track == 0x0)
  51   {
  52     Error("AliHBTClusterMap","Pointer to TPC track is NULL");
  53     return;
  54   }
  55  Int_t prevrow = -1;
  56  Int_t i = 0;
  57  for ( ; i < track->GetNumberOfClusters(); i++)
  58   {
  59     Int_t idx = track->GetClusterIndex(i);
  60     Int_t sect = (idx&0xff000000)>>24;
  61     Int_t row = (idx&0x00ff0000)>>16;
  62     if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
  63     if (AliHBTParticle::GetDebug() > 9)
  64       Info("AliHBTClusterMap","Cl.idx is %d, sect %d, row %d",idx,sect,row);
  65
  66     fPadRawMap.SetBitNumber(row,kTRUE);
  67
  68     //Fill the gap between previous row and this row with 0 bits
  69     if (prevrow < 0)
  70      {
  71        prevrow = row;//if previous bit was not assigned yet == this is the first one
  72      }
  73     else
  74      { //we don't know the order (inner to outer or reverse)
  75        //just to be save in case it is going change
  76        Int_t n = 0, m = 0;
  77        if (prevrow < row)
  78         {
  79           n = prevrow;
  80           m = row;
  81         }
  82        else
  83         {
  84           n = row;
  85           m = prevrow;
  86         }
  87        for (Int_t j = n+1; j < m; j++)
  88         {
  89           fPadRawMap.SetBitNumber(j,kFALSE);
  90         }
  91        prevrow = row;
  92      }
  93   }
  94   if (AliHBTParticle::GetDebug() > 2)
  95    {
  96      Print();
  97    }
  98 }
  99 /***********************************************************************/
 100
 101 void AliHBTClusterMap::Print() const
 102 {
 103 //Prints the bit map
 104   TString msg;
 105   for ( Int_t i = 0; i < fNPadRows; i++)
 106    {
 107      if ( fPadRawMap.TestBitNumber(i) )
 108       {
 109         msg+="1";
 110       }
 111      else
 112       {
 113         msg+="0";
 114       }
 115    }
 116   Info("AliHBTClusterMap","BitMap is\n  %s",msg.Data());
 117
 118 }
 119
 120 /***********************************************************************/
 121
 122 Float_t AliHBTClusterMap::GetOverlapFactor(const AliHBTClusterMap& clmap) const
 123 {
 124   //Returns quality factor FQ = Sum(An)/Sum(clusters)
 125   //      | -1; if both tracks have a cluster on padrow n
 126   //An = <  0; if neither track has a cluster on padrow n
 127   //     |  1; if only one trackhas a cluster on padrow n
 128   // Returned value ranges between
 129   //  -0.5 (low probability that these tracks are a split track)
 130   //  and
 131   //   1.0 (high probability that these tracks are a split track)
 132   TString msg1;
 133   TString msg2;
 134
 135   Int_t nh = 0;
 136   Int_t an = 0;
 137   for ( Int_t i = 0; i < fNPadRows; i++)
 138    {
 139      Bool_t x = HasClAtPadRow(i);
 140      Bool_t y = clmap.HasClAtPadRow(i);
 141
 142      if (x) msg1+="1";else msg1+="0";
 143      if (y) msg2+="1";else msg2+="0";
 144
 145      if (x && y)//both have clasters
 146       {
 147        an--;
 148        nh+=2;
 149       }
 150      else
 151       {
 152
 153        if (x || y)//only one have cluters
 154         {
 155           an++;
 156           nh++;
 157         }
 158       }
 159    }
 160
 161
 162   Float_t retval = 0.0;
 163   if (nh > 0) retval = ((Float_t)an)/((Float_t)nh);
 164   else Warning("GetOverlapFactor","Number of counted cluters is 0.");
 165
 166   if (AliHBTParticle::GetDebug() > 2)
 167    {
 168      Info("GetOverlapFactor","Splitting Quality Factor is %f. SumAn = %d, SumClusters %d",retval,an,nh);
 169      if (retval == 1.0)
 170       {
 171         Print();
 172         Info("AliHBTClusterMap","BitMap is\n  %s\n",msg1.Data());
 173         clmap.Print();
 174         Info("AliHBTClusterMap","BitMap is\n  %s\n\n\n\n",msg2.Data());
 175       }
 176      if (retval == -.5)
 177       {
 178         Print();
 179         Info("AliHBTClusterMap","BitMap is\n  %s\n",msg1.Data());
 180         clmap.Print();
 181         Info("AliHBTClusterMap","BitMap is\n  %s\n\n\n\n",msg2.Data());
 182       }
 183    }
 184
 185   return retval;
 186 }