TRD/AliTRDclusterizerV0.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
  16 /* $Id$ */
  17
  18 ///////////////////////////////////////////////////////////////////////////////
  19 //                                                                           //
  20 // TRD cluster finder for the fast simulator. It takes the hits from the     //
  21 // fast simulator (one hit per plane) and transforms them                    //
  22 // into cluster, by applying position smearing and merging                   //
  23 // of nearby cluster. The searing is done uniformly in z-direction           //
  24 // over the length of a readout pad. In rphi-direction a Gaussian            //
  25 // smearing is applied with a sigma given by fRphiSigma.                     //
  26 // Clusters are considered as overlapping when they are closer in            //
  27 // rphi-direction than the value defined in fRphiDist.                       //
  28 // Use the macro fastClusterCreate.C to create the cluster.                  //
  29 //                                                                           //
  30 ///////////////////////////////////////////////////////////////////////////////
  31
  32 #include <TRandom.h>
  33 #include <TTree.h>
  34
  35 #include "AliRun.h"
  36
  37 #include "AliTRD.h"
  38 #include "AliTRDclusterizerV0.h"
  39 #include "AliTRDhit.h"
  40 #include "AliTRDgeometry.h"
  41 #include "AliTRDrecPoint.h"
  42 #include "AliTRDparameter.h"
  43
  44 ClassImp(AliTRDclusterizerV0)
  45
  46 //_____________________________________________________________________________
  47   AliTRDclusterizerV0::AliTRDclusterizerV0():AliTRDclusterizer(), fTRD(NULL)
  48 {
  49   //
  50   // AliTRDclusterizerV0 default constructor
  51   //
  52
  53 }
  54
  55 //_____________________________________________________________________________
  56 AliTRDclusterizerV0::AliTRDclusterizerV0(const Text_t* name, const Text_t* title)
  57                     :AliTRDclusterizer(name,title), fTRD(NULL)
  58 {
  59   //
  60   // AliTRDclusterizerV0 default constructor
  61   //
  62
  63   Init();
  64
  65 }
  66
  67 //_____________________________________________________________________________
  68 AliTRDclusterizerV0::~AliTRDclusterizerV0()
  69 {
  70   //
  71   // AliTRDclusterizerV0 destructor
  72   //
  73
  74 }
  75
  76 //_____________________________________________________________________________
  77 void AliTRDclusterizerV0::Init()
  78 {
  79   //
  80   // Initializes the cluster finder
  81   //
  82
  83   // Position resolution in rphi-direction
  84   fRphiSigma  = 0.02;
  85   // Minimum distance of non-overlapping cluster
  86   fRphiDist   = 1.0;
  87
  88 }
  89
  90 //_____________________________________________________________________________
  91 Bool_t AliTRDclusterizerV0::MakeClusters()
  92 {
  93   //
  94   // Generates the cluster
  95   //
  96
  97   // Get the TRD object
  98   fTRD = (AliTRD*) gAlice->GetDetector("TRD");
  99   if (!fTRD) {
 100     printf("AliTRDclusterizerV0::makClusters -- ");
 101     printf("No TRD detector object found\n");
 102     return kFALSE;
 103   }
 104
 105   if (fTRD->IsVersion() != 0) {
 106     printf("AliTRDclusterizerV0::MakeCluster -- ");
 107     printf("TRD must be version 0 (fast simulator).\n");
 108     return kFALSE;
 109   }
 110
 111   // Get the geometry
 112   AliTRDgeometry *geo = fTRD->GetGeometry();
 113
 114   // Create a default parameter class if none is defined
 115   if (!fPar) {
 116     fPar = new AliTRDparameter("TRDparameter","Standard TRD parameter");
 117     if (fVerbose > 0) {
 118       printf("<AliTRDclusterizerV0::MakeCluster> ");
 119       printf("Create the default parameter object.\n");
 120     }
 121   }
 122
 123   printf("AliTRDclusterizerV0::MakeCluster -- ");
 124   printf("Start creating cluster.\n");
 125
 126   Int_t nBytes = 0;
 127
 128   AliTRDhit *hit;
 129
 130   // Get the pointer to the hit tree
 131   TTree     *hitTree      = fTRD->TreeH();
 132   // Get the pointer to the reconstruction tree
 133   TTree     *clusterTree  = gAlice->TreeR();
 134
 135   TObjArray *chamberArray = new TObjArray();
 136
 137   // Get the number of entries in the hit tree
 138   // (Number of primary particles creating a hit somewhere)
 139   Int_t nTrack = (Int_t) hitTree->GetEntries();
 140
 141   // Loop through all the chambers
 142   for (Int_t icham = 0; icham < AliTRDgeometry::Ncham(); icham++) {
 143     for (Int_t iplan = 0; iplan < AliTRDgeometry::Nplan(); iplan++) {
 144       for (Int_t isect = 0; isect < AliTRDgeometry::Nsect(); isect++) {
 145
 146         Int_t   nColMax     = fPar->GetColMax(iplan);
 147         Float_t row0        = fPar->GetRow0(iplan,icham,isect);
 148         Float_t col0        = fPar->GetCol0(iplan);
 149         Float_t time0       = fPar->GetTime0(iplan);
 150
 151         Float_t rowPadSize  = fPar->GetRowPadSize(iplan,icham,isect);
 152         Float_t colPadSize  = fPar->GetColPadSize(iplan);
 153         Float_t timeBinSize = fPar->GetDriftVelocity()
 154                             / fPar->GetSamplingFrequency();
 155
 156         // Loop through all entries in the tree
 157         for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
 158
 159           gAlice->ResetHits();
 160           nBytes += hitTree->GetEvent(iTrack);
 161
 162           // Get the number of hits in the TRD created by this particle
 163           Int_t nHit = fTRD->Hits()->GetEntriesFast();
 164
 165           // Loop through the TRD hits
 166           for (Int_t iHit = 0; iHit < nHit; iHit++) {
 167
 168             if (!(hit = (AliTRDhit *) fTRD->Hits()->UncheckedAt(iHit)))
 169               continue;
 170
 171             Float_t pos[3];
 172                     pos[0]   = hit->X();
 173                     pos[1]   = hit->Y();
 174                     pos[2]   = hit->Z();
 175             Int_t   track    = hit->Track();
 176             Int_t   detector = hit->GetDetector();
 177             Int_t   plane    = geo->GetPlane(detector);
 178             Int_t   sector   = geo->GetSector(detector);
 179             Int_t   chamber  = geo->GetChamber(detector);
 180
 181             if ((sector  != isect) ||
 182                 (plane   != iplan) ||
 183                 (chamber != icham))
 184               continue;
 185
 186             // Rotate the sectors on top of each other
 187             Float_t rot[3];
 188             geo->Rotate(detector,pos,rot);
 189
 190             // Add this recPoint to the temporary array for this chamber
 191             AliTRDrecPoint *recPoint = new AliTRDrecPoint("");
 192             recPoint->SetLocalRow(rot[2]);
 193             recPoint->SetLocalCol(rot[1]);
 194             recPoint->SetLocalTime(rot[0]);
 195             recPoint->SetEnergy(0);
 196             recPoint->SetDetector(detector);
 197             recPoint->AddDigit(track);
 198             chamberArray->Add(recPoint);
 199
 200           }
 201
 202         }
 203
 204         // Loop through the temporary cluster-array
 205         for (Int_t iClus1 = 0; iClus1 < chamberArray->GetEntries(); iClus1++) {
 206
 207           AliTRDrecPoint *recPoint1 = (AliTRDrecPoint *)
 208                                       chamberArray->UncheckedAt(iClus1);
 209           Float_t row1  = recPoint1->GetLocalRow();
 210           Float_t col1  = recPoint1->GetLocalCol();
 211           Float_t time1 = recPoint1->GetLocalTime();
 212
 213           if (recPoint1->GetEnergy() < 0) continue;        // Skip marked cluster
 214
 215           const Int_t kNsave  = 5;
 216           Int_t idxSave[kNsave];
 217           Int_t iSave = 0;
 218
 219           const Int_t kNsaveTrack = 3;
 220           Int_t tracks[kNsaveTrack];
 221           tracks[0] = recPoint1->GetDigit(0);
 222
 223           // Check the other cluster to see, whether there are close ones
 224           for (Int_t iClus2 = iClus1 + 1; iClus2 < chamberArray->GetEntries(); iClus2++) {
 225
 226             AliTRDrecPoint *recPoint2 = (AliTRDrecPoint *)
 227                                         chamberArray->UncheckedAt(iClus2);
 228             Float_t row2 = recPoint2->GetLocalRow();
 229             Float_t col2 = recPoint2->GetLocalCol();
 230
 231             if ((TMath::Abs(row1 - row2) < rowPadSize) ||
 232                 (TMath::Abs(col1 - col2) <  fRphiDist)) {
 233               if (iSave == kNsave) {
 234                 printf("AliTRDclusterizerV0::MakeCluster -- ");
 235                 printf("Boundary error: iSave = %d, kNsave = %d.\n"
 236                       ,iSave,kNsave);
 237               }
 238               else {
 239                 idxSave[iSave]  = iClus2;
 240                 iSave++;
 241                 if (iSave < kNsaveTrack) tracks[iSave] = recPoint2->GetDigit(0);
 242               }
 243             }
 244           }
 245
 246           // Merge close cluster
 247           Float_t rowMerge = row1;
 248           Float_t colMerge = col1;
 249           if (iSave) {
 250             for (Int_t iMerge = 0; iMerge < iSave; iMerge++) {
 251               AliTRDrecPoint *recPoint2 =
 252                 (AliTRDrecPoint *) chamberArray->UncheckedAt(idxSave[iMerge]);
 253               rowMerge += recPoint2->GetLocalRow();
 254               colMerge += recPoint2->GetLocalCol();
 255               recPoint2->SetEnergy(-1);     // Mark merged cluster
 256             }
 257             rowMerge /= (iSave + 1);
 258             colMerge /= (iSave + 1);
 259           }
 260
 261           Float_t smear[3];
 262
 263           // The position smearing in row-direction (uniform over pad width)
 264           Int_t row = (Int_t) ((rowMerge - row0) / rowPadSize);
 265           smear[0]  = (row + gRandom->Rndm()) * rowPadSize + row0;
 266
 267           // The position smearing in rphi-direction (Gaussian)
 268           smear[1] = 0;
 269           do
 270             smear[1] = gRandom->Gaus(colMerge,fRphiSigma);
 271           while ((smear[1] < col0                        ) ||
 272                  (smear[1] > col0 + nColMax * colPadSize));
 273
 274           // Time direction stays unchanged
 275           smear[2] = time1;
 276
 277           // Transform into local coordinates
 278           smear[0] = (Int_t) ((smear[0] -  row0) /  rowPadSize);
 279           smear[1] = (Int_t) ((smear[1] -  col0) /  colPadSize);
 280           smear[2] = (Int_t) ((time0 - smear[2]) / timeBinSize);
 281
 282           // Add the smeared cluster to the output array
 283           Int_t   detector  = recPoint1->GetDetector();
 284           Int_t   tr[9]     = { -1   };
 285           Float_t pos[3];
 286           Float_t sigma[2]  = {  0.0 };
 287           pos[0] = smear[1];
 288           pos[1] = smear[0];
 289           pos[2] = (time0 - smear[2]) / timeBinSize;
 290           AddCluster(pos,detector,0.0,tr,sigma,0);
 291
 292         }
 293
 294         // Clear the temporary cluster-array and delete the cluster
 295         chamberArray->Delete();
 296
 297       }
 298     }
 299   }
 300
 301   printf("AliTRDclusterizerV0::MakeCluster -- ");
 302   printf("Found %d points.\n",RecPoints()->GetEntries());
 303   printf("AliTRDclusterizerV0::MakeCluster -- ");
 304   printf("Fill the cluster tree.\n");
 305   clusterTree->Fill();
 306
 307   return kTRUE;
 308
 309 }