]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Major update for the TRD tracking code
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 25 Feb 2009 17:22:09 +0000 (17:22 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 25 Feb 2009 17:22:09 +0000 (17:22 +0000)
1. Detach offline tracklet from the old AliTRDseed. Prepare the junction
with AliTRDtrackletBase and remove obsolete data members
2. Store track position at lower (TPC) and upper (HMPID) entrances for
monitoring systematic effects
3. Restrict chi2 definition for tracklets to the y direction to cover up
the TPC z shift (temporary)
4. Rename and extend tracking resolution task.

Alex

17 files changed:
EVE/EveDet/AliEveTRDData.cxx
HLT/TRD/AliHLTTRDTracklet.cxx
HLT/TRD/AliHLTTRDTracklet.h
TRD/AliTRDCalibraFillHisto.cxx
TRD/AliTRDseedV1.cxx
TRD/AliTRDseedV1.h
TRD/AliTRDtrackV1.cxx
TRD/AliTRDtrackV1.h
TRD/AliTRDtrackerDebug.cxx
TRD/AliTRDtrackerV1.cxx
TRD/TRDqaRecLinkDef.h
TRD/qaRec/AliTRDcalibration.cxx
TRD/qaRec/AliTRDcheckDetector.cxx
TRD/qaRec/AliTRDresolution.cxx [moved from TRD/qaRec/AliTRDtrackingResolution.cxx with 84% similarity]
TRD/qaRec/AliTRDresolution.h [moved from TRD/qaRec/AliTRDtrackingResolution.h with 58% similarity]
TRD/qaRec/run.C
TRD/qaRec/run.h

index 198ab482c7b5c8fd5016641bfbbbda67850de632..f19127267d432dbc908d5697c0a4269b8480d395 100644 (file)
@@ -352,7 +352,7 @@ AliEveTRDTracklet::AliEveTRDTracklet(AliTRDseedV1 *trklt):TEveLine()
   Float_t tilt = trklt->GetTilt();
   Float_t g[3];
   AliTRDcluster *c = 0x0;
-  for(Int_t ic=0; ic<AliTRDseed::knTimebins; ic++){
+  for(Int_t ic=0; ic<AliTRDseedV1::kNTimeBins; ic++){
     if(!(c = trklt->GetClusters(ic))) continue;
     if(!fClusters) AddElement(fClusters = new AliEveTRDClusters());
     dx = x0 - c->GetX();
index ddf3b307e77d14a32ae9019e359b6198e0719c53..3696ecb0e3e9e4723b32f3205e8e2d57a0271727 100644 (file)
@@ -8,23 +8,24 @@
 AliHLTTRDTracklet::AliHLTTRDTracklet():
   fTRDtracklet(NULL),
   fSize(sizeof(AliHLTTRDTracklet)),
-  fSigmaY(-1),
+  //fSigmaY(-1),
   fSigmaY2(-1),
   fTilt(-1),
   fPadLength(-1),
   fX0(-1),
-  fMeanz(-1),
-  fZProb(-1),
-  fN(-1),
+  fUsable(0),
+  //fMeanz(-1),
+  //fZProb(-1),
+  //fN(-1),
   fN2(-1),
   fNUsed(-1),
-  fFreq(-1),
-  fNChange(-1),
-  fMPads(-1),
+  //fFreq(-1),
+  //fNChange(-1),
+  //fMPads(-1),
   fC(-1),
-  fCC(-1),
+  //fCC(-1),
   fChi2(-1),
-  fChi2Z(-1),
+  //fChi2Z(-1),
   fDet(-1),
   fMom(-1),
   fdX(-1)
@@ -39,23 +40,24 @@ AliHLTTRDTracklet::AliHLTTRDTracklet():
 AliHLTTRDTracklet::AliHLTTRDTracklet(AliTRDseedV1 * inTracklet):
   fTRDtracklet(NULL),
   fSize(sizeof(AliHLTTRDTracklet)),
-  fSigmaY(-1),
+  //fSigmaY(-1),
   fSigmaY2(-1),
   fTilt(-1),
   fPadLength(-1),
   fX0(-1),
-  fMeanz(-1),
-  fZProb(-1),
-  fN(-1),
+  fUsable(0),
+  //fMeanz(-1),
+  //fZProb(-1),
+  //fN(-1),
   fN2(-1),
   fNUsed(-1),
-  fFreq(-1),
-  fNChange(-1),
-  fMPads(-1),
+  //fFreq(-1),
+  //fNChange(-1),
+  //fMPads(-1),
   fC(-1),
-  fCC(-1),
+  //fCC(-1),
   fChi2(-1),
-  fChi2Z(-1),
+  //fChi2Z(-1),
   fDet(-1),
   fMom(-1),
   fdX(-1)
@@ -74,7 +76,7 @@ AliHLTTRDTracklet::AliHLTTRDTracklet(AliTRDseedV1 * inTracklet):
 //============================================================================
 void AliHLTTRDTracklet::AddClusters()
 {
-  for (Int_t iTimeBin = 0; iTimeBin < knTimebins; iTimeBin++)
+  for (Int_t iTimeBin = 0; iTimeBin < AliTRDseedV1::kNTimeBins; iTimeBin++)
     {
 //       if (fClusters[iTimeBin])
 //     HLTWarning("Trying to rewrite cluster in tracklet. Not good.");
@@ -101,46 +103,36 @@ void AliHLTTRDTracklet::CopyDataMembers()
     fYref[i] = fTRDtracklet->GetYref(i);
     fZref[i] = fTRDtracklet->GetZref(i);
   }
-  fSigmaY = fTRDtracklet->GetSigmaY();
-  fSigmaY2 = fTRDtracklet->GetSigmaY2();
+  //fSigmaY = fTRDtracklet->GetSigmaY();
+  fSigmaY2 = fTRDtracklet->GetS2Y();
   fTilt = fTRDtracklet->GetTilt();
   fPadLength = fTRDtracklet->GetPadLength();
   
   fX0 = fTRDtracklet->GetX0();
-  for (Int_t i = 0; i < knTimebins; i++){
-    fX[i] = fTRDtracklet->GetX(i);
-    fY[i] = fTRDtracklet->GetY(i);
-    fZ[i] = fTRDtracklet->GetZ(i);
+  for (Int_t i = 0; i < AliTRDseedV1::kNTimeBins; i++){
+//     fX[i] = fTRDtracklet->GetX(i);
+//     fY[i] = fTRDtracklet->GetY(i);
+//     fZ[i] = fTRDtracklet->GetZ(i);
     fIndexes[i] = fTRDtracklet->GetIndexes(i);
-    fUsable[i] = fTRDtracklet->IsUsable(i);
   }
+  fUsable = fTRDtracklet->GetUsabilityMap();
 
   for (Int_t i=0; i < 2; i++){
     fYfit[i] = fTRDtracklet->GetYfit(i);
     fZfit[i] = fTRDtracklet->GetZfit(i);
-    fYfitR[i] = fTRDtracklet->GetYfitR(i);
-    fZfitR[i] = fTRDtracklet->GetZfitR(i);
     fLabels[i] = fTRDtracklet->GetLabels(i);
   }
-  
-  fMeanz   = fTRDtracklet->GetMeanz();
-  fZProb   = fTRDtracklet->GetZProb();
-  fN       = fTRDtracklet->GetNbClusters();
+  fLabels[2] = fTRDtracklet->GetLabels(2);
   fN2      = fTRDtracklet->GetN2();
   fNUsed   = fTRDtracklet->GetNUsed();
-  fFreq    = fTRDtracklet->GetFreq();
-  fNChange = fTRDtracklet->GetNChange();
-  fMPads = fTRDtracklet->GetMPads();
    
   fC = fTRDtracklet->GetC();
-  fCC = fTRDtracklet->GetCC();
   fChi2 = fTRDtracklet->GetChi2();
-  fChi2Z = fTRDtracklet->GetChi2Z();
+  //fChi2Z = fTRDtracklet->GetChi2Z();
   
   fDet = fTRDtracklet->GetDetector();
   fMom = fTRDtracklet->GetMomentum();
   fdX = fTRDtracklet->GetdX();
-  
 }
 
 /**
@@ -161,48 +153,38 @@ void AliHLTTRDTracklet::ExportTRDTracklet(AliTRDseedV1 *outTracklet)
     outTracklet->SetZref(i, fZref[i]);
   }
   
-  outTracklet->SetSigmaY(fSigmaY);
-  outTracklet->SetSigmaY2(fSigmaY2);
+  //outTracklet->SetSigmaY(fSigmaY);
+  //outTracklet->SetSigmaY2(fSigmaY2);
   outTracklet->SetTilt(fTilt);
   outTracklet->SetPadLength(fPadLength);
   outTracklet->SetX0(fX0);
 
-  for (Int_t i=0; i < knTimebins; i++){
-    outTracklet->SetX(i,fX[i]);
-    outTracklet->SetX(i,fY[i]);
-    outTracklet->SetX(i,fZ[i]);
+  for (Int_t i=0; i < AliTRDseedV1::kNTimeBins; i++){
+//     outTracklet->SetX(i,fX[i]);
+//     outTracklet->SetX(i,fY[i]);
+//     outTracklet->SetX(i,fZ[i]);
     outTracklet->SetIndexes(i, fIndexes[i]);
-    outTracklet->SetUsable(i, fUsable[i]);
-  }
-  for (Int_t i=0; i < 2; i++){
-    outTracklet->SetYfit(i,fYfit[i]);
-    outTracklet->SetYfitR(i,fYfitR[i]);
-    outTracklet->SetZfit(i,fZfit[i]);
-    outTracklet->SetZfitR(i,fZfitR[i]);
-    outTracklet->SetLabels(i,fLabels[i]);
   }
+  outTracklet->SetUsabilityMap(fUsable);
+
+//   for (Int_t i=0; i < 2; i++){
+//     outTracklet->SetYfit(i,fYfit[i]);
+//     outTracklet->SetZfit(i,fZfit[i]);
+//   }
+  outTracklet->SetLabels(fLabels);
   
-  outTracklet->SetMeanz(fMeanz);
-  outTracklet->SetZProb(fZProb);
-  outTracklet->SetN(fN);
-  outTracklet->SetN2(fN2);
-  outTracklet->SetNUsed(fNUsed);
-  outTracklet->SetFreq(fFreq);
-  outTracklet->SetNChange(fNChange);
-  outTracklet->SetMPads(fMPads);
+  //outTracklet->SetN2(fN2);
+  //outTracklet->SetNUsed(fNUsed);
   outTracklet->SetC(fC);
-  outTracklet->SetCC(fCC);
   outTracklet->SetChi2(fChi2);
-  outTracklet->SetChi2Z(fChi2Z);
 
-  for (Int_t iCluster = 0; iCluster < knTimebins; iCluster++){
+  for (Int_t iCluster = 0; iCluster < AliTRDseedV1::kNTimeBins; iCluster++){
     if (fClusters[iCluster]){
       AliTRDcluster *trdCluster = new AliTRDcluster();
       fClusters[iCluster]->ExportTRDCluster(trdCluster);
-      outTracklet->SetClusters(iCluster, trdCluster);
+      //outTracklet->SetClusters(iCluster, trdCluster);
     }
   }
-  
 }
 
 
@@ -212,7 +194,7 @@ void AliHLTTRDTracklet::ExportTRDTracklet(AliTRDseedV1 *outTracklet)
 //============================================================================
 void AliHLTTRDTracklet::InitArrays()
 {
-  for (Int_t i=0; i < knTimebins; i++){
+  for (Int_t i=0; i < AliTRDseedV1::kNTimeBins; i++){
     fClusters[i] = NULL;
   }
 
@@ -220,22 +202,22 @@ void AliHLTTRDTracklet::InitArrays()
     fYref[i] = -1;
     fZref[i] = -1;
   }
-  for (Int_t i = 0; i < knTimebins; i++){
-    fX[i] = -1;
-    fY[i] = -1;
-    fZ[i] = -1;
+  for (Int_t i = 0; i < AliTRDseedV1::kNTimeBins; i++){
+//     fX[i] = -1;
+//     fY[i] = -1;
+//     fZ[i] = -1;
     fIndexes[i] = -1;
-    fUsable[i] = -1;
   }
+  fUsable = 0;
 
   for (Int_t i=0; i < 2; i++){
     fYfit[i] = -1;
     fZfit[i] = -1;
-    fYfitR[i] = -1;
-    fZfitR[i] = -1;
+//     fYfitR[i] = -1;
+//     fZfitR[i] = -1;
     fLabels[i] = -1;
   }
-  
+  fLabels[2] = 0;
 }
 
 /**
@@ -245,20 +227,14 @@ void AliHLTTRDTracklet::InitArrays()
 void AliHLTTRDTracklet::Print(Bool_t printClusters)
 {
   //printf("--hltTracklet-- addr 0x%p(%i); fSize %i\n", this, (int)this, fSize);
-  printf("      fDet %i; dMom %f; fdX %f fN %i\n", fDet, fMom, fdX, fN);
+  printf("      fDet %i; dMom %f; fdX %f fN %i\n", fDet, fMom, fdX, fN2);
 
-  if (printClusters){
-    
-    for (Int_t iCluster = 0; iCluster < knTimebins; iCluster++){
-      printf(" [%i] ",iCluster);
-      if (fClusters[iCluster]){
-       fClusters[iCluster]->Print();
-      }
-      else 
-       printf("      NULL\n");
-    }
+  if(!printClusters) return;
+  for (Int_t iCluster = 0; iCluster < AliTRDseedV1::kNTimeBins; iCluster++){
+    printf(" [%i] ",iCluster);
+    if (fClusters[iCluster]) fClusters[iCluster]->Print();
+    else printf("      NULL\n");
   }
-  
 }
 
 /**
@@ -270,7 +246,7 @@ void AliHLTTRDTracklet::ReadClustersFromMemory(void *input)
   AliHLTUInt8_t *iterPtr = (AliHLTUInt8_t*) input;
   AliHLTTRDCluster* hltCluster = NULL;
   
-  for (Int_t iCluster = 0; iCluster < knTimebins; iCluster++){
+  for (Int_t iCluster = 0; iCluster < AliTRDseedV1::kNTimeBins; iCluster++){
     // if we had something in the fClusters[iCluster] before copying,
     // then this entry in the array should not be empty. Fill it.
     if (fClusters[iCluster]){
index 133c5a03dac06ecbf4bcfc595ab8aba37b94c48f..07177eb4900b6f71d199b89ad474306537424346 100644 (file)
@@ -29,44 +29,40 @@ class AliHLTTRDTracklet
   AliTRDseedV1* fTRDtracklet;
   AliHLTUInt32_t fSize; // Size of the tracklet with clusters in the memory
   
-
-  /* ======= From AliTRDseed ======== */
-  enum { knTimebins = 35 };
-
   /* Defenitely need */
-  AliHLTTRDCluster *fClusters[knTimebins]; // Clusters
+  AliHLTTRDCluster *fClusters[AliTRDseedV1::kNTimeBins]; // Clusters
   Float_t        fYref[2];              //  Reference y
   Float_t        fZref[2];              //  Reference z
-  Float_t        fSigmaY;               //  "Robust" sigma in Y - constant fit
+  //Float_t        fSigmaY;               //  "Robust" sigma in Y - constant fit
   Float_t        fSigmaY2;              //  "Robust" sigma in Y - line fit
 
   /* Probably need */
   Float_t        fTilt;                 //  Tilting angle
   Float_t        fPadLength;            //  Pad length
   Float_t        fX0;                   //  X0 position
-  Float_t        fX[knTimebins];        //! X position
-  Float_t        fY[knTimebins];        //! Y position
-  Float_t        fZ[knTimebins];        //! Z position
-  Int_t          fIndexes[knTimebins];  //! Indexes
-  Bool_t         fUsable[knTimebins];   //  Indication  - usable cluster
+//  Float_t        fX[knTimebins];        //! X position
+//  Float_t        fY[knTimebins];        //! Y position
+//  Float_t        fZ[knTimebins];        //! Z position
+  Int_t          fIndexes[AliTRDseedV1::kNTimeBins];  //! Indexes
+  Long_t         fUsable;                //  Indication  - usable cluster
   Float_t        fYfit[2];              //  Y fit position +derivation
-  Float_t        fYfitR[2];             //  Y fit position +derivation
+  //Float_t        fYfitR[2];             //  Y fit position +derivation
   Float_t        fZfit[2];              //  Z fit position
-  Float_t        fZfitR[2];             //  Z fit position
-  Float_t        fMeanz;                //  Mean vaue of z
-  Float_t        fZProb;                //  Max probbable z
-  Int_t          fLabels[2];            //  Labels
-  Int_t          fN;                    //  Number of associated clusters
+  //Float_t        fZfitR[2];             //  Z fit position
+  //Float_t        fMeanz;                //  Mean vaue of z
+  //Float_t        fZProb;                //  Max probbable z
+  Int_t          fLabels[3];            //  Labels
+  //Int_t          fN;                    //  Number of associated clusters
   Int_t          fN2;                   //  Number of not crossed
   Int_t          fNUsed;                //  Number of used clusters
-  Int_t          fFreq;                 //  Frequency
-  Int_t          fNChange;              //  Change z counter
-  Float_t        fMPads;                //  Mean number of pads per cluster
+  //Int_t          fFreq;                 //  Frequency
+  //Int_t          fNChange;              //  Change z counter
+  //Float_t        fMPads;                //  Mean number of pads per cluster
 
   Float_t        fC;                    //  Curvature
-  Float_t        fCC;                   //  Curvature with constrain
+  //Float_t        fCC;                   //  Curvature with constrain
   Float_t        fChi2;                 //  Global chi2
-  Float_t        fChi2Z;                //  Global chi2
+  //Float_t        fChi2Z;                //  Global chi2
 
   /* ======= From AliTRDseedV1 ======== */
 
index f78b97f547f927a4e2784c9fed553ae7002984f3..140c65277b2329fbda56f1f92165bd0e3896a082 100644 (file)
@@ -691,7 +691,7 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(AliTRDtrackV1 *t)
     // loop over the clusters
     ////////////////////////////
     Int_t nbclusters = 0;
-    for(int jc=0; jc<AliTRDseed::kNtb; jc++){
+    for(int jc=0; jc<AliTRDseedV1::kNtb; jc++){
       if(!(cl = tracklet->GetClusters(jc))) continue;
       nbclusters++;
       
@@ -927,7 +927,7 @@ Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *track
   ////////////////////////////
   Int_t  nbli = 0;
   AliTRDcluster *cl                   = 0x0;
-  for(int ic=0; ic<AliTRDseed::kNtb; ic++){
+  for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
     if(!(cl = tracklet->GetClusters(ic))) continue;
     if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue;
     
@@ -942,16 +942,13 @@ Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *track
     fLinearFitterTracklet->AddPoint(&timeis,ycluster,1);
     nbli++;  
 
+    //////////////////////////////
+    // Check no shared clusters
+    //////////////////////////////
+    for(int ic=AliTRDseedV1::kNtb; ic<AliTRDseedV1::kNTimeBins; ic++){
+      if((cl = tracklet->GetClusters(ic)))  crossrow = 1;
+    }
   }
-
-  //////////////////////////////
-  // Check no shared clusters
-  //////////////////////////////
-
-  for(int ic=AliTRDseed::kNtb; ic<AliTRDseed::knTimebins; ic++){
-    if((cl = tracklet->GetClusters(ic)))  crossrow = 1;
-  }
-
   
   ////////////////////////////////////
   // Do the straight line fit now
@@ -1397,11 +1394,9 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet,
   //////////////////////
 
   Int_t nb3pc    = 0;              // number of three pads clusters used for fit 
-  Double_t *padPositions;          // to see the difference between the fit and the 3 pad clusters position
-  padPositions = new Double_t[AliTRDseed::kNtb];
-  for(Int_t k = 0; k < AliTRDseed::kNtb; k++){
-    padPositions[k] = 0.0;
-  } 
+  // to see the difference between the fit and the 3 pad clusters position
+  Double_t padPositions[AliTRDseedV1::kNtb];
+  memset(padPositions, 0, AliTRDseedV1::kNtb*sizeof(Double_t)); 
   fLinearFitterTracklet->ClearPoints();  
   
   //printf("loop clusters \n");
@@ -1409,13 +1404,10 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet,
   // loop over the clusters
   ////////////////////////////
   AliTRDcluster *cl                   = 0x0;
-  for(int ic=0; ic<AliTRDseed::kNtb; ic++){
+  for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
     // reject shared clusters on pad row
-    if(((ic+AliTRDseed::kNtb) < AliTRDseed::knTimebins) && (cl = tracklet->GetClusters(ic+AliTRDseed::kNtb))) continue;
-    //    
-    if(!(cl = tracklet->GetClusters(ic))) continue;
+    if(((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNTimeBins) && (cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
    
-    
     Double_t     time  = cl->GetLocalTimeBin();
     if((time<=7) || (time>=21)) continue;
     Short_t  *signals  = cl->GetSignals(); 
@@ -1461,9 +1453,7 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet,
   // fit with a straight line
   /////////////////////////////
   if(nb3pc < 3){ 
-    delete [] padPositions;
     fLinearFitterTracklet->ClearPoints();  
-    delete [] padPositions;
     return kFALSE;
   }
   fLinearFitterTracklet->Eval();
@@ -1479,9 +1469,9 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet,
   ////////////////////////////////////////////////
   // Fill the PRF: Second loop over clusters
   //////////////////////////////////////////////
-  for(int ic=0; ic<AliTRDseed::kNtb; ic++){
+  for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
     // reject shared clusters on pad row
-    if(((ic+AliTRDseed::kNtb) < AliTRDseed::knTimebins) && (cl = tracklet->GetClusters(ic+AliTRDseed::kNtb))) continue;
+    if(((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNTimeBins) && (cl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) continue;
     //
     if(!(cl = tracklet->GetClusters(ic))) continue;
 
@@ -1679,9 +1669,7 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet,
   } // second loop over clusters
 
 
-  delete [] padPositions;
   return kTRUE;
-  
 }
 ///////////////////////////////////////////////////////////////////////////////////////
 // Pad row col stuff: see if masked or not
index caa076758c81bbe04fd693e9a65d19a1b5c29626..f028c49e9ac75c6827db911577da6a9b73b24112 100644 (file)
@@ -55,23 +55,46 @@ ClassImp(AliTRDseedV1)
 
 //____________________________________________________________________
 AliTRDseedV1::AliTRDseedV1(Int_t det) 
-  :AliTRDseed()
+  :TObject()
   ,fReconstructor(0x0)
   ,fClusterIter(0x0)
+  ,fExB(0.)
+  ,fVD(0.)
+  ,fT0(0.)
+  ,fS2PRF(0.)
+  ,fDiffL(0.)
+  ,fDiffT(0.)
   ,fClusterIdx(0)
+  ,fUsable(0)
+  ,fN2(0)
+  ,fNUsed(0)
   ,fDet(det)
+  ,fTilt(0.)
+  ,fPadLength(0.)
   ,fMom(0.)
-  ,fSnp(0.)
-  ,fTgl(0.)
   ,fdX(0.)
-  ,fXref(0.)
-  ,fExB(0.)
+  ,fX0(0.)
+  ,fX(0.)
+  ,fY(0.)
+  ,fZ(0.)
+  ,fS2Y(0.)
+  ,fS2Z(0.)
+  ,fC(0.)
+  ,fChi2(0.)
 {
   //
   // Constructor
   //
-  for(int islice=0; islice < knSlices; islice++) fdEdx[islice] = 0.;
+  for(Int_t ic=kNTimeBins; ic--;) fIndexes[ic] = -1;
+  memset(fClusters, 0, kNTimeBins*sizeof(AliTRDcluster*));
+  fYref[0] = 0.; fYref[1] = 0.; 
+  fZref[0] = 0.; fZref[1] = 0.; 
+  fYfit[0] = 0.; fYfit[1] = 0.; 
+  fZfit[0] = 0.; fZfit[1] = 0.; 
+  memset(fdEdx, 0, kNSlices*sizeof(Float_t)); 
   for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec]  = -1.;
+  fLabels[0]=-1; fLabels[1]=-1; // most freq MC labels
+  fLabels[2]=0;  // number of different labels for tracklet
   fRefCov[0] = 1.; fRefCov[1] = 0.; fRefCov[2] = 1.;
   // covariance matrix [diagonal]
   // default sy = 200um and sz = 2.3 cm 
@@ -80,28 +103,40 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
 
 //____________________________________________________________________
 AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &ref)
-  :AliTRDseed((AliTRDseed&)ref)
-  ,fReconstructor(ref.fReconstructor)
+  :TObject((TObject&)ref)
+  ,fReconstructor(0x0)
   ,fClusterIter(0x0)
+  ,fExB(0.)
+  ,fVD(0.)
+  ,fT0(0.)
+  ,fS2PRF(0.)
+  ,fDiffL(0.)
+  ,fDiffT(0.)
   ,fClusterIdx(0)
-  ,fDet(ref.fDet)
-  ,fMom(ref.fMom)
-  ,fSnp(ref.fSnp)
-  ,fTgl(ref.fTgl)
-  ,fdX(ref.fdX)
-  ,fXref(ref.fXref)
-  ,fExB(ref.fExB)
+  ,fUsable(0)
+  ,fN2(0)
+  ,fNUsed(0)
+  ,fDet(-1)
+  ,fTilt(0.)
+  ,fPadLength(0.)
+  ,fMom(0.)
+  ,fdX(0.)
+  ,fX0(0.)
+  ,fX(0.)
+  ,fY(0.)
+  ,fZ(0.)
+  ,fS2Y(0.)
+  ,fS2Z(0.)
+  ,fC(0.)
+  ,fChi2(0.)
 {
   //
   // Copy Constructor performing a deep copy
   //
-
-  //printf("AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &)\n");
+  if(this != &ref){
+    ref.Copy(*this);
+  }
   SetBit(kOwner, kFALSE);
-  for(int islice=0; islice < knSlices; islice++) fdEdx[islice] = ref.fdEdx[islice];
-  for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec] = ref.fProb[ispec];
-  memcpy(fRefCov, ref.fRefCov, 3*sizeof(Double_t));
-  memcpy(fCov, ref.fCov, 3*sizeof(Double_t));
 }
 
 
@@ -118,7 +153,6 @@ AliTRDseedV1& AliTRDseedV1::operator=(const AliTRDseedV1 &ref)
   SetBit(kOwner, kFALSE);
 
   return *this;
-
 }
 
 //____________________________________________________________________
@@ -130,13 +164,14 @@ AliTRDseedV1::~AliTRDseedV1()
 
   //printf("I-AliTRDseedV1::~AliTRDseedV1() : Owner[%s]\n", IsOwner()?"YES":"NO");
 
-  if(IsOwner()) 
-    for(int itb=0; itb<knTimebins; itb++){
+  if(IsOwner()) {
+    for(int itb=0; itb<kNTimeBins; itb++){
       if(!fClusters[itb]) continue; 
       //AliInfo(Form("deleting c %p @ %d", fClusters[itb], itb));
       delete fClusters[itb];
       fClusters[itb] = 0x0;
     }
+  }
 }
 
 //____________________________________________________________________
@@ -149,23 +184,45 @@ void AliTRDseedV1::Copy(TObject &ref) const
   //AliInfo("");
   AliTRDseedV1 &target = (AliTRDseedV1 &)ref; 
 
+  target.fReconstructor = fReconstructor;
   target.fClusterIter   = 0x0;
+  target.fExB           = fExB;
+  target.fVD            = fVD;
+  target.fT0            = fT0;
+  target.fS2PRF         = fS2PRF;
+  target.fDiffL         = fDiffL;
+  target.fDiffT         = fDiffT;
   target.fClusterIdx    = 0;
+  target.fUsable        = fUsable;
+  target.fN2            = fN2;
+  target.fNUsed         = fNUsed;
   target.fDet           = fDet;
+  target.fTilt          = fTilt;
+  target.fPadLength     = fPadLength;
   target.fMom           = fMom;
-  target.fSnp           = fSnp;
-  target.fTgl           = fTgl;
   target.fdX            = fdX;
-  target.fXref          = fXref;
-  target.fExB           = fExB;
-  target.fReconstructor = fReconstructor;
+  target.fX0            = fX0;
+  target.fX             = fX;
+  target.fY             = fY;
+  target.fZ             = fZ;
+  target.fS2Y           = fS2Y;
+  target.fS2Z           = fS2Z;
+  target.fC             = fC;
+  target.fChi2          = fChi2;
   
-  for(int islice=0; islice < knSlices; islice++) target.fdEdx[islice] = fdEdx[islice];
-  for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) target.fProb[ispec] = fProb[ispec];
-  memcpy(target.fRefCov, fRefCov, 3*sizeof(Double_t));
-  memcpy(target.fCov, fCov, 3*sizeof(Double_t));
+  memcpy(target.fIndexes, fIndexes, kNTimeBins*sizeof(Int_t));
+  memcpy(target.fClusters, fClusters, kNTimeBins*sizeof(AliTRDcluster*));
+  target.fYref[0] = fYref[0]; target.fYref[1] = fYref[1]; 
+  target.fZref[0] = fZref[0]; target.fZref[1] = fZref[1]; 
+  target.fYfit[0] = fYfit[0]; target.fYfit[1] = fYfit[1]; 
+  target.fZfit[0] = fZfit[0]; target.fZfit[1] = fZfit[1]; 
+  memcpy(target.fdEdx, fdEdx, kNSlices*sizeof(Float_t)); 
+  memcpy(target.fProb, fProb, AliPID::kSPECIES*sizeof(Float_t)); 
+  memcpy(target.fLabels, fLabels, 3*sizeof(Int_t)); 
+  memcpy(target.fRefCov, fRefCov, 3*sizeof(Double_t)); 
+  memcpy(target.fCov, fCov, 3*sizeof(Double_t)); 
   
-  AliTRDseed::Copy(target);
+  TObject::Copy(ref);
 }
 
 
@@ -193,14 +250,46 @@ Bool_t AliTRDseedV1::Init(AliTRDtrackV1 *track)
 }
 
 
+//_____________________________________________________________________________
+void AliTRDseedV1::Reset()
+{
+  //
+  // Reset seed
+  //
+  fExB=0.;fVD=0.;fT0=0.;fS2PRF=0.;
+  fDiffL=0.;fDiffT=0.;
+  fClusterIdx=0;fUsable=0;
+  fN2=0;fNUsed=0;
+  fDet=-1;fTilt=0.;fPadLength=0.;
+  fMom=0.;
+  fdX=0.;fX0=0.; fX=0.; fY=0.; fZ=0.;
+  fS2Y=0.; fS2Z=0.;
+  fC=0.; fChi2 = 0.;
+
+  for(Int_t ic=kNTimeBins; ic--;) fIndexes[ic] = -1;
+  memset(fClusters, 0, kNTimeBins*sizeof(AliTRDcluster*));
+  fYref[0] = 0.; fYref[1] = 0.; 
+  fZref[0] = 0.; fZref[1] = 0.; 
+  fYfit[0] = 0.; fYfit[1] = 0.; 
+  fZfit[0] = 0.; fZfit[1] = 0.; 
+  memset(fdEdx, 0, kNSlices*sizeof(Float_t)); 
+  for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec]  = -1.;
+  fLabels[0]=-1; fLabels[1]=-1; // most freq MC labels
+  fLabels[2]=0;  // number of different labels for tracklet
+  fRefCov[0] = 1.; fRefCov[1] = 0.; fRefCov[2] = 1.;
+  // covariance matrix [diagonal]
+  // default sy = 200um and sz = 2.3 cm 
+  fCov[0] = 4.e-4; fCov[1] = 0.; fCov[2] = 5.3; 
+}
+
 //____________________________________________________________________
 void AliTRDseedV1::UpDate(const AliTRDtrackV1 *trk)
 { 
   // update tracklet reference position from the TRD track
   // Funny name to avoid the clash with the function AliTRDseed::Update() (has to be made obselete)
 
-  fSnp = trk->GetSnp();
-  fTgl = trk->GetTgl();
+  Double_t fSnp = trk->GetSnp();
+  Double_t fTgl = trk->GetTgl();
   fMom = trk->GetP();
   fYref[1] = fSnp/(1. - fSnp*fSnp);
   fZref[1] = fTgl;
@@ -211,6 +300,35 @@ void AliTRDseedV1::UpDate(const AliTRDtrackV1 *trk)
   fZref[0] = trk->GetZ() - dx*fZref[1];
 }
 
+//_____________________________________________________________________________
+void AliTRDseedV1::UpdateUsed()
+{
+  //
+  // Update used seed
+  //
+
+  fNUsed = 0;
+  for (Int_t i = kNTimeBins; i--; ) {
+    if (!fClusters[i]) continue;
+    if(!TESTBIT(fUsable, i)) continue;
+    if((fClusters[i]->IsUsed())) fNUsed++;
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDseedV1::UseClusters()
+{
+  //
+  // Use clusters
+  //
+  AliTRDcluster **c = &fClusters[0];
+  for (Int_t ic=kNTimeBins; ic--; c++) {
+    if(!(*c)) continue;
+    if(!((*c)->IsUsed())) (*c)->Use();
+  }
+}
+
+
 //____________________________________________________________________
 void AliTRDseedV1::CookdEdx(Int_t nslices)
 {
@@ -233,11 +351,10 @@ void AliTRDseedV1::CookdEdx(Int_t nslices)
 // 3. cluster size
 //
 
-  Int_t nclusters[knSlices];
-  for(int i=0; i<knSlices; i++){ 
-    fdEdx[i]     = 0.;
-    nclusters[i] = 0;
-  }
+  Int_t nclusters[kNSlices]; 
+  memset(nclusters, 0, kNSlices*sizeof(Int_t));
+  memset(fdEdx, 0, kNSlices*sizeof(Float_t));
+
   const Double_t kDriftLength = (.5 * AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
 
   AliTRDcluster *c = 0x0;
@@ -274,6 +391,32 @@ void AliTRDseedV1::CookdEdx(Int_t nslices)
   }
 }
 
+//_____________________________________________________________________________
+void AliTRDseedV1::CookLabels()
+{
+  //
+  // Cook 2 labels for seed
+  //
+
+  Int_t labels[200];
+  Int_t out[200];
+  Int_t nlab = 0;
+  for (Int_t i = 0; i < kNTimeBins; i++) {
+    if (!fClusters[i]) continue;
+    for (Int_t ilab = 0; ilab < 3; ilab++) {
+      if (fClusters[i]->GetLabel(ilab) >= 0) {
+        labels[nlab] = fClusters[i]->GetLabel(ilab);
+        nlab++;
+      }
+    }
+  }
+
+  fLabels[2] = AliTRDtrackerV1::Freq(nlab,labels,out,kTRUE);
+  fLabels[0] = out[0];
+  if ((fLabels[2]  > 1) && (out[3] > 1)) fLabels[1] = out[2];
+}
+
+
 //____________________________________________________________________
 void AliTRDseedV1::GetClusterXY(const AliTRDcluster *c, Double_t &x, Double_t &y)
 {
@@ -331,7 +474,7 @@ Float_t AliTRDseedV1::GetdQdl(Int_t ic) const
 }
 
 //____________________________________________________________________
-Double_t* AliTRDseedV1::GetProbability()
+Float_t* AliTRDseedV1::GetProbability()
 {      
 // Fill probability array for tracklet from the DB.
 //
@@ -385,12 +528,12 @@ Float_t AliTRDseedV1::GetQuality(Bool_t kZcorr) const
   // Returns a quality measurement of the current seed
   //
 
-  Float_t zcorr = kZcorr ? fTilt * (fZProb - fZref[0]) : 0.;
+  Float_t zcorr = kZcorr ? fTilt * (fZfit[0] - fZref[0]) : 0.;
   return 
       .5 * TMath::Abs(18.0 - fN2)
     + 10.* TMath::Abs(fYfit[1] - fYref[1])
     + 5. * TMath::Abs(fYfit[0] - fYref[0] + zcorr)
-    + 2. * TMath::Abs(fMeanz - fZref[0]) / fPadLength;
+    + 2. * TMath::Abs(fZfit[0] - fZref[0]) / fPadLength;
 }
 
 //____________________________________________________________________
@@ -461,14 +604,20 @@ void AliTRDseedV1::GetCovAt(Double_t x, Double_t *cov) const
 
 
 //____________________________________________________________________
-void AliTRDseedV1::SetExB()
+void AliTRDseedV1::Calibrate()
 {
-// Retrive the tg(a_L) from OCDB. The following information are used
+// Retrieve calibration and position parameters from OCDB. 
+// The following information are used
 //  - detector index
-//  - column and row position of first attached cluster. 
-// 
-// If no clusters are attached to the tracklet a random central chamber 
-// position (c=70, r=7) will be used to retrieve the drift velocity.
+//  - column and row position of first attached cluster. If no clusters are attached 
+// to the tracklet a random central chamber position (c=70, r=7) will be used.
+//
+// The following information is cached in the tracklet
+//   t0 (trigger delay)
+//   drift velocity
+//   PRF width
+//   omega*tau = tg(a_L)
+//   diffusion coefficients (longitudinal and transversal)
 //
 // Author :
 // Alex Bercuci <A.Bercuci@gsi.de> 
@@ -481,23 +630,30 @@ void AliTRDseedV1::SetExB()
     return;
   }
 
-  AliTRDcalibDB *fCalib = AliTRDcalibDB::Instance();
-  AliTRDCalROC  *fCalVdriftROC = fCalib->GetVdriftROC(fDet);
-  const AliTRDCalDet  *fCalVdriftDet = fCalib->GetVdriftDet();
+  AliTRDcalibDB *calib = AliTRDcalibDB::Instance();
+  AliTRDCalROC  *vdROC = calib->GetVdriftROC(fDet),
+                *t0ROC = calib->GetT0ROC(fDet);;
+  const AliTRDCalDet *vdDet = calib->GetVdriftDet();
+  const AliTRDCalDet *t0Det = calib->GetT0Det();
 
   Int_t col = 70, row = 7;
   AliTRDcluster **c = &fClusters[0];
-  if(fN){ 
+  if(fN2){ 
     Int_t ic = 0;
-    while (ic<AliTRDseed::knTimebins && !(*c)){ic++; c++;} 
+    while (ic<kNTimeBins && !(*c)){ic++; c++;} 
     if(*c){
       col = (*c)->GetPadCol();
       row = (*c)->GetPadRow();
     }
   }
 
-  Double_t vd = fCalVdriftDet->GetValue(fDet) * fCalVdriftROC->GetValue(col, row);
-  fExB   = AliTRDCommonParam::Instance()->GetOmegaTau(vd);
+  fT0    = t0Det->GetValue(fDet) + t0ROC->GetValue(col,row);
+  fVD    = vdDet->GetValue(fDet) * vdROC->GetValue(col, row);
+  fS2PRF = calib->GetPRFWidth(fDet, col, row); fS2PRF *= fS2PRF;
+  fExB   = AliTRDCommonParam::Instance()->GetOmegaTau(fVD);
+  AliTRDCommonParam::Instance()->GetDiffCoeff(fDiffL,
+  fDiffT, fVD);
+  SetBit(kCalib, kTRUE);
 }
 
 //____________________________________________________________________
@@ -506,7 +662,7 @@ void AliTRDseedV1::SetOwner()
   //AliInfo(Form("own [%s] fOwner[%s]", own?"YES":"NO", fOwner?"YES":"NO"));
   
   if(TestBit(kOwner)) return;
-  for(int ic=0; ic<knTimebins; ic++){
+  for(int ic=0; ic<kNTimeBins; ic++){
     if(!fClusters[ic]) continue;
     fClusters[ic] = new AliTRDcluster(*fClusters[ic]);
   }
@@ -545,7 +701,7 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
   Double_t kroadz = fPadLength * .5 + 1.;
   
   // initialize configuration parameters
-  Float_t zcorr = kZcorr ? fTilt * (fZProb - fZref[0]) : 0.;
+  Float_t zcorr = kZcorr ? fTilt * (fZfit[0] - fZref[0]) : 0.;
   Int_t   niter = kZcorr ? 1 : 2;
   
   Double_t yexp, zexp;
@@ -577,15 +733,15 @@ Bool_t    AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
       
       fIndexes[iTime]  = layer->GetGlobalIndex(index);
       fClusters[iTime] = cl;
-      fY[iTime]        = cl->GetY();
-      fZ[iTime]        = cl->GetZ();
+//       fY[iTime]        = cl->GetY();
+//       fZ[iTime]        = cl->GetZ();
       ncl++;
     }
     if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d ncl [%d] = %d", iter, fDet, ncl));
     
     if(ncl>1){ 
       // calculate length of the time bin (calibration aware)
-      Int_t irp = 0; Float_t x[2]; Int_t tb[2];
+      Int_t irp = 0; Float_t x[2]={0., 0.}; Int_t tb[2] = {0,0};
       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
         if(!fClusters[iTime]) continue;
         x[irp]  = fClusters[iTime]->GetX();
@@ -593,8 +749,9 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
         irp++;
         if(irp==2) break;
       } 
-      fdX = (x[1] - x[0]) / (tb[0] - tb[1]);
-  
+      Int_t dtb = tb[1] - tb[0];
+      fdX = dtb ? (x[0] - x[1]) / dtb : 0.15;
+
       // update X0 from the clusters (calibration/alignment aware)
       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
         if(!(layer = chamber->GetTB(iTime))) continue;
@@ -615,12 +772,12 @@ Bool_t    AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
       // TODO
       
       // update x reference positions (calibration/alignment aware)
-      for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
-        if(!fClusters[iTime]) continue;
-        fX[iTime] = fX0 - fClusters[iTime]->GetX();
-      } 
+//       for (Int_t iTime = 0; iTime < AliTRDtrackerV1::GetNTimeBins(); iTime++) {
+//         if(!fClusters[iTime]) continue;
+//         fX[iTime] = fX0 - fClusters[iTime]->GetX();
+//       } 
       
-      AliTRDseed::Update();
+      FitMI();
     }
     if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=7) AliInfo(Form("iter = %d nclFit [%d] = %d", iter, fDet, fN2));
     
@@ -635,8 +792,8 @@ Bool_t      AliTRDseedV1::AttachClustersIter(AliTRDtrackingChamber *chamber, Float_t
 
   if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker)>=1) CookLabels();
 
-  // set ExB angle
-  SetExB();
+  // load calibration params
+  Calibrate();
   UpdateUsed();
   return kTRUE;        
 }
@@ -678,12 +835,12 @@ Bool_t    AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
 
   // working variables
   const Int_t kNrows = 16;
-  AliTRDcluster *clst[kNrows][knTimebins];
+  AliTRDcluster *clst[kNrows][kNTimeBins];
   Double_t cond[4], dx, dy, yt, zt,
-    yres[kNrows][knTimebins];
-  Int_t idxs[kNrows][knTimebins], ncl[kNrows], ncls = 0;
+    yres[kNrows][kNTimeBins];
+  Int_t idxs[kNrows][kNTimeBins], ncl[kNrows], ncls = 0;
   memset(ncl, 0, kNrows*sizeof(Int_t));
-  memset(clst, 0, kNrows*knTimebins*sizeof(AliTRDcluster*));
+  memset(clst, 0, kNrows*kNTimeBins*sizeof(AliTRDcluster*));
 
   // Do cluster projection
   AliTRDcluster *c = 0x0;
@@ -719,8 +876,8 @@ Bool_t      AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
       yres[r][ncl[r]] = dy;
       ncl[r]++; ncls++;
 
-      if(ncl[r] >= knTimebins) {
-        AliWarning(Form("Cluster candidates reached limit %d. Some may be lost.", knTimebins));
+      if(ncl[r] >= kNTimeBins) {
+        AliWarning(Form("Cluster candidates reached limit %d. Some may be lost.", kNTimeBins));
         kBUFFER = kTRUE;
         break;
       }
@@ -796,8 +953,8 @@ Bool_t      AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
       if(!(c = clst[jr][ic])) continue;
       Int_t it = c->GetPadTime();
       // TODO proper indexing of clusters !!
-      fIndexes[it+35*ir]  = chamber->GetTB(it)->GetGlobalIndex(idxs[jr][ic]);
-      fClusters[it+35*ir] = c;
+      fIndexes[it+kNtb*ir]  = chamber->GetTB(it)->GetGlobalIndex(idxs[jr][ic]);
+      fClusters[it+kNtb*ir] = c;
   
       //printf("\tid[%2d] it[%d] idx[%d]\n", ic, it, fIndexes[it]);
   
@@ -816,7 +973,7 @@ Bool_t      AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
   fNUsed = 0;
   for (Int_t it = 0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
     if(fClusters[it] && fClusters[it]->IsUsed()) fNUsed++;
-    if(fClusters[it+35] && fClusters[it+35]->IsUsed()) fNUsed++;
+    if(fClusters[it+kNtb] && fClusters[it+kNtb]->IsUsed()) fNUsed++;
   }
   if (fN2-fNUsed < kClmin){
     //AliWarning(Form("Too many clusters already in use %d (from %d).", fNUsed, fN2));
@@ -824,18 +981,18 @@ Bool_t    AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *chamber, Bool_t tilt)
     return kFALSE;
   }
 
-  // set the Lorentz angle for this tracklet  
-  SetExB();
+  // Load calibration parameters for this tracklet  
+  Calibrate();
 
   // calculate dx for time bins in the drift region (calibration aware)
-  Int_t irp = 0; Float_t x[2]={0., 0.}; Int_t tb[2] = {0, 0};
+  Int_t irp = 0; Float_t x[2] = {0.,0.}; Int_t tb[2]={0,0};
   for (Int_t it = t0; it < AliTRDtrackerV1::GetNTimeBins(); it++) {
     if(!fClusters[it]) continue;
     x[irp]  = fClusters[it]->GetX();
     tb[irp] = it;
     irp++;
     if(irp==2) break;
-  }
+  }  
   Int_t dtb = tb[1] - tb[0];
   fdX = dtb ? (x[0] - x[1]) / dtb : 0.15;
 
@@ -877,18 +1034,19 @@ void AliTRDseedV1::Bootstrap(const AliTRDReconstructor *rec)
   AliTRDpadPlane *pp = g.GetPadPlane(fDet);
   fTilt      = TMath::Tan(TMath::DegToRad()*pp->GetTiltingAngle());
   fPadLength = pp->GetLengthIPad();
-  fSnp = fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);
-  fTgl = fZref[1];
-  fN = 0; fN2 = 0; fMPads = 0.;
+  //fSnp = fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);
+  //fTgl = fZref[1];
+  fN2 = 0;// fMPads = 0.;
   AliTRDcluster **cit = &fClusters[0];
-  for(Int_t ic = knTimebins; ic--; cit++){
+  for(Int_t ic = kNTimeBins; ic--; cit++){
     if(!(*cit)) return;
-    fN++; fN2++;
-    fX[ic] = (*cit)->GetX() - fX0;
+    fN2++;
+/*    fX[ic] = (*cit)->GetX() - fX0;
     fY[ic] = (*cit)->GetY();
-    fZ[ic] = (*cit)->GetZ();
+    fZ[ic] = (*cit)->GetZ();*/
   }
-  Update(); // Fit();
+  //Update(); // 
+  Fit();
   CookLabels();
   GetProbability();
 }
@@ -911,6 +1069,11 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
   // 3. Do a Least Square Fit to the data
   //
 
+  if(!IsCalibrated()){
+    AliWarning("Tracklet fit failed. Call Calibrate().");
+    return kFALSE;
+  }
+
   const Int_t kClmin = 8;
 
 
@@ -952,17 +1115,20 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
   Double_t yt, zt;
 
   const Int_t kNtb = AliTRDtrackerV1::GetNTimeBins();
+  // calculation of tg^2(phi - a_L) and tg^2(a_L)
+  Double_t tgg = (dydx-fExB)/(1.+dydx*fExB); tgg *= tgg;
+  //Double_t exb2= fExB*fExB;
+
   //AliTRDtrackerV1::AliTRDLeastSquare fitterZ;
   TLinearFitter  fitterY(1, "pol1");
   // convertion factor from square to gauss distribution for sigma
   //Double_t convert = 1./TMath::Sqrt(12.);
   
   // book cluster information
-  Double_t q, xc[knTimebins], yc[knTimebins], zc[knTimebins], sy[knTimebins]/*, sz[knTimebins]*/;
-//   Int_t zRow[knTimebins];
-  
+  Double_t qc[kNTimeBins], xc[kNTimeBins], yc[kNTimeBins], zc[kNTimeBins], sy[kNTimeBins];
+
   Int_t ily = AliTRDgeometry::GetLayer(fDet);
-  fN = 0; //fXref = 0.; Double_t ssx = 0.;
+  Int_t fN = 0;
   AliTRDcluster *c=0x0, **jc = &fClusters[0];
   for (Int_t ic=0; ic<kNtb; ic++, ++jc) {
     //zRow[ic] = -1;
@@ -979,7 +1145,15 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     if(c->GetNPads()>5) w = .2;
 
     //zRow[fN] = c->GetPadRow();
+    qc[fN]   = TMath::Abs(c->GetQ());
     // correct cluster position for PRF and v drift
+    //Int_t jc = TMath::Max(fN-3, 0);
+    //xc[fN]   = c->GetXloc(fT0, fVD, &qc[jc], &xc[jc]/*, z0 - c->GetX()*dzdx*/);
+    //Double_t s2 = fS2PRF + fDiffL*fDiffL*xc[fN]/(1.+2.*exb2)+tgg*xc[fN]*xc[fN]*exb2/12.;
+    //yc[fN]   = c->GetYloc(s2, fPadLength, xc[fN], fExB);
+    
+    // uncalibrated cluster correction 
+    // TODO remove
     Double_t x, y; GetClusterXY(c, x, y);
     xc[fN]   = fX0 - x;
     yc[fN]   = y;
@@ -994,13 +1168,11 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
 
     // ELABORATE CLUSTER ERROR
     // TODO to be moved to AliTRDcluster
-    q = TMath::Abs(c->GetQ());
-    Double_t tgg = (dydx-fExB)/(1.+dydx*fExB); tgg *= tgg;
     // basic y error (|| to track).
     sy[fN]  = xc[fN] < AliTRDgeometry::CamHght() ? 2. : sy0 + sya*TMath::Exp(1./(xc[fN]+syb));
     //printf("cluster[%d]\n\tsy[0] = %5.3e [um]\n", fN,  sy[fN]*1.e4);
     // y error due to total charge
-    sy[fN] += sqb*(1./q - sq0inv);
+    sy[fN] += sqb*(1./qc[fN] - sq0inv);
     //printf("\tsy[1] = %5.3e [um]\n", sy[fN]*1.e4);
     // y error due to PRF
     sy[fN] += scy[ily][0]*TMath::Gaus(c->GetCenter(), scy[ily][1], scy[ily][2]) - scy[ily][3];
@@ -1015,8 +1187,6 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     // error of drift length perpendicular to the track
     //sx += sxd0 + sxd1*d + sxd2*d*d;
     sx *= sx; // square sx
-    // update xref
-    //fXref += xc[fN]/sx; ssx+=1./sx;
 
     // add error from ExB 
     if(errors>0) sy[fN] += fExB*fExB*sx;
@@ -1031,7 +1201,7 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
     c->SetSigmaY2(sy[fN]);
 
     sy[fN]  = TMath::Sqrt(sy[fN]);
-    fitterY.AddPoint(&xc[fN], yc[fN]/*-yt*/, sy[fN]);
+    fitterY.AddPoint(&xc[fN], yc[fN], sy[fN]);
     fN++;
   }
   // to few clusters
@@ -1048,7 +1218,7 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
   fCov[2] = p[3]; // variance of dydx
   // the ref radial position is set at the minimum of 
   // the y variance of the tracklet
-  fXref = -fCov[1]/fCov[2]; //fXref = fX0 - fXref;
+  fX = -fCov[1]/fCov[2]; //fXref = fX0 - fXref;
 
   // fit XZ
   if(IsRowCross()){ 
@@ -1105,6 +1275,262 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Int_t errors)
 }
 
 
+
+
+
+
+//_____________________________________________________________________________
+void AliTRDseedV1::FitMI()
+{
+//
+// Fit the seed.
+// Marian Ivanov's version 
+//
+// linear fit on the y direction with respect to the reference direction. 
+// The residuals for each x (x = xc - x0) are deduced from:
+// dy = y - yt             (1)
+// the tilting correction is written :
+// y = yc + h*(zc-zt)      (2)
+// yt = y0+dy/dx*x         (3)
+// zt = z0+dz/dx*x         (4)
+// from (1),(2),(3) and (4)
+// dy = yc - y0 - (dy/dx + h*dz/dx)*x + h*(zc-z0)
+// the last term introduces the correction on y direction due to tilting pads. There are 2 ways to account for this:
+// 1. use tilting correction for calculating the y
+// 2. neglect tilting correction here and account for it in the error parametrization of the tracklet.
+  const Float_t kRatio  = 0.8;
+  const Int_t   kClmin  = 5;
+  const Float_t kmaxtan = 2;
+
+  if (TMath::Abs(fYref[1]) > kmaxtan){
+               //printf("Exit: Abs(fYref[1]) = %3.3f, kmaxtan = %3.3f\n", TMath::Abs(fYref[1]), kmaxtan);
+               return;              // Track inclined too much
+       }
+
+  Float_t  sigmaexp  = 0.05 + TMath::Abs(fYref[1] * 0.25); // Expected r.m.s in y direction
+  Float_t  ycrosscor = fPadLength * fTilt * 0.5;           // Y correction for crossing 
+  Int_t fNChange = 0;
+
+  Double_t sumw;
+  Double_t sumwx;
+  Double_t sumwx2;
+  Double_t sumwy;
+  Double_t sumwxy;
+  Double_t sumwz;
+  Double_t sumwxz;
+
+       // Buffering: Leave it constant fot Performance issues
+  Int_t    zints[kNtb];            // Histograming of the z coordinate 
+                                         // Get 1 and second max probable coodinates in z
+  Int_t    zouts[2*kNtb];       
+  Float_t  allowedz[kNtb];         // Allowed z for given time bin
+  Float_t  yres[kNtb];             // Residuals from reference
+  //Float_t  anglecor = fTilt * fZref[1];  // Correction to the angle
+  
+  Float_t pos[3*kNtb]; memset(pos, 0, 3*kNtb*sizeof(Float_t));
+  Float_t *fX = &pos[0], *fY = &pos[kNtb], *fZ = &pos[2*kNtb];
+  
+  Int_t fN  = 0; AliTRDcluster *c = 0x0; 
+  fN2 = 0;
+  for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
+    yres[i] = 10000.0;
+    if (!(c = fClusters[i])) continue;
+    if(!c->IsInChamber()) continue;
+    // Residual y
+    //yres[i] = fY[i] - fYref[0] - (fYref[1] + anglecor) * fX[i] + fTilt*(fZ[i] - fZref[0]);
+    fX[i] = fX0 - c->GetX();
+    fY[i] = c->GetY();
+    fZ[i] = c->GetZ();
+    yres[i] = fY[i] - fTilt*(fZ[i] - (fZref[0] - fX[i]*fZref[1]));
+    zints[fN] = Int_t(fZ[i]);
+    fN++;
+  }
+
+  if (fN < kClmin){
+    //printf("Exit fN < kClmin: fN = %d\n", fN);
+    return; 
+  }
+  Int_t nz = AliTRDtrackerV1::Freq(fN, zints, zouts, kFALSE);
+  Float_t fZProb   = zouts[0];
+  if (nz <= 1) zouts[3] = 0;
+  if (zouts[1] + zouts[3] < kClmin) {
+    //printf("Exit zouts[1] = %d, zouts[3] = %d\n",zouts[1],zouts[3]);
+    return;
+  }
+  
+  // Z distance bigger than pad - length
+  if (TMath::Abs(zouts[0]-zouts[2]) > 12.0) zouts[3] = 0;
+  
+  Int_t  breaktime = -1;
+  Bool_t mbefore   = kFALSE;
+  Int_t  cumul[kNtb][2];
+  Int_t  counts[2] = { 0, 0 };
+  
+  if (zouts[3] >= 3) {
+
+    //
+    // Find the break time allowing one chage on pad-rows
+    // with maximal number of accepted clusters
+    //
+    fNChange = 1;
+    for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
+      cumul[i][0] = counts[0];
+      cumul[i][1] = counts[1];
+      if (TMath::Abs(fZ[i]-zouts[0]) < 2) counts[0]++;
+      if (TMath::Abs(fZ[i]-zouts[2]) < 2) counts[1]++;
+    }
+    Int_t  maxcount = 0;
+    for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins(); i++) {
+      Int_t after  = cumul[AliTRDtrackerV1::GetNTimeBins()][0] - cumul[i][0];
+      Int_t before = cumul[i][1];
+      if (after + before > maxcount) { 
+        maxcount  = after + before; 
+        breaktime = i;
+        mbefore   = kFALSE;
+      }
+      after  = cumul[AliTRDtrackerV1::GetNTimeBins()-1][1] - cumul[i][1];
+      before = cumul[i][0];
+      if (after + before > maxcount) { 
+        maxcount  = after + before; 
+        breaktime = i;
+        mbefore   = kTRUE;
+      }
+    }
+    breaktime -= 1;
+  }
+
+  for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+    if (i >  breaktime) allowedz[i] =   mbefore  ? zouts[2] : zouts[0];
+    if (i <= breaktime) allowedz[i] = (!mbefore) ? zouts[2] : zouts[0];
+  }  
+
+  if (((allowedz[0] > allowedz[AliTRDtrackerV1::GetNTimeBins()]) && (fZref[1] < 0)) ||
+      ((allowedz[0] < allowedz[AliTRDtrackerV1::GetNTimeBins()]) && (fZref[1] > 0))) {
+    //
+    // Tracklet z-direction not in correspondance with track z direction 
+    //
+    fNChange = 0;
+    for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+      allowedz[i] = zouts[0];  // Only longest taken
+    } 
+  }
+  
+  if (fNChange > 0) {
+    //
+    // Cross pad -row tracklet  - take the step change into account
+    //
+    for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+      if (!fClusters[i]) continue; 
+      if(!fClusters[i]->IsInChamber()) continue;
+      if (TMath::Abs(fZ[i] - allowedz[i]) > 2) continue;
+      // Residual y
+      //yres[i] = fY[i] - fYref[0] - (fYref[1] + anglecor) * fX[i] /*+ fTilt*(fZ[i] - fZref[0])*/;   
+      yres[i] = fY[i] - fTilt*(fZ[i] - (fZref[0] - fX[i]*fZref[1]));
+/*      if (TMath::Abs(fZ[i] - fZProb) > 2) {
+        if (fZ[i] > fZProb) yres[i] += fTilt * fPadLength;
+        if (fZ[i] < fZProb) yres[i] -= fTilt * fPadLength;
+      }*/
+    }
+  }
+  
+  Double_t yres2[kNtb];
+  Double_t mean;
+  Double_t sigma;
+  for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+    if (!fClusters[i]) continue;
+    if(!fClusters[i]->IsInChamber()) continue;
+    if (TMath::Abs(fZ[i] - allowedz[i]) > 2) continue;
+    yres2[fN2] = yres[i];
+    fN2++;
+  }
+  if (fN2 < kClmin) {
+               //printf("Exit fN2 < kClmin: fN2 = %d\n", fN2);
+    fN2 = 0;
+    return;
+  }
+  AliMathBase::EvaluateUni(fN2,yres2,mean,sigma, Int_t(fN2*kRatio-2.));
+  if (sigma < sigmaexp * 0.8) {
+    sigma = sigmaexp;
+  }
+  //Float_t fSigmaY = sigma;
+
+  // Reset sums
+  sumw   = 0; 
+  sumwx  = 0; 
+  sumwx2 = 0;
+  sumwy  = 0; 
+  sumwxy = 0; 
+  sumwz  = 0;
+  sumwxz = 0;
+
+  fN2    = 0;
+  Float_t fMeanz = 0;
+  Float_t fMPads = 0;
+  fUsable = 0;
+  for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+    if (!fClusters[i]) continue;
+    if (!fClusters[i]->IsInChamber()) continue;
+    if (TMath::Abs(fZ[i] - allowedz[i]) > 2){fClusters[i] = 0x0; continue;}
+    if (TMath::Abs(yres[i] - mean) > 4.0 * sigma){fClusters[i] = 0x0;  continue;}
+    SETBIT(fUsable,i);
+    fN2++;
+    fMPads += fClusters[i]->GetNPads();
+    Float_t weight = 1.0;
+    if (fClusters[i]->GetNPads() > 4) weight = 0.5;
+    if (fClusters[i]->GetNPads() > 5) weight = 0.2;
+   
+       
+    Double_t x = fX[i];
+    //printf("x = %7.3f dy = %7.3f fit %7.3f\n", x, yres[i], fY[i]-yres[i]);
+    
+    sumw   += weight; 
+    sumwx  += x * weight; 
+    sumwx2 += x*x * weight;
+    sumwy  += weight * yres[i];  
+    sumwxy += weight * (yres[i]) * x;
+    sumwz  += weight * fZ[i];    
+    sumwxz += weight * fZ[i] * x;
+
+  }
+
+  if (fN2 < kClmin){
+               //printf("Exit fN2 < kClmin(2): fN2 = %d\n",fN2);
+    fN2 = 0;
+    return;
+  }
+  fMeanz = sumwz / sumw;
+  Float_t correction = 0;
+  if (fNChange > 0) {
+    // Tracklet on boundary
+    if (fMeanz < fZProb) correction =  ycrosscor;
+    if (fMeanz > fZProb) correction = -ycrosscor;
+  }
+
+  Double_t det = sumw * sumwx2 - sumwx * sumwx;
+  fYfit[0]    = (sumwx2 * sumwy  - sumwx * sumwxy) / det;
+  fYfit[1]    = (sumw   * sumwxy - sumwx * sumwy)  / det;
+  
+  fS2Y = 0;
+  for (Int_t i = 0; i < AliTRDtrackerV1::GetNTimeBins()+1; i++) {
+    if (!TESTBIT(fUsable,i)) continue;
+    Float_t delta = yres[i] - fYfit[0] - fYfit[1] * fX[i];
+    fS2Y += delta*delta;
+  }
+  fS2Y = TMath::Sqrt(fS2Y / Float_t(fN2-2));
+       // TEMPORARY UNTIL covariance properly calculated
+       fS2Y = TMath::Max(fS2Y, Float_t(.1));
+  
+  fZfit[0]   = (sumwx2 * sumwz  - sumwx * sumwxz) / det;
+  fZfit[1]   = (sumw   * sumwxz - sumwx * sumwz)  / det;
+//   fYfitR[0] += fYref[0] + correction;
+//   fYfitR[1] += fYref[1];
+//  fYfit[0]   = fYfitR[0];
+  fYfit[1]   = -fYfit[1];
+
+  UpdateUsed();
+}
+
+
 //___________________________________________________________________
 void AliTRDseedV1::Print(Option_t *o) const
 {
@@ -1113,7 +1539,7 @@ void AliTRDseedV1::Print(Option_t *o) const
   //
 
   AliInfo(Form("Det[%3d] Tilt[%+6.2f] Pad[%5.2f]", fDet, fTilt, fPadLength));
-  AliInfo(Form("Nattach[%2d] Nfit[%2d] Nuse[%2d] pads[%f]", fN, fN2, fNUsed, fMPads));
+  AliInfo(Form("N[%2d] Nuse[%2d]", fN2, fNUsed));
   AliInfo(Form("x[%7.2f] y[%7.2f] z[%7.2f] dydx[%5.2f] dzdx[%5.2f]", fX0, fYfit[0], fZfit[0], fYfit[1], fZfit[1]));
   AliInfo(Form("Ref        y[%7.2f] z[%7.2f] dydx[%5.2f] dzdx[%5.2f]", fYref[0], fZref[0], fYref[1], fZref[1]))
 
@@ -1121,7 +1547,7 @@ void AliTRDseedV1::Print(Option_t *o) const
   if(strcmp(o, "a")!=0) return;
 
   AliTRDcluster* const* jc = &fClusters[0];
-  for(int ic=0; ic<AliTRDtrackerV1::GetNTimeBins(); ic++, jc++) {
+  for(int ic=0; ic<kNTimeBins; ic++, jc++) {
     if(!(*jc)) continue;
     (*jc)->Print(o);
   }
@@ -1139,51 +1565,47 @@ Bool_t AliTRDseedV1::IsEqual(const TObject *o) const
   if(!inTracklet) return kFALSE;
 
   for (Int_t i = 0; i < 2; i++){
-    if ( fYref[i] != inTracklet->GetYref(i) ) return kFALSE;
-    if ( fZref[i] != inTracklet->GetZref(i) ) return kFALSE;
+    if ( fYref[i] != inTracklet->fYref[i] ) return kFALSE;
+    if ( fZref[i] != inTracklet->fZref[i] ) return kFALSE;
   }
   
-  if ( fSigmaY != inTracklet->GetSigmaY() ) return kFALSE;
-  if ( fSigmaY2 != inTracklet->GetSigmaY2() ) return kFALSE;
-  if ( fTilt != inTracklet->GetTilt() ) return kFALSE;
-  if ( fPadLength != inTracklet->GetPadLength() ) return kFALSE;
+  if ( fS2Y != inTracklet->fS2Y ) return kFALSE;
+  if ( fTilt != inTracklet->fTilt ) return kFALSE;
+  if ( fPadLength != inTracklet->fPadLength ) return kFALSE;
   
-  for (Int_t i = 0; i < knTimebins; i++){
-    if ( fX[i] != inTracklet->GetX(i) ) return kFALSE;
-    if ( fY[i] != inTracklet->GetY(i) ) return kFALSE;
-    if ( fZ[i] != inTracklet->GetZ(i) ) return kFALSE;
-    if ( fIndexes[i] != inTracklet->GetIndexes(i) ) return kFALSE;
-    if ( fUsable[i] != inTracklet->IsUsable(i) ) return kFALSE;
+  for (Int_t i = 0; i < kNTimeBins; i++){
+//     if ( fX[i] != inTracklet->GetX(i) ) return kFALSE;
+//     if ( fY[i] != inTracklet->GetY(i) ) return kFALSE;
+//     if ( fZ[i] != inTracklet->GetZ(i) ) return kFALSE;
+    if ( fIndexes[i] != inTracklet->fIndexes[i] ) return kFALSE;
   }
+  if ( fUsable != inTracklet->fUsable ) return kFALSE;
 
   for (Int_t i=0; i < 2; i++){
-    if ( fYfit[i] != inTracklet->GetYfit(i) ) return kFALSE;
-    if ( fZfit[i] != inTracklet->GetZfit(i) ) return kFALSE;
-    if ( fYfitR[i] != inTracklet->GetYfitR(i) ) return kFALSE;
-    if ( fZfitR[i] != inTracklet->GetZfitR(i) ) return kFALSE;
-    if ( fLabels[i] != inTracklet->GetLabels(i) ) return kFALSE;
+    if ( fYfit[i] != inTracklet->fYfit[i] ) return kFALSE;
+    if ( fZfit[i] != inTracklet->fZfit[i] ) return kFALSE;
+    if ( fLabels[i] != inTracklet->fLabels[i] ) return kFALSE;
   }
   
-  if ( fMeanz != inTracklet->GetMeanz() ) return kFALSE;
-  if ( fZProb != inTracklet->GetZProb() ) return kFALSE;
-  if ( fN2 != inTracklet->GetN2() ) return kFALSE;
-  if ( fNUsed != inTracklet->GetNUsed() ) return kFALSE;
-  if ( fFreq != inTracklet->GetFreq() ) return kFALSE;
-  if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
-  if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
+/*  if ( fMeanz != inTracklet->GetMeanz() ) return kFALSE;
+  if ( fZProb != inTracklet->GetZProb() ) return kFALSE;*/
+  if ( fN2 != inTracklet->fN2 ) return kFALSE;
+  if ( fNUsed != inTracklet->fNUsed ) return kFALSE;
+  //if ( fFreq != inTracklet->GetFreq() ) return kFALSE;
+  //if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
    
-  if ( fC != inTracklet->GetC() ) return kFALSE;
-  if ( fCC != inTracklet->GetCC() ) return kFALSE;
-  if ( fChi2 != inTracklet->GetChi2() ) return kFALSE;
+  if ( fC != inTracklet->fC ) return kFALSE;
+  //if ( fCC != inTracklet->GetCC() ) return kFALSE;
+  if ( fChi2 != inTracklet->fChi2 ) return kFALSE;
   //  if ( fChi2Z != inTracklet->GetChi2Z() ) return kFALSE;
 
-  if ( fDet != inTracklet->GetDetector() ) return kFALSE;
-  if ( fMom != inTracklet->GetMomentum() ) return kFALSE;
-  if ( fdX != inTracklet->GetdX() ) return kFALSE;
+  if ( fDet != inTracklet->fDet ) return kFALSE;
+  if ( fMom != inTracklet->fMom ) return kFALSE;
+  if ( fdX != inTracklet->fdX ) return kFALSE;
   
-  for (Int_t iCluster = 0; iCluster < knTimebins; iCluster++){
+  for (Int_t iCluster = 0; iCluster < kNTimeBins; iCluster++){
     AliTRDcluster *curCluster = fClusters[iCluster];
-    AliTRDcluster *inCluster = inTracklet->GetClusters(iCluster);
+    AliTRDcluster *inCluster = inTracklet->fClusters[iCluster];
     if (curCluster && inCluster){
       if (! curCluster->IsEqual(inCluster) ) {
         curCluster->Print();
index a9fa408b3bf98f45beab259709cf2ad0ed53e346..d56ca089b09d0ab4c6d8e7642b311acc5cbdea85 100644 (file)
 //                                                                        //
 ////////////////////////////////////////////////////////////////////////////
 
-#ifndef ALITRDSEED_H
-#include "AliTRDseed.h"
+#ifndef ROOT_TObject
+#include "TObject.h"
+#endif
+
+#ifndef ROOT_TMath
+#include "TMath.h"
 #endif
 
 #ifndef ALITRDGEOMETRY_H
@@ -35,18 +39,20 @@ class AliTRDtrackingChamber;
 class AliTRDcluster;
 class AliTRDtrackV1;
 class AliTRDReconstructor;
-class AliTRDseedV1 : public AliTRDseed
-{
-
-  public:
-
-  enum {
-    knSlices = 10
+class AliTRDseedV1 : public TObject //TODO we should inherit 
+{                                   // AliTRDtrackletBase
+public:
+  enum ETRDtrackletBuffers {    
+    kNtb = 32           // max clusters/pad row
+   ,kNTimeBins = 2*kNtb // max number of clusters/tracklet
+   ,kNSlices = 10       // max dEdx slices
   };
+
   // bits from 0-13 are reserved by ROOT (see TObject.h)
-  enum AliTRDtrackletStatus {
-    kOwner    = BIT(14)
-  , kRowCross = BIT(15) 
+  enum ETRDtrackletStatus {
+    kOwner    = BIT(14) // owner of its clusters
+   ,kRowCross = BIT(15) // pad row cross tracklet
+   ,kCalib    = BIT(16) // calibrated tracklet
   };
 
   AliTRDseedV1(Int_t det = -1);
@@ -60,94 +66,150 @@ class AliTRDseedV1 : public AliTRDseed
   Bool_t         AttachClusters(
               AliTRDtrackingChamber *chamber, Bool_t tilt = kFALSE);
   void      Bootstrap(const AliTRDReconstructor *rec);
+  void      Calibrate();
   void      CookdEdx(Int_t nslices);
+  void      CookLabels();
   Bool_t    Fit(Bool_t tilt=kTRUE, Int_t errors = 2);
-
+  void      FitMI();
   Bool_t    Init(AliTRDtrackV1 *track);
   inline void      Init(const AliRieman *fit);
   Bool_t    IsEqual(const TObject *inTracklet) const;
+  Bool_t    IsCalibrated() const     { return TestBit(kCalib);}
   Bool_t    IsOwner() const          { return TestBit(kOwner);}
+  Bool_t    IsOK() const             { return fN2 > 4;}
   Bool_t    IsRowCross() const       { return TestBit(kRowCross);}
+  Bool_t    IsUsable(Int_t i) const  { return TESTBIT(fUsable, i);}
 
-  inline Float_t   GetChi2Z(const Float_t z = 999.) const;
-  inline Float_t   GetChi2Y(const Float_t y = 999.) const;
+  Float_t   GetC() const             { return fC; }
+  Float_t   GetChi2() const          { return fChi2; }
+  inline Float_t   GetChi2Z() const;
+  inline Float_t   GetChi2Y() const;
   static void      GetClusterXY(const AliTRDcluster *c, Double_t &x, Double_t &y);
   void      GetCovAt(Double_t x, Double_t *cov) const;
   void      GetCovXY(Double_t *cov) const { memcpy(cov, &fCov[0], 3*sizeof(Double_t));}
   void      GetCovRef(Double_t *cov) const { memcpy(cov, &fRefCov[0], 3*sizeof(Double_t));}
-  Double_t* GetCrossXYZ()            { return &fCross[0];}
-  Double_t  GetCrossSz2() const      { return fCross[3];}
   Float_t   GetdX() const            { return fdX;}
   Float_t*  GetdEdx()                { return &fdEdx[0];}
   Float_t   GetdQdl(Int_t ic) const;
   Int_t     GetDetector() const      { return fDet;}
-  Float_t   GetExB() const           { return fExB;}
+  void      GetCalibParam(Float_t &exb, Float_t &vd, Float_t &t0, Float_t &s2, Float_t &dl, Float_t &dt) const    { 
+              exb = fExB; vd = fVD; t0 = fT0; s2 = fS2PRF; dl = fDiffL; dt = fDiffT;}
+  AliTRDcluster*  GetClusters(Int_t i) const               { return i<0 || i>=kNTimeBins ? 0x0 : fClusters[i];}
+  Int_t     GetIndexes(Int_t i) const{ return i<0 || i>=kNTimeBins ? -1 : fIndexes[i];}
+  Int_t     GetLabels(Int_t i) const { return fLabels[i];}  
   Double_t  GetMomentum() const      { return fMom;}
   Int_t     GetN() const             { return fN2;}
+  Int_t     GetN2() const            { return fN2;}
+  Int_t     GetNUsed() const         { return fNUsed;}
   Float_t   GetQuality(Bool_t kZcorr) const;
+  Float_t   GetPadLength() const     { return fPadLength;}
   Int_t     GetPlane() const         { return AliTRDgeometry::GetLayer(fDet);    }
 
-  Double_t* GetProbability();
-  Double_t  GetSnp() const           { return fSnp;}
-  Double_t  GetTgl() const           { return fTgl;}
-  Float_t   GetXref() const          { return fX0 - fXref;}
+  Float_t*  GetProbability();
+  Float_t   GetS2Y() const           { return fS2Y;}
+  Float_t   GetS2Z() const           { return fS2Z;}
+  Float_t   GetSigmaY() const        { return fS2Y > 0. ? TMath::Sqrt(fS2Y) : 0.2;}
+  Float_t   GetSnp() const           { return fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);}
+  Float_t   GetTgl() const           { return fZref[1];}
+  Float_t   GetTilt() const          { return fTilt;}
+  Float_t   GetX0() const            { return fX0;}
+  Float_t   GetX() const             { return fX0 - fX;}
+  Float_t   GetY() const             { return GetYat(fX); }
   Double_t  GetYat(Double_t x) const { return fYfit[0] - fYfit[1] * (fX0-x);}
+  Float_t   GetYfit(Int_t id) const { return fYfit[id];}
+  Float_t   GetYref(Int_t id) const { return fYref[id];}
+  Float_t   GetZ() const             { return GetZat(fX); }
   Double_t  GetZat(Double_t x) const { return fZfit[0] - fZfit[1] * (fX0-x);}
-  
+  Float_t   GetZfit(Int_t id) const { return fZfit[id];}
+  Float_t   GetZref(Int_t id) const { return fZref[id];}
+  Long_t    GetUsabilityMap() const { return fUsable; }
+
   inline AliTRDcluster* NextCluster();
   inline AliTRDcluster* PrevCluster();
   void      Print(Option_t *o = "") const;
   inline void ResetClusterIter(Bool_t forward = kTRUE);
+  void      Reset();
 
+  void      SetC(Float_t c)         { fC = c;}
+  void      SetChi2(Float_t chi2)   { fChi2 = chi2;}
   void      SetCovRef(const Double_t *cov) { memcpy(&fRefCov[0], cov, 3*sizeof(Double_t));}
+  void      SetIndexes(Int_t i, Int_t idx) { fIndexes[i]  = idx; }
+  void      SetLabels(Int_t *lbls)   { memcpy(fLabels, lbls, 3*sizeof(Int_t)); }
   void      SetMomentum(Double_t mom){ fMom = mom;}
   void      SetOwner();
+  void      SetTilt(Float_t tilt)    { fTilt = tilt; }
+  void      SetPadLength(Float_t len){ fPadLength = len;}
   void      SetDetector(Int_t d)     { fDet = d;  }
   void      SetDX(Float_t inDX)      { fdX = inDX;}
-  void      SetSnp(Double_t snp)     { fSnp = snp;}
-  void      SetTgl(Double_t tgl)     { fTgl = tgl;}
   void      SetReconstructor(const AliTRDReconstructor *rec) {fReconstructor = rec;}
+  void      SetX0(Float_t x0)        { fX0 = x0; }
+  void      SetYref(Int_t i, Float_t y) { fYref[i]     = y;}
+  void      SetZref(Int_t i, Float_t z) { fZref[i]     = z;}
+  void      SetUsabilityMap(Long_t um)  { fUsable = um; }
   void      UpDate(const AliTRDtrackV1* trk);
+  void      UpdateUsed();
+  void      UseClusters();
 
 protected:
   void Copy(TObject &ref) const;
-  void SetExB();
 
 private:
   const AliTRDReconstructor *fReconstructor;//! local reconstructor
-  AliTRDcluster    **fClusterIter;          //! clusters iterator
+  AliTRDcluster  **fClusterIter;            //! clusters iterator
+  Int_t            fIndexes[kNTimeBins];    //! Indexes
+  Float_t          fExB;                    //! tg(a_L) @ tracklet location
+  Float_t          fVD;                     //! drift velocity @ tracklet location
+  Float_t          fT0;                     //! time 0 @ tracklet location
+  Float_t          fS2PRF;                  //! sigma^2 PRF for xd->0 and phi=a_L 
+  Float_t          fDiffL;                  //! longitudinal diffusion coefficient
+  Float_t          fDiffT;                  //! transversal diffusion coefficient
   Char_t           fClusterIdx;             //! clusters iterator
-  Int_t            fDet;                    //  TRD detector
-  Float_t          fMom;                    //  Momentum estimate for tracklet [GeV/c]
-  Float_t          fSnp;                    // sin of track with respect to x direction in XY plane    
-  Float_t          fTgl;                    // tg of track with respect to x direction in XZ plane     
+  Long_t           fUsable;                 //! bit map of usable clusters
+  UChar_t          fN2;                     // number of clusters attached
+  UChar_t          fNUsed;                  // number of used usable clusters
+  Short_t          fDet;                    // TRD detector
+  Float_t          fTilt;                   // local tg of the tilt angle 
+  Float_t          fPadLength;              // local pad length 
+  AliTRDcluster   *fClusters[kNTimeBins];   // Clusters
+  Float_t          fYref[2];                //  Reference y
+  Float_t          fZref[2];                //  Reference z
+  Float_t          fYfit[2];                //  Y fit position +derivation
+  Float_t          fZfit[2];                //  Z fit position
+  Float_t          fMom;                    //  Momentum estimate @ tracklet [GeV/c]
   Float_t          fdX;                     // length of time bin
-  Float_t          fXref;                   // average radial position of clusters
-  Float_t          fExB;                    // tg(a_L) for the tracklet reagion
-  Float_t          fdEdx[knSlices];         // dE/dx measurements for tracklet
-  Double_t         fCross[4];               // spatial parameters of the pad row crossing
+  Float_t          fX0;                     // anode wire position
+  Float_t          fX;                      // radial position of the tracklet
+  Float_t          fY;                      // r-phi position of the tracklet
+  Float_t          fZ;                      // z position of the tracklet
+  Float_t          fS2Y;                    // estimated resolution in the r-phi direction 
+  Float_t          fS2Z;                    // estimated resolution in the z direction 
+  Float_t          fC;                      // Curvature
+  Float_t          fChi2;                   // Global chi2  
+  Float_t          fdEdx[kNSlices];         // dE/dx measurements for tracklet
+  Float_t          fProb[AliPID::kSPECIES]; //  PID probabilities
+  Int_t            fLabels[3];              // most frequent MC labels and total number of different labels
   Double_t         fRefCov[3];              // covariance matrix of the track in the yz plane
   Double_t         fCov[3];                 // covariance matrix of the tracklet in the xy plane
-  Double_t         fProb[AliPID::kSPECIES]; //  PID probabilities
-
-  ClassDef(AliTRDseedV1, 4)                 //  New TRD seed 
 
+  ClassDef(AliTRDseedV1, 5)                 // The offline TRD tracklet 
 };
 
 //____________________________________________________________
-inline Float_t AliTRDseedV1::GetChi2Z(const Float_t z) const
+inline Float_t AliTRDseedV1::GetChi2Z() const
 {
-  Float_t z1  = (z == 999.) ? fMeanz : z;
-  Float_t chi = fZref[0] - z1;
-  return chi*chi;
+  Double_t dz = fZref[0]-fZfit[0]; dz*=dz;
+  Double_t cov[3]; GetCovAt(fX, cov);
+  Double_t s2 = fRefCov[2]+cov[2];
+  return s2 > 0. ? dz/s2 : 0.; 
 }
 
 //____________________________________________________________
-inline Float_t AliTRDseedV1::GetChi2Y(const Float_t y) const
+inline Float_t AliTRDseedV1::GetChi2Y() const
 {
-  Float_t y1  = (y == 999.) ? fYfitR[0] : y;
-  Float_t chi = fYref[0] - y1;
-  return chi*chi;
+  Double_t dy = fYref[0]-fYfit[0]; dy*=dy;
+  Double_t cov[3]; GetCovAt(fX, cov);
+  Double_t s2 = fRefCov[0]+cov[0];
+  return s2 > 0. ? dy/s2 : 0.; 
 }
 
 //____________________________________________________________
@@ -168,7 +230,7 @@ inline AliTRDcluster* AliTRDseedV1::NextCluster()
 // Forward iterator
 
   fClusterIdx++; fClusterIter++;
-  while(fClusterIdx < AliTRDseed::knTimebins){
+  while(fClusterIdx < kNTimeBins){
     if(!(*fClusterIter)){ 
       fClusterIdx++; 
       fClusterIter++;
@@ -208,8 +270,8 @@ inline void AliTRDseedV1::ResetClusterIter(Bool_t forward)
     fClusterIter = &fClusters[0]; fClusterIter--; 
     fClusterIdx=-1;
   } else {
-    fClusterIter = &fClusters[AliTRDseed::knTimebins-1]; fClusterIter++; 
-    fClusterIdx=AliTRDseed::knTimebins;
+    fClusterIter = &fClusters[kNTimeBins-1]; fClusterIter++; 
+    fClusterIdx=kNTimeBins;
   }
 }
 
index 8be1deeb24066a56a74ab568bcd1a5daacd3be48..ce696ffd0704f10dd133d28314625c9be597cd88 100644 (file)
@@ -44,6 +44,8 @@ AliTRDtrackV1::AliTRDtrackV1() : AliKalmanTrack()
   ,fDE(0.)
   ,fReconstructor(0x0)
   ,fBackupTrack(0x0)
+  ,fTrackLow(0x0)
+  ,fTrackHigh(0x0)
 {
   //
   // Default constructor
@@ -67,6 +69,8 @@ AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref) : AliKalmanTrack(ref)
   ,fDE(ref.fDE)
   ,fReconstructor(ref.fReconstructor)
   ,fBackupTrack(0x0)
+  ,fTrackLow(0x0)
+  ,fTrackHigh(0x0)
 {
   //
   // Copy constructor
@@ -78,7 +82,9 @@ AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref) : AliKalmanTrack(ref)
     fTrackletIndex[ip] = ref.fTrackletIndex[ip];
     fTracklet[ip]      = ref.fTracklet[ip];
   }
-
+  if(ref.fTrackLow) fTrackLow = new AliExternalTrackParam(*ref.fTrackLow);
+  if(ref.fTrackHigh) fTrackHigh = new AliExternalTrackParam(*ref.fTrackHigh);
   for (Int_t i = 0; i < 3;i++) fBudget[i]      = ref.fBudget[i];
 
        for(Int_t is = 0; is<AliPID::kSPECIES; is++) fPID[is] = ref.fPID[is];
@@ -92,6 +98,8 @@ AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t) : AliKalmanTrack()
   ,fDE(0.)
   ,fReconstructor(0x0)
   ,fBackupTrack(0x0)
+  ,fTrackLow(0x0)
+  ,fTrackHigh(0x0)
 {
   //
   // Constructor from AliESDtrack
@@ -140,6 +148,8 @@ AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5], const Do
   ,fDE(0.)
   ,fReconstructor(0x0)
   ,fBackupTrack(0x0)
+  ,fTrackLow(0x0)
+  ,fTrackHigh(0x0)
 {
   //
   // The stand alone tracking constructor
@@ -200,8 +210,10 @@ AliTRDtrackV1::~AliTRDtrackV1()
   //AliInfo("");
   //printf("I-AliTRDtrackV1::~AliTRDtrackV1() : Owner[%s]\n", TestBit(kOwner)?"YES":"NO");
 
-  if(fBackupTrack) delete fBackupTrack;
-  fBackupTrack = 0x0;
+  if(fBackupTrack) delete fBackupTrack; fBackupTrack = 0x0;
+
+  if(fTrackLow) delete fTrackLow; fTrackLow = 0x0;
+  if(fTrackHigh) delete fTrackHigh; fTrackHigh = 0x0;
 
   for(Int_t ip=0; ip<kNplane; ip++){
     if(TestBit(kOwner) && fTracklet[ip]) delete fTracklet[ip];
@@ -226,7 +238,7 @@ Bool_t AliTRDtrackV1::CookLabel(Float_t wrong)
   AliTRDcluster *c    = 0x0;
   for (Int_t ip = 0; ip < kNplane; ip++) {
     if(fTrackletIndex[ip] == 0xffff) continue;
-    for (Int_t ic = 0; ic < AliTRDseed::knTimebins; ic++) {
+    for (Int_t ic = 0; ic < AliTRDseedV1::kNTimeBins; ic++) {
       if(!(c = fTracklet[ip]->GetClusters(ic))) continue;
       for (Int_t k = 0; k < 3; k++) { 
         label      = c->GetLabel(k);
@@ -277,7 +289,7 @@ Bool_t AliTRDtrackV1::CookPID()
   fPIDquality = 0;
   
   // steer PID calculation @ tracklet level
-  Double_t *prob = 0x0;
+  Float_t *prob = 0x0;
   for(int ip=0; ip<kNplane; ip++){
     if(fTrackletIndex[ip] == 0xffff) continue;
     if(!fTracklet[ip]->IsOK()) continue;
@@ -326,7 +338,7 @@ AliTRDcluster* AliTRDtrackV1::GetCluster(Int_t id)
       continue;
     }
     AliTRDcluster *c = 0x0;
-    for(Int_t ic=AliTRDseed::knTimebins-1; ic>=0; ic--){
+    for(Int_t ic=AliTRDseedV1::kNTimeBins; ic--;){
       if(!(c = fTracklet[ip]->GetClusters(ic))) continue;
 
       if(n<id){n++; continue;}
@@ -347,7 +359,7 @@ Int_t  AliTRDtrackV1::GetClusterIndex(Int_t id) const
       continue;
     }
     AliTRDcluster *c = 0x0;
-    for(Int_t ic=AliTRDseed::knTimebins-1; ic>=0; ic--){
+    for(Int_t ic=AliTRDseedV1::kNTimeBins; ic--;){
       if(!(c = fTracklet[ip]->GetClusters(ic))) continue;
 
       if(n<id){n++; continue;}
@@ -370,7 +382,11 @@ Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt) const
   Double_t cov[3];
   trklt->GetCovAt(x, cov);
   
-  return AliExternalTrackParam::GetPredictedChi2(p, cov);
+  const Double_t *cc = GetCovariance();
+  Double_t dy = p[0]-GetY(); dy*=dy;
+  Double_t s2 = cov[0]+cc[0];
+  return s2 > 0. ? dy/s2 : 0.; 
+  //return AliExternalTrackParam::GetPredictedChi2(p, cov);
 }
 
 //_______________________________________________________________
@@ -699,6 +715,20 @@ void AliTRDtrackV1::SetTracklet(AliTRDseedV1 *trklt, Int_t index)
   fTrackletIndex[plane] = index;
 }
 
+//_______________________________________________________________
+void AliTRDtrackV1::SetTrackLow()
+{
+  const AliExternalTrackParam *op = dynamic_cast<const AliExternalTrackParam*>(this);
+  fTrackLow = fTrackLow ? new(fTrackLow) AliExternalTrackParam(*op) : new AliExternalTrackParam(*op);
+}
+
+//_______________________________________________________________
+void AliTRDtrackV1::SetTrackHigh(const AliExternalTrackParam *op)
+{
+  if(!op) op = dynamic_cast<const AliExternalTrackParam*>(this);
+  fTrackHigh = fTrackHigh ? new(fTrackHigh) AliExternalTrackParam(*op) : new AliExternalTrackParam(*op);
+}
+
 //_______________________________________________________________
 void AliTRDtrackV1::UnsetTracklet(Int_t plane)
 {
index 3a36075d09528269f1259f92e9e40de843d45424..524773dc7dcb0a94451a7cd2cb0dc66847097ffe 100644 (file)
@@ -29,20 +29,21 @@ class AliTRDReconstructor;
 class AliTRDtrackV1 : public AliKalmanTrack
 {
 public:
-  enum { kMAXCLUSTERSPERTRACK = 210 };
-    
-  enum { kNdet      = 540
-        , kNstacks   =  90
-        , kNplane    =   AliESDtrack::kTRDnPlanes
-        , kNcham     =   5
-        , kNsect     =  18
-        , kNslice    =   3
-        , kNMLPslice =   8 };
+  enum ETRDtrackV1Size { 
+    kNdet      = AliTRDgeometry::kNdet
+   ,kNstacks   = AliTRDgeometry::kNstack*AliTRDgeometry::kNsector
+   ,kNplane    = AliTRDgeometry::kNlayer
+   ,kNcham     = AliTRDgeometry::kNstack
+   ,kNsect     = AliTRDgeometry::kNsector
+   ,kNslice    =   3
+   ,kNMLPslice =   8 
+   ,kMAXCLUSTERSPERTRACK = 210
+  };
   
   // bits from 0-13 are reserved by ROOT (see TObject.h)
-  enum AliTRDtrackStatus {
-          kOwner   = BIT(14)
-        , kStopped = BIT(15) 
+  enum ETRDtrackV1Status {
+    kOwner   = BIT(14)
+   ,kStopped = BIT(15) 
   };
 
   AliTRDtrackV1();
@@ -69,8 +70,12 @@ public:
   Double_t       GetPredictedChi2(const AliCluster* /*c*/) const                   { return 0.0; }
   Int_t          GetProlongation(Double_t xk, Double_t &y, Double_t &z);
   AliTRDseedV1*  GetTracklet(Int_t plane) const {return plane >=0 && plane <kNplane ? fTracklet[plane] : 0x0;}
-  Int_t          GetTrackletIndex(Int_t plane) const          { return (plane>=0 && plane<kNplane) ? fTrackletIndex[plane] : -1;} 
-  UShort_t*      GetTrackletIndexes() {return &fTrackletIndex[0];}
+  Int_t          GetTrackletIndex(Int_t plane) const          { return (plane>=0 && plane<kNplane) ? fTrackletIndex[plane] : -1;}
+  AliExternalTrackParam*
+                 GetTrackLow() const  { return fTrackLow;} 
+  AliExternalTrackParam*
+                 GetTrackHigh() const  { return fTrackHigh;} 
+  UShort_t*      GetTrackletIndexes() { return &fTrackletIndex[0];}
   
   Bool_t         IsEqual(const TObject *inTrack) const;
   Bool_t         IsOwner() const {return TestBit(kOwner);};
@@ -91,6 +96,8 @@ public:
   void           SetPIDquality(UChar_t inPIDquality){fPIDquality = inPIDquality;};
   void           SetStopped(Bool_t stop) {SetBit(kStopped, stop);}
   void           SetTracklet(AliTRDseedV1 *trklt,  Int_t index);
+  void           SetTrackLow();
+  void           SetTrackHigh(const AliExternalTrackParam *op=0x0);
   inline void    SetReconstructor(const AliTRDReconstructor *rec);
   inline Float_t StatusForTOF();
   void           UnsetTracklet(Int_t plane);
@@ -108,9 +115,10 @@ private:
   const AliTRDReconstructor *fReconstructor;//! reconstructor link 
   AliTRDseedV1 *fTracklet[kNplane];   //  Tracklets array defining the track
   AliTRDtrackV1 *fBackupTrack;        // Backup track
+  AliExternalTrackParam *fTrackLow;   // parameters of the track which enter TRD from below (TPC) 
+  AliExternalTrackParam *fTrackHigh;  // parameters of the track which enter TRD from above (HMPID, PHOS) 
 
-
-  ClassDef(AliTRDtrackV1, 3)          // new TRD track
+  ClassDef(AliTRDtrackV1, 4)          // new TRD track
 };
 
 //____________________________________________________
index 7aa821489723a7823444923dc3bb95fe1ed54e9f..4e704a48e9a1a2b00dfca3851472f9fb92f58f32 100644 (file)
@@ -1,17 +1,17 @@
 /**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+*                                                                        *
+* Author: The ALICE Off-line Project.                                    *
+* Contributors are mentioned in the code where appropriate.              *
+*                                                                        *
+* Permission to use, copy, modify and distribute this software and its   *
+* documentation strictly for non-commercial purposes is hereby granted   *
+* without fee, provided that the above copyright notice appears in all   *
+* copies and that both the copyright notice and this permission notice   *
+* appear in the supporting documentation. The authors make no claims     *
+* about the suitability of this software for any purpose. It is          *
+* provided "as is" without express or implied warranty.                  *
+**************************************************************************/
 
 /* $Id: AliTRDtrackerDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */
 
@@ -51,25 +51,25 @@ Int_t AliTRDtrackerDebug::fgCandidateNumber = 0;
 
 //____________________________________________________
 AliTRDtrackerDebug::AliTRDtrackerDebug() : AliTRDtrackerV1()
-       ,fOutputStreamer(0x0)
-       ,fTree(0x0)
-       ,fTracklet(0x0)
-       ,fTrack(0x0)
-       ,fNClusters(0)
-       ,fAlpha(0.)
+  ,fOutputStreamer(0x0)
+  ,fTree(0x0)
+  ,fTracklet(0x0)
+  ,fTrack(0x0)
+  ,fNClusters(0)
+  ,fAlpha(0.)
 {
         //
-       // Default constructor
-       //
-       fOutputStreamer = new TTreeSRedirector("TRD.Debug.root");
+  // Default constructor
+  //
+  fOutputStreamer = new TTreeSRedirector("TRD.Debug.root");
 }
 
 //____________________________________________________
 AliTRDtrackerDebug::~AliTRDtrackerDebug()
 {
-       // destructor
-       
-       delete fOutputStreamer;
+  // destructor
+  
+  delete fOutputStreamer;
 }
 
 
@@ -84,52 +84,52 @@ void AliTRDtrackerDebug::Draw(Option_t *)
 Bool_t AliTRDtrackerDebug::Init()
 {
 // steer linking data for various debug streams        
-       fTrack = new AliTRDtrackV1();
-       fTree->SetBranchAddress("ncl", &fNClusters);
-       fTree->SetBranchAddress("track.", &fTrack);
-       return kTRUE;
+  fTrack = new AliTRDtrackV1();
+  fTree->SetBranchAddress("ncl", &fNClusters);
+  fTree->SetBranchAddress("track.", &fTrack);
+  return kTRUE;
 }
 
 //____________________________________________________
 Bool_t AliTRDtrackerDebug::Open(const char *method)
 {
-       // Connect to the tracker debug file
-       
-       TDirectory *savedir = gDirectory; 
-       TFile::Open("TRD.TrackerDebugger.root");
-       fTree = (TTree*)gFile->Get(method);
-       if(!fTree){
-               AliInfo(Form("Can not find debug stream for the %s method.\n", method));
-               savedir->cd();
-               return kFALSE;
-       }
-       savedir->cd();
-       return kTRUE;
+  // Connect to the tracker debug file
+  
+  TDirectory *savedir = gDirectory; 
+  TFile::Open("TRD.TrackerDebugger.root");
+  fTree = (TTree*)gFile->Get(method);
+  if(!fTree){
+    AliInfo(Form("Can not find debug stream for the %s method.\n", method));
+    savedir->cd();
+    return kFALSE;
+  }
+  savedir->cd();
+  return kTRUE;
 }
 
 //____________________________________________________
 Int_t AliTRDtrackerDebug::Process()
 {
 // steer debug process threads
-       
-       for(int it = 0; it<fTree->GetEntries(); it++){
-               if(!fTree->GetEntry(it)) continue;
-               if(!fNClusters) continue;
-               fAlpha = fTrack->GetAlpha();
-               //printf("Processing track %d [%d] ...\n", it, fNClusters);
-               ResidualsTrackletsTrack();
-
-               const AliTRDseedV1 *tracklet = 0x0;
-               for(int ip = 5; ip>=0; ip--){
-                       if(!(tracklet = fTrack->GetTracklet(ip))) continue;
-                       if(!tracklet->GetN()) continue;
-                       
-                       ResidualsClustersTrack(tracklet);
-                       ResidualsClustersTracklet(tracklet);
-                       ResidualsClustersParametrisation(tracklet);
-               }
-       }
-       return kTRUE;
+  
+  for(int it = 0; it<fTree->GetEntries(); it++){
+    if(!fTree->GetEntry(it)) continue;
+    if(!fNClusters) continue;
+    fAlpha = fTrack->GetAlpha();
+    //printf("Processing track %d [%d] ...\n", it, fNClusters);
+    ResidualsTrackletsTrack();
+
+    const AliTRDseedV1 *tracklet = 0x0;
+    for(int ip = 5; ip>=0; ip--){
+      if(!(tracklet = fTrack->GetTracklet(ip))) continue;
+      if(!tracklet->GetN()) continue;
+      
+      ResidualsClustersTrack(tracklet);
+      ResidualsClustersTracklet(tracklet);
+      ResidualsClustersParametrisation(tracklet);
+    }
+  }
+  return kTRUE;
 }
 
 
@@ -137,99 +137,99 @@ Int_t AliTRDtrackerDebug::Process()
 void AliTRDtrackerDebug::ResidualsClustersTrack(const AliTRDseedV1 *tracklet)
 {
 // Calculate averange distances from clusters to the TRD track 
-       
-       Double_t x[3]; 
-       AliTRDcluster *c = 0x0;
-       for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
-               if(!(c = tracklet->GetClusters(ic))) continue;
-               Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
-
-               // propagate track to cluster 
-               PropagateToX(*fTrack, xc, 2.); 
-               fTrack->GetXYZ(x);
-               
-               // transform to local tracking coordinates
-               //Double_t xg =  x[0] * TMath::Cos(fAlpha) + x[1] * TMath::Sin(fAlpha); 
-               Double_t yg = -x[0] * TMath::Sin(fAlpha) + x[1] * TMath::Cos(fAlpha);
-
-               // apply tilt pad correction
-               yc+= (zc - x[2]) * tracklet->GetTilt();
-               
-               Double_t dy = yc-yg;
-
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "ResidualsClustersTrack"
-                       << "c.="   << c
-                       << "dy="   << dy
-                       << "\n";
-       }
+  
+  Double_t x[3]; 
+  AliTRDcluster *c = 0x0;
+  for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
+    if(!(c = tracklet->GetClusters(ic))) continue;
+    Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
+
+    // propagate track to cluster 
+    PropagateToX(*fTrack, xc, 2.); 
+    fTrack->GetXYZ(x);
+    
+    // transform to local tracking coordinates
+    //Double_t xg =  x[0] * TMath::Cos(fAlpha) + x[1] * TMath::Sin(fAlpha); 
+    Double_t yg = -x[0] * TMath::Sin(fAlpha) + x[1] * TMath::Cos(fAlpha);
+
+    // apply tilt pad correction
+    yc+= (zc - x[2]) * tracklet->GetTilt();
+    
+    Double_t dy = yc-yg;
+
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "ResidualsClustersTrack"
+      << "c.="   << c
+      << "dy="   << dy
+      << "\n";
+  }
 }
 
 //____________________________________________________
 void AliTRDtrackerDebug::ResidualsClustersTracklet(const AliTRDseedV1 *tracklet) const
 {
 // Calculates distances from clusters to tracklets
-       
-       Double_t x0 = tracklet->GetX0(), 
-                y0 = tracklet->GetYfit(0), 
-                ys = tracklet->GetYfit(1);
-                //z0 = tracklet->GetZfit(0), 
-                //zs = tracklet->GetZfit(1);
-       
-       AliTRDcluster *c = 0x0;
-       for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
-               if(!(c = tracklet->GetClusters(ic))) continue;
-               Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
-               Double_t dy = yc- (y0-(x0-xc)*ys);
-
-               //To draw  use : 
+  
+  Double_t x0 = tracklet->GetX0(), 
+          y0 = tracklet->GetYfit(0), 
+          ys = tracklet->GetYfit(1);
+          //z0 = tracklet->GetZfit(0), 
+          //zs = tracklet->GetZfit(1);
+  
+  AliTRDcluster *c = 0x0;
+  for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
+    if(!(c = tracklet->GetClusters(ic))) continue;
+    Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
+    Double_t dy = yc- (y0-(x0-xc)*ys);
+
+    //To draw  use : 
     //ResidualsClustersTracklet->Draw("TMath::Abs(10.*dy):TMath::ATan(ys)*TMath::RadToDeg()>>h(20, -40, 40)", "", "prof");
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "ResidualsClustersTracklet"
-                       << "c.="   << c
-                       << "ys="   << ys
-                       << "dy="   << dy
-                       << "\n";
-       }
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "ResidualsClustersTracklet"
+      << "c.="   << c
+      << "ys="   << ys
+      << "dy="   << dy
+      << "\n";
+  }
 }
 
 //____________________________________________________
 void AliTRDtrackerDebug::ResidualsClustersParametrisation(const AliTRDseedV1 *tracklet) const
 {
 // Calculates distances from clusters to Rieman fit.
-       
-       // store cluster positions
-       Double_t x0 = tracklet->GetX0();
-       AliTRDcluster *c = 0x0;
-       
-       Double_t x[2]; Int_t ncl, mcl, jc;
-       TLinearFitter fitter(3, "hyp2");
-       for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
-               if(!(c = tracklet->GetClusters(ic))) continue;
-               Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
-               
-               jc = ic; ncl = 0; mcl=0; fitter.ClearPoints();
-               while(ncl<6){
-                       // update index
-                       mcl++;
-                       jc = ic + ((mcl&1)?-1:1)*(mcl>>1);
-
-                       if(jc<0 || jc>=35) continue;
-                       if(!(c = tracklet->GetClusters(jc))) continue;
-
-                       x[0] = c->GetX()-x0;
-                       x[1] = x[0]*x[0];
-                       fitter.AddPoint(x, c->GetY(), c->GetSigmaY2());
-                       ncl++;
-               }
-               fitter.Eval();
-               Double_t dy = yc - fitter.GetParameter(0) -fitter.GetParameter(1) * (xc-x0) - fitter.GetParameter(2)* (xc-x0)*(xc-x0); 
-       
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "ResidualsClustersParametrisation"
-                       << "dy="   << dy
-                       << "\n";
-       }
+  
+  // store cluster positions
+  Double_t x0 = tracklet->GetX0();
+  AliTRDcluster *c = 0x0;
+  
+  Double_t x[2]; Int_t ncl, mcl, jc;
+  TLinearFitter fitter(3, "hyp2");
+  for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
+    if(!(c = tracklet->GetClusters(ic))) continue;
+    Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
+    
+    jc = ic; ncl = 0; mcl=0; fitter.ClearPoints();
+    while(ncl<6){
+      // update index
+      mcl++;
+      jc = ic + ((mcl&1)?-1:1)*(mcl>>1);
+
+      if(jc<0 || jc>=35) continue;
+      if(!(c = tracklet->GetClusters(jc))) continue;
+
+      x[0] = c->GetX()-x0;
+      x[1] = x[0]*x[0];
+      fitter.AddPoint(x, c->GetY(), c->GetSigmaY2());
+      ncl++;
+    }
+    fitter.Eval();
+    Double_t dy = yc - fitter.GetParameter(0) -fitter.GetParameter(1) * (xc-x0) - fitter.GetParameter(2)* (xc-x0)*(xc-x0); 
+  
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "ResidualsClustersParametrisation"
+      << "dy="   << dy
+      << "\n";
+  }
 }
 
 
@@ -237,83 +237,82 @@ void AliTRDtrackerDebug::ResidualsClustersParametrisation(const AliTRDseedV1 *tr
 void AliTRDtrackerDebug::ResidualsTrackletsTrack() const
 {
 // Calculates distances from tracklets to the TRD track.
-       
-       if(fTrack->GetNumberOfTracklets() < 6) return;
-
-       // build a working copy of the tracklets attached to the track 
-       // and initialize working variables fX, fY and fZ
-       AliTRDcluster *c = 0x0;
-       AliTRDseedV1 tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
-       const AliTRDseedV1 *ctracklet = 0x0;
-       for(int ip = 0; ip<6; ip++){
-               if(!(ctracklet = fTrack->GetTracklet(ip))) continue;
-               tracklet[ip] = (*ctracklet); 
-               Double_t x0 = tracklet[ip].GetX0();
-               for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
-                       if(!(c = tracklet[ip].GetClusters(ic))) continue;
-                       Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
-                       tracklet[ip].SetX(ic, xc-x0);
-                       tracklet[ip].SetY(ic, yc);
-                       tracklet[ip].SetZ(ic, zc);
-               }
-       }
-       
-       // Do a Rieman fit (with tilt correction) for all tracklets 
-       // except the one which is tested. 
-       // (Based on AliTRDseed::IsOK() return false)
-       for(int ip=0; ip<6; ip++){
-               // reset tracklet to be tested
-               Double_t x0 = tracklet[ip].GetX0();
-               new(&tracklet[ip]) AliTRDseedV1();
-               tracklet[ip].SetX0(x0);
-
-               // fit Rieman with tilt correction
-               AliTRDseedV1::FitRiemanTilt(&tracklet[0], kTRUE);
-
-               // make a copy of the fit result
-               Double_t 
-                       y0   = tracklet[ip].GetYref(0),
-                       dydx = tracklet[ip].GetYref(1),
-                       z0   = tracklet[ip].GetZref(0),
-                       dzdx = tracklet[ip].GetZref(1);
-
-               // restore tracklet
-               tracklet[ip] = (*fTrack->GetTracklet(ip)); 
-               for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
-                       if(!(c = tracklet[ip].GetClusters(ic))) continue;
-                       Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
-                       tracklet[ip].SetX(ic, xc-x0);
-                       tracklet[ip].SetY(ic, yc);
-                       tracklet[ip].SetZ(ic, zc);
-               }               
-               
-               // fit clusters
-               AliTRDseedV1 ts(tracklet[ip]);
-               ts.SetYref(0, y0); ts.SetYref(1, dydx);
-               ts.SetZref(0, z0); ts.SetZref(1, dzdx);
-               ts.Update();
-
-               // save results for plotting
-               Int_t plane   = tracklet[ip].GetPlane();
-               Double_t dy   = tracklet[ip].GetYfit(0) - ts.GetYfit(0);
-               Double_t tgy  = tracklet[ip].GetYfit(1);
-               Double_t dtgy = (tracklet[ip].GetYfit(1) - ts.GetYfit(1))/(1. + tracklet[ip].GetYfit(1) * ts.GetYfit(1));
-               Double_t dz   = tracklet[ip].GetZfit(0) - ts.GetZfit(0);
-               Double_t tgz  = tracklet[ip].GetZfit(1);
-               Double_t dtgz = (tracklet[ip].GetZfit(1) - ts.GetZfit(1))/(1. + tracklet[ip].GetZfit(1) * ts.GetZfit(1));
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "ResidualsTrackletsTrack"
-                       << "ref.="   << &tracklet[ip]
-                       << "fit.="   << &ts
-                       << "plane="  << plane
-                       << "dy="     << dy
-                       << "tgy="    << tgy
-                       << "dtgy="   << dtgy
-                       << "dz="     << dz
-                       << "tgz="    << tgz
-                       << "dtgz="   << dtgz
-                       << "\n";
-       }
+  
+  if(fTrack->GetNumberOfTracklets() < 6) return;
+
+  // build a working copy of the tracklets attached to the track 
+  // and initialize working variables fX, fY and fZ
+  AliTRDseedV1 tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
+  const AliTRDseedV1 *ctracklet = 0x0;
+  for(int ip = 0; ip<6; ip++){
+    if(!(ctracklet = fTrack->GetTracklet(ip))) continue;
+    tracklet[ip] = (*ctracklet); 
+//             Double_t x0 = tracklet[ip].GetX0();
+//             for(int ic=0; ic<AliTRDseedV1:knTimebins; ic++){
+//                     if(!(c = tracklet[ip].GetClusters(ic))) continue;
+//                     Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
+//                     tracklet[ip].SetX(ic, xc-x0);
+//                     tracklet[ip].SetY(ic, yc);
+//                     tracklet[ip].SetZ(ic, zc);
+//             }
+  }
+  
+  // Do a Rieman fit (with tilt correction) for all tracklets 
+  // except the one which is tested. 
+  // (Based on AliTRDseed::IsOK() return false)
+  for(int ip=0; ip<6; ip++){
+    // reset tracklet to be tested
+    Double_t x0 = tracklet[ip].GetX0();
+    new(&tracklet[ip]) AliTRDseedV1();
+    tracklet[ip].SetX0(x0);
+
+    // fit Rieman with tilt correction
+    AliTRDtrackerV1::FitRiemanTilt(0x0, &tracklet[0], kTRUE);
+
+    // make a copy of the fit result
+    Double_t 
+      y0   = tracklet[ip].GetYref(0),
+      dydx = tracklet[ip].GetYref(1),
+      z0   = tracklet[ip].GetZref(0),
+      dzdx = tracklet[ip].GetZref(1);
+
+    // restore tracklet
+    tracklet[ip] = (*fTrack->GetTracklet(ip)); 
+//             for(int ic=0; ic<AliTRDseedV1:knTimebins; ic++){
+//                     if(!(c = tracklet[ip].GetClusters(ic))) continue;
+//                     Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
+//                     tracklet[ip].SetX(ic, xc-x0);
+//                     tracklet[ip].SetY(ic, yc);
+//                     tracklet[ip].SetZ(ic, zc);
+//             }               
+    
+    // fit clusters
+    AliTRDseedV1 ts(tracklet[ip]);
+    ts.SetYref(0, y0); ts.SetYref(1, dydx);
+    ts.SetZref(0, z0); ts.SetZref(1, dzdx);
+    ts.FitMI();
+
+    // save results for plotting
+    Int_t plane   = tracklet[ip].GetPlane();
+    Double_t dy   = tracklet[ip].GetYfit(0) - ts.GetYfit(0);
+    Double_t tgy  = tracklet[ip].GetYfit(1);
+    Double_t dtgy = (tracklet[ip].GetYfit(1) - ts.GetYfit(1))/(1. + tracklet[ip].GetYfit(1) * ts.GetYfit(1));
+    Double_t dz   = tracklet[ip].GetZfit(0) - ts.GetZfit(0);
+    Double_t tgz  = tracklet[ip].GetZfit(1);
+    Double_t dtgz = (tracklet[ip].GetZfit(1) - ts.GetZfit(1))/(1. + tracklet[ip].GetZfit(1) * ts.GetZfit(1));
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "ResidualsTrackletsTrack"
+      << "ref.="   << &tracklet[ip]
+      << "fit.="   << &ts
+      << "plane="  << plane
+      << "dy="     << dy
+      << "tgy="    << tgy
+      << "dtgy="   << dtgy
+      << "dz="     << dz
+      << "tgz="    << tgz
+      << "dtgz="   << dtgz
+      << "\n";
+  }
 }
 
 //____________________________________________________
@@ -330,46 +329,46 @@ void AliTRDtrackerDebug::AnalyseFindable(Char_t *treename){
 // A new TTree containing the number of findable tracklets and the number of clusters
 // attached to the full track is stored to disk
 //
-       // Link the File
-       TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
-       fTree = (TTree *)(debfile->Get(treename));
-       if(!fTree){
-               AliError(Form("Tree %s not found in file TRDdebug.root. Abborting!", treename));
-               debfile->Close();
-               return;
-       }
-       
-       AliTRDseedV1 *tracklets[kNPlanes];
-       for(Int_t iPlane = 0; iPlane < AliTRDtrackerV1::kNPlanes; iPlane++)
-               tracklets[iPlane] = 0x0;
-       for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++)
-               fTree->SetBranchAddress(Form("S%d.", iPlane), &tracklets[iPlane]);
-       fTree->SetBranchAddress("EventNumber", &fgEventNumber);
-       fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
-       
-       Int_t findable = 0, nClusters = 0;
-       Int_t nEntries = fTree->GetEntriesFast();
-       for(Int_t iEntry = 0; iEntry < nEntries; iEntry++){
-               printf("Entry %d\n", iEntry);
-               fTree->GetEntry(iEntry);
-               findable = 0;
-               nClusters = 0;
-               // Calculate Findable
-               for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++){
-                       if (TMath::Abs(tracklets[iPlane]->GetYref(0) / tracklets[iPlane]->GetX0()) < 0.15) findable++;
-                       if (!tracklets[iPlane]->IsOK()) continue;
-                       nClusters += tracklets[iPlane]->GetN2();
-               }
-               
-               // Fill Histogramms
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "AnalyseFindable"
-                       << "EventNumber="               << fgEventNumber
-                       << "CandidateNumber="   << fgCandidateNumber
-                       << "Findable="                  << findable
-                       << "NClusters="                 << nClusters
-                       << "\n";
-       }
+  // Link the File
+  TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
+  fTree = (TTree *)(debfile->Get(treename));
+  if(!fTree){
+    AliError(Form("Tree %s not found in file TRDdebug.root. Abborting!", treename));
+    debfile->Close();
+    return;
+  }
+  
+  AliTRDseedV1 *tracklets[kNPlanes];
+  for(Int_t iPlane = 0; iPlane < AliTRDtrackerV1::kNPlanes; iPlane++)
+    tracklets[iPlane] = 0x0;
+  for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++)
+    fTree->SetBranchAddress(Form("S%d.", iPlane), &tracklets[iPlane]);
+  fTree->SetBranchAddress("EventNumber", &fgEventNumber);
+  fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
+  
+  Int_t findable = 0, nClusters = 0;
+  Int_t nEntries = fTree->GetEntriesFast();
+  for(Int_t iEntry = 0; iEntry < nEntries; iEntry++){
+    printf("Entry %d\n", iEntry);
+    fTree->GetEntry(iEntry);
+    findable = 0;
+    nClusters = 0;
+    // Calculate Findable
+    for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++){
+      if (TMath::Abs(tracklets[iPlane]->GetYref(0) / tracklets[iPlane]->GetX0()) < 0.15) findable++;
+      if (!tracklets[iPlane]->IsOK()) continue;
+      nClusters += tracklets[iPlane]->GetN2();
+    }
+    
+    // Fill Histogramms
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "AnalyseFindable"
+      << "EventNumber="                << fgEventNumber
+      << "CandidateNumber="    << fgCandidateNumber
+      << "Findable="                   << findable
+      << "NClusters="                  << nClusters
+      << "\n";
+  }
 }
 //____________________________________________________
 void AliTRDtrackerDebug::AnalyseTiltedRiemanFit(){
@@ -382,35 +381,35 @@ void AliTRDtrackerDebug::AnalyseTiltedRiemanFit(){
 //
 // TODO: Match everything with Init and Process
 //
-       TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
-       fTree = (TTree *)(debfile->Get("MakeSeeds2"));
-       if(!fTree) return;
-       Int_t nEntries = fTree->GetEntries();
-       TLinearFitter *tiltedRiemanFitter = 0x0;
-       fTree->SetBranchAddress("FitterT.", &tiltedRiemanFitter);
-       fTree->SetBranchAddress("EventNumber", &fgEventNumber);
-       fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
-       for(Int_t entry = 0; entry < nEntries; entry++){
-               fTree->GetEntry(entry);
-               Double_t a = tiltedRiemanFitter->GetParameter(0);
-               Double_t b = tiltedRiemanFitter->GetParameter(1);
-               Double_t c = tiltedRiemanFitter->GetParameter(2);
-               Double_t offset = tiltedRiemanFitter->GetParameter(3);
-               Double_t slope  = tiltedRiemanFitter->GetParameter(4);
-               Float_t radius = GetTrackRadius(a, b, c);
-               Float_t curvature = GetTrackCurvature(a, b, c);
-               Float_t dca = GetDCA(a, b, c);
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "AnalyseTiltedRiemanFit"
-               << "EventNumber="               << fgEventNumber
-               << "CandidateNumber=" << fgCandidateNumber
-               << "Radius="                                    << radius
-               << "Curvature="                         << curvature
-               << "DCA="                                                       << dca
-               << "Offset="                                    << offset
-               << "Slope="                                             << slope
-               << "\n";
-       }
+  TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
+  fTree = (TTree *)(debfile->Get("MakeSeeds2"));
+  if(!fTree) return;
+  Int_t nEntries = fTree->GetEntries();
+  TLinearFitter *tiltedRiemanFitter = 0x0;
+  fTree->SetBranchAddress("FitterT.", &tiltedRiemanFitter);
+  fTree->SetBranchAddress("EventNumber", &fgEventNumber);
+  fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
+  for(Int_t entry = 0; entry < nEntries; entry++){
+    fTree->GetEntry(entry);
+    Double_t a = tiltedRiemanFitter->GetParameter(0);
+    Double_t b = tiltedRiemanFitter->GetParameter(1);
+    Double_t c = tiltedRiemanFitter->GetParameter(2);
+    Double_t offset = tiltedRiemanFitter->GetParameter(3);
+    Double_t slope  = tiltedRiemanFitter->GetParameter(4);
+    Float_t radius = GetTrackRadius(a, b, c);
+    Float_t curvature = GetTrackCurvature(a, b, c);
+    Float_t dca = GetDCA(a, b, c);
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "AnalyseTiltedRiemanFit"
+    << "EventNumber="          << fgEventNumber
+    << "CandidateNumber=" << fgCandidateNumber
+    << "Radius="                                       << radius
+    << "Curvature="                            << curvature
+    << "DCA="                                                  << dca
+    << "Offset="                                       << offset
+    << "Slope="                                                << slope
+    << "\n";
+  }
 }
 
 //____________________________________________________
@@ -421,10 +420,10 @@ Float_t AliTRDtrackerDebug::GetTrackRadius(Float_t a, Float_t b, Float_t c) cons
 // Parameters: The three parameters from the Rieman fit
 // Output:     The track radius
 //
-       Float_t radius = 0;
-       if(1.0 + b*b - c*a > 0.0)
-               radius = TMath::Sqrt(1.0 + b*b - c*a )/a;
-       return radius;
+  Float_t radius = 0;
+  if(1.0 + b*b - c*a > 0.0)
+    radius = TMath::Sqrt(1.0 + b*b - c*a )/a;
+  return radius;
 }
 
 //____________________________________________________
@@ -435,10 +434,10 @@ Float_t AliTRDtrackerDebug::GetTrackCurvature(Float_t a, Float_t b, Float_t c) c
 // Parameters: the three parameters from the tilted Rieman fitter
 // Output:             the full track curvature
 //
-       Float_t curvature =  1.0 + b*b - c*a;
-       if (curvature > 0.0) 
-               curvature  =  a / TMath::Sqrt(curvature);
-       return curvature;
+  Float_t curvature =  1.0 + b*b - c*a;
+  if (curvature > 0.0) 
+    curvature  =  a / TMath::Sqrt(curvature);
+  return curvature;
 }
 
 //____________________________________________________
@@ -450,53 +449,53 @@ Float_t AliTRDtrackerDebug::GetDCA(Float_t a, Float_t b, Float_t c) const {
 // Parameters: the three parameters from the tilted Rieman fitter
 // Output:     the Distance to Closest Approach
 //
-       Float_t dca  =  0.0;
-       if (1.0 + b*b - c*a > 0.0) {
-               dca = -c / (TMath::Sqrt(1.0 + b*b - c*a) + TMath::Sqrt(1.0 + b*b));
-       }
-       return dca;
+  Float_t dca  =  0.0;
+  if (1.0 + b*b - c*a > 0.0) {
+    dca = -c / (TMath::Sqrt(1.0 + b*b - c*a) + TMath::Sqrt(1.0 + b*b));
+  }
+  return dca;
 }
 
 //____________________________________________________
 void AliTRDtrackerDebug::AnalyseMinMax()
 {
 //
-       TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
-       if(!debfile){
-               AliError("File TRD.TrackerDebug.root not found!");
-               return; 
-       }
-       fTree = (TTree *)(debfile->Get("MakeSeeds0"));
-       if(!fTree){
-               AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
-               return;
-       }
-       AliTRDseedV1 *cseed[4] = {0x0, 0x0, 0x0, 0x0};
-       AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
-       for(Int_t il = 0; il < 4; il++){
-               fTree->SetBranchAddress(Form("Seed%d.", il),    &cseed[il]);
-               fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
-       }
-       fTree->SetBranchAddress("CandidateNumber",      &fgCandidateNumber);
-       fTree->SetBranchAddress("EventNumber",  &fgEventNumber);
-       Int_t entries = fTree->GetEntries();
-       for(Int_t ientry = 0; ientry < entries; ientry++){
-               fTree->GetEntry(ientry);
-               Float_t minmax[2] = { -100.0,  100.0 };
-               for (Int_t iLayer = 0; iLayer < 4; iLayer++) {
-                       Float_t max = c[iLayer]->GetZ() + cseed[iLayer]->GetPadLength() * 0.5 + 1.0 - cseed[iLayer]->GetZref(0);
-                       if (max < minmax[1]) minmax[1] = max;
-                       Float_t min = c[iLayer]->GetZ()-cseed[iLayer]->GetPadLength() * 0.5 - 1.0 - cseed[iLayer]->GetZref(0);
-                       if (min > minmax[0]) minmax[0] = min;
-               }
-               TTreeSRedirector &cstreamer = *fOutputStreamer;
-               cstreamer << "AnalyseMinMaxLayer"
-               << "EventNumber="                               << fgEventNumber
-               << "CandidateNumber="           << fgCandidateNumber
-               << "Min="                                                               << minmax[0]
-               << "Max="                                                               << minmax[1]
-               << "\n";
-       }
+  TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
+  if(!debfile){
+    AliError("File TRD.TrackerDebug.root not found!");
+    return; 
+  }
+  fTree = (TTree *)(debfile->Get("MakeSeeds0"));
+  if(!fTree){
+    AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
+    return;
+  }
+  AliTRDseedV1 *cseed[4] = {0x0, 0x0, 0x0, 0x0};
+  AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
+  for(Int_t il = 0; il < 4; il++){
+    fTree->SetBranchAddress(Form("Seed%d.", il),       &cseed[il]);
+    fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
+  }
+  fTree->SetBranchAddress("CandidateNumber",   &fgCandidateNumber);
+  fTree->SetBranchAddress("EventNumber",       &fgEventNumber);
+  Int_t entries = fTree->GetEntries();
+  for(Int_t ientry = 0; ientry < entries; ientry++){
+    fTree->GetEntry(ientry);
+    Float_t minmax[2] = { -100.0,  100.0 };
+    for (Int_t iLayer = 0; iLayer < 4; iLayer++) {
+      Float_t max = c[iLayer]->GetZ() + cseed[iLayer]->GetPadLength() * 0.5 + 1.0 - cseed[iLayer]->GetZref(0);
+      if (max < minmax[1]) minmax[1] = max;
+      Float_t min = c[iLayer]->GetZ()-cseed[iLayer]->GetPadLength() * 0.5 - 1.0 - cseed[iLayer]->GetZref(0);
+      if (min > minmax[0]) minmax[0] = min;
+    }
+    TTreeSRedirector &cstreamer = *fOutputStreamer;
+    cstreamer << "AnalyseMinMaxLayer"
+    << "EventNumber="                          << fgEventNumber
+    << "CandidateNumber="              << fgCandidateNumber
+    << "Min="                                                          << minmax[0]
+    << "Max="                                                          << minmax[1]
+    << "\n";
+  }
 }
 
 //____________________________________________________
@@ -509,95 +508,95 @@ TCanvas* AliTRDtrackerDebug::PlotSeedingConfiguration(const Char_t *direction, I
 //                             -Event Nr
 //             -Candidate that has to be plotted
 //
-       const Float_t kxmin = 280;
-       const Float_t kxmax = 380;
-       const Float_t kxdelta = (kxmax - kxmin)/1000;
-       
-       if((strcmp(direction, "y") != 0) && (strcmp(direction, "z") != 0)){
-               AliError(Form("Direction %s does not exist. Abborting!", direction));
-               return 0x0;
-       }
-
-       TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
-       if(!debfile){
-               AliError("File TRD.TrackerDebug.root not found!");
-               return 0x0; 
-       }
-       fTree = (TTree *)(debfile->Get("MakeSeeds0"));
-       if(!fTree){
-               AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
-               return 0x0;
-       }
-       
-       TGraph *seedcl = new TGraph(4);
-       TGraph *seedRef = new TGraph(4);
-       TGraph *riemanFit = new TGraph(1000);
-       seedcl->SetMarkerStyle(20);
-       seedcl->SetMarkerColor(kRed);
-       seedRef->SetMarkerStyle(2);
-
-       AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
-       AliRieman *rim = 0x0;
-       Bool_t found = kFALSE;
-       for(Int_t il = 0; il < 4; il++) fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
-       fTree->SetBranchAddress("EventNumber", &fgEventNumber);
-       fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
-       fTree->SetBranchAddress("RiemanFitter.", &rim);
-       Int_t entries = fTree->GetEntries();
-       for(Int_t entry = 0; entry < entries; entry++){
-               fTree->GetEntry(entry);
-               if(fgEventNumber < event) continue;
-               if(fgEventNumber > event) break;
-               // EventNumber matching: Do the same for the candidate number
-               if(fgCandidateNumber < candidate) continue;
-               if(fgCandidateNumber > candidate) break;
-               found = kTRUE;
-               Int_t nPoints = 0;
-               for(Int_t il = 0; il < 4; il++){
-                       Float_t cluster = 0.0, reference = 0.0;
-                       if(!strcmp(direction, "y")){
-                               cluster = c[il]->GetY();
-                               reference = rim->GetYat(c[il]->GetX());
-                       }
-                       else{
-                               cluster = c[il]->GetZ();
-                               reference = rim->GetZat(c[il]->GetX());
-                       }
-                       seedcl->SetPoint(nPoints, cluster, c[il]->GetX());
-                       seedRef->SetPoint(nPoints, reference , c[il]->GetX());
-                       nPoints++;
-               }
-               // evaluate the fitter Function numerically
-               nPoints = 0;
-               for(Int_t ipt = 0; ipt < 1000; ipt++){
-                       Float_t x = kxmin + ipt * kxdelta;
-                       Float_t point = 0.0;
-                       if(!strcmp(direction, "y"))
-                               point = rim->GetYat(x);
-                       else
-                               point = rim->GetZat(x);
-                       riemanFit->SetPoint(nPoints++, point, x);
-               }
-               // We reached the End: break
-               break;
-       }
-       if(found){
-               seedcl->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               seedRef->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               riemanFit->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               TCanvas *c1 = new TCanvas();
-               seedcl->Draw("ap");
-               seedRef->Draw("psame");
-               riemanFit->Draw("lpsame");
-               return c1;
-       }
-       else{
-               AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
-               delete seedcl;
-               delete seedRef;
-               delete riemanFit;
-               return 0x0;
-       }
+  const Float_t kxmin = 280;
+  const Float_t kxmax = 380;
+  const Float_t kxdelta = (kxmax - kxmin)/1000;
+  
+  if((strcmp(direction, "y") != 0) && (strcmp(direction, "z") != 0)){
+    AliError(Form("Direction %s does not exist. Abborting!", direction));
+    return 0x0;
+  }
+
+  TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
+  if(!debfile){
+    AliError("File TRD.TrackerDebug.root not found!");
+    return 0x0; 
+  }
+  fTree = (TTree *)(debfile->Get("MakeSeeds0"));
+  if(!fTree){
+    AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
+    return 0x0;
+  }
+  
+  TGraph *seedcl = new TGraph(4);
+  TGraph *seedRef = new TGraph(4);
+  TGraph *riemanFit = new TGraph(1000);
+  seedcl->SetMarkerStyle(20);
+  seedcl->SetMarkerColor(kRed);
+  seedRef->SetMarkerStyle(2);
+
+  AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
+  AliRieman *rim = 0x0;
+  Bool_t found = kFALSE;
+  for(Int_t il = 0; il < 4; il++) fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
+  fTree->SetBranchAddress("EventNumber", &fgEventNumber);
+  fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
+  fTree->SetBranchAddress("RiemanFitter.", &rim);
+  Int_t entries = fTree->GetEntries();
+  for(Int_t entry = 0; entry < entries; entry++){
+    fTree->GetEntry(entry);
+    if(fgEventNumber < event) continue;
+    if(fgEventNumber > event) break;
+    // EventNumber matching: Do the same for the candidate number
+    if(fgCandidateNumber < candidate) continue;
+    if(fgCandidateNumber > candidate) break;
+    found = kTRUE;
+    Int_t nPoints = 0;
+    for(Int_t il = 0; il < 4; il++){
+      Float_t cluster = 0.0, reference = 0.0;
+      if(!strcmp(direction, "y")){
+        cluster = c[il]->GetY();
+        reference = rim->GetYat(c[il]->GetX());
+      }
+      else{
+        cluster = c[il]->GetZ();
+        reference = rim->GetZat(c[il]->GetX());
+      }
+      seedcl->SetPoint(nPoints, cluster, c[il]->GetX());
+      seedRef->SetPoint(nPoints, reference , c[il]->GetX());
+      nPoints++;
+    }
+    // evaluate the fitter Function numerically
+    nPoints = 0;
+    for(Int_t ipt = 0; ipt < 1000; ipt++){
+      Float_t x = kxmin + ipt * kxdelta;
+      Float_t point = 0.0;
+      if(!strcmp(direction, "y"))
+        point = rim->GetYat(x);
+      else
+        point = rim->GetZat(x);
+      riemanFit->SetPoint(nPoints++, point, x);
+    }
+    // We reached the End: break
+    break;
+  }
+  if(found){
+    seedcl->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    seedRef->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    riemanFit->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    TCanvas *c1 = new TCanvas();
+    seedcl->Draw("ap");
+    seedRef->Draw("psame");
+    riemanFit->Draw("lpsame");
+    return c1;
+  }
+  else{
+    AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
+    delete seedcl;
+    delete seedRef;
+    delete riemanFit;
+    return 0x0;
+  }
 }
 
 //____________________________________________________
@@ -612,170 +611,170 @@ TCanvas* AliTRDtrackerDebug::PlotFullTrackFit(Int_t event, Int_t candidate, Int_
 //                        -direction (default: y)
 // Output:     -TCanvas (containing the Picture);
 //
-       const Float_t kxmin = 280;
-       const Float_t kxmax = 380;
-       const Float_t kxdelta = (kxmax - kxmin)/1000;
-       
-       if(strcmp(direction, "y") && strcmp(direction, "z")){
-               AliError(Form("Direction %s does not exist. Abborting!", direction));
-               return 0x0;
-       }
-
-       TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
-       if(!debfile){
-               AliError("File TRD.TrackerDebug.root not found.");
-               return 0x0;
-       }
-       TString *treename = 0x0;
-       if(iteration > -1)
-               treename = new TString("ImproveSeedQuality");
-       else
-               treename = new TString("MakeSeeds1");
-       fTree = (TTree *)(debfile->Get(treename->Data()));
-       if(!fTree){
-               AliError(Form("Tree %s not found in File TRD.TrackerDebug.root.", treename->Data()));
-               delete treename;
-               return 0x0;
-       }
-       delete treename;
-
-       TGraph *fitfun = new TGraph(1000);
-       // Prepare containers
-       Float_t x0[AliTRDtrackerV1::kNPlanes],
-                       refP[AliTRDtrackerV1::kNPlanes],
-                       clx[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins],
-                       clp[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins];
-       Int_t nLayers = 0, ncls = 0;
-       
-       TLinearFitter *fitter = 0x0;
-       AliTRDseedV1 *tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
-       for(Int_t iLayer = 0; iLayer < 6; iLayer++)
-               fTree->SetBranchAddress(Form("S%d.", iLayer), &tracklet[iLayer]);
-       fTree->SetBranchAddress("FitterT.", &fitter);
-       fTree->SetBranchAddress("EventNumber", &fgEventNumber);
-       fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
-       
-       Int_t nEntries = fTree->GetEntriesFast();
-       Bool_t found = kFALSE;
-       for(Int_t entry = 0; entry < nEntries; entry++){
-               fTree->GetEntry(entry);
-               if(fgEventNumber < event) continue;
-               if(fgEventNumber > event) break;
-               // EventNumber matching: Do the same for the candidate number
-               if(fgCandidateNumber < candidate) continue;
-               if(fgCandidateNumber > candidate) break;
-               found = kTRUE;
-               
-               for(Int_t iLayer = 0; iLayer < 6; iLayer++){
-                       if(!tracklet[iLayer]->IsOK()) continue;
-                       x0[nLayers] = tracklet[iLayer]->GetX0();
-                       if(!strcmp(direction, "y"))
-                               refP[nLayers] = tracklet[iLayer]->GetYref(0);
-                       else
-                               refP[nLayers] = tracklet[iLayer]->GetZref(0);
-                       nLayers++;
-                       for(Int_t itb = 0; itb < 30; itb++){
-                               if(!tracklet[iLayer]->IsUsable(itb)) continue;
-                               AliTRDcluster *cl = tracklet[iLayer]->GetClusters(itb);
-                               if(!strcmp(direction, "y"))
-                                       clp[ncls] = cl->GetY();
-                               else
-                                       clp[ncls] = cl->GetZ();
-                               clx[ncls] = cl->GetX();
-                               ncls++;
-                       }
-               }
-               // Add function derived by the tilted Rieman fit (Defined by the curvature)
-               Int_t nPoints = 0;
-               if(!strcmp(direction, "y")){
-                       Double_t a = fitter->GetParameter(0);
-                       Double_t b = fitter->GetParameter(1);
-                       Double_t c = fitter->GetParameter(2);
-                       Double_t curvature =  1.0 + b*b - c*a;
-                       if (curvature > 0.0) {
-                               curvature  =  a / TMath::Sqrt(curvature);
-                       }
-                       // Numerical evaluation of the function:
-                       for(Int_t ipt = 0; ipt < 1000; ipt++){
-                               Float_t x = kxmin + ipt * kxdelta;
-                               Double_t res = (x * a + b);                                                             // = (x - x0)/y0
-                               res *= res;
-                               res  = 1.0 - c * a + b * b - res;                                       // = (R^2 - (x - x0)^2)/y0^2
-                               Double_t y = 0.;
-                               if (res >= 0) {
-                                       res = TMath::Sqrt(res);
-                                       y    = (1.0 - res) / a;
-                               }
-                               fitfun->SetPoint(nPoints++, y, x);
-                       }
-               }
-               else{
-                       Double_t offset = fitter->GetParameter(3);
-                       Double_t slope  = fitter->GetParameter(4);       
-                       // calculate the reference x (defined as medium between layer 2 and layer 3)
-                       // same procedure as in the tracker code
-                       Float_t medx = 0, xref = 0;
-                       Int_t startIndex = 5, nDistances = 0;
-                       for(Int_t iLayer = 5; iLayer > 0; iLayer--){
-                               if(tracklet[iLayer]->IsOK() && tracklet[iLayer - 1]->IsOK()){
-                                       medx += tracklet[iLayer]->GetX0() - tracklet[iLayer - 1]->GetX0();
-                                       startIndex = iLayer - 1;
-                                       nDistances++;
-                               }
-                       }
-                       if(nDistances){
-                               medx /= nDistances;
-                       }
-                       else{
-                               Float_t xpos[2];        memset(xpos, 0, sizeof(Float_t) * 2);
-                               Int_t ien = 0, idiff = 0;
-                               for(Int_t iLayer = 5; iLayer > 0; iLayer--){
-                                       if(tracklet[iLayer]->IsOK()){
-                                               xpos[ien++] = tracklet[iLayer]->GetX0();
-                                               startIndex = iLayer;
-                                       }
-                                       if(ien)
-                                               idiff++;
-                                       if(ien >=2)
-                                               break;
-                               }
-                               medx = (xpos[0] - xpos[1])/idiff;
-                       }
-                       xref = tracklet[startIndex]->GetX0() + medx * (2.5 - startIndex) - 0.5 * (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
-
-                       for(Int_t ipt = 0; ipt < 1000; ipt++){
-                               Float_t x = kxmin + ipt * kxdelta;
-                               Float_t z = offset + slope * (x - xref);
-                               fitfun->SetPoint(nPoints++, z, x);
-                       }
-               }
-               break;
-       }
-       if(found){
-               TGraph *trGraph         = new TGraph(ncls);
-               TGraph *refPoints       = new TGraph(nLayers);
-               trGraph->SetMarkerStyle(20);
-               trGraph->SetMarkerColor(kRed);
-               refPoints->SetMarkerStyle(21);
-               refPoints->SetMarkerColor(kBlue);
-               // fill the graphs
-               for(Int_t iLayer = 0; iLayer < nLayers; iLayer++)
-                       refPoints->SetPoint(iLayer, refP[iLayer], x0[iLayer]);
-               for(Int_t icls = 0; icls < ncls; icls++)
-                       trGraph->SetPoint(icls, clp[icls], clx[icls]);
-               TCanvas *c1 = new TCanvas();
-               trGraph->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               refPoints->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               fitfun->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
-               trGraph->Draw("ap");
-               refPoints->Draw("psame");
-               fitfun->Draw("lpsame");
-               return c1;
-       }
-       else{
-               AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
-               delete fitfun;
-               return 0x0;
-       }
+  const Float_t kxmin = 280;
+  const Float_t kxmax = 380;
+  const Float_t kxdelta = (kxmax - kxmin)/1000;
+  
+  if(strcmp(direction, "y") && strcmp(direction, "z")){
+    AliError(Form("Direction %s does not exist. Abborting!", direction));
+    return 0x0;
+  }
+
+  TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
+  if(!debfile){
+    AliError("File TRD.TrackerDebug.root not found.");
+    return 0x0;
+  }
+  TString *treename = 0x0;
+  if(iteration > -1)
+    treename = new TString("ImproveSeedQuality");
+  else
+    treename = new TString("MakeSeeds1");
+  fTree = (TTree *)(debfile->Get(treename->Data()));
+  if(!fTree){
+    AliError(Form("Tree %s not found in File TRD.TrackerDebug.root.", treename->Data()));
+    delete treename;
+    return 0x0;
+  }
+  delete treename;
+
+  TGraph *fitfun = new TGraph(1000);
+  // Prepare containers
+  Float_t x0[AliTRDtrackerV1::kNPlanes],
+      refP[AliTRDtrackerV1::kNPlanes],
+      clx[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins],
+      clp[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins];
+  Int_t nLayers = 0, ncls = 0;
+  
+  TLinearFitter *fitter = 0x0;
+  AliTRDseedV1 *tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
+  for(Int_t iLayer = 0; iLayer < 6; iLayer++)
+    fTree->SetBranchAddress(Form("S%d.", iLayer), &tracklet[iLayer]);
+  fTree->SetBranchAddress("FitterT.", &fitter);
+  fTree->SetBranchAddress("EventNumber", &fgEventNumber);
+  fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
+  
+  Int_t nEntries = fTree->GetEntriesFast();
+  Bool_t found = kFALSE;
+  for(Int_t entry = 0; entry < nEntries; entry++){
+    fTree->GetEntry(entry);
+    if(fgEventNumber < event) continue;
+    if(fgEventNumber > event) break;
+    // EventNumber matching: Do the same for the candidate number
+    if(fgCandidateNumber < candidate) continue;
+    if(fgCandidateNumber > candidate) break;
+    found = kTRUE;
+    
+    for(Int_t iLayer = 0; iLayer < 6; iLayer++){
+      if(!tracklet[iLayer]->IsOK()) continue;
+      x0[nLayers] = tracklet[iLayer]->GetX0();
+      if(!strcmp(direction, "y"))
+        refP[nLayers] = tracklet[iLayer]->GetYref(0);
+      else
+        refP[nLayers] = tracklet[iLayer]->GetZref(0);
+      nLayers++;
+      for(Int_t itb = 0; itb < 30; itb++){
+        if(!tracklet[iLayer]->IsUsable(itb)) continue;
+        AliTRDcluster *cl = tracklet[iLayer]->GetClusters(itb);
+        if(!strcmp(direction, "y"))
+          clp[ncls] = cl->GetY();
+        else
+          clp[ncls] = cl->GetZ();
+        clx[ncls] = cl->GetX();
+        ncls++;
+      }
+    }
+    // Add function derived by the tilted Rieman fit (Defined by the curvature)
+    Int_t nPoints = 0;
+    if(!strcmp(direction, "y")){
+      Double_t a = fitter->GetParameter(0);
+      Double_t b = fitter->GetParameter(1);
+      Double_t c = fitter->GetParameter(2);
+      Double_t curvature =  1.0 + b*b - c*a;
+      if (curvature > 0.0) {
+        curvature  =  a / TMath::Sqrt(curvature);
+      }
+      // Numerical evaluation of the function:
+      for(Int_t ipt = 0; ipt < 1000; ipt++){
+        Float_t x = kxmin + ipt * kxdelta;
+        Double_t res = (x * a + b);                                                            // = (x - x0)/y0
+        res *= res;
+        res  = 1.0 - c * a + b * b - res;                                      // = (R^2 - (x - x0)^2)/y0^2
+        Double_t y = 0.;
+        if (res >= 0) {
+          res = TMath::Sqrt(res);
+          y    = (1.0 - res) / a;
+        }
+        fitfun->SetPoint(nPoints++, y, x);
+      }
+    }
+    else{
+      Double_t offset  = fitter->GetParameter(3);
+      Double_t slope   = fitter->GetParameter(4);       
+      // calculate the reference x (defined as medium between layer 2 and layer 3)
+      // same procedure as in the tracker code
+      Float_t medx = 0, xref = 0;
+      Int_t startIndex = 5, nDistances = 0;
+      for(Int_t iLayer = 5; iLayer > 0; iLayer--){
+        if(tracklet[iLayer]->IsOK() && tracklet[iLayer - 1]->IsOK()){
+          medx += tracklet[iLayer]->GetX0() - tracklet[iLayer - 1]->GetX0();
+          startIndex = iLayer - 1;
+          nDistances++;
+        }
+      }
+      if(nDistances){
+        medx /= nDistances;
+      }
+      else{
+        Float_t xpos[2];       memset(xpos, 0, sizeof(Float_t) * 2);
+        Int_t ien = 0, idiff = 0;
+        for(Int_t iLayer = 5; iLayer > 0; iLayer--){
+          if(tracklet[iLayer]->IsOK()){
+            xpos[ien++] = tracklet[iLayer]->GetX0();
+            startIndex = iLayer;
+          }
+          if(ien)
+            idiff++;
+          if(ien >=2)
+            break;
+        }
+        medx = (xpos[0] - xpos[1])/idiff;
+      }
+      xref = tracklet[startIndex]->GetX0() + medx * (2.5 - startIndex) - 0.5 * (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
+
+      for(Int_t ipt = 0; ipt < 1000; ipt++){
+        Float_t x = kxmin + ipt * kxdelta;
+        Float_t z = offset + slope * (x - xref);
+        fitfun->SetPoint(nPoints++, z, x);
+      }
+    }
+    break;
+  }
+  if(found){
+    TGraph *trGraph            = new TGraph(ncls);
+    TGraph *refPoints  = new TGraph(nLayers);
+    trGraph->SetMarkerStyle(20);
+    trGraph->SetMarkerColor(kRed);
+    refPoints->SetMarkerStyle(21);
+    refPoints->SetMarkerColor(kBlue);
+    // fill the graphs
+    for(Int_t iLayer = 0; iLayer < nLayers; iLayer++)
+      refPoints->SetPoint(iLayer, refP[iLayer], x0[iLayer]);
+    for(Int_t icls = 0; icls < ncls; icls++)
+      trGraph->SetPoint(icls, clp[icls], clx[icls]);
+    TCanvas *c1 = new TCanvas();
+    trGraph->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    refPoints->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    fitfun->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
+    trGraph->Draw("ap");
+    refPoints->Draw("psame");
+    fitfun->Draw("lpsame");
+    return c1;
+  }
+  else{
+    AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
+    delete fitfun;
+    return 0x0;
+  }
 }
 
index d7f854f02c7891e002e848c0ef32f1fe69713f68..042f4d81e09335a6294db991b30eb5a209f59986 100644 (file)
@@ -314,14 +314,10 @@ Int_t AliTRDtrackerV1::PropagateBack(AliESDEvent *event)
     }
 
     if ((TMath::Abs(track.GetC(track.GetBz()) - p4) / TMath::Abs(p4) < 0.2) ||(track.Pt() > 0.8)) {
-      //
+
       // Make backup for back propagation
-      //
       Int_t foundClr = track.GetNumberOfClusters();
       if (foundClr >= foundMin) {
-        //AliInfo(Form("Making backup track ncls [%d]...", foundClr));
-        //track.CookdEdx();
-        //track.CookdEdxTimBin(seed->GetID());
         track.CookLabel(1. - fgkLabelFraction);
         if(track.GetBackupTrack()) UseClusters(track.GetBackupTrack());
 
@@ -444,12 +440,11 @@ Int_t AliTRDtrackerV1::RefitInward(AliESDEvent *event)
       continue;
     }
 
+    // reject tracks which failed propagation in the TRD or
+    // are produced by the TRD stand alone tracker
     ULong_t status = seed->GetStatus();
-    // reject tracks which failed propagation in the TRD
-    if((status & AliESDtrack::kTRDout) == 0) continue;
-
-    // reject tracks which are produced by the TRD stand alone track finder.
-    if((status & AliESDtrack::kTRDin)  == 0) continue;
+    if(!(status & AliESDtrack::kTRDout)) continue;
+    if(!(status & AliESDtrack::kTRDin)) continue;
     nseed++; 
 
     track.ResetCovariance(50.0);
@@ -460,16 +455,26 @@ Int_t AliTRDtrackerV1::RefitInward(AliESDEvent *event)
     if(FollowProlongation(track)){     
       // Prolongate to TPC
       if (PropagateToX(track, xTPC, fgkMaxStep)) { //  -with update
-  seed->UpdateTrackParams(&track, AliESDtrack::kTRDrefit);
-  found++;
-  kUPDATE = kTRUE;
+        seed->UpdateTrackParams(&track, AliESDtrack::kTRDrefit);
+        found++;
+        kUPDATE = kTRUE;
+      }
+
+      // Update the friend track
+      if (fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker) > 0){ 
+        TObject *o = 0x0; Int_t ic = 0;
+        AliTRDtrackV1 *calibTrack = 0x0; 
+        while((o = seed->GetCalibObject(ic++))){
+          if(!(calibTrack = dynamic_cast<AliTRDtrackV1*>(o))) continue;
+          calibTrack->SetTrackHigh(track.GetTrackHigh());
+        }
       }
-    }   
+    }
     
     // Prolongate to TPC without update
     if(!kUPDATE) {
       AliTRDtrackV1 tt(*seed);
-      if (PropagateToX(tt, xTPC, fgkMaxStep)) seed->UpdateTrackParams(&tt, AliESDtrack::kTRDrefit);
+      if (PropagateToX(tt, xTPC, fgkMaxStep)) seed->UpdateTrackParams(&tt, AliESDtrack::kTRDbackup);
     }
   }
   AliInfo(Form("Number of loaded seeds: %d",nseed));
@@ -505,15 +510,15 @@ Int_t AliTRDtrackerV1::FollowProlongation(AliTRDtrackV1 &t)
   // Debug level 2
   //
   
+  Bool_t kStoreIn = kTRUE;
   Int_t    nClustersExpected = 0;
-  Int_t lastplane = 5; //GetLastPlane(&t);
-  for (Int_t iplane = lastplane; iplane >= 0; iplane--) {
+  for (Int_t iplane = kNPlanes; iplane--;) {
     Int_t   index   = 0;
     AliTRDseedV1 *tracklet = GetTracklet(&t, iplane, index);
     if(!tracklet) continue;
     if(!tracklet->IsOK()) AliWarning("tracklet not OK");
     
-    Double_t x  = tracklet->GetXref();//GetX0();
+    Double_t x  = tracklet->GetX();//GetX0();
     // reject tracklets which are not considered for inward refit
     if(x > t.GetX()+fgkMaxStep) continue;
 
@@ -551,7 +556,11 @@ Int_t AliTRDtrackerV1::FollowProlongation(AliTRDtrackV1 &t)
       t.PropagateTo(x, xx0, xrho);
       if (!AdjustSector(&t)) break;
     }
-    
+    if(kStoreIn){
+      t.SetTrackHigh(); 
+      kStoreIn = kFALSE;
+    }
+
     Double_t maxChi2 = t.GetPredictedChi2(tracklet);
     if (maxChi2 < 1e+10 && t.Update(tracklet, maxChi2)){ 
       nClustersExpected += tracklet->GetN();
@@ -618,9 +627,11 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     tracklets[ip] = t.GetTracklet(ip);
     t.UnsetTracklet(ip);
   } 
+  Bool_t kStoreIn = kTRUE;
+
 
   // Loop through the TRD layers
-  for (Int_t ilayer = 0; ilayer < AliTRDgeometry::Nlayer(); ilayer++) {
+  for (Int_t ilayer = 0; ilayer < kNPlanes; ilayer++) {
     // BUILD TRACKLET IF NOT ALREADY BUILT
     Double_t x = 0., x0, y, z, alpha;
     ptrTracklet  = tracklets[ilayer];
@@ -699,7 +710,7 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     Double_t xyz0[3]; // entry point 
     t.GetXYZ(xyz0);
     alpha = t.GetAlpha();
-    x = ptrTracklet->GetXref(); //GetX0();
+    x = ptrTracklet->GetX(); //GetX0();
     if (!t.GetProlongation(x, y, z)) return -1/*nClustersExpected*/;
     Double_t xyz1[3]; // exit point
     xyz1[0] =  x * TMath::Cos(alpha) - y * TMath::Sin(alpha); 
@@ -714,6 +725,11 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
     // Propagate and update track
     if (!t.PropagateTo(x, xx0, xrho)) return -1/*nClustersExpected*/;
     if (!AdjustSector(&t)) return -1/*nClustersExpected*/;
+
+    if(kStoreIn){
+      t.SetTrackLow(); 
+      kStoreIn = kFALSE;
+    }
     Double_t maxChi2 = t.GetPredictedChi2(ptrTracklet);
     if (!t.Update(ptrTracklet, maxChi2)) return -1/*nClustersExpected*/;
     ptrTracklet->UpDate(&t);
@@ -744,8 +760,9 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
         //(ratio0+ratio1           >   1.5) && 
         (t.GetNCross()           ==    0) && 
         (TMath::Abs(t.GetSnp())  <  0.85) &&
-        (t.GetNumberOfClusters() >    20)) t.MakeBackupTrack();
-    
+        (t.GetNumberOfClusters() >    20)){
+      t.MakeBackupTrack();
+    }
   } // end layers loop
 
   if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker) > 1){
@@ -787,7 +804,7 @@ Float_t AliTRDtrackerV1::FitRieman(AliTRDseedV1 *tracklets, Double_t *chi2, Int_
   }
   for(Int_t il = 0; il < maxLayers; il++){
     if(!tracklets[ppl[il]].IsOK()) continue;
-    fitter->AddPoint(tracklets[ppl[il]].GetX0(), tracklets[ppl[il]].GetYfitR(0), tracklets[ppl[il]].GetZProb(),1,10);
+    fitter->AddPoint(tracklets[ppl[il]].GetX0(), tracklets[ppl[il]].GetYfit(0), tracklets[ppl[il]].GetZfit(0),1,10);
   }
   fitter->Update();
   // Set the reference position of the fit and calculate the chi2 values
@@ -900,7 +917,7 @@ Float_t AliTRDtrackerV1::FitTiltedRiemanConstraint(AliTRDseedV1 *tracklets, Doub
 
   Float_t chi2track = fitter->GetChisquare()/Double_t(nPoints);
   for(Int_t ip = 0; ip < AliTRDtrackerV1::kNPlanes; ip++)
-    tracklets[ip].SetCC(curvature);
+    tracklets[ip].SetC(curvature);
 
 /*  if(fReconstructor->GetStreamLevel() >= 5){
     //Linear Model on z-direction
@@ -1004,7 +1021,7 @@ Float_t AliTRDtrackerV1::FitTiltedRieman(AliTRDseedV1 *tracklets, Bool_t sigErro
   for (Int_t iLayer = 0; iLayer < kNPlanes; iLayer++) {
     if(!tracklets[iLayer].IsOK()) continue;
     zref = offset + slope * (tracklets[iLayer].GetX0() - xref);
-    if (TMath::Abs(tracklets[iLayer].GetZProb() - zref) > tracklets[iLayer].GetPadLength() * 0.5 + 1.0) 
+    if (TMath::Abs(tracklets[iLayer].GetZfit(0) - zref) > tracklets[iLayer].GetPadLength() * 0.5 + 1.0) 
       acceptablez = kFALSE;
   }
   if (!acceptablez) {
@@ -1249,7 +1266,7 @@ Double_t AliTRDtrackerV1::FitRiemanTilt(const AliTRDtrackV1 *track, AliTRDseedV1
   for (Int_t iLayer = 0; iLayer < kNPlanes; iLayer++) {
     if(!tracklets[iLayer].IsOK()) continue;
     zref = z0 + dzdx * (tracklets[iLayer].GetX0() - xref);
-    if (TMath::Abs(tracklets[iLayer].GetZProb() - zref) > tracklets[iLayer].GetPadLength() * 0.5 + 1.0) 
+    if (TMath::Abs(tracklets[iLayer].GetZfit(0) - zref) > tracklets[iLayer].GetPadLength() * 0.5 + 1.0) 
       accept = kFALSE;
   }
   if (!accept) {
@@ -1448,7 +1465,7 @@ Float_t AliTRDtrackerV1::CalculateChi2Z(AliTRDseedV1 *tracklets, Double_t offset
   for (Int_t iLayer = 0; iLayer < AliTRDgeometry::kNlayer; iLayer++) {
     if(!tracklets[iLayer].IsOK()) continue;
     Double_t z = offset + slope * (tracklets[iLayer].GetX0() - xref);
-    chi2Z += TMath::Abs(tracklets[iLayer].GetMeanz() - z);
+    chi2Z += TMath::Abs(tracklets[iLayer].GetZfit(0) - z);
     nLayers++;
   }
   chi2Z /= TMath::Max((nLayers - 3.0),1.0);
@@ -2558,10 +2575,10 @@ Int_t AliTRDtrackerV1::MakeSeeds(AliTRDtrackingChamber **stack, AliTRDseedV1 *ss
         Freq(nlab,labels,outlab,kFALSE);
         Int_t label     = outlab[0];
         Int_t frequency = outlab[1];
-        for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
-          cseed[iLayer].SetFreq(frequency);
-          cseed[iLayer].SetChi2Z(chi2[1]);
-        }
+//         for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
+//           cseed[iLayer].SetFreq(frequency);
+//           cseed[iLayer].SetChi2Z(chi2[1]);
+//         }
             
         if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker) >= 2){
           TTreeSRedirector &cstreamer = *fReconstructor->GetDebugStream(AliTRDReconstructor::kTracker);
@@ -2801,7 +2818,7 @@ Double_t AliTRDtrackerV1::CalculateTrackLikelihood(AliTRDseedV1 *tracklets, Doub
   Double_t sumdaf = 0, nLayers = 0;
   for (Int_t iLayer = 0; iLayer < kNPlanes; iLayer++) {
     if(!tracklets[iLayer].IsOK()) continue;
-    sumdaf += TMath::Abs((tracklets[iLayer].GetYfit(1) - tracklets[iLayer].GetYref(1))/ tracklets[iLayer].GetSigmaY2());
+    sumdaf += TMath::Abs((tracklets[iLayer].GetYfit(1) - tracklets[iLayer].GetYref(1))/ tracklets[iLayer].GetS2Y());
     nLayers++;
   }
   sumdaf /= Float_t (nLayers - 2.0);
@@ -2869,7 +2886,7 @@ Double_t AliTRDtrackerV1::CookLikelihood(AliTRDseedV1 *cseed, Int_t planes[4])
   for(UChar_t ilayer = 0; ilayer < 4; ilayer++){
     Int_t jlayer = planes[ilayer];
     nclusters += cseed[jlayer].GetN2();
-    sumda += TMath::Abs(cseed[jlayer].GetYfitR(1) - cseed[jlayer].GetYref(1));
+    sumda += TMath::Abs(cseed[jlayer].GetYfit(1) - cseed[jlayer].GetYref(1));
   }
   nclusters *= .25;
 
@@ -3273,22 +3290,6 @@ Int_t AliTRDtrackerV1::Freq(Int_t n, const Int_t *inlist
 }
 
 
-//____________________________________________________________________
-
-//_____________________________________________________________________________
-Float_t AliTRDtrackerV1::GetChi2Y(AliTRDseedV1 *tracklets) const
-{
-  //   Chi2 definition on y-direction
-
-  Float_t chi2 = 0;
-  for(Int_t ipl = 0; ipl < kNPlanes; ipl++){
-    if(!tracklets[ipl].IsOK()) continue;
-    Double_t distLayer = (tracklets[ipl].GetYfit(0) - tracklets[ipl].GetYref(0));// /tracklets[ipl].GetSigmaY(); 
-    chi2 += distLayer * distLayer;
-  }
-  return chi2;
-}
-
 //____________________________________________________________________
 void AliTRDtrackerV1::ResetSeedTB()
 {
@@ -3301,20 +3302,35 @@ void AliTRDtrackerV1::ResetSeedTB()
   }
 }
 
+
+//_____________________________________________________________________________
+Float_t AliTRDtrackerV1::GetChi2Y(AliTRDseedV1 *tracklets) const
+{
+  //   Calculates normalized chi2 in y-direction
+  // chi2 = Sum chi2 / n_tracklets
+
+  Double_t chi2 = 0.; Int_t n = 0;
+  for(Int_t ipl = kNPlanes; ipl--;){
+    if(!tracklets[ipl].IsOK()) continue;
+    chi2 += tracklets[ipl].GetChi2Y();
+    n++;
+  }
+  return n ? chi2/n : 0.;
+}
+
 //_____________________________________________________________________________
 Float_t AliTRDtrackerV1::GetChi2Z(AliTRDseedV1 *tracklets) const 
 {
   //   Calculates normalized chi2 in z-direction
+  // chi2 = Sum chi2 / n_tracklets
 
-  Float_t chi2 = 0;
-  // chi2 = Sum ((z - zmu)/sigma)^2
-  // Sigma for the z direction is defined as half of the padlength
-  for(Int_t ipl = 0; ipl < kNPlanes; ipl++){
+  Double_t chi2 = 0; Int_t n = 0;
+  for(Int_t ipl = kNPlanes; ipl--;){
     if(!tracklets[ipl].IsOK()) continue;
-    Double_t distLayer = (tracklets[ipl].GetMeanz() - tracklets[ipl].GetZref(0)); // /(tracklets[ipl].GetPadLength()/2); 
-    chi2 += distLayer * distLayer;
+    chi2 += tracklets[ipl].GetChi2Z();
+    n++;
   }
-  return chi2;
+  return n ? chi2/n : 0.;
 }
 
 ///////////////////////////////////////////////////////
index ecdb004b1774056a96098db3932886b69f41a094..67921276f3ec757428a740a802b5efa6a309d653 100644 (file)
@@ -15,7 +15,7 @@
 #pragma link C++ class  AliTRDtrackInfoGen+;
 #pragma link C++ class  AliTRDrecoTask+;
 #pragma link C++ class  AliTRDcheckDetector+;
-#pragma link C++ class  AliTRDtrackingResolution+;
+#pragma link C++ class  AliTRDresolution+;
 #pragma link C++ class  AliTRDtrackingEfficiency+;
 #pragma link C++ class  AliTRDtrackingEfficiencyCombined+;
 #pragma link C++ class  AliTRDpidChecker+;
index e7884edb4768853c4f9576a7ac7e78fe4a312a93..2906ae9f08b68834f975834b70d2909d5d3df8ac 100644 (file)
@@ -147,21 +147,20 @@ void AliTRDcalibration::Exec(Option_t *)
     fTRDCalibraFillHisto->UpdateHistogramsV1(ftrdTrack);
     
     const AliTRDseedV1 *tracklet = 0x0;
-    for(Int_t itr = 0; itr < 6; itr++){
+    for(Int_t itr = 0; itr < AliTRDgeometry::kNlayer; itr++){
       if(!(tracklet = ftrdTrack->GetTracklet(itr))) continue;
       if(!tracklet->IsOK()) continue;
       Int_t nbclusters = 0;
       // For PH
-      Double_t phtb[35];
-      for(Int_t k=0; k < 35; k++){
-  phtb[k] = 0.0;
-      }
+      Double_t phtb[AliTRDseedV1::kNtb];
+      memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
+
       // For CH
       Double_t sum = 0.0;
       // normalisation
       Float_t normalisation = 6.67;
       Int_t detector = 0;
-      for(int ic=0; ic<AliTRDseed::knTimebins; ic++){
+      for(int ic=0; ic<AliTRDseedV1::kNTimeBins; ic++){
         if(!(fcl = tracklet->GetClusters(ic))) continue;
         nbclusters++;
         Int_t time = fcl->GetPadTime();
index 290c946f10b11444caf13fdbcc6d0c6a505692a9..a8a4c29cd7c78c1e057aa06e0d28f761cc27f80a 100644 (file)
@@ -661,7 +661,7 @@ TH1 *AliTRDcheckDetector::PlotFindableTracklets(const AliTRDtrackV1 *track){
     // ignore y-crossing (material)
     if((z + delta_z > zmin && z - delta_z < zmax) && (y + delta_y > ymin && y - delta_y < ymax)) nFindable++;
       if(fDebugLevel > 3){
-        Double_t pos_tracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetMeanz() : 0};
+        Double_t pos_tracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetZfit(0) : 0};
         Int_t hasTracklet = tracklet ? 1 : 0;
         (*fDebugStream)   << "FindableTracklets"
           << "layer="     << il
@@ -865,7 +865,7 @@ TH1 *AliTRDcheckDetector::PlotChargeTracklet(const AliTRDtrackV1 *track){
   for(Int_t itl = 0x0; itl < AliTRDgeometry::kNlayer; itl++){
     if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
     Qtot = 0.;
-    for(Int_t ic = AliTRDseed::knTimebins; ic--;){
+    for(Int_t ic = AliTRDseedV1::kNTimeBins; ic--;){
       if(!(c = tracklet->GetClusters(ic))) continue;
       Qtot += TMath::Abs(c->GetQ());
     }
similarity index 84%
rename from TRD/qaRec/AliTRDtrackingResolution.cxx
rename to TRD/qaRec/AliTRDresolution.cxx
index 90c98d780953fde13b65ce83903b3cdcc70ff82c..dda142fc20ff189fd92ed07218db54766bc37301 100644 (file)
@@ -13,7 +13,7 @@
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
 
-/* $Id: AliTRDtrackingResolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
+/* $Id: AliTRDresolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
 
 ////////////////////////////////////////////////////////////////////////////
 //                                                                        //
@@ -33,7 +33,7 @@
 // {  
 //   gSystem->Load("libANALYSIS.so");
 //   gSystem->Load("libTRDqaRec.so");
-//   AliTRDtrackingResolution *res = new AliTRDtrackingResolution();
+//   AliTRDresolution *res = new AliTRDresolution();
 //   //res->SetMCdata();
 //   //res->SetVerbose();
 //   //res->SetVisual();
 
 #include "AliTRDtrackInfo/AliTRDclusterInfo.h"
 #include "AliTRDtrackInfo/AliTRDtrackInfo.h"
-#include "AliTRDtrackingResolution.h"
+#include "AliTRDresolution.h"
 
-ClassImp(AliTRDtrackingResolution)
+ClassImp(AliTRDresolution)
 
 //________________________________________________________
-AliTRDtrackingResolution::AliTRDtrackingResolution()
-  :AliTRDrecoTask("Resolution", "Tracking Resolution")
+AliTRDresolution::AliTRDresolution()
+  :AliTRDrecoTask("Resolution", "Spatial and Momentum Resolution")
   ,fStatus(0)
   ,fReconstructor(0x0)
   ,fGeo(0x0)
@@ -112,7 +112,7 @@ AliTRDtrackingResolution::AliTRDtrackingResolution()
 }
 
 //________________________________________________________
-AliTRDtrackingResolution::~AliTRDtrackingResolution()
+AliTRDresolution::~AliTRDresolution()
 {
   if(fGraphS){fGraphS->Delete(); delete fGraphS;}
   if(fGraphM){fGraphM->Delete(); delete fGraphM;}
@@ -127,7 +127,7 @@ AliTRDtrackingResolution::~AliTRDtrackingResolution()
 
 
 //________________________________________________________
-void AliTRDtrackingResolution::CreateOutputObjects()
+void AliTRDresolution::CreateOutputObjects()
 {
   // spatial resolution
   OpenFile(0, "RECREATE");
@@ -145,7 +145,7 @@ void AliTRDtrackingResolution::CreateOutputObjects()
 }
 
 //________________________________________________________
-void AliTRDtrackingResolution::Exec(Option_t *opt)
+void AliTRDresolution::Exec(Option_t *opt)
 {
   fCl->Delete();
   fTrklt->Delete();
@@ -161,7 +161,7 @@ void AliTRDtrackingResolution::Exec(Option_t *opt)
 }
 
 //________________________________________________________
-TH1* AliTRDtrackingResolution::PlotCluster(const AliTRDtrackV1 *track)
+TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track)
 {
   if(track) fTrack = track;
   if(!fTrack){
@@ -231,7 +231,7 @@ TH1* AliTRDtrackingResolution::PlotCluster(const AliTRDtrackV1 *track)
 
 
 //________________________________________________________
-TH1* AliTRDtrackingResolution::PlotTracklet(const AliTRDtrackV1 *track)
+TH1* AliTRDresolution::PlotTracklet(const AliTRDtrackV1 *track)
 {
   if(track) fTrack = track;
   if(!fTrack){
@@ -245,17 +245,17 @@ TH1* AliTRDtrackingResolution::PlotTracklet(const AliTRDtrackV1 *track)
   }
 
   Double_t cov[3], covR[3];
-  Float_t xref, y0, dx, dy, dydx;
+  Float_t x, y0, dx, dy, dydx;
   AliTRDseedV1 *fTracklet = 0x0;  
   for(Int_t il=AliTRDgeometry::kNlayer; il--;){
     if(!(fTracklet = fTrack->GetTracklet(il))) continue;
     if(!fTracklet->IsOK()) continue;
     y0   = fTracklet->GetYref(0);
     dydx = fTracklet->GetYref(1);
-    xref = fTracklet->GetXref();
-    dx   = fTracklet->GetX0() - xref;
-    dy   = y0-dx*dydx - fTracklet->GetYat(xref);
-    fTracklet->GetCovAt(xref, cov);
+    x    = fTracklet->GetX();
+    dx   = fTracklet->GetX0() - x;
+    dy   = y0-dx*dydx - fTracklet->GetY();
+    fTracklet->GetCovAt(x, cov);
     fTracklet->GetCovRef(covR);
     h->Fill(dydx, dy/TMath::Sqrt(cov[0] + covR[0]));
   }
@@ -263,7 +263,7 @@ TH1* AliTRDtrackingResolution::PlotTracklet(const AliTRDtrackV1 *track)
 }
 
 //________________________________________________________
-TH1* AliTRDtrackingResolution::PlotTrackletPhi(const AliTRDtrackV1 *track)
+TH1* AliTRDresolution::PlotTrackletPhi(const AliTRDtrackV1 *track)
 {
   if(track) fTrack = track;
   if(!fTrack){
@@ -287,7 +287,7 @@ TH1* AliTRDtrackingResolution::PlotTrackletPhi(const AliTRDtrackV1 *track)
 
 
 //________________________________________________________
-TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
+TH1* AliTRDresolution::PlotResolution(const AliTRDtrackV1 *track)
 {
   if(!HasMCdata()){ 
     AliWarning("No MC defined. Results will not be available.");
@@ -302,7 +302,9 @@ TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
   UChar_t s;
   Int_t pdg = fMC->GetPDG(), det=-1;
   Int_t label = fMC->GetLabel();
-  Float_t p, pt, xref, x0, y0, z0, dx, dy, dz, dydx, dzdx;
+  Double_t x, y, z;
+  Float_t p, pt, x0, y0, z0, dx, dy, dz, dydx, dzdx;
+  Double_t covR[3];
 
   if(fDebugLevel>=1){
     Double_t DX[12], DY[12], DZ[12], DPt[12], COV[12][15];
@@ -332,41 +334,57 @@ TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
     det = fTracklet->GetDetector();
     x0  = fTracklet->GetX0();
     //radial shift with respect to the MC reference (radial position of the pad plane)
-    xref= fTracklet->GetXref();
-    dx  = x0 - xref;
+    x= fTracklet->GetX();
+    dx  = x0 - x;
     if(!fMC->GetDirections(x0, y0, z0, dydx, dzdx, pt, s)) continue;
+    if(fDebugLevel>=1){
+      (*fDebugStream) << "MC"
+        << "det="     << det
+        << "pdg="     << pdg
+        << "pt="      << pt
+        << "x0="      << x0
+        << "y0="      << y0
+        << "z0="      << z0
+        << "dydx="    << dydx
+        << "dzdx="    << dzdx
+        << "\n";
+    }
     // MC track position at reference radial position
-    Float_t yt = y0 - (x0-xref)*dydx;
-    Float_t zt = z0 - (x0-xref)*dzdx;
+    Float_t yt = y0 - dx*dydx;
+    Float_t zt = z0 - dx*dzdx;
     p = pt*(1.+dzdx*dzdx); // pt -> p
 
     // add Kalman residuals for y, z and pt
-    Float_t yr = fTracklet->GetYref(0) - dx*fTracklet->GetYref(1) + .05;
+    Float_t yr = fTracklet->GetYref(0) - dx*fTracklet->GetYref(1);
     dy = yt - yr;
     Float_t zr = fTracklet->GetZref(0) - dx*fTracklet->GetZref(1);
     dz = zt - zr;
     Float_t tgl = fTracklet->GetTgl();
-    Float_t ptr = fTracklet->GetMomentum()/(1.+tgl*tgl);
+    Float_t dpt = pt - fTracklet->GetMomentum()/(1.+tgl*tgl);
+    fTracklet->GetCovRef(covR);
 
     ((TH2I*)fContainer->At(kMCtrackY))->Fill(dydx, dy);
-    ((TH2I*)fContainer->At(kMCtrackZ))->Fill(dzdx, dz);
-    if(pdg!=kElectron && pdg!=kPositron) ((TH2I*)fContainer->At(kMCtrackPt))->Fill(1./pt, ptr-pt);
+    ((TH2I*)fContainer->At(kMCtrackYPull))->Fill(dydx, dy/TMath::Sqrt(covR[0]));
+    if(ily==0){
+      ((TH2I*)fContainer->At(kMCtrackZIn))->Fill(dzdx, dz);
+      ((TH2I*)fContainer->At(kMCtrackZInPull))->Fill(dzdx, dz/TMath::Sqrt(covR[2]));
+    } else if(ily==AliTRDgeometry::kNlayer-1) {
+      ((TH2I*)fContainer->At(kMCtrackZOut))->Fill(dzdx, dz);
+      ((TH2I*)fContainer->At(kMCtrackZOutPull))->Fill(dzdx, dz/TMath::Sqrt(covR[2]));
+    }
+    if(pdg!=kElectron && pdg!=kPositron) ((TH2I*)fContainer->At(kMCtrackPt))->Fill(1./pt, dpt);
     // Fill Debug stream for Kalman track
     if(fDebugLevel>=1){
       Float_t dydxr = fTracklet->GetYref(1);
       (*fDebugStream) << "MCtrack"
-        << "det="     << det
-        << "pdg="     << pdg
-        << "pt="      << pt
-        << "yt="      << yt
-        << "zt="      << zt
-        << "dydx="    << dydx
-        << "dzdx="    << dzdx
-        << "ptr="     << ptr
+        << "x="       << x
+        << "dpt="     << dpt
         << "dy="      << dy
         << "dz="      << dz
         << "dydxr="   << dydxr
         << "dzdxr="   << tgl
+        << "s2y="     << covR[0]
+        << "s2z="     << covR[2]
         << "\n";
     }
 
@@ -374,41 +392,39 @@ TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
     AliTRDseedV1 tt(*fTracklet);
     tt.SetZref(0, z0);
     tt.SetZref(1, dzdx); 
-    if(!tt.Fit(kTRUE)) continue;
-    xref= fTracklet->GetXref(); // the true one 
-    yt = y0 - (x0-xref)*dydx;
-
+    tt.Fit(kTRUE);
+    x= tt.GetX(); // the true one 
+    dx = x0 - x;
+    yt = y0 - dx*dydx;
+    Bool_t rc = tt.IsRowCross(); 
+    
     // add tracklet residuals for y and dydx
-    Float_t yf = tt.GetYat(xref) + .05;
-    dy = yt - yf;
+    Float_t yf = tt.GetY();
+    dy = yt - yf; dz = 100.;
     Float_t dphi   = (tt.GetYfit(1) - dydx);
     dphi /= 1.- tt.GetYfit(1)*dydx;
-    ((TH2I*)fContainer->At(kMCtrackletY))->Fill(dydx, dy);
-    ((TH2I*)fContainer->At(kMCtrackletPhi))->Fill(dydx, dphi*TMath::RadToDeg());
-
-    Float_t dz = 100.;
-    Bool_t rc = fTracklet->IsRowCross(); 
-    if(rc){
+    Double_t s2y = tt.GetS2Y(), s2z = tt.GetS2Z();
+    if(!rc){
+      ((TH2I*)fContainer->At(kMCtrackletY))->Fill(dydx, dy);
+      ((TH2I*)fContainer->At(kMCtrackletYPull))->Fill(dydx, dy/TMath::Sqrt(s2y));
+      ((TH2I*)fContainer->At(kMCtrackletPhi))->Fill(dydx, dphi*TMath::RadToDeg());
+    } else {
       // add tracklet residuals for z
-      Double_t *xyz = tt.GetCrossXYZ();
-      dz = xyz[2] - (z0 - (x0 - xyz[0])*dzdx) ;
+      dz = tt.GetZ() - (z0 - dx*dzdx) ;
       ((TH2I*)fContainer->At(kMCtrackletZ))->Fill(dzdx, dz);
+      ((TH2I*)fContainer->At(kMCtrackletZPull))->Fill(dzdx, dz/TMath::Sqrt(s2z));
     }
   
     // Fill Debug stream for tracklet
     if(fDebugLevel>=1){
       (*fDebugStream) << "MCtracklet"
-        << "det="     << det
-        << "pdg="     << pdg
-        << "p="       << p
-        << "yt="      << yt
-        << "zt="      << zt
-        << "dydx="    << dydx
-        << "dzdx="    << dzdx
-        << "rowc="    << rc
-        << "dy="      << dy
-        << "dz="      << dz
-        << "dphi="    << dphi
+        << "rc="    << rc
+        << "x="     << x
+        << "dy="    << dy
+        << "dz="    << dz
+        << "dphi="  << dphi
+        << "s2y="   << s2y
+        << "s2z="   << s2z
         << "\n";
     }
 
@@ -417,16 +433,15 @@ TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
     Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset();
     Float_t tilt = fTracklet->GetTilt();
 
-    Double_t x,y;
     AliTRDcluster *c = 0x0;
     fTracklet->ResetClusterIter(kFALSE);
     while((c = fTracklet->PrevCluster())){
       Float_t  q = TMath::Abs(c->GetQ());
       //AliTRDseedV1::GetClusterXY(c,x,y);
-      x = c->GetX(); y = c->GetY();
+      x = c->GetX(); y = c->GetY(); z = c->GetZ();
       Float_t xc = x;
       Float_t yc = y;
-      Float_t zc = c->GetZ();
+      Float_t zc = z;
       dx = x0 - xc; 
       yt = y0 - dx*dydx;
       zt = z0 - dx*dzdx;
@@ -463,7 +478,7 @@ TH1* AliTRDtrackingResolution::PlotResolution(const AliTRDtrackV1 *track)
 
 
 //________________________________________________________
-Bool_t AliTRDtrackingResolution::GetRefFigure(Int_t ifig)
+Bool_t AliTRDresolution::GetRefFigure(Int_t ifig)
 {
   Float_t y[2] = {0., 0.};
   TBox *b = 0x0;
@@ -584,7 +599,7 @@ Bool_t AliTRDtrackingResolution::GetRefFigure(Int_t ifig)
     b->SetFillStyle(3002);b->SetFillColor(kBlue);
     b->SetLineColor(0); b->Draw();
     return kTRUE;
-  case kMCtrackZ:
+  case kMCtrackZIn:
     if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
     ax = g->GetHistogram()->GetYaxis();
     ax->SetRangeUser(-.5, 2.);
@@ -613,7 +628,7 @@ Bool_t AliTRDtrackingResolution::GetRefFigure(Int_t ifig)
 
 
 //________________________________________________________
-Bool_t AliTRDtrackingResolution::PostProcess()
+Bool_t AliTRDresolution::PostProcess()
 {
   //fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
   if (!fContainer) {
@@ -843,9 +858,9 @@ Bool_t AliTRDtrackingResolution::PostProcess()
   }
 
   // track z resolution
-  h2 = (TH2I*)fContainer->At(kMCtrackZ);
-  gm = (TGraphErrors*)fGraphM->At(kMCtrackZ);
-  gs = (TGraphErrors*)fGraphS->At(kMCtrackZ);
+  h2 = (TH2I*)fContainer->At(kMCtrackZIn);
+  gm = (TGraphErrors*)fGraphM->At(kMCtrackZIn);
+  gs = (TGraphErrors*)fGraphS->At(kMCtrackZIn);
   for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
     h = h2->ProjectionY("pz", iphi, iphi);
     if(h->GetEntries()<70) continue;
@@ -905,7 +920,7 @@ Bool_t AliTRDtrackingResolution::PostProcess()
 
 
 //________________________________________________________
-void AliTRDtrackingResolution::Terminate(Option_t *)
+void AliTRDresolution::Terminate(Option_t *)
 {
   if(fDebugStream){ 
     delete fDebugStream;
@@ -916,7 +931,7 @@ void AliTRDtrackingResolution::Terminate(Option_t *)
 }
 
 //________________________________________________________
-void AliTRDtrackingResolution::AdjustF1(TH1 *h, TF1 *f)
+void AliTRDresolution::AdjustF1(TH1 *h, TF1 *f)
 {
 // Helper function to avoid duplication of code
 // Make first guesses on the fit parameters
@@ -946,11 +961,11 @@ void AliTRDtrackingResolution::AdjustF1(TH1 *h, TF1 *f)
 }
 
 //________________________________________________________
-TObjArray* AliTRDtrackingResolution::Histos()
+TObjArray* AliTRDresolution::Histos()
 {
   if(fContainer) return fContainer;
 
-  fContainer  = new TObjArray(10);
+  fContainer  = new TObjArray(kNhistos);
   //fContainer->SetOwner(kTRUE);
 
   TH1 *h = 0x0;
@@ -1003,6 +1018,15 @@ TObjArray* AliTRDtrackingResolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackletY);
 
+  // tracklet y resolution [0]
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltYPull"))){
+    h = new TH2I("hMCtrkltYPull", "Tracklet Pulls (Y)", 31, -.48, .48, 100, -3.5, 3.5);
+    h->GetXaxis()->SetTitle("tg(#phi)");
+    h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackletYPull);
+
   // tracklet y resolution [0]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZ"))){
     h = new TH2I("hMCtrkltZ", "Tracklet Resolution (Z)", 31, -.48, .48, 100, -.5, .5);
@@ -1012,6 +1036,15 @@ TObjArray* AliTRDtrackingResolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackletZ);
 
+  // tracklet y resolution [0]
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZPull"))){
+    h = new TH2I("hMCtrkltZ", "Tracklet Pulls (Z)", 31, -.48, .48, 100, -3.5, 3.5);
+    h->GetXaxis()->SetTitle("tg(#theta)");
+    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackletZPull);
+
   // tracklet angular resolution [1]
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltPhi"))){
     h = new TH2I("hMCtrkltPhi", "Tracklet Resolution (#Phi)", 31, -.48, .48, 100, -10., 10.);
@@ -1030,14 +1063,50 @@ TObjArray* AliTRDtrackingResolution::Histos()
   } else h->Reset();
   fContainer->AddAt(h, kMCtrackY);
 
+  // Kalman track y resolution
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkYPull"))){
+    h = new TH2I("hMCtrkYPull", "Kalman Track Pulls (Y)", 31, -.48, .48, 100, -3.5, 3.5);
+    h->GetXaxis()->SetTitle("tg(#phi)");
+    h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackYPull);
+
+  // Kalman track Z resolution
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZIn"))){
+    h = new TH2I("hMCtrkZIn", "Kalman Track Resolution (Zin)", 20, -1., 1., 100, -1.5, 1.5);
+    h->GetXaxis()->SetTitle("tg(#theta)");
+    h->GetYaxis()->SetTitle("#Delta z [cm]");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackZIn);
+
   // Kalman track Z resolution
-  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZ"))){
-    h = new TH2I("hMCtrkZ", "Kalman Track Resolution (Z)", 20, -1., 1., 100, -1.5, 1.5);
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOut"))){
+    h = new TH2I("hMCtrkZOut", "Kalman Track Resolution (Zout)", 20, -1., 1., 100, -1.5, 1.5);
     h->GetXaxis()->SetTitle("tg(#theta)");
     h->GetYaxis()->SetTitle("#Delta z [cm]");
     h->GetZaxis()->SetTitle("entries");
   } else h->Reset();
-  fContainer->AddAt(h, kMCtrackZ);
+  fContainer->AddAt(h, kMCtrackZOut);
+
+  // Kalman track Z resolution
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZInPull"))){
+    h = new TH2I("hMCtrkZInPull", "Kalman Track Pulls (Zin)", 20, -1., 1., 100, -4.5, 4.5);
+    h->GetXaxis()->SetTitle("tg(#theta)");
+    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackZInPull);
+
+  // Kalman track Z resolution
+  if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOutPull"))){
+    h = new TH2I("hMCtrkZOutPull", "Kalman Track Pulls (Zout)", 20, -1., 1., 100, -4.5, 4.5);
+    h->GetXaxis()->SetTitle("tg(#theta)");
+    h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
+    h->GetZaxis()->SetTitle("entries");
+  } else h->Reset();
+  fContainer->AddAt(h, kMCtrackZOutPull);
 
   // Kalman track Pt resolution
   if(!(h = (TH2I*)gROOT->FindObject("hMCtrkPt"))){
@@ -1053,7 +1122,7 @@ TObjArray* AliTRDtrackingResolution::Histos()
 
 
 //________________________________________________________
-void AliTRDtrackingResolution::SetRecoParam(AliTRDrecoParam *r)
+void AliTRDresolution::SetRecoParam(AliTRDrecoParam *r)
 {
 
   fReconstructor->SetRecoParam(r);
similarity index 58%
rename from TRD/qaRec/AliTRDtrackingResolution.h
rename to TRD/qaRec/AliTRDresolution.h
index 63b412e8d27740e713d322555608ef31327d8327..2d3720e1a288d59528faeeaf27c573572127c8e3 100644 (file)
@@ -1,13 +1,13 @@
-#ifndef ALITRDTRACKINGRESOLUTION_H
-#define ALITRDTRACKINGRESOLUTION_H
+#ifndef ALITRDRESOLUTION_H
+#define ALITRDRESOLUTION_H
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
-/* $Id: AliTRDtrackingResolution.h 27496 2008-07-22 08:35:45Z cblume $ */
+/* $Id: AliTRDresolution.h 27496 2008-07-22 08:35:45Z cblume $ */
 
 ////////////////////////////////////////////////////////////////////////////
 //                                                                        //
-//  Reconstruction QA                                                     //
+//  Resolution QA                                                     //
 //                                                                        //
 ////////////////////////////////////////////////////////////////////////////
 
@@ -23,28 +23,35 @@ class AliTRDgeometry;
 class AliTRDrecoParam;
 class AliTRDseedV1;
 class AliTRDtrackInfo;
-class AliTRDtrackingResolution : public AliTRDrecoTask
+class AliTRDresolution : public AliTRDrecoTask
 {
 public:
-  enum{
-    kCluster        = 0 // cluster - track
-    ,kTrackletY     = 1 // tracklet - track y pools
-    ,kTrackletPhi   = 2 // tracklet - track angular pools residuals
-    ,kMCcluster     = 3 // cluster - mc residuals/systematics
-    ,kMCtrackletY   = 4 // tracklet - mc y resolution/systematics
-    ,kMCtrackletZ   = 5 // tracklet - mc z resolution/systematics (pad row cross)
-    ,kMCtrackletPhi = 6 // tracklet - mc phi resolution/systematics
-    ,kMCtrackY      = 7 // Kalman Y resolution
-    ,kMCtrackZ      = 8 // Kalman Z resolution
-    ,kMCtrackPt     = 9 // Kalman Pt resolution
+  enum ETRDresolutionPlots{
+    kCluster          =  0 // cluster - track
+    ,kTrackletY       =  1 // tracklet - track y pulls
+    ,kTrackletPhi     =  2 // tracklet - track angular pulls residuals
+    ,kMCcluster       =  3 // cluster - mc residuals/systematics
+    ,kMCtrackletY     =  4 // tracklet - mc y resolution/systematics
+    ,kMCtrackletYPull =  5 // tracklet - mc y resolution/systematics
+    ,kMCtrackletZ     =  6 // tracklet - mc z resolution/systematics (pad row cross)
+    ,kMCtrackletZPull =  7 // tracklet - mc z resolution/systematics (pad row cross)
+    ,kMCtrackletPhi   =  8 // tracklet - mc phi resolution/systematics
+    ,kMCtrackY        =  9 // Kalman Y resolution
+    ,kMCtrackYPull    = 10 // Kalman Y resolution
+    ,kMCtrackZIn      = 11 // Kalman Z resolution
+    ,kMCtrackZOut     = 12 // Kalman Z resolution
+    ,kMCtrackZInPull  = 13 // Kalman Z resolution
+    ,kMCtrackZOutPull = 14 // Kalman Z resolution
+    ,kMCtrackPt       = 15 // Kalman Pt resolution
+    ,kNhistos         = 16
   };
-  enum{
+  enum ETRDresolutionSteer {
     kVerbose  = 0
     ,kVisual  = 1
   };
 
-  AliTRDtrackingResolution();
-  virtual ~AliTRDtrackingResolution();
+  AliTRDresolution();
+  virtual ~AliTRDresolution();
   
   void    CreateOutputObjects();
   Bool_t  GetRefFigure(Int_t ifig);
@@ -66,8 +73,8 @@ public:
   void    Terminate(Option_t *);
   
 private:
-  AliTRDtrackingResolution(const AliTRDtrackingResolution&);
-  AliTRDtrackingResolution& operator=(const AliTRDtrackingResolution&);
+  AliTRDresolution(const AliTRDresolution&);
+  AliTRDresolution& operator=(const AliTRDresolution&);
   void        AdjustF1(TH1 *h, TF1 *f);
 
 private:
@@ -83,6 +90,6 @@ private:
   TObjArray           *fMCcl;   //! cluster2mc calib
   TObjArray           *fMCtrklt;//! tracklet2mc calib
   
-  ClassDef(AliTRDtrackingResolution, 1) // TRD tracking resolution task
+  ClassDef(AliTRDresolution, 1) // TRD tracking resolution task
 };
 #endif
index cc57343a60e0f8639c4855ef4f639b84400cb7f5..a13cc6cea27a0cb44f7b1d6e63c648b83368c6a5 100644 (file)
@@ -63,7 +63,7 @@
 #include "TRD/qaRec/AliTRDtrackInfoGen.h"
 #include "TRD/qaRec/AliTRDtrackingEfficiency.h"
 #include "TRD/qaRec/AliTRDtrackingEfficiencyCombined.h"
-#include "TRD/qaRec/AliTRDtrackingResolution.h"
+#include "TRD/qaRec/AliTRDresolution.h"
 #include "TRD/qaRec/AliTRDcalibration.h"
 #include "TRD/qaRec/AliTRDalignmentTask.h"
 #include "TRD/qaRec/AliTRDpidChecker.h"
@@ -186,14 +186,15 @@ void run(Char_t *tasks="ALL", const Char_t *files=0x0)
 
   //____________________________________________
   // TRD track summary generator
+  AliAnalysisDataContainer *coutput1 = 0x0, *coutput1a = 0x0;
        if(TSTBIT(fSteerTask, kInfoGen)){
     mgr->AddTask(task = new AliTRDtrackInfoGen());
     taskPtr[(Int_t)kInfoGen] = task;
     task->SetDebugLevel(0);
     task->SetMCdata(fHasMCdata);
     mgr->ConnectInput( task, 0, mgr->GetCommonInputContainer());
-    AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("trackInfo", TObjArray::Class(), AliAnalysisManager::kExchangeContainer);
-    AliAnalysisDataContainer *coutput1a = mgr->CreateContainer("eventInfo", AliTRDeventInfo::Class(), AliAnalysisManager::kExchangeContainer);
+    coutput1 = mgr->CreateContainer("trackInfo", TObjArray::Class(), AliAnalysisManager::kExchangeContainer);
+    coutput1a = mgr->CreateContainer("eventInfo", AliTRDeventInfo::Class(), AliAnalysisManager::kExchangeContainer);
     mgr->ConnectOutput(task, 0, coutput1);
     mgr->ConnectOutput(task, 1, coutput1a);
   }
@@ -239,7 +240,7 @@ void run(Char_t *tasks="ALL", const Char_t *files=0x0)
   //____________________________________________
   // TRD tracking resolution
   if(TSTBIT(fSteerTask, kResolution)){
-    mgr->AddTask(task = new AliTRDtrackingResolution());
+    mgr->AddTask(task = new AliTRDresolution());
     taskPtr[(Int_t)kResolution] = task;
     task->SetMCdata(fHasMCdata);
     task->SetPostProcess(kFALSE);
index 542a2faf4b3f98292f60b75849448d199205cf9e..00d0025d29fb131581f2ec59309c3a06504f9bac 100644 (file)
@@ -28,7 +28,7 @@ const Char_t* fgkTRDtaskClassName[NQATASKS+NCALIBTASKS] = {
   ,"AliTRDcheckDetector"
   ,"AliTRDtrackingEfficiency"
   ,"AliTRDtrackingEfficiencyCombined"
-  ,"AliTRDtrackingResolution"
+  ,"AliTRDresolution"
   ,"AliTRDpidChecker"
   ,"AliTRDcalibration"
   ,"AliTRDalignmentTask"