]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDseedV1.cxx
- fix cluster ownership problems in time bin container
[u/mrichter/AliRoot.git] / TRD / AliTRDseedV1.cxx
index 8a040a011363e67f163ec5d53588ef5983d3e809..43792bc2278996ce31e336f708ab6c64eb71a164 100644 (file)
@@ -375,18 +375,19 @@ Bool_t    AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
       << "\n"; 
   }
 
-  Float_t  tquality=0., zcorr=0.;
+  Float_t  tquality;
   Double_t kroady = fReconstructor->GetRecoParam() ->GetRoad1y();
   Double_t kroadz = fPadLength * .5 + 1.;
+  
+  // initialize configuration parameters
+  Float_t zcorr = kZcorr ? fTilt * (fZProb - fZref[0]) : 0.;
+  Int_t   niter = kZcorr ? 1 : 2;
+  
   Double_t yexp, zexp;
   Int_t ncl = 0;
-
-  // looking for clusters
-  for (Int_t iter = 0; iter < 2; iter++) {
+  // start seed update
+  for (Int_t iter = 0; iter < niter; iter++) {
     ncl = 0;
-    // recalculate correction for tilt pad
-    zcorr = (kZcorr&&fN) ? fTilt * (fZProb - fZref[0]) : 0.;
-    
     for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
       if(!(layer = chamber->GetTB(iTime))) continue;
       if(!Int_t(*layer)) continue;
@@ -402,13 +403,13 @@ Bool_t    AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
           if (zexp < c->GetZ()) zexp = c->GetZ() - fPadLength*0.5;
         }
       } else zexp = fZref[0] + (kZcorr ? fZref[1] * dxlayer : 0.);
-      yexp  = fYref[0] + fYref[1]* dxlayer + (kZcorr ? fZref[1] * dxlayer : 0.) + zcorr;
+      yexp  = fYref[0] + fYref[1] * dxlayer - zcorr;
       
       // Get and register cluster
       Int_t    index = layer->SearchNearestCluster(yexp, zexp, kroady, kroadz);
       if (index < 0) continue;
       AliTRDcluster *cl = (*layer)[index];
-
+      
       fIndexes[iTime]  = layer->GetGlobalIndex(index);
       fClusters[iTime] = cl;
       fY[iTime]        = cl->GetY();
@@ -455,7 +456,6 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
       } 
       
       AliTRDseed::Update();
-      //Fit(kZcorr);
     }
     if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fDet, fN2));
     
@@ -638,6 +638,12 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt)
   const Float_t clSigma0 = 2.E-2;    //[cm]
   const Float_t clSlopeQ = -1.19E-2; //[1/cm]
 
+  // get track direction
+  Double_t y0   = fYref[0];
+  Double_t dydx = fYref[1]; 
+  Double_t z0   = fZref[0];
+  Double_t dzdx = fZref[1];
+  Double_t yt, zt;
 
   const Int_t kNtb = AliTRDtrackerV1::GetNTimeBins();
   AliTRDtrackerV1::AliTRDLeastSquare fitterY, fitterZ;
@@ -646,54 +652,69 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt)
   Double_t convert = 1./TMath::Sqrt(12.);
   
   // book cluster information
-  Double_t xc[knTimebins+1], yc[knTimebins], zc[knTimebins+1], sy[knTimebins], sz[knTimebins+1];
+  Double_t xc[knTimebins], yc[knTimebins], zc[knTimebins], sy[knTimebins], sz[knTimebins];
   Int_t zRow[knTimebins];
-  AliTRDcluster *c = 0x0;
-  Int_t nc = 0;
-  for (Int_t ic=0; ic<kNtb; ic++) {
+
+
+  fN = 0;
+  AliTRDcluster *c=0x0, **jc = &fClusters[0];
+  for (Int_t ic=0; ic<kNtb; ic++, ++jc) {
     zRow[ic] = -1;
     xc[ic]  = -1.;
     yc[ic]  = 999.;
     zc[ic]  = 999.;
     sy[ic]  = 0.;
     sz[ic]  = 0.;
-    if(!(c = fClusters[ic])) continue;
+    if(!(c = (*jc))) continue;
     if(!c->IsInChamber()) continue;
     Float_t w = 1.;
     if(c->GetNPads()>4) w = .5;
     if(c->GetNPads()>5) w = .2;
-    zRow[nc] = c->GetPadRow();
-    xc[nc]   = c->GetX() - fX0;
-    yc[nc]   = c->GetY();
-    zc[nc]   = c->GetZ();
+    zRow[fN] = c->GetPadRow();
+    xc[fN]   = fX0 - c->GetX();
+    yc[fN]   = c->GetY();
+    zc[fN]   = c->GetZ();
+
+    // extrapolated y value for the track
+    yt = y0 - xc[fN]*dydx; 
+    // extrapolated z value for the track
+    zt = z0 - xc[fN]*dzdx; 
+    // tilt correction
+    if(tilt) yc[fN] -= fTilt*(zc[fN] - zt); 
+
+    // elaborate cluster error
     Float_t qr = c->GetQ() - q0;
-    sy[nc]   = qr < 0. ? clSigma0*TMath::Exp(clSlopeQ*qr) : clSigma0;
-    sz[nc]   = fPadLength*convert;
-    fitterZ.AddPoint(&xc[nc], zc[nc], sz[nc]);
-    nc++;
+    sy[fN]   = qr < 0. ? clSigma0*TMath::Exp(clSlopeQ*qr) : clSigma0;
+
+    fitterY.AddPoint(&xc[fN], yc[fN]-yt, sy[fN]);
+
+    sz[fN]   = fPadLength*convert;
+    fitterZ.AddPoint(&xc[fN], zc[fN], sz[fN]);
+    fN++;
   }
   // to few clusters
-  if (nc < kClmin) return kFALSE; 
-  
+  if (fN < kClmin) return kFALSE; 
 
-  Int_t zN[2*35];
-  Int_t nz = AliTRDtrackerV1::Freq(nc, zRow, zN, kFALSE);
+  // fit XY plane
+  fitterY.Eval();
+  fYfit[0] = y0+fitterY.GetFunctionParameter(0);
+  fYfit[1] = dydx-fitterY.GetFunctionParameter(1);
+
+  // check par row crossing
+  Int_t zN[2*AliTRDseed::knTimebins];
+  Int_t nz = AliTRDtrackerV1::Freq(fN, zRow, zN, kFALSE);
   // more than one pad row crossing
   if(nz>2) return kFALSE; 
-  
-  // estimate reference parameter at average x
-  Double_t y0 = fYref[0];
-  Double_t dydx = fYref[1]; 
-  Double_t dzdx = fZref[1];
-  zc[nc]  = fZref[0];
+
 
   // determine z offset of the fit
+  Float_t zslope = 0.;
   Int_t nchanges = 0, nCross = 0;
   if(nz==2){ // tracklet is crossing pad row
     // Find the break time allowing one chage on pad-rows
     // with maximal number of accepted clusters
     Int_t padRef = zRow[0];
-    for (Int_t ic=1; ic<nc; ic++) {
+    for (Int_t ic=1; ic<fN; ic++) {
       if(zRow[ic] == padRef) continue;
       
       // debug
@@ -703,12 +724,12 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt)
     
       // evaluate parameters of the crossing point
       Float_t sx = (xc[ic-1] - xc[ic])*convert;
-      xc[nc] = .5 * (xc[ic-1] + xc[ic]);
-      zc[nc] = .5 * (zc[ic-1] + zc[ic]);
-      sz[nc] = TMath::Max(dzdx * sx, .01);
-      dzdx   = zc[ic-1] > zc[ic] ? 1. : -1.;
-      padRef = zRow[ic];
-      nCross = ic;
+      fCross[0] = .5 * (xc[ic-1] + xc[ic]);
+      fCross[2] = .5 * (zc[ic-1] + zc[ic]);
+      fCross[3] = TMath::Max(dzdx * sx, .01);
+      zslope    = zc[ic-1] > zc[ic] ? 1. : -1.;
+      padRef    = zRow[ic];
+      nCross    = ic;
       nchanges++;
     }
   }
@@ -716,38 +737,21 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt)
   // condition on nCross and reset nchanges TODO
 
   if(nchanges==1){
-    if(dzdx * fZref[1] < 0.){
+    if(dzdx * zslope < 0.){
       AliInfo("tracklet direction does not correspond to the track direction. TODO.");
     }
     SetBit(kRowCross, kTRUE); // mark pad row crossing
-    fCross[0] = xc[nc]; fCross[2] = zc[nc]; fCross[3] = sz[nc]; 
-    fitterZ.AddPoint(&xc[nc], zc[nc], sz[nc]);
+    fitterZ.AddPoint(&fCross[0], fCross[2], fCross[3]);
     fitterZ.Eval();
-    dzdx = fZref[1]; // we don't trust Parameter[1] ??;
-    zc[nc] = fitterZ.GetFunctionParameter(0); 
+    //zc[nc] = fitterZ.GetFunctionParameter(0); 
+    fCross[1] = fYfit[0] - fCross[0] * fYfit[1];
+    fCross[0] = fX0 - fCross[0];
   } else if(nchanges > 1){ // debug
-    AliInfo("ERROR in n changes!!!");
+    AliError("N pad row crossing > 1.");
     return kFALSE;
   }
 
-  
-  // estimate deviation from reference direction
-  dzdx *= fTilt;
-  for (Int_t ic=0; ic<nc; ic++) {
-    yc[ic] -= (y0 + xc[ic]*dydx);
-    if(tilt) yc[ic] -= fTilt*(xc[ic]*dzdx + (zc[ic] - zc[nc])); 
-    fitterY.AddPoint(&xc[ic], yc[ic], sy[ic]);
-  }
-  fitterY.Eval();
-  fYfit[0] = y0+fitterY.GetFunctionParameter(0);
-  fYfit[1] = dydx+fitterY.GetFunctionParameter(1);
-
-  if(nchanges) fCross[1] = fYfit[0] + fCross[0] * fYfit[1];
-
-//     printf("\nnz = %d\n", nz);
-//     for(int ic=0; ic<35; ic++) printf("%d row[%d]\n", ic, zRow[ic]);        
-// 
-//     for(int ic=0; ic<nz; ic++) printf("%d n[%d]\n", ic, zN[ic]);    
+  UpdateUsed();
 
   return kTRUE;
 }