]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCtrackerMI.cxx
Added a cut on PtHard at 2.76 GeV/c (Nicole)
[u/mrichter/AliRoot.git] / TPC / AliTPCtrackerMI.cxx
index 695ec8274ff39af108ce8c50be017962dfc948fe..ebb763c45246021957ba642adfb38bf32e4d4adb 100644 (file)
 //
 //   Origin: Marian Ivanov   Marian.Ivanov@cern.ch
 // 
-//  AliTPC parallel tracker - 
-//  How to use?  - 
-//  run AliTPCFindClusters.C macro - clusters neccessary for tracker are founded
-//  run AliTPCFindTracksMI.C macro - to find tracks
-//  tracks are written to AliTPCtracks.root file
-//  for comparison also seeds are written to the same file - to special branch
+//  AliTPC parallel tracker
 //-------------------------------------------------------
 
 
@@ -128,13 +123,13 @@ Int_t AliTPCtrackerMI::UpdateTrack(AliTPCseed * track, Int_t accept){
   track->fClusterPointer[track->fRow] = c;  
   //
 
-  Float_t angle2 = track->GetSnp()*track->GetSnp();
-  angle2 = TMath::Sqrt(angle2/(1-angle2)); 
+  Double_t angle2 = track->GetSnp()*track->GetSnp();
   //
   //SET NEW Track Point
   //
-  //  if (debug)
+  if (angle2<1) //PH sometimes angle2 is very big. To be investigated...
   {
+    angle2 = TMath::Sqrt(angle2/(1-angle2)); 
     AliTPCTrackerPoint   &point =*(track->GetTrackPoint(track->fRow));
     //
     point.SetSigmaY(c->GetSigmaY2()/track->fCurrentSigmaY2);
@@ -645,7 +640,7 @@ Double_t AliTPCtrackerMI::ErrY2(AliTPCseed* seed, AliTPCclusterMI * cl){
   Float_t z = TMath::Abs(fParam->GetZLength()-TMath::Abs(seed->GetZ()));
   Int_t ctype = cl->GetType();  
   Float_t padlength= GetPadPitchLength(seed->fRow);
-  Float_t angle2 = seed->GetSnp()*seed->GetSnp();
+  Double_t angle2 = seed->GetSnp()*seed->GetSnp();
   angle2 = angle2/(1-angle2); 
   //
   //cluster "quality"
@@ -780,7 +775,7 @@ Double_t AliTPCtrackerMI::ErrZ2(AliTPCseed* seed, AliTPCclusterMI * cl){
   Int_t ctype = cl->GetType();  
   Float_t padlength= GetPadPitchLength(seed->fRow);
   //
-  Float_t angle2 = seed->GetSnp()*seed->GetSnp();
+  Double_t angle2 = seed->GetSnp()*seed->GetSnp();
   //  if (angle2<0.6) angle2 = 0.6;
   angle2 = seed->GetTgl()*seed->GetTgl()*(1+angle2/(1-angle2)); 
   //
@@ -1131,6 +1126,9 @@ Int_t  AliTPCtrackerMI::LoadClusters()
     //  
     Int_t sec,row;
     fParam->AdjustSectorRow(clrow->GetID(),sec,row);
+    for (Int_t icl=0; icl<clrow->GetArray()->GetEntriesFast(); icl++){
+      Transform((AliCluster*)(clrow->GetArray()->At(icl)));
+    }
     //
     AliTPCRow * tpcrow=0;
     Int_t left=0;
@@ -1193,6 +1191,19 @@ void AliTPCtrackerMI::UnloadClusters()
   return ;
 }
 
+void   AliTPCtrackerMI::Transform(AliCluster * cluster){
+  //
+  // 
+  //
+  //if (!fParam->IsGeoRead()) fParam->ReadGeoMatrices();
+  TGeoHMatrix  *mat = fParam->GetClusterMatrix(cluster->GetDetector());
+  Double_t pos[3]= {cluster->GetX(),cluster->GetY(),cluster->GetZ()};
+  Double_t posC[3];
+  mat->LocalToMaster(pos,posC);
+  cluster->SetX(posC[0]);
+  cluster->SetY(posC[1]);
+  cluster->SetZ(posC[2]);
+}
 
 //_____________________________________________________________________________
 Int_t AliTPCtrackerMI::LoadOuterSectors() {
@@ -1305,6 +1316,11 @@ AliTPCclusterMI *AliTPCtrackerMI::GetClusterMI(Int_t index) const {
   const AliTPCRow * tpcrow=0;
   AliTPCclusterMI * clrow =0;
 
+  if (sec<0 || sec>=fkNIS*4) {
+    AliWarning(Form("Wrong sector %d",sec));
+    return 0x0;
+  }
+
   if (sec<fkNIS*2){
     tpcrow = &(fInnerSec[sec%fkNIS][row]);
     if (tpcrow==0) return 0;
@@ -1370,9 +1386,9 @@ Int_t AliTPCtrackerMI::FollowToNext(AliTPCseed& t, Int_t nr) {
       if (TMath::Abs(angle-t.GetAlpha())>0.001){
        Double_t rotation = angle-t.GetAlpha();
        t.fRelativeSector= relativesector;
-       t.Rotate(rotation);     
+       if (!t.Rotate(rotation)) return 0;      
       }
-      t.PropagateTo(x);
+      if (!t.PropagateTo(x)) return 0;
       //
       t.fCurrentCluster = cl; 
       t.fRow = nr;
@@ -1392,10 +1408,15 @@ Int_t AliTPCtrackerMI::FollowToNext(AliTPCseed& t, Int_t nr) {
       }    
     }
   }
-  if (fIteration>1) return 0;  // not look for new cluster during refitting
+  if (TMath::Abs(t.GetSnp())>AliTPCReconstructor::GetMaxSnpTracker()) return 0;  // cut on angle
+  if (fIteration>1){
+    // not look for new cluster during refitting    
+    t.fNFoundable++; 
+    return 0;
+  }
   //
   UInt_t index=0;
-  if (TMath::Abs(t.GetSnp())>0.95 || TMath::Abs(x*t.GetC()-t.GetEta())>0.95) return 0;
+  //  if (TMath::Abs(t.GetSnp())>0.95 || TMath::Abs(x*t.GetC()-t.GetEta())>0.95) return 0;// patch 28 fev 06
   Double_t  y=t.GetYat(x);
   if (TMath::Abs(y)>ymax){
     if (y > ymax) {
@@ -1429,7 +1450,8 @@ Int_t AliTPCtrackerMI::FollowToNext(AliTPCseed& t, Int_t nr) {
   } 
   else
     {
-      if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*x+10) && TMath::Abs(z)<fParam->GetZLength() ) t.fNFoundable++;
+      if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*x+10) && TMath::Abs(z)<fParam->GetZLength() && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker())) 
+       t.fNFoundable++;
       else
        return 0;
     }   
@@ -1558,9 +1580,10 @@ Bool_t AliTPCtrackerMI::GetTrackPoint(Int_t index, AliTrackPoint &p ) const
   AliTPCclusterMI *cl = GetClusterMI(index);
   if (!cl) return kFALSE;
   Int_t sector = (index&0xff000000)>>24;
-  Int_t row = (index&0x00ff0000)>>16;
+  //  Int_t row = (index&0x00ff0000)>>16;
   Float_t xyz[3];
-  xyz[0] = fParam->GetPadRowRadii(sector,row);
+  //  xyz[0] = fParam->GetPadRowRadii(sector,row);
+  xyz[0] = cl->GetX();
   xyz[1] = cl->GetY();
   xyz[2] = cl->GetZ();
   Float_t sin,cos;
@@ -1608,14 +1631,14 @@ Int_t AliTPCtrackerMI::UpdateClusters(AliTPCseed& t,  Int_t nr) {
   Double_t  ymax= GetMaxY(nr);
 
   if (row < nr) return 1; // don't prolongate if not information until now -
-  if (TMath::Abs(t.GetSnp())>0.9 && t.GetNumberOfClusters()>40. && fIteration!=2) {
-    t.fRemoval =10;
-    return 0;  // not prolongate strongly inclined tracks
-  } 
-  if (TMath::Abs(t.GetSnp())>0.95) {
-    t.fRemoval =10;
-    return 0;  // not prolongate strongly inclined tracks
-  }
+//   if (TMath::Abs(t.GetSnp())>0.9 && t.GetNumberOfClusters()>40. && fIteration!=2) {
+//     t.fRemoval =10;
+//     return 0;  // not prolongate strongly inclined tracks
+//   } 
+//   if (TMath::Abs(t.GetSnp())>0.95) {
+//     t.fRemoval =10;
+//     return 0;  // not prolongate strongly inclined tracks
+//   }// patch 28 fev 06
 
   Double_t x= GetXrow(nr);
   Double_t y,z;
@@ -1645,7 +1668,7 @@ Int_t AliTPCtrackerMI::UpdateClusters(AliTPCseed& t,  Int_t nr) {
     //y = t.GetY();    
   }
   //
-
+  if (TMath::Abs(t.GetSnp())>AliTPCReconstructor::GetMaxSnpTracker()) return 0;
   AliTPCRow &krow=GetRow(t.fRelativeSector,nr);
 
   if (TMath::Abs(TMath::Abs(y)-ymax)<krow.fDeadZone){
@@ -1655,7 +1678,7 @@ Int_t AliTPCtrackerMI::UpdateClusters(AliTPCseed& t,  Int_t nr) {
   } 
   else
     {
-      if (TMath::Abs(t.GetZ())<(AliTPCReconstructor::GetCtgRange()*t.GetX()+10)) t.fNFoundable++;
+      if (TMath::Abs(t.GetZ())<(AliTPCReconstructor::GetCtgRange()*t.GetX()+10) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker())) t.fNFoundable++;
       else
        return 0;      
     }
@@ -1833,7 +1856,7 @@ Int_t AliTPCtrackerMI::FollowBackProlongation(AliTPCseed& t, Int_t rf) {
   //
   if (first<0) first=0;
   for (Int_t nr=first; nr<=rf; nr++) {
-    if ( (TMath::Abs(t.GetSnp())>0.95)) break;
+    //    if ( (TMath::Abs(t.GetSnp())>0.95)) break;//patch 28 fev 06
     if (t.GetKinkIndexes()[0]<0){
       for (Int_t i=0;i<3;i++){
        Int_t index = t.GetKinkIndexes()[i];
@@ -2235,21 +2258,28 @@ void AliTPCtrackerMI::RemoveUsed2(TObjArray * arr, Float_t factor1,  Float_t fac
     Int_t found,foundable,shared;
     pt->GetClusterStatistic(first,last, found, foundable,shared,kFALSE);
     Float_t sharedfactor = Float_t(shared+1)/Float_t(found+1);
-    //
-    if (Float_t(shared+1)/Float_t(found+1)>factor){
-      if (pt->GetKinkIndexes()[0]!=0) continue;  //don't remove tracks  - part of the kinks
-      delete arr->RemoveAt(trackindex);
-      continue;
+    Bool_t itsgold =kFALSE;
+    if (pt->fEsd){
+      Int_t dummy[12];
+      if (pt->fEsd->GetITSclusters(dummy)>4) itsgold= kTRUE;
     }
-
-    if (pt->GetNumberOfClusters()<50&&(found-0.5*shared)<minimal){  //remove short tracks
-      if (pt->GetKinkIndexes()[0]!=0) continue;  //don't remove tracks  - part of the kinks
-      delete arr->RemoveAt(trackindex);
-      continue;
+    if (!itsgold){
+      //
+      if (Float_t(shared+1)/Float_t(found+1)>factor){
+       if (pt->GetKinkIndexes()[0]!=0) continue;  //don't remove tracks  - part of the kinks
+       delete arr->RemoveAt(trackindex);
+       continue;
+      }      
+      if (pt->GetNumberOfClusters()<50&&(found-0.5*shared)<minimal){  //remove short tracks
+       if (pt->GetKinkIndexes()[0]!=0) continue;  //don't remove tracks  - part of the kinks
+       delete arr->RemoveAt(trackindex);
+       continue;
+      }
     }
 
     good++;
     if (sharedfactor>0.4) continue;
+    if (pt->GetKinkIndexes()[0]>0) continue;
     for (Int_t i=first; i<last; i++) {
       Int_t index=pt->GetClusterIndex2(i);
       // if (index<0 || index&0x8000 ) continue;
@@ -2539,7 +2569,7 @@ Int_t AliTPCtrackerMI::RefitInward(AliESD *event)
     seed->CookdEdx(0.02,0.6);
     CookLabel(seed,0.1); //For comparison only
     //
-    if (0 && seed!=0&&esd!=0) {
+    if (AliTPCReconstructor::StreamLevel()>0 && seed!=0&&esd!=0) {
       TTreeSRedirector &cstream = *fDebugStreamer;
       cstream<<"Crefit"<<
        "Esd.="<<esd<<
@@ -2549,6 +2579,18 @@ Int_t AliTPCtrackerMI::RefitInward(AliESD *event)
     if (seed->GetNumberOfClusters()>15){
       esd->UpdateTrackParams(seed,AliESDtrack::kTPCrefit); 
       esd->SetTPCPoints(seed->GetPoints());
+      esd->SetTPCPointsF(seed->fNFoundable);
+      Int_t ndedx   = seed->fNCDEDX[0]+seed->fNCDEDX[1]+seed->fNCDEDX[2]+seed->fNCDEDX[3];
+      Float_t sdedx = (seed->fSDEDX[0]+seed->fSDEDX[1]+seed->fSDEDX[2]+seed->fSDEDX[3])*0.25;
+      Float_t dedx  = seed->GetdEdx();
+      esd->SetTPCsignal(dedx, sdedx, ndedx);
+      //
+      // add seed to the esd track in Calib level
+      //
+      if (AliTPCReconstructor::StreamLevel()>0){
+       AliTPCseed * seedCopy = new AliTPCseed(*seed, kTRUE); 
+       esd->AddCalibObject(seedCopy);
+      }
       ntracks++;
     }
     else{
@@ -2572,7 +2614,9 @@ Int_t AliTPCtrackerMI::PropagateBack(AliESD *event)
   fEvent = event;
   fIteration = 1;
   ReadSeeds(event,1);
-  PropagateBack(fSeeds);
+  PropagateBack(fSeeds); 
+  RemoveUsed2(fSeeds,0.4,0.4,20);
+  //
   Int_t nseed = fSeeds->GetEntriesFast();
   Int_t ntracks=0;
   for (Int_t i=0;i<nseed;i++){
@@ -2586,7 +2630,17 @@ Int_t AliTPCtrackerMI::PropagateBack(AliESD *event)
     if (seed->GetNumberOfClusters()>15){
       esd->UpdateTrackParams(seed,AliESDtrack::kTPCout);
       esd->SetTPCPoints(seed->GetPoints());
+      esd->SetTPCPointsF(seed->fNFoundable);
+      Int_t ndedx   = seed->fNCDEDX[0]+seed->fNCDEDX[1]+seed->fNCDEDX[2]+seed->fNCDEDX[3];
+      Float_t sdedx = (seed->fSDEDX[0]+seed->fSDEDX[1]+seed->fSDEDX[2]+seed->fSDEDX[3])*0.25;
+      Float_t dedx  = seed->GetdEdx();
+      esd->SetTPCsignal(dedx, sdedx, ndedx);
       ntracks++;
+      Int_t eventnumber = event->GetEventNumber();// patch 28 fev 06
+      (*fDebugStreamer)<<"Cback"<<
+       "Tr0.="<<seed<<
+       "EventNr="<<eventnumber<<
+       "\n"; // patch 28 fev 06   
     }
   }
   //FindKinks(fSeeds,event);
@@ -2689,20 +2743,22 @@ void AliTPCtrackerMI::ReadSeeds(AliESD *event, Int_t direction)
     }
     seed->fEsd = esd;
     // sign clusters
-    for (Int_t irow=0;irow<160;irow++){
-      Int_t index = seed->GetClusterIndex2(irow);    
-      if (index>0){ 
-       //
-       AliTPCclusterMI * cl = GetClusterMI(index);
-       seed->fClusterPointer[irow] = cl;
-       if (cl){
-         if ((index & 0x8000)==0){
-           cl->Use(10);  // accepted cluster     
+    if (esd->GetKinkIndex(0)<=0){
+      for (Int_t irow=0;irow<160;irow++){
+       Int_t index = seed->GetClusterIndex2(irow);    
+       if (index>0){ 
+         //
+         AliTPCclusterMI * cl = GetClusterMI(index);
+         seed->fClusterPointer[irow] = cl;
+         if (cl){
+           if ((index & 0x8000)==0){
+             cl->Use(10);  // accepted cluster   
+           }else{
+             cl->Use(6);   // close cluster not accepted
+           }   
          }else{
-           cl->Use(6);   // close cluster not accepted
-         }     
-       }else{
-          Info("ReadSeeds","Not found cluster");
+           Info("ReadSeeds","Not found cluster");
+         }
        }
       }
     }
@@ -4103,8 +4159,8 @@ void  AliTPCtrackerMI::FindKinks(TObjArray * array, AliESD *esd)
            track0->fCircular += 2;
          }
        }               
-       if (sign&&0){     
-         //debug stream
+       if (sign&&AliTPCReconstructor::StreamLevel()>1){          
+         //debug stream          
          cstream<<"Curling"<<
            "lab0="<<track0->fLab<<
            "lab1="<<track1->fLab<<   
@@ -4718,18 +4774,20 @@ void  AliTPCtrackerMI::FindV0s(TObjArray * array, AliESD *esd)
     if (sign[i]==0) continue;
     AliTPCseed * track0 = (AliTPCseed*)array->At(i);
     if (!track0) continue;
-    cstream<<"Tracks"<<
-      "Tr0.="<<track0<<
-      "dca="<<dca[i]<<
-      "z0="<<z0[i]<<
-      "zvertex="<<zvertex<<
-      "sdcar0="<<sdcar[i]<<
-      "cdcar0="<<cdcar[i]<<
-      "pulldcar0="<<pulldcar[i]<<
-      "pulldcaz0="<<pulldcaz[i]<<
-      "pulldca0="<<pulldca[i]<<
-      "isPrim="<<isPrim[i]<<
-      "\n";
+    if (AliTPCReconstructor::StreamLevel()>0){
+      cstream<<"Tracks"<<
+       "Tr0.="<<track0<<
+       "dca="<<dca[i]<<
+       "z0="<<z0[i]<<
+       "zvertex="<<zvertex<<
+       "sdcar0="<<sdcar[i]<<
+       "cdcar0="<<cdcar[i]<<
+       "pulldcar0="<<pulldcar[i]<<
+       "pulldcaz0="<<pulldcaz[i]<<
+       "pulldca0="<<pulldca[i]<<
+       "isPrim="<<isPrim[i]<<
+       "\n";
+    }
     //
     if (track0->fP4<0) continue;
     if (track0->GetKinkIndex(0)>0||isPrim[i]) continue;   //daughter kink
@@ -4821,7 +4879,8 @@ void  AliTPCtrackerMI::FindV0s(TObjArray * array, AliESD *esd)
        Int_t eventNr = esd->GetEventNumber();
        Double_t radiusm= (delta[0]<delta[1])? TMath::Sqrt(radius[0]):TMath::Sqrt(radius[1]);  
        Double_t deltam= (delta[0]<delta[1])? TMath::Sqrt(delta[0]):TMath::Sqrt(delta[1]);  
-       cstream<<"V0"<<
+       if (AliTPCReconstructor::StreamLevel()>0)
+         cstream<<"V0"<<
          "Event="<<eventNr<<
          "vertex.="<<&vertex<<
          "Tr0.="<<track0<<
@@ -4922,7 +4981,8 @@ void  AliTPCtrackerMI::FindV0s(TObjArray * array, AliESD *esd)
     }
     {
       Int_t eventNr = esd->GetEventNumber();
-      cstream<<"V02"<<
+      if (AliTPCReconstructor::StreamLevel()>0)
+       cstream<<"V02"<<
        "Event="<<eventNr<<
        "vertex.="<<v0<<        
        "vertex2.="<<v02<<
@@ -5734,7 +5794,11 @@ void AliTPCtrackerMI::SumTracks(TObjArray *arr1,TObjArray *arr2) const
   for (Int_t i=0;i<nseed;i++){
     AliTPCseed *pt=(AliTPCseed*)arr2->UncheckedAt(i);    
     if (pt){
-      
+       // REMOVE VERY SHORT  TRACKS
+      if (pt->GetNumberOfClusters()<20){ 
+       delete arr2->RemoveAt(i);
+       continue;
+      }// patch 28 fev06
       // NORMAL ACTIVE TRACK
       if (pt->IsActive()){
        arr1->AddLast(arr2->RemoveAt(i));
@@ -5745,11 +5809,7 @@ void AliTPCtrackerMI::SumTracks(TObjArray *arr1,TObjArray *arr2) const
        delete arr2->RemoveAt(i);
        continue;
       }
-      // REMOVE VERY SHORT  TRACKS
-      if (pt->GetNumberOfClusters()<20){ 
-       delete arr2->RemoveAt(i);
-       continue;
-      }
+     
       // ENABLE ONLY ENOUGH GOOD STOPPED TRACKS
       if (pt->GetDensityFirst(20)>0.8 || pt->GetDensityFirst(30)>0.8 || pt->GetDensityFirst(40)>0.7)
        arr1->AddLast(arr2->RemoveAt(i));
@@ -5988,7 +6048,7 @@ void  AliTPCtrackerMI::GetShape(AliTPCseed * seed, Int_t row)
   Float_t padlength =  GetPadPitchLength(row);
   //
   Float_t sresy = (seed->fSector < fParam->GetNSector()/2) ? 0.2 :0.3;
-  Float_t angulary  = seed->GetSnp();
+  Double_t angulary  = seed->GetSnp();
   angulary = angulary*angulary/(1-angulary*angulary);
   seed->fCurrentSigmaY2 = sd2+padlength*padlength*angulary/12.+sresy*sresy;  
   //