ITS/UPGRADE/AliITSURecoLayer.cxx

   1 #include <TClonesArray.h>
   2 #include "AliITSURecoLayer.h"
   3 #include "AliITSUGeomTGeo.h"
   4 #include "AliITSsegmentation.h"
   5 #include "AliITSUAux.h"
   6 #include "AliITSUClusterPix.h"
   7
   8 using namespace AliITSUAux;
   9 using namespace TMath;
  10
  11 ClassImp(AliITSURecoLayer)
  12
  13
  14 //______________________________________________________
  15 AliITSURecoLayer::AliITSURecoLayer(const char* name)
  16   :fActiveID(-1)
  17   ,fNSensors(0)
  18   ,fNSensInLadder(0)
  19   ,fNLadders(0)
  20   ,fR(0)
  21   ,fRMax(0)
  22   ,fRMin(0)
  23   ,fZMax(0)
  24   ,fZMin(0)
  25   ,fPhiLadMax(0)
  26   ,fPhiLadMin(0)
  27   ,fPhiOffs(0)
  28   ,fSensDZInv(0)
  29   ,fDPhiLadInv(0)
  30   ,fMaxStep(0.5)
  31   ,fSensors(0)
  32   ,fITSGeom(0)
  33   ,fClusters(0)
  34 {
  35   // def. c-tor
  36   SetNameTitle(name,name);
  37 }
  38
  39 //______________________________________________________
  40 AliITSURecoLayer::AliITSURecoLayer(const char* name, Int_t activeID, AliITSUGeomTGeo* gm)
  41   :fActiveID(activeID)
  42   ,fNSensors(0)
  43   ,fNSensInLadder(0)
  44   ,fNLadders(0)
  45   ,fR(0)
  46   ,fRMax(0)
  47   ,fRMin(0)
  48   ,fZMax(0)
  49   ,fZMin(0)
  50   ,fPhiLadMax(0)
  51   ,fPhiLadMin(0)
  52   ,fPhiOffs(0)
  53   ,fSensDZInv(0)
  54   ,fDPhiLadInv(0)
  55   ,fMaxStep(0.5)
  56   ,fSensors(0)
  57   ,fITSGeom(gm)
  58   ,fClusters(0)
  59 {
  60   // def. c-tor
  61   SetNameTitle(name,name);
  62   Build();
  63 }
  64
  65 //______________________________________________________
  66 AliITSURecoLayer::~AliITSURecoLayer()
  67 {
  68   // def. d-tor
  69   delete[] fSensors;
  70   delete[] fPhiLadMax;
  71   delete[] fPhiLadMin;
  72 }
  73
  74 //______________________________________________________
  75 void AliITSURecoLayer::Print(Option_t* opt) const
  76 {
  77   //print
  78   printf("Lr %-15s %d (act:%+d), NSens: %4d | MaxStep:%.2f ",GetName(),GetID(),GetActiveID(),GetNSensors(),fMaxStep);
  79   printf("%6.3f<R<%6.3f | %+8.3f<Z<%+8.3f dZ:%6.3f\n",fRMin,fRMax,fZMin,fZMax,fSensDZInv>0 ? 1/fSensDZInv : 0);
  80   TString opts = opt; opts.ToLower();
  81   if (opts.Contains("sn")) for (int i=0;i<GetNSensors();i++) GetSensor(i)->Print(opt);
  82 }
  83
  84 //______________________________________________________
  85 void AliITSURecoLayer::Build()
  86 {
  87   // build internal structures
  88   if (fActiveID<0) return;
  89   fNLadders = fITSGeom->GetNLadders(fActiveID);
  90   fNSensInLadder = fITSGeom->GetNDetectors(fActiveID);
  91   fNSensors = fNLadders*fNSensInLadder;
  92   fSensors = new AliITSURecoSens*[fNSensors];
  93   const AliITSsegmentation* kSegm = fITSGeom->GetSegmentation(fActiveID);
  94   //
  95   // name layer according its active id, detector type and segmentation tyoe
  96   TGeoHMatrix mmod;
  97   const TGeoHMatrix* mt2l;
  98   fRMin=fZMin=1e9;
  99   fRMax=fZMax=-1e9;
 100   double phiTF,rTF, loc[3]={0,0,0},glo[3];
 101   fNSensors = 0;
 102   fPhiLadMin = new Double_t[fNLadders];
 103   fPhiLadMax = new Double_t[fNLadders];
 104   fSensDZInv = 0;
 105   fDPhiLadInv = fNLadders/TwoPi();
 106   //
 107   for (int ild=0;ild<fNLadders;ild++) {
 108     fPhiLadMin[ild] = 1e9;
 109     fPhiLadMax[ild] = -1e9;
 110     //
 111     for (int idt=0;idt<fNSensInLadder;idt++) {
 112       AliITSURecoSens* sens = new AliITSURecoSens(fNSensors++);
 113       fSensors[ild*fNSensInLadder+idt] = sens;
 114       //
 115       double phiMin=1e9,phiMax=-1e9,zMin=1e9,zMax=-1e9;
 116       mmod = *fITSGeom->GetMatrix(fActiveID,ild,idt);
 117       for (int ix=0;ix<2;ix++) {
 118         loc[0] = (ix-0.5)*kSegm->Dx(); // +-DX/2
 119         for (int iy=0;iy<2;iy++) {
 120           loc[1] = (iy-0.5)*kSegm->Dy(); // +-DY/2
 121           for (int iz=0;iz<2;iz++) {
 122             loc[2] = (iz-0.5)*kSegm->Dz(); // +-DZ/2
 123             //
 124             mmod.LocalToMaster(loc,glo);
 125             double phi = ATan2(glo[1],glo[0]);
 126             double r   = glo[0]*glo[0] + glo[1]*glo[1];
 127             if (fRMin>r) fRMin = r;
 128             if (fRMax<r) fRMax = r;
 129             BringTo02Pi(phi);
 130             if      (phiMin>1e8) phiMin=phi;
 131             else if (!OKforPhiMin(phiMin,phi)) phiMin=phi;
 132             if      (phiMax<-1e8) phiMax=phi;
 133             else if (!OKforPhiMax(phiMax,phi)) phiMax=phi;
 134             if (glo[2]>zMax) zMax=glo[2];
 135             if (glo[2]<zMin) zMin=glo[2];
 136           }
 137         }
 138       }
 139       sens->SetBoundaries(phiMin,phiMax,zMin,zMax);
 140       mt2l = fITSGeom->GetMatrixT2L(fActiveID,ild,idt);
 141       mmod.Multiply(mt2l);
 142       loc[0]=loc[1]=loc[2]=0;
 143       mmod.LocalToMaster(loc,glo);
 144       rTF   = Sqrt(glo[0]*glo[0] + glo[1]*glo[1]);  //  tracking params (misaligned)
 145       phiTF = ATan2(glo[1],glo[0]);
 146       BringTo02Pi(phiTF);
 147       //
 148       sens->SetXTF(rTF);
 149       sens->SetPhiTF(phiTF);
 150       //
 151       if      (fPhiLadMin[ild]>1e8)  fPhiLadMin[ild] = phiMin;
 152       else if (!OKforPhiMin(fPhiLadMin[ild],phiMin)) fPhiLadMin[ild] = phiMin;
 153       if      (fPhiLadMax[ild]<-1e8) fPhiLadMax[ild] = phiMax;
 154       else if (!OKforPhiMax(fPhiLadMax[ild],phiMax)) fPhiLadMax[ild] = phiMax;
 155       if (fZMin>zMin) fZMin = zMin;
 156       if (fZMax<zMax) fZMax = zMax;
 157       //
 158       if (idt>0) fSensDZInv += zMax - GetSensor(ild,idt-1)->GetZMax(); // z interval to previous
 159     }
 160   }
 161   //
 162   fRMin = Sqrt(fRMin);
 163   fRMax = Sqrt(fRMax);
 164   fR = 0.5*(fRMin+fRMax);
 165   double dz = fNSensInLadder>0 ? fSensDZInv/(fNSensInLadder-1)/fNLadders : fZMax-fZMin;
 166   fSensDZInv = 1./dz;
 167
 168   const int kNBId[3][3] = {
 169     {AliITSURecoSens::kNghbBL,AliITSURecoSens::kNghbB,AliITSURecoSens::kNghbBR},
 170     {AliITSURecoSens::kNghbL,          -1            ,AliITSURecoSens::kNghbR },
 171     {AliITSURecoSens::kNghbTL,AliITSURecoSens::kNghbT,AliITSURecoSens::kNghbTR}
 172   };
 173
 174   // add neighbours info
 175   double zTol = 0.45*dz, phiTol = 0.45*TwoPi()/fNLadders;
 176   for (int ild=0;ild<fNLadders;ild++) {
 177     for (int idt=0;idt<fNSensInLadder;idt++) {
 178       AliITSURecoSens* sens = GetSensor(ild,idt);
 179       //
 180       for (int ils=-1;ils<=1;ils++) {
 181         int ildN = ild+ils;  // ladders of neighbouring sensors
 182         if (ildN<0) ildN = fNLadders-1; else if (ildN==fNLadders) ildN = 0;
 183         for (int ids=-1;ids<=1;ids++) {
 184           int idtN = idt+ids;
 185           if (idtN<0 || idtN==fNSensInLadder || (ids==0&&ils==0)) continue;
 186           AliITSURecoSens* sensN = GetSensor(ildN,idtN); // potential neighbor
 187           int neighbID = ildN*fNSensInLadder+idtN;
 188           //
 189           int zType = 1;  // side
 190           if (sens->GetZMin()-zTol  > sensN->GetZMax()) continue; // too large distance
 191           if (sensN->GetZMin()-zTol > sens->GetZMax() ) continue; // too large distance
 192           if      (sens->GetZMin()-zTol>sensN->GetZMin()) zType =  0;     // bottom
 193           else if (sensN->GetZMin()-zTol>sens->GetZMin()) zType =  2;     // top
 194           //
 195           int phiType = 1;
 196
 197           double phiTstMn = sensN->GetPhiMin()-phiTol;
 198           BringTo02Pi(phiTstMn);
 199           if (!OKforPhiMax(sens->GetPhiMax(),phiTstMn)) continue; // too large angle
 200           double phiTstMx = sensN->GetPhiMax()+phiTol;
 201           BringTo02Pi(phiTstMx);
 202           if (!OKforPhiMin(sens->GetPhiMin(),phiTstMx)) continue; // too large angle
 203           //
 204           phiTstMn = sensN->GetPhiMin()+phiTol;
 205           BringTo02Pi(phiTstMn);
 206           phiTstMx = sensN->GetPhiMax()-phiTol;
 207           BringTo02Pi(phiTstMx);
 208           if      (!OKforPhiMax(sens->GetPhiMax(),phiTstMx)) phiType = 0; // left
 209           else if (!OKforPhiMin(sens->GetPhiMin(),phiTstMn)) phiType = 2; // right
 210           //
 211           sens->SetNeighborID(kNBId[zType][phiType], neighbID);
 212         } // phi scan
 213       } // z scan
 214     } // sensors
 215   } // ladders
 216   //
 217 }
 218
 219 //______________________________________________________
 220 Int_t AliITSURecoLayer::FindSensors(const double* impPar, AliITSURecoSens *sensors[AliITSURecoSens::kNNeighbors])
 221 {
 222   // find sensors having intersection with track
 223   // impPar contains: lab phi of track, dphi, labZ, dz
 224   double z = impPar[2];
 225   if (z>fZMax+impPar[3]) return 0; // outside of Z coverage
 226   z -= fZMin;
 227   if (z<-impPar[3]) return 0; // outside of Z coverage
 228   int sensInLad = int(z*fSensDZInv);
 229   if      (sensInLad<0) sensInLad = 0;
 230   else if (sensInLad>=fNSensInLadder) sensInLad = fNSensInLadder-1;
 231   //
 232   double phi = impPar[0] - fPhiOffs;
 233   BringTo02Pi(phi);
 234   int ladID = int(phi*fDPhiLadInv);  // ladder id
 235   int nsens = 0;
 236   //
 237   AliITSURecoSens* sensN,*sens = GetSensor(ladID*fNSensInLadder+sensInLad);
 238   sensors[nsens++] = sens;
 239   // check neighbours
 240   double zMn=impPar[2]-impPar[3], zMx=impPar[2]+impPar[3], phiMn=impPar[0]-impPar[1], phiMx=impPar[0]+impPar[1];
 241   BringTo02Pi(phiMn);
 242   BringTo02Pi(phiMx);
 243   //
 244   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbR)); // neighbor on the right (smaller phi)
 245   if (sensN && OKforPhiMin(phiMn,sensN->GetPhiMax())) sensors[nsens++] = sensN;
 246   //
 247   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbTR)); // neighbor on the top right (smaller phi, larger Z)
 248   if (sensN && OKforPhiMin(phiMn,sensN->GetPhiMax()) && sensN->GetZMin()<zMx) sensors[nsens++] = sensN;
 249   //
 250   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbT)); // neighbor on the top (larger Z)
 251   if (sensN && sensN->GetZMin()<zMx) sensors[nsens++] = sensN;
 252   //
 253   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbTL)); // neighbor on the top left (larger Z, larger phi)
 254   if (sensN && OKforPhiMax(phiMx,sensN->GetPhiMin()) && sensN->GetZMin()<zMx) sensors[nsens++] = sensN;
 255   //
 256   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbL)); // neighbor on the left (larger phi)
 257   if (sensN && OKforPhiMax(phiMx,sensN->GetPhiMin())) sensors[nsens++] = sensN;
 258   //
 259   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbBL)); // neighbor on the bottom left (smaller Z, larger phi)
 260   if (sensN && OKforPhiMax(phiMx,sensN->GetPhiMin()) && sensN->GetZMax()>zMn) sensors[nsens++] = sensN;
 261   //
 262   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbB));  // neighbor on the bottom (smaller Z)
 263   if (sensN && sensN->GetZMax()>zMn) sensors[nsens++] = sensN;
 264   //
 265   sensN = GetSensor(sens->GetNeighborID(AliITSURecoSens::kNghbBR)); // neighbor on the bottom right (smaller Z, smaller phi)
 266   if (sensN && OKforPhiMin(phiMn,sensN->GetPhiMax()) && sensN->GetZMax()>zMn) sensors[nsens++] = sensN;
 267   //
 268   return nsens;
 269 }
 270
 271 //______________________________________________________
 272 void AliITSURecoLayer::ProcessClusters(Int_t mode)
 273 {
 274   // register in each sensor of the layer its cluster.
 275   // the clusters of the layer must be sorted per sensor
 276   int ncl = fClusters->GetEntriesFast();
 277   int curSensID = -1;
 278   for (int i=fNSensors;i--;) GetSensor(i)->SetNClusters(0);
 279   AliITSURecoSens* curSens = 0;
 280   for (int icl=0;icl<ncl;icl++) {
 281     AliITSUClusterPix* cl = (AliITSUClusterPix*) fClusters->UncheckedAt(icl);
 282     cl->GoToFrameTrk();
 283     int vID = cl->GetVolumeId();
 284     if (vID<curSensID) {AliFatal("Clusters are not sorted in increasing sensorID");}
 285     if (vID>curSensID) {
 286       if (curSens) curSens->ProcessClusters(mode);    // prepare clusters for reconstruction
 287       curSens   = GetSensor(vID - fITSGeom->GetFirstModIndex(fActiveID));
 288       curSensID = vID;
 289       curSens->SetFirstClusterId(icl);
 290     }
 291     curSens->IncNClusters();
 292   }
 293   if (curSens) curSens->ProcessClusters(mode); // last sensor was not processed yet
 294   //
 295 }
 296
 297 //______________________________________________________
 298 Bool_t AliITSURecoLayer::IsEqual(const TObject* obj) const
 299 {
 300   // check if layers are equal in R
 301   const AliITSURecoLayer* lr = (const AliITSURecoLayer*)obj;
 302   return Abs(lr->GetR()-GetR())<1e-6 ? kTRUE : kFALSE;
 303 }
 304
 305 //______________________________________________________
 306 Int_t  AliITSURecoLayer::Compare(const TObject* obj) const
 307 {
 308   // compare two layers
 309   const AliITSURecoLayer* lr = (const AliITSURecoLayer*)obj;
 310   double dr = GetR() - lr->GetR();
 311   if (Abs(dr)<1e-6) return 0;
 312   return dr>0 ? 1:-1;
 313   //
 314 }