]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Possibility to calculate the DCA between two ESD track. The V0 and cascade vertexes...
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 30 Mar 2006 11:50:54 +0000 (11:50 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 30 Mar 2006 11:50:54 +0000 (11:50 +0000)
30 files changed:
ANALYSIS/AliD0toKpiAnalysis.cxx
ANALYSIS/AliD0toKpiAnalysis.h
ITS/AliCascadeVertex.cxx [deleted file]
ITS/AliCascadeVertex.h [deleted file]
ITS/AliCascadeVertexer.cxx [deleted file]
ITS/AliITStrackV2.cxx
ITS/AliITStrackV2.h
ITS/AliITStrackerMI.h
ITS/AliV0FindVertices.C [deleted file]
ITS/AliV0vertex.cxx [deleted file]
ITS/AliV0vertex.h [deleted file]
ITS/AliV0vertexer.cxx [deleted file]
ITS/ITSrecLinkDef.h
ITS/libITSrec.pkg
STEER/AliCascadeComparison.C [moved from ITS/AliCascadeComparison.C with 100% similarity]
STEER/AliCascadeVertexer.cxx [new file with mode: 0644]
STEER/AliCascadeVertexer.h [moved from ITS/AliCascadeVertexer.h with 93% similarity]
STEER/AliESDcascade.cxx
STEER/AliESDcascade.h
STEER/AliESDv0.cxx
STEER/AliESDv0.h
STEER/AliExternalTrackParam.cxx
STEER/AliExternalTrackParam.h
STEER/AliKalmanTrack.cxx
STEER/AliKalmanTrack.h
STEER/AliV0Comparison.C [moved from ITS/AliV0Comparison.C with 100% similarity]
STEER/AliV0vertexer.cxx [new file with mode: 0644]
STEER/AliV0vertexer.h [moved from ITS/AliV0vertexer.h with 94% similarity]
STEER/STEERLinkDef.h
STEER/libSTEER.pkg

index 4da29ecd1367fa1f47258165ce2b1afd9be0c701..00f6e914d16f4b7ffbf7783aae2af22c68283cb3 100644 (file)
 #include "AliESD.h"
 #include "AliStack.h"
 #include "AliRunLoader.h"
-#include "AliITStrackV2.h"
 #include "AliITSVertexerTracks.h"
 #include "AliESDVertex.h"
-#include "AliV0vertexer.h"
-#include "AliV0vertex.h"
+#include "AliESDv0.h"
 #include "AliD0toKpi.h"
 #include "AliD0toKpiAnalysis.h"
 #include "AliLog.h"
@@ -135,6 +133,7 @@ Double_t AliD0toKpiAnalysis::CalculateTOFmass(Double_t mom,Double_t length,
 
   return mom*a;
 }
+/*
 //----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::FindCandidates(Int_t evFirst,Int_t evLast,
                                        const Char_t *outName) {
@@ -365,6 +364,7 @@ void AliD0toKpiAnalysis::FindCandidates(Int_t evFirst,Int_t evLast,
 
   return;
 }
+*/
 //----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
                                           const Char_t *outName) {
@@ -387,14 +387,14 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
   Int_t    nTotEv=0,nD0rec=0,nD0rec1ev=0;
   Double_t dca;
   Double_t v2[3],mom[6],d0[2];
-  Double_t alphaP,alphaN,ptP,ptN,phiP,phiN;
+  //Double_t alphaP,alphaN,ptP,ptN,phiP,phiN;
   Int_t    iTrkP,iTrkN,trkEntries;
   Int_t    nTrksP=0,nTrksN=0;
   Int_t    trkNum[2];
   Double_t tofmass[2];
   Int_t    okD0=0,okD0bar=0;
-  AliITStrackV2 *postrack = 0;
-  AliITStrackV2 *negtrack = 0;
+  AliESDtrack *postrack = 0;
+  AliESDtrack *negtrack = 0;
 
   // create the AliITSVertexerTracks object
   // (it will be used only if fVertexOnTheFly=kTrue)
@@ -403,7 +403,7 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
   Int_t  skipped[2];
   Bool_t goodVtx1;
   
-
+  /*
   // define the cuts for vertexing
   Double_t vtxcuts[]={50., // max. allowed chi2
                      0.0, // min. allowed negative daughter's impact param 
@@ -415,6 +415,7 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
   
   // create the AliV0vertexer object
   AliV0vertexer *vertexer2 = new AliV0vertexer(vtxcuts);
+  */
 
   // create tree for reconstructed D0s
   AliD0toKpi *ioD0toKpi=0;
@@ -475,39 +476,47 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
       if(iTrkP%10==0) AliDebugClass(1,Form("  Processing positive track number %d of %d",iTrkP,nTrksP));
          
       // get track from track array
-      postrack = (AliITStrackV2*)trksP.UncheckedAt(iTrkP);
+      postrack = (AliESDtrack*)trksP.UncheckedAt(iTrkP);
       trkNum[0] = trkEntryP[iTrkP];      
 
       // loop on negative tracks 
       for(iTrkN=0; iTrkN<nTrksN; iTrkN++) {
 
        // get track from tracks array
-       negtrack = (AliITStrackV2*)trksN.UncheckedAt(iTrkN);
+       negtrack = (AliESDtrack*)trksN.UncheckedAt(iTrkN);
        trkNum[1] = trkEntryN[iTrkN];      
 
-       AliITStrackV2 nt(*negtrack), pt(*postrack), *pnt=&nt, *ppt=&pt;
-
+        {
        //
        // ----------- DCA MINIMIZATION ------------------
        //
-       // find the DCA and propagate the tracks to the DCA 
-       dca = vertexer2->PropagateToDCA(pnt,ppt);
+       // find the DCA and propagate the tracks to the DCA
+       Double_t b=event->GetMagneticField(); 
+       AliESDtrack nt(*negtrack), pt(*postrack);
+       dca = nt.PropagateToDCA(&pt,b);
 
-       // define the AliV0vertex object
-       AliV0vertex *vertex2 = new AliV0vertex(*pnt,*ppt);
+       // define the AliESDv0 object
+       AliESDv0 vertex2(nt,trkNum[0],pt,trkNum[1]);
          
        // get position of the secondary vertex
-       vertex2->GetXYZ(v2[0],v2[1],v2[2]);
-       
-       delete vertex2;
-  
+       vertex2.GetXYZ(v2[0],v2[1],v2[2]);
+        vertex2.GetPPxPyPz(mom[0],mom[1],mom[2]);
+        vertex2.GetNPxPyPz(mom[3],mom[4],mom[5]);
+       // impact parameters of the tracks w.r.t. the primary vertex
+
+       d0[0] =  10000.*pt.GetD(fV1[0],fV1[1],b);
+       d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b);
+       }
+        /*
        // momenta of the tracks at the vertex
-       ptP = 1./TMath::Abs(ppt->Get1Pt());
-       alphaP = ppt->GetAlpha();
-       phiP = alphaP+TMath::ASin(ppt->GetSnp());
+        //Double_t x,par[5]; postrack->GetExternalParameters(x,par); 
+       //ptP = 1./TMath::Abs(par[4]);
+       //alphaP = postrack->GetAlpha();
+       //phiP = alphaP+TMath::ASin(par[2]);
+         postrack->GetPxPyPz();
        mom[0] = ptP*TMath::Cos(phiP); 
        mom[1] = ptP*TMath::Sin(phiP);
-       mom[2] = ptP*ppt->GetTgl();
+       mom[2] = ptP*par[3];
          
        ptN = 1./TMath::Abs(pnt->Get1Pt());
        alphaN = pnt->GetAlpha();
@@ -515,7 +524,7 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
        mom[3] = ptN*TMath::Cos(phiN); 
        mom[4] = ptN*TMath::Sin(phiN);
        mom[5] = ptN*pnt->GetTgl();
-         
+       */
        goodVtx1 = kTRUE;
        // no vertexing if DeltaMass > fMassCut 
        if(fVertexOnTheFly) {
@@ -535,10 +544,6 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
          }
        }         
 
-       // impact parameters of the tracks w.r.t. the primary vertex
-       d0[0] =  10000.*ppt->GetD(fV1[0],fV1[1]);
-       d0[1] = -10000.*pnt->GetD(fV1[0],fV1[1]);
-
        // create the object AliD0toKpi
        AliD0toKpi theD0(ev,trkNum,fV1,v2,dca,mom,d0);
 
@@ -598,7 +603,7 @@ void AliD0toKpiAnalysis::FindCandidatesESD(Int_t evFirst,Int_t evLast,
   printf("\n+++\n+++ Total number of D0 candidates: %d\n+++\n",nD0rec);
 
   delete vertexer1;
-  delete vertexer2;
+  //delete vertexer2;
 
   esdFile->Close();
 
@@ -694,6 +699,7 @@ Bool_t AliD0toKpiAnalysis::SelectInvMass(const Double_t p[6]) const {
   if(TMath::Abs(minvD0bar-mD0) < fMassCut) return kTRUE;
   return kFALSE;
 }
+/*
 //-----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::SelectTracks(TTree &trkTree,
                  TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
@@ -730,6 +736,7 @@ void AliD0toKpiAnalysis::SelectTracks(TTree &trkTree,
 
   return;
 }
+*/
 //-----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::SelectTracksESD(AliESD &event,
         TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
@@ -744,23 +751,20 @@ void AliD0toKpiAnalysis::SelectTracksESD(AliESD &event,
   // transfer ITS tracks from ESD to arrays and to a tree
   for(Int_t i=0; i<entr; i++) {
 
-    AliESDtrack *esdtrack = (AliESDtrack*)event.GetTrack(i);
+    AliESDtrack *esdtrack = event.GetTrack(i);
     UInt_t status = esdtrack->GetStatus();
 
     if(!(status&AliESDtrack::kITSrefit)) continue;
 
-    AliITStrackV2 *itstrack = new AliITStrackV2(*esdtrack); 
-
     // single track selection
-    if(!SingleTrkCuts(*itstrack)) 
-      { delete itstrack; continue; }
+    if(!SingleTrkCuts(*esdtrack,event.GetMagneticField())) continue;
 
-    if(itstrack->Get1Pt()>0.) { // negative track
-      trksN.AddLast(itstrack);
+    if(esdtrack->GetSign()<0) { // negative track
+      trksN.AddLast(esdtrack);
       trkEntryN[nTrksN] = i;
       nTrksN++;
     } else {                 // positive track
-      trksP.AddLast(itstrack);
+      trksP.AddLast(esdtrack);
       trkEntryP[nTrksP] = i;
       nTrksP++;
     }
@@ -780,42 +784,40 @@ void AliD0toKpiAnalysis::SelectTracksESDvtx(AliESD &event,TTree *trkTree,
 
   Int_t entr = event.GetNumberOfTracks();
  
-  AliITStrackV2 *itstrackfortree = 0;
-  trkTree->Branch("tracks","AliITStrackV2",&itstrackfortree,entr,0);
+  AliESDtrack *esdtrackfortree = 0;
+  trkTree->Branch("tracks","AliESDtrack",&esdtrackfortree,entr,0);
 
 
-  // transfer ITS tracks from ESD to arrays and to a tree
+  // transfer the tracks from ESD to arrays and to a tree
   for(Int_t i=0; i<entr; i++) {
 
-    AliESDtrack *esdtrack = (AliESDtrack*)event.GetTrack(i);
+    AliESDtrack *esdtrack = event.GetTrack(i);
     UInt_t status = esdtrack->GetStatus();
 
     if(!(status&AliESDtrack::kITSrefit)) continue;
 
-    AliITStrackV2 *itstrack = new AliITStrackV2(*esdtrack);
-
     // store track in the tree to be used for primary vertex finding
-    itstrackfortree = new AliITStrackV2(*esdtrack);
+    esdtrackfortree = new AliESDtrack(*esdtrack);
     trkTree->Fill();
-    itstrackfortree = 0;
+    //itstrackfortree = 0;
+    delete esdtrackfortree;
 
     // single track selection
-    if(!SingleTrkCuts(*itstrack)) 
-      { delete itstrack; continue; }
+    if(!SingleTrkCuts(*esdtrack,event.GetMagneticField())) continue;
 
-    if(itstrack->Get1Pt()>0.) { // negative track
-      trksN.AddLast(itstrack);
+    if(esdtrack->GetSign()<0) { // negative track
+      trksN.AddLast(esdtrack);
       trkEntryN[nTrksN] = i;
       nTrksN++;
     } else {                 // positive track
-      trksP.AddLast(itstrack);
+      trksP.AddLast(esdtrack);
       trkEntryP[nTrksP] = i;
       nTrksP++;
     }
 
   } // loop on esd tracks
 
-  delete itstrackfortree;
+  //delete itstrackfortree;
 
   return;
 }
@@ -846,16 +848,19 @@ void AliD0toKpiAnalysis::SetD0Cuts(const Double_t cuts[9]) {
   return;
 }
 //-----------------------------------------------------------------------------
-Bool_t AliD0toKpiAnalysis::SingleTrkCuts(const AliITStrackV2& trk) const {
+Bool_t 
+AliD0toKpiAnalysis::SingleTrkCuts(const AliESDtrack& trk, Double_t b) const {
   // Check if track passes some kinematical cuts  
+  // Magnetic field "b" (kG)
 
-  if(TMath::Abs(1./trk.Get1Pt()) < fPtCut) 
+  if(TMath::Abs(1./trk.GetParameter()[4]) < fPtCut) 
     return kFALSE;
-  if(TMath::Abs(10000.*trk.GetD(fV1[0],fV1[1])) < fd0Cut) 
+  if(TMath::Abs(10000.*trk.GetD(fV1[0],fV1[1],b)) < fd0Cut) 
     return kFALSE;
 
   return kTRUE;
 }
+/*
 //----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::MakeTracksRefFile(Int_t evFirst,Int_t evLast) 
   const {
@@ -933,6 +938,7 @@ void AliD0toKpiAnalysis::MakeTracksRefFile(Int_t evFirst,Int_t evLast)
 
   return;
 }
+*/
 //----------------------------------------------------------------------------
 void AliD0toKpiAnalysis::MakeTracksRefFileESD() const {
   // Create a file with simulation info for the reconstructed tracks
index f906d901b67bd4d40e27ef5a1cb652ed6baf8a01..d92f5e01766464da033885befcdf261e3ea60032 100644 (file)
@@ -13,9 +13,7 @@
 #include <TString.h>
 #include <TNamed.h>
 #include "AliESD.h"
-#include "AliMagF.h"
 #include "AliTracker.h"
-#include "AliITStrackV2.h"
 
 //-----------------------------------------------------------------------------
 class AliD0toKpiAnalysis : public TNamed {
@@ -26,8 +24,8 @@ class AliD0toKpiAnalysis : public TNamed {
 
   void ApplySelection(const Char_t *inName="AliD0toKpi.root",
                      const Char_t *outName="AliD0toKpi_sele.root") const;
-  void FindCandidates(Int_t evFirst=0,Int_t evLast=0,
-                     const Char_t *outName="AliD0toKpi.root");
+  //  void FindCandidates(Int_t evFirst=0,Int_t evLast=0,
+  //                 const Char_t *outName="AliD0toKpi.root");
   void FindCandidatesESD(Int_t evFirst=0,Int_t evLast=0,
                         const Char_t *outName="AliD0toKpi.root");
   void PrintStatus() const;
@@ -73,12 +71,12 @@ class AliD0toKpiAnalysis : public TNamed {
 
   //
   Double_t CalculateTOFmass(Double_t mom,Double_t length,Double_t time) const;
-  void     MakeTracksRefFile(Int_t evFirst=0,Int_t evLast=0) const;
+  //void     MakeTracksRefFile(Int_t evFirst=0,Int_t evLast=0) const;
   void     MakeTracksRefFileESD() const;
   Bool_t   SelectInvMass(const Double_t p[6]) const;
-  void     SelectTracks(TTree &trkTree,
-                       TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
-                       TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const;
+  //void     SelectTracks(TTree &trkTree,
+  //                   TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
+  //                   TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const;
   void     SelectTracksESD(AliESD &event,
                           TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP,
                           TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const;
@@ -88,7 +86,7 @@ class AliD0toKpiAnalysis : public TNamed {
   void     SetVertex1(Double_t x=0.,Double_t y=0.,Double_t z=0.) 
     { fV1[0]=x;fV1[1]=y;fV1[2]=z; }
   void     SimulationInfo(TTree *treeD0in,TTree *treeD0out) const;
-  Bool_t   SingleTrkCuts(const AliITStrackV2& trk) const;
+  Bool_t   SingleTrkCuts(const AliESDtrack& trk, Double_t b) const;
   //
   ClassDef(AliD0toKpiAnalysis,2)  // Reconstruction of D0 candidates class
 };
diff --git a/ITS/AliCascadeVertex.cxx b/ITS/AliCascadeVertex.cxx
deleted file mode 100644 (file)
index 21d252b..0000000
+++ /dev/null
@@ -1,96 +0,0 @@
-/**************************************************************************
- * 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.                  *
- **************************************************************************/
-
-////////////////////////////////////////////////////////////////////////////
-//               Implementation of the cascade vertex class               //
-//                                                                        //
-//  Origin: Christian Kuhn, IReS, Strasbourg, christian.kuhn@ires.in2p3.fr//
-////////////////////////////////////////////////////////////////////////////
-#include <TMath.h>
-
-#include "AliCascadeVertex.h"
-#include "AliITStrackV2.h"
-#include "AliV0vertex.h"
-
-ClassImp(AliCascadeVertex)
-
-
-AliCascadeVertex::AliCascadeVertex(const AliV0vertex &v,const AliITStrackV2 &t) {
-  //--------------------------------------------------------------------
-  // Main constructor
-  //--------------------------------------------------------------------
-  fPdgCode=kXiMinus;
-
-  fV0idx[0]=v.GetNindex(); fV0idx[1]=v.GetPindex();
-  fBachIdx=t.GetLabel(); 
-
-  //Trivial estimation of the vertex parameters
-  Double_t pt, phi, x, par[5];
-  Double_t alpha, cs, sn;
-
-  t.GetExternalParameters(x,par); alpha=t.GetAlpha();
-  pt=1./TMath::Abs(par[4]);
-  phi=TMath::ASin(par[2]) + alpha;  
-
-  // momentum of the bachelor track
-
-  Double_t px1=pt*TMath::Cos(phi), py1=pt*TMath::Sin(phi), pz1=pt*par[3];
-
-  cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
-
-  Double_t x1=x*cs - par[0]*sn; // position of the bachelor at dca (bachelor,V0)
-  Double_t y1=x*sn + par[0]*cs;
-  Double_t z1=par[1];
-
-  Double_t x2,y2,z2;          // position of the V0 
-  v.GetXYZ(x2,y2,z2);
-    
-  Double_t px2,py2,pz2;       // momentum of V0
-  v.GetPxPyPz(px2,py2,pz2);
-
-  Double_t a2=((x1-x2)*px2+(y1-y2)*py2+(z1-z2)*pz2)/(px2*px2+py2*py2+pz2*pz2);
-
-  Double_t xm=x2+a2*px2;
-  Double_t ym=y2+a2*py2;
-  Double_t zm=z2+a2*pz2;
-
-  // position of the cascade decay
-  
-  fPos[0]=0.5*(x1+xm); fPos[1]=0.5*(y1+ym); fPos[2]=0.5*(z1+zm);
-    
-
-  // invariant mass of the cascade (default is Ximinus)
-  
-  Double_t e1=TMath::Sqrt(0.13957*0.13957 + px1*px1 + py1*py1 + pz1*pz1);
-  Double_t e2=TMath::Sqrt(1.11568*1.11568 + px2*px2 + py2*py2 + pz2*pz2);
-  
-  fEffMass=TMath::Sqrt((e1+e2)*(e1+e2)-
-    (px1+px2)*(px1+px2)-(py1+py2)*(py1+py2)-(pz1+pz2)*(pz1+pz2));
-
-
-  // momenta of the bachelor and the V0
-  
-  fBachMom[0]=px1; fBachMom[1]=py1; fBachMom[2]=pz1; 
-  v.GetNPxPyPz(px2,py2,pz2);
-  fV0mom[0][0]=px2; fV0mom[0][1]=py2; fV0mom[0][2]=pz2;
-  v.GetPPxPyPz(px2,py2,pz2);
-  fV0mom[1][0]=px2; fV0mom[1][1]=py2; fV0mom[1][2]=pz2;
-
-
-  fChi2=7.;   
-
-}
-
-
diff --git a/ITS/AliCascadeVertex.h b/ITS/AliCascadeVertex.h
deleted file mode 100644 (file)
index f78ed43..0000000
+++ /dev/null
@@ -1,33 +0,0 @@
-#ifndef ALICASCADEVERTEX_H
-#define ALICASCADEVERTEX_H
-
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-////////////////////////////////////////////////////////////////////////////
-//                          Cascade Vertex Class                          //
-//                                                                        //
-// Origin: Christian Kuhn, IReS, Strasbourg, christian.kuhn@ires.in2p3.fr //
-////////////////////////////////////////////////////////////////////////////
-
-#include "AliESDcascade.h"
-
-class AliITStrackV2;
-class AliV0vertex;
-
-#define kXiMinus       3312
-#define kXiPlusBar    -3312
-#define kOmegaMinus    3334
-#define kOmegaPlusBar -3334
-
-class AliCascadeVertex : public AliESDcascade {
-public:
-  AliCascadeVertex():AliESDcascade(){;}
-  AliCascadeVertex(const AliV0vertex &vtx, const AliITStrackV2 &trk);
-
-  ClassDef(AliCascadeVertex,1)   // reconstructed cascade vertex
-};
-
-#endif
-
-
diff --git a/ITS/AliCascadeVertexer.cxx b/ITS/AliCascadeVertexer.cxx
deleted file mode 100644 (file)
index 59af374..0000000
+++ /dev/null
@@ -1,393 +0,0 @@
-/**************************************************************************
- * 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.                  *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-//               Implementation of the cascade vertexer class
-//          Reads V0s and tracks, writes out cascade vertices
-//                     Fills the ESD with the cascades 
-//    Origin: Christian Kuhn, IReS, Strasbourg, christian.kuhn@ires.in2p3.fr
-//-------------------------------------------------------------------------
-#include <TObjArray.h>
-#include <TTree.h>
-
-#include "AliESD.h"
-#include "AliESDv0.h"
-#include "AliCascadeVertex.h"
-#include "AliCascadeVertexer.h"
-#include "AliITStrackV2.h"
-#include "AliV0vertex.h"
-
-ClassImp(AliCascadeVertexer)
-
-Int_t AliCascadeVertexer::V0sTracks2CascadeVertices(AliESD *event) {
-  //--------------------------------------------------------------------
-  // This function reconstructs cascade vertices
-  //      Adapted to the ESD by I.Belikov (Jouri.Belikov@cern.ch)
-  //--------------------------------------------------------------------
-
-   Int_t nV0=(Int_t)event->GetNumberOfV0s();
-   TObjArray vtcs(nV0);
-   Int_t i;
-   for (i=0; i<nV0; i++) {
-       const AliESDv0 *esdV0=event->GetV0(i);
-       vtcs.AddLast(new AliV0vertex(*esdV0));
-   }
-
-
-   Int_t ntr=(Int_t)event->GetNumberOfTracks();
-   TObjArray trks(ntr);
-   for (i=0; i<ntr; i++) {
-       AliESDtrack *esdtr=event->GetTrack(i);
-       UInt_t status=esdtr->GetStatus();
-       UInt_t flags=AliESDtrack::kITSin|AliESDtrack::kTPCin|
-                    AliESDtrack::kTPCpid|AliESDtrack::kESDpid;
-
-       if ((status&AliESDtrack::kITSrefit)==0)
-          if (flags!=status) continue;
-
-       AliITStrackV2 *iotrack=new AliITStrackV2(*esdtr);
-       iotrack->SetLabel(i);  // now it is the index in array of ESD tracks
-       if ((status&AliESDtrack::kITSrefit)==0)   //correction for the beam pipe
-          if (!iotrack->PropagateTo(3.,0.0023,65.19)) continue; 
-       if (!iotrack->PropagateTo(2.5,0.,0.)) continue;
-       trks.AddLast(iotrack);
-   }   
-   ntr=trks.GetEntriesFast();
-
-   Double_t massLambda=1.11568;
-   Int_t ncasc=0;
-
-   // Looking for the cascades...
-   for (i=0; i<nV0; i++) {
-      AliV0vertex *v=(AliV0vertex*)vtcs.UncheckedAt(i);
-      v->ChangeMassHypothesis(kLambda0); // the v0 must be Lambda 
-      if (TMath::Abs(v->GetEffMass()-massLambda)>fMassWin) continue; 
-      if (v->GetD(fX,fY,fZ)<fDV0min) continue;
-      for (Int_t j=0; j<ntr; j++) {
-        AliITStrackV2 *b=(AliITStrackV2*)trks.UncheckedAt(j);
-
-         if (b->Get1Pt()>0.) continue;  // bachelor's charge 
-         if (TMath::Abs(b->GetD(fX,fY))<fDBachMin) continue;
-          
-        AliV0vertex v0(*v), *pv0=&v0;
-         AliITStrackV2 bt(*b), *pbt=&bt;
-
-         Double_t dca=PropagateToDCA(pv0,pbt);
-         if (dca > fDCAmax) continue;
-
-         AliCascadeVertex cascade(*pv0,*pbt);
-         if (cascade.GetChi2() > fChi2max) continue;
-
-        Double_t x,y,z; cascade.GetXYZ(x,y,z); 
-         Double_t r2=x*x + y*y; 
-         if (r2 > fRmax*fRmax) continue;   // condition on fiducial zone
-         if (r2 < fRmin*fRmin) continue;
-
-         {
-         Double_t x1,y1,z1; pv0->GetXYZ(x1,y1,z1);
-         if (r2 > (x1*x1+y1*y1)) continue;
-         //if ((z-fZ)*(z-fZ) > (z1-fZ)*(z1-fZ)) continue;
-         }
-
-        Double_t px,py,pz; cascade.GetPxPyPz(px,py,pz);
-         Double_t p2=px*px+py*py+pz*pz;
-         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
-               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
-        if (cost<fCPAmax) continue; //condition on the cascade pointing angle 
-         //cascade.ChangeMassHypothesis(); //default is Xi
-
-         event->AddCascade(&cascade);
-
-         ncasc++;
-
-      }
-   }
-
-   // Looking for the anti-cascades...
-   for (i=0; i<nV0; i++) {
-      AliV0vertex *v=(AliV0vertex*)vtcs.UncheckedAt(i);
-      v->ChangeMassHypothesis(kLambda0Bar); //the v0 must be anti-Lambda 
-      if (TMath::Abs(v->GetEffMass()-massLambda)>fMassWin) continue; 
-      if (v->GetD(fX,fY,fZ)<fDV0min) continue;
-      for (Int_t j=0; j<ntr; j++) {
-        AliITStrackV2 *b=(AliITStrackV2*)trks.UncheckedAt(j);
-
-         if (TMath::Abs(b->GetD(fX,fY))<fDBachMin) continue;
-         if (b->Get1Pt()<0.) continue;  // bachelor's charge 
-          
-        AliV0vertex v0(*v), *pv0=&v0;
-         AliITStrackV2 bt(*b), *pbt=&bt;
-
-         Double_t dca=PropagateToDCA(pv0,pbt);
-         if (dca > fDCAmax) continue;
-
-         AliCascadeVertex cascade(*pv0,*pbt);
-         if (cascade.GetChi2() > fChi2max) continue;
-
-        Double_t x,y,z; cascade.GetXYZ(x,y,z); 
-         Double_t r2=x*x + y*y; 
-         if (r2 > fRmax*fRmax) continue;   // condition on fiducial zone
-         if (r2 < fRmin*fRmin) continue;
-
-         {
-         Double_t x1,y1,z1; pv0->GetXYZ(x1,y1,z1);
-         if (r2 > (x1*x1+y1*y1)) continue;
-         if (z*z > z1*z1) continue;
-         }
-
-        Double_t px,py,pz; cascade.GetPxPyPz(px,py,pz);
-         Double_t p2=px*px+py*py+pz*pz;
-         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
-               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
-
-         if (cost<fCPAmax) continue; //condition on the cascade pointing angle 
-         //cascade.ChangeMassHypothesis(); //default is Xi
-
-         event->AddCascade(&cascade);
-
-         ncasc++;
-
-      }
-   }
-
-Info("V0sTracks2CascadeVertices","Number of reconstructed cascades: %d",ncasc);
-
-   trks.Delete();
-   vtcs.Delete();
-
-   return 0;
-}
-
-Int_t AliCascadeVertexer::
-V0sTracks2CascadeVertices(TTree *vTree,TTree *tTree, TTree *xTree) {
-  //--------------------------------------------------------------------
-  // This function reconstructs cascade vertices
-  //--------------------------------------------------------------------
-  Warning("V0sTracks2CascadeVertices(TTree*,TTree*,TTree*)",
-  "Will be removed soon !  Use V0sTracks2CascadeVertices(AliESD*) instead");
-
-   TBranch *branch=vTree->GetBranch("vertices");
-   if (!branch) {
-      Error("V0sTracks2CascadeVertices","Can't get the V0 branch !");
-      return 1;
-   }   
-   Int_t nentrV0=(Int_t)vTree->GetEntries();
-   
-   TObjArray vtxV0(nentrV0);
-
-   // fill TObjArray vtxV0 with vertices
-
-   Int_t i, nV0=0;
-   for (i=0; i<nentrV0; i++) {
-
-       AliV0vertex *ioVertex=new AliV0vertex;
-       branch->SetAddress(&ioVertex);
-       vTree->GetEvent(i);
-       nV0++; 
-       vtxV0.AddLast(ioVertex);
-       
-   }
-
-   branch=tTree->GetBranch("tracks");
-   if (!branch) {
-      Error("V0sTracks2CascadeVertices","Can't get the track branch !");
-      return 2;
-   }
-   Int_t nentr=(Int_t)tTree->GetEntries();
-
-   TObjArray trks(nentr);
-
-   // fill TObjArray trks with tracks
-
-   Int_t ntrack=0;
-
-   for (i=0; i<nentr; i++) {
-
-       AliITStrackV2 *iotrack=new AliITStrackV2;
-       branch->SetAddress(&iotrack);
-       tTree->GetEvent(i);
-
-       if (!iotrack->PropagateTo(3.,0.0023,65.19)) continue; 
-       if (!iotrack->PropagateTo(2.5,0.,0.)) continue;
-
-       ntrack++; trks.AddLast(iotrack);
-       
-   }   
-
-  AliCascadeVertex *ioCascade=0;
-  branch=xTree->GetBranch("cascades");
-  if (!branch) xTree->Branch("cascades","AliCascadeVertex",&ioCascade,32000,3);
-  else branch->SetAddress(&ioCascade); 
-
-   // loop on all vertices
-
-   Double_t massLambda=1.11568;
-   Int_t ncasc=0;
-
-   for (i=0; i<nV0; i++) {
-
-       AliV0vertex *lV0ver=(AliV0vertex *)vtxV0.UncheckedAt(i);
-
-       lV0ver->ChangeMassHypothesis(kLambda0); //I.B.
-
-       if (lV0ver->GetEffMass()<massLambda-fMassWin ||       // condition of the V0 mass window (cut fMassWin)
-           lV0ver->GetEffMass()>massLambda+fMassWin) continue; 
-
-       if (lV0ver->GetD(0,0,0)<fDV0min) continue;          // condition of minimum impact parameter of the V0 (cut fDV0min) 
-                                                          // here why not cuting on pointing angle ???
-
-   // for each vertex in the good mass range, loop on all tracks (= bachelor candidates)
-
-       for (Int_t k=0; k<ntrack; k++) {
-
-         AliITStrackV2 *bachtrk=(AliITStrackV2 *)trks.UncheckedAt(k);
-
-          if (TMath::Abs(bachtrk->GetD())<fDBachMin) continue;        // eliminate to small impact parameters
-
-          if (lV0ver->GetPdgCode()==kLambda0 && bachtrk->Get1Pt()>0.) continue;     // condition on V0 label 
-          if (lV0ver->GetPdgCode()==kLambda0Bar && bachtrk->Get1Pt()<0.) continue;  // + good sign for bachelor
-          
-         AliV0vertex lV0(*lV0ver), *pV0=&lV0;
-          AliITStrackV2 bt(*bachtrk), *pbt=&bt;
-
-   // calculation of the distance of closest approach between the V0 and the bachelor
-
-          Double_t dca=PropagateToDCA(pV0,pbt);
-          if (dca > fDCAmax) continue;                         // cut on dca
-
-   // construction of a cascade object
-
-          AliCascadeVertex cascade(*pV0,*pbt);
-          if (cascade.GetChi2() > fChi2max) continue;
-
-   // get cascade decay position (V0, bachelor)
-          
-        Double_t x,y,z; cascade.GetXYZ(x,y,z); 
-         Double_t r2=x*x + y*y; 
-         if (r2 > fRmax*fRmax) continue;   // condition on fiducial zone
-         if (r2 < fRmin*fRmin) continue;
-
-         {
-   //I.B.
-         Double_t x1,y1,z1; lV0ver->GetXYZ(x1,y1,z1);
-         if (r2 > (x1*x1+y1*y1)) continue;
-         if (z*z > z1*z1) continue;
-         }
-
-   // get cascade momentum
-         
-        Double_t px,py,pz; cascade.GetPxPyPz(px,py,pz);
-         Double_t p2=px*px+py*py+pz*pz;
-         Double_t cost=(x*px+y*py+z*pz)/TMath::Sqrt(p2*(r2+z*z));
-
-         if (cost<fCPAmax) continue; // condition on the cascade pointing angle 
-
-         //cascade.ChangeMassHypothesis(); //default is Xi
-
-
-         ioCascade=&cascade; xTree->Fill();
-
-         ncasc++;
-
-      }
-   }
-
-Info("V0sTracks2CascadeVertices","Number of reconstructed cascades: %d",ncasc);
-
-   trks.Delete();
-   vtxV0.Delete();
-
-   return 0;
-}
-
-inline Double_t det(Double_t a00, Double_t a01, Double_t a10, Double_t a11){
-  // determinant 2x2
-  return a00*a11 - a01*a10;
-}
-
-inline Double_t det (Double_t a00,Double_t a01,Double_t a02,
-                     Double_t a10,Double_t a11,Double_t a12,
-                     Double_t a20,Double_t a21,Double_t a22) {
-  // determinant 3x3
-  return 
-  a00*det(a11,a12,a21,a22)-a01*det(a10,a12,a20,a22)+a02*det(a10,a11,a20,a21);
-}
-
-
-
-
-Double_t 
-AliCascadeVertexer::PropagateToDCA(AliV0vertex *v, AliITStrackV2 *t) {
-  //--------------------------------------------------------------------
-  // This function returns the DCA between the V0 and the track
-  //--------------------------------------------------------------------
-
-  Double_t phi, x, par[5];
-  Double_t alpha, cs1, sn1;
-
-  t->GetExternalParameters(x,par); alpha=t->GetAlpha();
-  phi=TMath::ASin(par[2]) + alpha;  
-  Double_t px1=TMath::Cos(phi), py1=TMath::Sin(phi), pz1=par[3];
-
-  cs1=TMath::Cos(alpha); sn1=TMath::Sin(alpha);
-  Double_t x1=x*cs1 - par[0]*sn1;
-  Double_t y1=x*sn1 + par[0]*cs1;
-  Double_t z1=par[1];
-
-  
-  Double_t x2,y2,z2;     // position and momentum of V0
-  Double_t px2,py2,pz2;
-  
-  v->GetXYZ(x2,y2,z2);
-  v->GetPxPyPz(px2,py2,pz2);
-// calculation dca
-   
-  Double_t dd= det(x2-x1,y2-y1,z2-z1,px1,py1,pz1,px2,py2,pz2);
-  Double_t ax= det(py1,pz1,py2,pz2);
-  Double_t ay=-det(px1,pz1,px2,pz2);
-  Double_t az= det(px1,py1,px2,py2);
-
-  Double_t dca=TMath::Abs(dd)/TMath::Sqrt(ax*ax + ay*ay + az*az);
-
-//points of the DCA
-  Double_t t1 = det(x2-x1,y2-y1,z2-z1,px2,py2,pz2,ax,ay,az)/
-                det(px1,py1,pz1,px2,py2,pz2,ax,ay,az);
-  
-  x1 += px1*t1; y1 += py1*t1; //z1 += pz1*t1;
-  
-
-  //propagate track to the points of DCA
-
-  x1=x1*cs1 + y1*sn1;
-  if (!t->PropagateTo(x1,0.,0.)) {
-    Error("PropagateToDCA","Propagation failed !");
-    return 1.e+33;
-  }  
-
-  return dca;
-}
-
-
-
-
-
-
-
-
-
-
-
index ae042babdd838b3ab2e4c087b1ab0415b047942f..36d4e7c9af792e9856d90f73fb7d799b44f7473d 100644 (file)
@@ -734,3 +734,30 @@ void AliITStrackV2::CookdEdx(Double_t low, Double_t up) {
 
   SetdEdx(dedx);
 }
+
+Double_t AliITStrackV2::
+PropagateToDCA(AliKalmanTrack *p, Double_t d, Double_t x0) {
+  //--------------------------------------------------------------
+  // Propagates this track and the argument track to the position of the
+  // distance of closest approach. 
+  // Returns the (weighed !) distance of closest approach.
+  //--------------------------------------------------------------
+  Double_t xthis, xp, dca;
+  {
+  //Temporary solution
+  Double_t b=1./GetLocalConvConst()/kB2C;
+  AliExternalTrackParam dummy1(*this), dummy2(*p); 
+  dca=dummy1.GetDCA(&dummy2,b,xthis,xp);
+  }
+  if (!PropagateTo(xthis,d,x0)) {
+    //AliWarning(" propagation failed !");
+    return 1e+33;
+  }  
+
+  if (!p->PropagateTo(xp,d,x0)) {
+    //AliWarning(" propagation failed !";
+    return 1e+33;
+  }  
+
+  return dca;
+} 
index 8a766ec46b7af117be7dc3fba87d929f164975cb..c8441d6da6b2583116d683cbdb67e1f76ea5eb06 100644 (file)
@@ -42,6 +42,7 @@ public:
   Int_t Propagate(Double_t alpha, Double_t xr);
   virtual Int_t CorrectForMaterial(Double_t d, Double_t x0=21.82);
   Int_t PropagateTo(Double_t xr, Double_t d, Double_t x0=21.82);
+  Double_t PropagateToDCA(AliKalmanTrack *p, Double_t d=0., Double_t x0=0.); 
   Int_t Update(const AliCluster* cl,Double_t chi2,UInt_t i);
   Int_t Improve(Double_t x0,Double_t xyz[3],Double_t ers[3]);
   void SetdEdx(Double_t dedx) {fdEdx=dedx;}
index dbdd8154822fb54db50b82bafb0c784f24a9c25c..48cf2a90f3e34f26aeec5bcc9ed86c27c04247e6 100644 (file)
@@ -16,13 +16,11 @@ class TTreeSRedirector;
 class AliESD;
 class AliHelix;
 class AliITSgeom;
-class AliV0vertex;
 
 #include <TObjArray.h>
 #include "AliITSRecPoint.h"
 #include "AliITStrackMI.h"
 #include "AliTracker.h"
-#include "AliV0vertex.h"
 
 //-------------------------------------------------------------------------
 class AliITStrackerMI : public AliTracker {
diff --git a/ITS/AliV0FindVertices.C b/ITS/AliV0FindVertices.C
deleted file mode 100644 (file)
index 5c917ec..0000000
+++ /dev/null
@@ -1,97 +0,0 @@
-/****************************************************************************
- *           Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch                 *
- ****************************************************************************/
-
-#if !defined(__CINT__) || defined(__MAKECINT__)
-  #include "Riostream.h"
-  #include "AliV0vertexer.h"
-  #include "TFile.h"
-  #include "TKey.h"
-  #include "TStopwatch.h"
-
-  #include "AliRun.h"
-  #include "AliTracker.h"
-  #include "AliMagF.h"
-  #include "AliESD.h"
-  #include "AliRunLoader.h"
-#endif
-
-extern AliRun *gAlice;
-
-Int_t AliV0FindVertices(Int_t nev=5) {
-   cerr<<"Looking for V0 vertices...\n";
-
-   if (gAlice) {
-      delete gAlice->GetRunLoader();
-      delete gAlice; 
-      gAlice=0;
-   } 
-
-   AliRunLoader* rl = AliRunLoader::Open("galice.root");
-   if (rl == 0x0) {
-      cerr<<"AliV0FindVertices.C : Can not open session RL=NULL"<< endl;
-      return 1;
-   }
-
-   if (rl->LoadgAlice()) {
-      cerr<<"AliV0FindVertices.C : LoadgAlice returned error"<<endl;
-      delete rl;
-      return 3;
-   }
-
-   // Magnetic field
-   AliTracker::SetFieldMap(gAlice->Field(),1); // 1 means uniform magnetic field
-       
-   Double_t cuts[]={33,  // max. allowed chi2
-                    0.16,// min. allowed negative daughter's impact parameter 
-                    0.05,// min. allowed positive daughter's impact parameter 
-                    0.080,// max. allowed DCA between the daughter tracks
-                    0.998,// max. allowed cosine of V0's pointing angle
-                    0.9,  // min. radius of the fiducial volume
-                    2.9   // max. radius of the fiducial volume
-                   };
-   TStopwatch timer;
-   AliV0vertexer vtxer(cuts);
-   Int_t rc=0;
-   if (nev>rl->GetNumberOfEvents()) nev=rl->GetNumberOfEvents();
-
-   TFile *v0f=TFile::Open("AliESDv0.root","RECREATE");
-   if ((!v0f)||(!v0f->IsOpen())) {
-      cerr<<"Can't AliESDv0.root !\n"; return 1;
-   }
-   TFile *itsf=TFile::Open("AliESDits.root");
-   if ((!itsf)||(!itsf->IsOpen())) {
-      cerr<<"Can't AliESDits.root !\n"; return 1;
-   }
-   TKey *key=0;
-   TIter next(itsf->GetListOfKeys());
-   for (Int_t i=0; i<nev; i++) {
-     itsf->cd();
-     if ((key=(TKey*)next())==0) break;
-     cerr<<"Processing event number: "<<i<<endl;
-     AliESD *event=(AliESD*)key->ReadObj();
-
-     //Double_t vtx[3]={0.,0.,0.}; vtxer.SetVertex(vtx); // primary vertex (cm)
-
-     rc=vtxer.Tracks2V0vertices(event);
-
-     if (rc==0) {
-        Char_t ename[100]; 
-        sprintf(ename,"%d",i);
-        v0f->cd();
-        if (!event->Write(ename)) rc++;
-     } 
-     if (rc) {
-        cerr<<"Something bad happened...\n";
-     }
-     delete event;
-   }
-   timer.Stop(); timer.Print();
-    
-   itsf->Close();
-   v0f->Close();
-
-   delete rl;
-
-   return rc;
-}
diff --git a/ITS/AliV0vertex.cxx b/ITS/AliV0vertex.cxx
deleted file mode 100644 (file)
index e81b580..0000000
+++ /dev/null
@@ -1,81 +0,0 @@
-/**************************************************************************
- * 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.                  *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-//               Implementation of the V0 vertex class
-//
-//     Origin: Iouri Belikov, IReS, Strasbourg, Jouri.Belikov@cern.ch
-//-------------------------------------------------------------------------
-#include <TMath.h>
-
-#include "AliV0vertex.h"
-#include "AliITStrackV2.h"
-
-ClassImp(AliV0vertex)
-
-AliV0vertex::AliV0vertex(const AliITStrackV2 &n, const AliITStrackV2 &p) {
-  //--------------------------------------------------------------------
-  // Main constructor
-  //--------------------------------------------------------------------
-  for (Int_t i=0; i<6; i++) {
-    fPosCov[i]= 0.;
-    fNmomCov[i] = 0.;
-    fPmomCov[i] = 0.;
-  }
-  fPdgCode=kK0Short;
-  fNidx=n.GetLabel(); fPidx=p.GetLabel(); //indices in the array of ESD tracks
-
-  //Trivial estimation of the vertex parameters
-  Double_t pt, phi, x, par[5];
-  Double_t alpha, cs, sn;
-
-  n.GetExternalParameters(x,par); alpha=n.GetAlpha();
-  pt=1./TMath::Abs(par[4]);
-  phi=TMath::ASin(par[2]) + alpha;  
-  Double_t px1=pt*TMath::Cos(phi), py1=pt*TMath::Sin(phi), pz1=pt*par[3];
-  cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
-  Double_t x1=x*cs - par[0]*sn;
-  Double_t y1=x*sn + par[0]*cs;
-  Double_t z1=par[1];
-  Double_t sx1=sn*sn*n.GetSigmaY2(), sy1=cs*cs*n.GetSigmaY2(); 
-
-  p.GetExternalParameters(x,par); alpha=p.GetAlpha();
-  pt=1./TMath::Abs(par[4]);
-  phi=TMath::ASin(par[2]) + alpha;  
-  Double_t px2=pt*TMath::Cos(phi), py2=pt*TMath::Sin(phi), pz2=pt*par[3];
-  cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
-  Double_t x2=x*cs - par[0]*sn;
-  Double_t y2=x*sn + par[0]*cs;
-  Double_t z2=par[1];
-  Double_t sx2=sn*sn*p.GetSigmaY2(), sy2=cs*cs*p.GetSigmaY2(); 
-    
-  Double_t sz1=n.GetSigmaZ2(), sz2=p.GetSigmaZ2();
-  Double_t wx1=sx2/(sx1+sx2), wx2=1.- wx1;
-  Double_t wy1=sy2/(sy1+sy2), wy2=1.- wy1;
-  Double_t wz1=sz2/(sz1+sz2), wz2=1.- wz1;
-  fPos[0]=wx1*x1 + wx2*x2; fPos[1]=wy1*y1 + wy2*y2; fPos[2]=wz1*z1 + wz2*z2;
-
-  //fPos[0]=0.5*(x1+x2); fPos[1]=0.5*(y1+y2); fPos[2]=0.5*(z1+z2);
-  fNmom[0]=px1; fNmom[1]=py1; fNmom[2]=pz1; 
-  fPmom[0]=px2; fPmom[1]=py2; fPmom[2]=pz2;
-
-  Double_t e1=TMath::Sqrt(0.13957*0.13957 + px1*px1 + py1*py1 + pz1*pz1);
-  Double_t e2=TMath::Sqrt(0.13957*0.13957 + px2*px2 + py2*py2 + pz2*pz2);
-  fEffMass=TMath::Sqrt((e1+e2)*(e1+e2)-
-    (px1+px2)*(px1+px2)-(py1+py2)*(py1+py2)-(pz1+pz2)*(pz1+pz2));
-
-  fChi2=7.;   
-}
-
diff --git a/ITS/AliV0vertex.h b/ITS/AliV0vertex.h
deleted file mode 100644 (file)
index 9391c6d..0000000
+++ /dev/null
@@ -1,28 +0,0 @@
-#ifndef ALIV0VERTEX_H
-#define ALIV0VERTEX_H
-
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-//-------------------------------------------------------------------------
-//                          V0 Vertex Class
-//
-//    Origin: Iouri Belikov, IReS, Strasbourg, Jouri.Belikov@cern.ch 
-//-------------------------------------------------------------------------
-
-#include <AliESDv0.h>
-
-class AliITStrackV2;
-
-class AliV0vertex : public AliESDv0 {
-public:
-  AliV0vertex() : AliESDv0() {;}
-  AliV0vertex(const AliESDv0 &v) : AliESDv0(v) {;}
-  AliV0vertex(const AliITStrackV2 &neg, const AliITStrackV2 &pos);
-
-  ClassDef(AliV0vertex,1) // reconstructed V0 vertex
-};
-
-#endif
-
-
diff --git a/ITS/AliV0vertexer.cxx b/ITS/AliV0vertexer.cxx
deleted file mode 100644 (file)
index 9e7c23d..0000000
+++ /dev/null
@@ -1,258 +0,0 @@
-/**************************************************************************
- * 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.                  *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-//               Implementation of the V0 vertexer class
-//                  reads tracks writes out V0 vertices
-//                      fills the ESD with the V0s       
-//     Origin: Iouri Belikov, IReS, Strasbourg, Jouri.Belikov@cern.ch
-//-------------------------------------------------------------------------
-#include <TObjArray.h>
-#include <TTree.h>
-
-#include "AliESD.h"
-#include "AliESDtrack.h"
-#include "AliITStrackV2.h"
-#include "AliV0vertex.h"
-#include "AliV0vertexer.h"
-
-ClassImp(AliV0vertexer)
-
-Int_t AliV0vertexer::Tracks2V0vertices(AliESD *event) {
-  //--------------------------------------------------------------------
-  //This function reconstructs V0 vertices
-  //--------------------------------------------------------------------
-
-   Int_t nentr=event->GetNumberOfTracks();
-
-   TObjArray negtrks(nentr/2);
-   TObjArray postrks(nentr/2);
-
-   Int_t nneg=0, npos=0, nvtx=0;
-
-   Int_t i;
-   for (i=0; i<nentr; i++) {
-     AliESDtrack *esd=event->GetTrack(i);
-     UInt_t status=esd->GetStatus();
-     UInt_t flags=AliESDtrack::kITSin|AliESDtrack::kTPCin|
-                  AliESDtrack::kTPCpid|AliESDtrack::kESDpid;
-
-     if ((status&AliESDtrack::kITSrefit)==0)
-        if (flags!=status) continue;
-
-     AliITStrackV2 *iotrack=new AliITStrackV2(*esd);
-     iotrack->SetLabel(i);  // now it is the index in array of ESD tracks
-     if ((status&AliESDtrack::kITSrefit)==0)   //correction for the beam pipe
-        if (!iotrack->PropagateTo(3.,0.0023,65.19)) {
-         delete iotrack;
-         continue;
-       }
-     if (!iotrack->PropagateTo(2.5,0.,0.)) {
-       delete iotrack;
-       continue;
-     }
-
-     if (iotrack->Get1Pt() < 0.) {nneg++; negtrks.AddLast(iotrack);}
-     else {npos++; postrks.AddLast(iotrack);}
-   }   
-
-
-   for (i=0; i<nneg; i++) {
-      //if (i%10==0) cerr<<nneg-i<<'\r';
-      AliITStrackV2 *ntrk=(AliITStrackV2 *)negtrks.UncheckedAt(i);
-
-      if (TMath::Abs(ntrk->GetD(fX,fY))<fDPmin) continue;
-      if (TMath::Abs(ntrk->GetD(fX,fY))>fRmax) continue;
-
-      for (Int_t k=0; k<npos; k++) {
-         AliITStrackV2 *ptrk=(AliITStrackV2 *)postrks.UncheckedAt(k);
-
-         if (TMath::Abs(ptrk->GetD(fX,fY))<fDPmin) continue;
-         if (TMath::Abs(ptrk->GetD(fX,fY))>fRmax) continue;
-
-         if (TMath::Abs(ntrk->GetD(fX,fY))<fDNmin)
-         if (TMath::Abs(ptrk->GetD(fX,fY))<fDNmin) continue;
-
-
-         AliITStrackV2 nt(*ntrk), pt(*ptrk), *pnt=&nt, *ppt=&pt;
-
-         Double_t dca=PropagateToDCA(pnt,ppt);
-         if (dca > fDCAmax) continue;
-
-         AliV0vertex vertex(*pnt,*ppt);
-         if (vertex.GetChi2() > fChi2max) continue;
-        
-         /*  Think of something better here ! 
-         nt.PropagateToVertex(); if (TMath::Abs(nt.GetZ())<0.04) continue;
-         pt.PropagateToVertex(); if (TMath::Abs(pt.GetZ())<0.04) continue;
-        */
-
-         Double_t x,y,z; vertex.GetXYZ(x,y,z); 
-         Double_t r2=x*x + y*y; 
-         if (r2 > fRmax*fRmax) continue;
-         if (r2 < fRmin*fRmin) continue;
-
-         Double_t px,py,pz; vertex.GetPxPyPz(px,py,pz);
-         Double_t p2=px*px+py*py+pz*pz;
-         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
-               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
-
-        //if (cost < (5*fCPAmax-0.9-TMath::Sqrt(r2)*(fCPAmax-1))/4.1) continue;
-         if (cost < fCPAmax) continue;
-        vertex.SetDcaDaughters(dca);
-         //vertex.ChangeMassHypothesis(); //default is Lambda0 
-
-         event->AddV0(&vertex);
-
-         nvtx++;
-      }
-   }
-
-   Info("Tracks2V0vertices","Number of reconstructed V0 vertices: %d",nvtx);
-
-   negtrks.Delete();
-   postrks.Delete();
-
-   return 0;
-}
-
-Int_t AliV0vertexer::Tracks2V0vertices(TTree *tTree, TTree *vTree) {
-  //--------------------------------------------------------------------
-  //This function reconstructs V0 vertices
-  //--------------------------------------------------------------------
-  Warning("Tracks2V0vertices(TTree*,TTree*)",
-         "Will be removed soon !  Use Tracks2V0vertices(AliESD*) instead");
-
-   TBranch *branch=tTree->GetBranch("tracks");
-   if (!branch) {
-      Error("Tracks2V0vertices","Can't get the branch !");
-      return 1;
-   }
-   Int_t nentr=(Int_t)tTree->GetEntries();
-
-   TObjArray negtrks(nentr/2);
-   TObjArray postrks(nentr/2);
-
-   Int_t nneg=0, npos=0, nvtx=0;
-
-   Int_t i;
-   for (i=0; i<nentr; i++) {
-       AliITStrackV2 *iotrack=new AliITStrackV2;
-       branch->SetAddress(&iotrack);
-       tTree->GetEvent(i);
-
-       if (!iotrack->PropagateTo(3.,0.0023,65.19)) {
-        delete iotrack;
-        continue; 
-       }
-       if (!iotrack->PropagateTo(2.5,0.,0.)) {
-        delete iotrack;
-        continue;
-       }
-
-       if (iotrack->Get1Pt() > 0.) {nneg++; negtrks.AddLast(iotrack);}
-       else {npos++; postrks.AddLast(iotrack);}
-   }   
-
-
-   AliV0vertex *ioVertex=0;
-   branch=vTree->GetBranch("vertices");
-   if (!branch) vTree->Branch("vertices","AliV0vertex",&ioVertex,32000,3);
-   else branch->SetAddress(&ioVertex);
-
-
-   for (i=0; i<nneg; i++) {
-      //if (i%10==0) cerr<<nneg-i<<'\r';
-      AliITStrackV2 *ntrk=(AliITStrackV2 *)negtrks.UncheckedAt(i);
-
-      if (TMath::Abs(ntrk->GetD(fX,fY))<fDPmin) continue;
-      if (TMath::Abs(ntrk->GetD(fX,fY))>fRmax) continue;
-
-      for (Int_t k=0; k<npos; k++) {
-         AliITStrackV2 *ptrk=(AliITStrackV2 *)postrks.UncheckedAt(k);
-
-         if (TMath::Abs(ptrk->GetD(fX,fY))<fDPmin) continue;
-         if (TMath::Abs(ptrk->GetD(fX,fY))>fRmax) continue;
-
-         if (TMath::Abs(ntrk->GetD(fX,fY))<fDNmin)
-         if (TMath::Abs(ptrk->GetD(fX,fY))<fDNmin) continue;
-
-
-         AliITStrackV2 nt(*ntrk), pt(*ptrk), *pnt=&nt, *ppt=&pt;
-
-         Double_t dca=PropagateToDCA(pnt,ppt);
-         if (dca > fDCAmax) continue;
-
-         AliV0vertex vertex(*pnt,*ppt);
-         if (vertex.GetChi2() > fChi2max) continue;
-        
-         /*  Think of something better here ! 
-         nt.PropagateToVertex(); if (TMath::Abs(nt.GetZ())<0.04) continue;
-         pt.PropagateToVertex(); if (TMath::Abs(pt.GetZ())<0.04) continue;
-        */
-
-         Double_t x,y,z; vertex.GetXYZ(x,y,z); 
-         Double_t r2=x*x + y*y; 
-         if (r2 > fRmax*fRmax) continue;
-         if (r2 < fRmin*fRmin) continue;
-
-         Double_t px,py,pz; vertex.GetPxPyPz(px,py,pz);
-         Double_t p2=px*px+py*py+pz*pz;
-         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
-               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
-
-        //if (cost < (5*fCPAmax-0.9-TMath::Sqrt(r2)*(fCPAmax-1))/4.1) continue;
-         if (cost < fCPAmax) continue;
-        vertex.SetDcaDaughters(dca);
-         //vertex.ChangeMassHypothesis(); //default is Lambda0 
-
-         ioVertex=&vertex; vTree->Fill();
-
-         nvtx++;
-      }
-   }
-
-   Info("Tracks2V0vertices","Number of reconstructed V0 vertices: %d",nvtx);
-
-   negtrks.Delete();
-   postrks.Delete();
-
-   return 0;
-}
-
-Double_t 
-AliV0vertexer::PropagateToDCA(AliITStrackV2 *n, AliITStrackV2 *p) const {
-  //--------------------------------------------------------------------
-  // This function returns the DCA between two tracks
-  // The tracks will be moved to the point of DCA ! 
-  //--------------------------------------------------------------------
-  return n->PropagateToDCA(p);
-}
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
index c0656ba2ed59142816696db501eca3eca1fcca61..a9d7fc7d979a065d162f9583cb3323e872a54c1f 100644 (file)
 #pragma link C++ class AliITStrackMI+;
 #pragma link C++ class AliITStrackerMI+;
 //#pragma link C++ class AliITSRecV0Info+;
-#pragma link C++ class  AliV0vertex+;
-#pragma link C++ class  AliV0vertexer+;
-#pragma link C++ class  AliCascadeVertex+;
-#pragma link C++ class  AliCascadeVertexer+;
 
 #pragma link C++ class  AliITSVertexer+;
 #pragma link C++ class  AliITSVertexerIons+;
index 884c0900e69cca93c60c769e08be96d07eff475b..fa1c5e0c3b6cf4f135a4d9ab162f8575901b20f8 100644 (file)
@@ -25,12 +25,8 @@ SRCS =       AliITSDetTypeRec.cxx \
                AliITSVertexerTracks.cxx \
                AliITSVertexerZ.cxx \
                AliITSVertexerFast.cxx \
-               AliV0vertex.cxx \
-               AliV0vertexer.cxx \
                AliITSPid.cxx \
                AliITStrackV2Pid.cxx \
-               AliCascadeVertex.cxx \
-               AliCascadeVertexer.cxx \
                AliITSreconstruction.cxx \
                AliITSFindClustersV2.cxx \
                AliITSRiemannFit.cxx \
@@ -63,5 +59,5 @@ DHDR=ITSrecLinkDef.h
 
 EINCLUDE:=$(ALICE)/geant3/TGeant3 TPC RAW
 
-EXPORT:=AliITStrackV2.h AliITSVertexerTracks.h AliV0vertexer.h \
-        AliV0vertex.h AliITSVertexer.h AliITSrecoV2.h
+EXPORT:=AliITStrackV2.h AliITSVertexerTracks.h \
+        AliITSVertexer.h AliITSrecoV2.h
diff --git a/STEER/AliCascadeVertexer.cxx b/STEER/AliCascadeVertexer.cxx
new file mode 100644 (file)
index 0000000..7451457
--- /dev/null
@@ -0,0 +1,233 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+
+//-------------------------------------------------------------------------
+//               Implementation of the cascade vertexer class
+//          Reads V0s and tracks, writes out cascade vertices
+//                     Fills the ESD with the cascades 
+//    Origin: Christian Kuhn, IReS, Strasbourg, christian.kuhn@ires.in2p3.fr
+//-------------------------------------------------------------------------
+#include <TObjArray.h>
+#include <TTree.h>
+
+#include "AliESD.h"
+#include "AliESDv0.h"
+#include "AliESDcascade.h"
+#include "AliCascadeVertexer.h"
+
+ClassImp(AliCascadeVertexer)
+
+Int_t AliCascadeVertexer::V0sTracks2CascadeVertices(AliESD *event) {
+  //--------------------------------------------------------------------
+  // This function reconstructs cascade vertices
+  //      Adapted to the ESD by I.Belikov (Jouri.Belikov@cern.ch)
+  //--------------------------------------------------------------------
+   Double_t b=event->GetMagneticField();
+   Int_t nV0=(Int_t)event->GetNumberOfV0s();
+   TObjArray vtcs(nV0);
+   Int_t i;
+   for (i=0; i<nV0; i++) {
+       AliESDv0 *v=event->GetV0(i);
+       if (v->GetD(fX,fY,fZ)<fDV0min) continue;
+       vtcs.AddLast(v);
+   }
+   nV0=vtcs.GetEntriesFast();
+
+   Int_t nentr=(Int_t)event->GetNumberOfTracks();
+   TArrayI trk(nentr); Int_t ntr=0;
+   for (i=0; i<nentr; i++) {
+       AliESDtrack *esdtr=event->GetTrack(i);
+       UInt_t status=esdtr->GetStatus();
+       UInt_t flags=AliESDtrack::kITSin|AliESDtrack::kTPCin|
+                    AliESDtrack::kTPCpid|AliESDtrack::kESDpid;
+
+       if ((status&AliESDtrack::kITSrefit)==0)
+          if (flags!=status) continue;
+
+       if (TMath::Abs(esdtr->GetD(fX,fY,b))<fDBachMin) continue;
+
+       trk[ntr++]=i;
+   }   
+
+   Double_t massLambda=1.11568;
+   Int_t ncasc=0;
+
+   // Looking for the cascades...
+   for (i=0; i<nV0; i++) {
+      AliESDv0 *v=(AliESDv0*)vtcs.UncheckedAt(i);
+      v->ChangeMassHypothesis(kLambda0); // the v0 must be Lambda 
+      if (TMath::Abs(v->GetEffMass()-massLambda)>fMassWin) continue; 
+      for (Int_t j=0; j<ntr; j++) {
+        Int_t bidx=trk[j];
+        AliESDtrack *btrk=event->GetTrack(bidx);
+
+         if (btrk->GetSign()>0) continue;  // bachelor's charge 
+          
+        AliESDv0 v0(*v), *pv0=&v0;
+         AliExternalTrackParam bt(*btrk), *pbt=&bt;
+
+         Double_t dca=PropagateToDCA(pv0,pbt,b);
+         if (dca > fDCAmax) continue;
+
+         AliESDcascade cascade(*pv0,*pbt,bidx);
+         if (cascade.GetChi2() > fChi2max) continue;
+
+        Double_t x,y,z; cascade.GetXYZ(x,y,z); 
+         Double_t r2=x*x + y*y; 
+         if (r2 > fRmax*fRmax) continue;   // condition on fiducial zone
+         if (r2 < fRmin*fRmin) continue;
+
+         {
+         Double_t x1,y1,z1; pv0->GetXYZ(x1,y1,z1);
+         if (r2 > (x1*x1+y1*y1)) continue;
+         //if ((z-fZ)*(z-fZ) > (z1-fZ)*(z1-fZ)) continue;
+         }
+
+        Double_t px,py,pz; cascade.GetPxPyPz(px,py,pz);
+         Double_t p2=px*px+py*py+pz*pz;
+         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
+               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
+        if (cost<fCPAmax) continue; //condition on the cascade pointing angle 
+         //cascade.ChangeMassHypothesis(); //default is Xi
+
+         event->AddCascade(&cascade);
+
+         ncasc++;
+
+      }
+   }
+
+   // Looking for the anti-cascades...
+   for (i=0; i<nV0; i++) {
+      AliESDv0 *v=(AliESDv0*)vtcs.UncheckedAt(i);
+      v->ChangeMassHypothesis(kLambda0Bar); //the v0 must be anti-Lambda 
+      if (TMath::Abs(v->GetEffMass()-massLambda)>fMassWin) continue; 
+      for (Int_t j=0; j<ntr; j++) {
+        Int_t bidx=trk[j];
+        AliESDtrack *btrk=event->GetTrack(bidx);
+
+         if (btrk->GetSign()<0) continue;  // bachelor's charge 
+          
+        AliESDv0 v0(*v), *pv0=&v0;
+         AliESDtrack bt(*btrk), *pbt=&bt;
+
+         Double_t dca=PropagateToDCA(pv0,pbt,b);
+         if (dca > fDCAmax) continue;
+
+         AliESDcascade cascade(*pv0,*pbt,bidx);
+         if (cascade.GetChi2() > fChi2max) continue;
+
+        Double_t x,y,z; cascade.GetXYZ(x,y,z); 
+         Double_t r2=x*x + y*y; 
+         if (r2 > fRmax*fRmax) continue;   // condition on fiducial zone
+         if (r2 < fRmin*fRmin) continue;
+
+         {
+         Double_t x1,y1,z1; pv0->GetXYZ(x1,y1,z1);
+         if (r2 > (x1*x1+y1*y1)) continue;
+         if (z*z > z1*z1) continue;
+         }
+
+        Double_t px,py,pz; cascade.GetPxPyPz(px,py,pz);
+         Double_t p2=px*px+py*py+pz*pz;
+         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
+               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
+
+         if (cost<fCPAmax) continue; //condition on the cascade pointing angle 
+         //cascade.ChangeMassHypothesis(); //default is Xi
+
+         event->AddCascade(&cascade);
+
+         ncasc++;
+
+      }
+   }
+
+Info("V0sTracks2CascadeVertices","Number of reconstructed cascades: %d",ncasc);
+
+   return 0;
+}
+
+
+inline Double_t det(Double_t a00, Double_t a01, Double_t a10, Double_t a11){
+  // determinant 2x2
+  return a00*a11 - a01*a10;
+}
+
+inline Double_t det (Double_t a00,Double_t a01,Double_t a02,
+                     Double_t a10,Double_t a11,Double_t a12,
+                     Double_t a20,Double_t a21,Double_t a22) {
+  // determinant 3x3
+  return 
+  a00*det(a11,a12,a21,a22)-a01*det(a10,a12,a20,a22)+a02*det(a10,a11,a20,a21);
+}
+
+
+
+
+Double_t AliCascadeVertexer::
+PropagateToDCA(AliESDv0 *v, AliExternalTrackParam *t, Double_t b) {
+  //--------------------------------------------------------------------
+  // This function returns the DCA between the V0 and the track
+  //--------------------------------------------------------------------
+  Double_t alpha=t->GetAlpha(), cs1=TMath::Cos(alpha), sn1=TMath::Sin(alpha);
+  Double_t r[3]; t->GetXYZ(r);
+  Double_t x1=r[0], y1=r[1], z1=r[2];
+  Double_t p[3]; t->GetPxPyPz(p);
+  Double_t px1=p[0], py1=p[1], pz1=p[2];
+  
+  Double_t x2,y2,z2;     // position and momentum of V0
+  Double_t px2,py2,pz2;
+  
+  v->GetXYZ(x2,y2,z2);
+  v->GetPxPyPz(px2,py2,pz2);
+// calculation dca
+   
+  Double_t dd= det(x2-x1,y2-y1,z2-z1,px1,py1,pz1,px2,py2,pz2);
+  Double_t ax= det(py1,pz1,py2,pz2);
+  Double_t ay=-det(px1,pz1,px2,pz2);
+  Double_t az= det(px1,py1,px2,py2);
+
+  Double_t dca=TMath::Abs(dd)/TMath::Sqrt(ax*ax + ay*ay + az*az);
+
+//points of the DCA
+  Double_t t1 = det(x2-x1,y2-y1,z2-z1,px2,py2,pz2,ax,ay,az)/
+                det(px1,py1,pz1,px2,py2,pz2,ax,ay,az);
+  
+  x1 += px1*t1; y1 += py1*t1; //z1 += pz1*t1;
+  
+
+  //propagate track to the points of DCA
+
+  x1=x1*cs1 + y1*sn1;
+  if (!t->PropagateTo(x1,b)) {
+    Error("PropagateToDCA","Propagation failed !");
+    return 1.e+33;
+  }  
+
+  return dca;
+}
+
+
+
+
+
+
+
+
+
+
+
similarity index 93%
rename from ITS/AliCascadeVertexer.h
rename to STEER/AliCascadeVertexer.h
index a88e392c820e642cd5b36877f734f8446a2daa51..c44d9a76e3a3082b6263fdd85d8b5db3499f2eef 100644 (file)
@@ -13,8 +13,8 @@
 
 class AliESD;
 class TTree;
-class AliITStrackV2;
-class AliV0vertex;
+class AliESDv0;
+class AliExternalTrackParam;
 
 //_____________________________________________________________________________
 class AliCascadeVertexer : public TObject {
@@ -25,8 +25,7 @@ public:
   void SetVertex(Double_t *vtx) { fX=vtx[0]; fY=vtx[1]; fZ=vtx[2]; }
 
   Int_t V0sTracks2CascadeVertices(AliESD *event);
-  Int_t V0sTracks2CascadeVertices(TTree *v, TTree *t, TTree *x);
-  Double_t PropagateToDCA(AliV0vertex *vtx, AliITStrackV2 *trk);
+  Double_t PropagateToDCA(AliESDv0 *vtx,AliExternalTrackParam *trk,Double_t b);
 
   void GetCuts(Double_t cuts[8]) const;
   void GetVertex(Double_t *vtx) const { vtx[0]=fX; vtx[1]=fY; vtx[2]=fZ; }
index 578beeaba68f7d99300f3ec584b4bd4af8c9e3b1..5a041992f3ca2587469fb2a1c653c81050b36a10 100644 (file)
@@ -27,6 +27,8 @@
 #include <TMath.h>
 
 #include "AliLog.h"
+#include "AliExternalTrackParam.h"
+#include "AliESDv0.h"
 #include "AliESDcascade.h"
 
 ClassImp(AliESDcascade)
@@ -45,6 +47,63 @@ AliESDcascade::AliESDcascade() :
   fPosCov[0]=fPosCov[1]=fPosCov[2]=fPosCov[3]=fPosCov[4]=fPosCov[5]=0.;
 }
 
+AliESDcascade::AliESDcascade(const AliESDv0 &v,
+                            const AliExternalTrackParam &t, Int_t i) : 
+  TObject(),
+  fPdgCode(kXiMinus),
+  fEffMass(TDatabasePDG::Instance()->GetParticle(kXiMinus)->Mass()),
+  fChi2(1.e+33),
+  fBachIdx(i)
+{
+  //--------------------------------------------------------------------
+  // Main constructor  (Xi-)
+  //--------------------------------------------------------------------
+
+  fV0idx[0]=v.GetNindex(); fV0idx[1]=v.GetPindex();
+
+  Double_t r[3]; t.GetXYZ(r);
+  Double_t x1=r[0], y1=r[1], z1=r[2]; // position of the bachelor
+  Double_t p[3]; t.GetPxPyPz(p);
+  Double_t px1=p[0], py1=p[1], pz1=p[2];// momentum of the bachelor track
+
+  Double_t x2,y2,z2;          // position of the V0 
+  v.GetXYZ(x2,y2,z2);    
+  Double_t px2,py2,pz2;       // momentum of V0
+  v.GetPxPyPz(px2,py2,pz2);
+
+  Double_t a2=((x1-x2)*px2+(y1-y2)*py2+(z1-z2)*pz2)/(px2*px2+py2*py2+pz2*pz2);
+
+  Double_t xm=x2+a2*px2;
+  Double_t ym=y2+a2*py2;
+  Double_t zm=z2+a2*pz2;
+
+  // position of the cascade decay
+  
+  fPos[0]=0.5*(x1+xm); fPos[1]=0.5*(y1+ym); fPos[2]=0.5*(z1+zm);
+    
+
+  // invariant mass of the cascade (default is Ximinus)
+  
+  Double_t e1=TMath::Sqrt(0.13957*0.13957 + px1*px1 + py1*py1 + pz1*pz1);
+  Double_t e2=TMath::Sqrt(1.11568*1.11568 + px2*px2 + py2*py2 + pz2*pz2);
+  
+  fEffMass=TMath::Sqrt((e1+e2)*(e1+e2)-
+    (px1+px2)*(px1+px2)-(py1+py2)*(py1+py2)-(pz1+pz2)*(pz1+pz2));
+
+
+  // momenta of the bachelor and the V0
+  
+  fBachMom[0]=px1; fBachMom[1]=py1; fBachMom[2]=pz1; 
+  v.GetNPxPyPz(px2,py2,pz2);
+  fV0mom[0][0]=px2; fV0mom[0][1]=py2; fV0mom[0][2]=pz2;
+  v.GetPPxPyPz(px2,py2,pz2);
+  fV0mom[1][0]=px2; fV0mom[1][1]=py2; fV0mom[1][2]=pz2;
+
+
+  fChi2=7.;   
+
+}
+
 Double_t AliESDcascade::ChangeMassHypothesis(Double_t &v0q, Int_t code) {
   //--------------------------------------------------------------------
   // This function changes the mass hypothesis for this cascade
index 6b2676b981692b83010c846e11e0b4804b121b6f..322e679add484646ff57fb43532a8f61747e9d69 100644 (file)
@@ -15,7 +15,7 @@
 #include <TObject.h>
 #include <TPDGCode.h>
 
-class AliESDtrack;
+class AliExternalTrackParam;
 class AliESDv0;
 
 #define kXiMinus       3312
@@ -26,6 +26,8 @@ class AliESDv0;
 class AliESDcascade : public TObject {
 public:
   AliESDcascade();
+  AliESDcascade(const AliESDv0 &v0, 
+                const AliExternalTrackParam &t, Int_t i);
 
   Double_t ChangeMassHypothesis(Double_t &v0q, Int_t code=kXiMinus); 
 
index e01c51e314f24f7caa16312962bee713eafb4beb..64a72255e24099e6d3d8eb71c4c5fccf250bbae4 100644 (file)
@@ -31,6 +31,7 @@
 
 #include "AliLog.h"
 #include "AliESDv0.h"
+#include "AliExternalTrackParam.h"
 
 ClassImp(AliESDv0)
 
@@ -60,6 +61,71 @@ AliESDv0::AliESDv0() :
   }
 }
 
+AliESDv0::AliESDv0(const AliExternalTrackParam &t1, Int_t i1,
+                   const AliExternalTrackParam &t2, Int_t i2) :
+  TObject(),
+  fPdgCode(kK0Short),
+  fEffMass(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass()),
+  fDcaDaughters(0),
+  fChi2(1.e+33),
+  fNidx(i1),
+  fPidx(i2)
+{
+  //--------------------------------------------------------------------
+  // Main constructor  (K0s)
+  //--------------------------------------------------------------------
+
+  for (Int_t i=0; i<6; i++) {
+    fPosCov[i]= 0.;
+    fNmomCov[i] = 0.;
+    fPmomCov[i] = 0.;
+  }
+
+  //Trivial estimation of the vertex parameters
+  Double_t x=t1.GetX(), alpha=t1.GetAlpha();
+  const Double_t *par=t1.GetParameter();
+  Double_t pt=1./TMath::Abs(par[4]), 
+           phi=TMath::ASin(par[2]) + alpha, 
+           cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
+
+  Double_t px1=pt*TMath::Cos(phi), py1=pt*TMath::Sin(phi), pz1=pt*par[3];
+  Double_t x1=x*cs - par[0]*sn;
+  Double_t y1=x*sn + par[0]*cs;
+  Double_t z1=par[1];
+  Double_t sx1=sn*sn*t1.GetSigmaY2(), sy1=cs*cs*t1.GetSigmaY2(); 
+
+
+
+  x=t2.GetX(); alpha=t2.GetAlpha(); par=t2.GetParameter();
+  pt=1./TMath::Abs(par[4]);
+  phi=TMath::ASin(par[2]) + alpha;  
+  cs=TMath::Cos(alpha); sn=TMath::Sin(alpha);
+
+  Double_t px2=pt*TMath::Cos(phi), py2=pt*TMath::Sin(phi), pz2=pt*par[3];
+  Double_t x2=x*cs - par[0]*sn;
+  Double_t y2=x*sn + par[0]*cs;
+  Double_t z2=par[1];
+  Double_t sx2=sn*sn*t2.GetSigmaY2(), sy2=cs*cs*t2.GetSigmaY2(); 
+    
+  Double_t sz1=t1.GetSigmaZ2(), sz2=t2.GetSigmaZ2();
+  Double_t wx1=sx2/(sx1+sx2), wx2=1.- wx1;
+  Double_t wy1=sy2/(sy1+sy2), wy2=1.- wy1;
+  Double_t wz1=sz2/(sz1+sz2), wz2=1.- wz1;
+  fPos[0]=wx1*x1 + wx2*x2; fPos[1]=wy1*y1 + wy2*y2; fPos[2]=wz1*z1 + wz2*z2;
+
+  //fPos[0]=0.5*(x1+x2); fPos[1]=0.5*(y1+y2); fPos[2]=0.5*(z1+z2);
+  fNmom[0]=px1; fNmom[1]=py1; fNmom[2]=pz1; 
+  fPmom[0]=px2; fPmom[1]=py2; fPmom[2]=pz2;
+
+  Double_t e1=TMath::Sqrt(0.13957*0.13957 + px1*px1 + py1*py1 + pz1*pz1);
+  Double_t e2=TMath::Sqrt(0.13957*0.13957 + px2*px2 + py2*py2 + pz2*pz2);
+  fEffMass=TMath::Sqrt((e1+e2)*(e1+e2)-
+    (px1+px2)*(px1+px2)-(py1+py2)*(py1+py2)-(pz1+pz2)*(pz1+pz2));
+
+  fChi2=7.;   
+
+}
+
 Double_t AliESDv0::ChangeMassHypothesis(Int_t code) {
   //--------------------------------------------------------------------
   // This function changes the mass hypothesis for this V0
index d9381c9295e72c4a12561bf2eca2e2bef5d4d55c..0599a455597063575bbdc4af9decfd1408600372 100644 (file)
 #include <TObject.h>
 #include <TPDGCode.h>
 
-class AliESDtrack;
+class AliExternalTrackParam;
 
 class AliESDv0 : public TObject {
 public:
   AliESDv0();
+  AliESDv0(const AliExternalTrackParam &t1, Int_t i1,
+           const AliExternalTrackParam &t2, Int_t i2);
 
   Double_t ChangeMassHypothesis(Int_t code=kK0Short); 
 
index abbb7c212f2aacf34445f937924f40d2fe108df8..8b94b8211cb074ffc28ea023d6940d7dd7fd2aab 100644 (file)
@@ -94,7 +94,7 @@ Double_t AliExternalTrackParam::GetP() const {
 }
 
 //_______________________________________________________________________
-Double_t AliExternalTrackParam::GetD(Double_t b,Double_t x,Double_t y) const {
+Double_t AliExternalTrackParam::GetD(Double_t x,Double_t y,Double_t b) const {
   //------------------------------------------------------------------
   // This function calculates the transverse impact parameter
   // with respect to a point with global coordinates (x,y)
@@ -131,6 +131,50 @@ Double_t AliExternalTrackParam::GetLinearD(Double_t xv,Double_t yv) const {
   return d;
 }
 
+Bool_t AliExternalTrackParam::
+CorrectForMaterial(Double_t d,  Double_t x0, Double_t mass) {
+  //------------------------------------------------------------------
+  // This function corrects the track parameters for the crossed material
+  // "d"    - the thickness (fraction of the radiation length)
+  // "x0"   - the radiation length (g/cm^2) 
+  // "mass" - the mass of this particle (GeV/c^2)
+  //------------------------------------------------------------------
+  Double_t &fP2=fP[2];
+  Double_t &fP3=fP[3];
+  Double_t &fP4=fP[4];
+
+  Double_t &fC22=fC[5];
+  Double_t &fC33=fC[9];
+  Double_t &fC43=fC[13];
+  Double_t &fC44=fC[14];
+
+  Double_t p2=(1.+ fP3*fP3)/(fP4*fP4);
+  Double_t beta2=p2/(p2 + mass*mass);
+  d*=TMath::Sqrt((1.+ fP3*fP3)/(1.- fP2*fP2));
+
+  //Multiple scattering******************
+  if (d!=0) {
+     Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(d);
+     //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*TMath::Abs(d)*9.36*2.33;
+     fC22 += theta2*(1.- fP2*fP2)*(1. + fP3*fP3);
+     fC33 += theta2*(1. + fP3*fP3)*(1. + fP3*fP3);
+     fC43 += theta2*fP3*fP4*(1. + fP3*fP3);
+     fC44 += theta2*fP3*fP4*fP3*fP4;
+  }
+
+  //Energy losses************************
+  if (x0!=0.) {
+     d*=x0;
+     Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d;
+     if (beta2/(1-beta2)>3.5*3.5)
+       dE=0.153e-3/beta2*(log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2)*d;
+
+     fP4*=(1.- TMath::Sqrt(p2 + mass*mass)/p2*dE);
+  }
+
+  return kTRUE;
+}
+
 Bool_t AliExternalTrackParam::Rotate(Double_t alpha) {
   //------------------------------------------------------------------
   // Transform this track to the local coord. system rotated
@@ -321,6 +365,168 @@ Bool_t AliExternalTrackParam::Update(Double_t p[2], Double_t cov[3]) {
   return kTRUE;
 }
 
+void 
+AliExternalTrackParam::GetHelixParameters(Double_t hlx[6], Double_t b) const {
+  //--------------------------------------------------------------------
+  // External track parameters -> helix parameters 
+  // "b" - magnetic field (kG)
+  //--------------------------------------------------------------------
+  Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
+  
+  hlx[0]=fP[0]; hlx[1]=fP[1]; hlx[2]=fP[2]; hlx[3]=fP[3]; hlx[4]=fP[4];
+
+  hlx[5]=fX*cs - hlx[0]*sn;               // x0
+  hlx[0]=fX*sn + hlx[0]*cs;               // y0
+//hlx[1]=                                 // z0
+  hlx[2]=TMath::ASin(hlx[2]) + fAlpha;    // phi0
+//hlx[3]=                                 // tgl
+  hlx[4]=hlx[4]*kB2C*b;                   // C
+}
+
+
+static void Evaluate(const Double_t *h, Double_t t,
+                     Double_t r[3],  //radius vector
+                     Double_t g[3],  //first defivatives
+                     Double_t gg[3]) //second derivatives
+{
+  //--------------------------------------------------------------------
+  // Calculate position of a point on a track and some derivatives
+  //--------------------------------------------------------------------
+  Double_t phase=h[4]*t+h[2];
+  Double_t sn=TMath::Sin(phase), cs=TMath::Cos(phase);
+
+  r[0] = h[5] + (sn - h[6])/h[4];
+  r[1] = h[0] - (cs - h[7])/h[4];  
+  r[2] = h[1] + h[3]*t;
+
+  g[0] = cs; g[1]=sn; g[2]=h[3];
+  
+  gg[0]=-h[4]*sn; gg[1]=h[4]*cs; gg[2]=0.;
+}
+
+Double_t AliExternalTrackParam::GetDCA(const AliExternalTrackParam *p, 
+Double_t b, Double_t &xthis, Double_t &xp) const {
+  //------------------------------------------------------------
+  // Returns the (weighed !) distance of closest approach between 
+  // this track and the track "p".
+  // Other returned values:
+  //   xthis, xt - coordinates of tracks' reference planes at the DCA 
+  //-----------------------------------------------------------
+  Double_t dy2=GetSigmaY2() + p->GetSigmaY2();
+  Double_t dz2=GetSigmaZ2() + p->GetSigmaZ2();
+  Double_t dx2=dy2; 
+
+  //dx2=dy2=dz2=1.;
+
+  Double_t p1[8]; GetHelixParameters(p1,b);
+  p1[6]=TMath::Sin(p1[2]); p1[7]=TMath::Cos(p1[2]);
+  Double_t p2[8]; p->GetHelixParameters(p2,b);
+  p2[6]=TMath::Sin(p2[2]); p2[7]=TMath::Cos(p2[2]);
+
+
+  Double_t r1[3],g1[3],gg1[3]; Double_t t1=0.;
+  Evaluate(p1,t1,r1,g1,gg1);
+  Double_t r2[3],g2[3],gg2[3]; Double_t t2=0.;
+  Evaluate(p2,t2,r2,g2,gg2);
+
+  Double_t dx=r2[0]-r1[0], dy=r2[1]-r1[1], dz=r2[2]-r1[2];
+  Double_t dm=dx*dx/dx2 + dy*dy/dy2 + dz*dz/dz2;
+
+  Int_t max=27;
+  while (max--) {
+     Double_t gt1=-(dx*g1[0]/dx2 + dy*g1[1]/dy2 + dz*g1[2]/dz2);
+     Double_t gt2=+(dx*g2[0]/dx2 + dy*g2[1]/dy2 + dz*g2[2]/dz2);
+     Double_t h11=(g1[0]*g1[0] - dx*gg1[0])/dx2 + 
+                  (g1[1]*g1[1] - dy*gg1[1])/dy2 +
+                  (g1[2]*g1[2] - dz*gg1[2])/dz2;
+     Double_t h22=(g2[0]*g2[0] + dx*gg2[0])/dx2 + 
+                  (g2[1]*g2[1] + dy*gg2[1])/dy2 +
+                  (g2[2]*g2[2] + dz*gg2[2])/dz2;
+     Double_t h12=-(g1[0]*g2[0]/dx2 + g1[1]*g2[1]/dy2 + g1[2]*g2[2]/dz2);
+
+     Double_t det=h11*h22-h12*h12;
+
+     Double_t dt1,dt2;
+     if (TMath::Abs(det)<1.e-33) {
+        //(quasi)singular Hessian
+        dt1=-gt1; dt2=-gt2;
+     } else {
+        dt1=-(gt1*h22 - gt2*h12)/det; 
+        dt2=-(h11*gt2 - h12*gt1)/det;
+     }
+
+     if ((dt1*gt1+dt2*gt2)>0) {dt1=-dt1; dt2=-dt2;}
+
+     //check delta(phase1) ?
+     //check delta(phase2) ?
+
+     if (TMath::Abs(dt1)/(TMath::Abs(t1)+1.e-3) < 1.e-4)
+     if (TMath::Abs(dt2)/(TMath::Abs(t2)+1.e-3) < 1.e-4) {
+        if ((gt1*gt1+gt2*gt2) > 1.e-4/dy2/dy2) 
+         AliWarning(" stopped at not a stationary point !");
+        Double_t lmb=h11+h22; lmb=lmb-TMath::Sqrt(lmb*lmb-4*det);
+        if (lmb < 0.) 
+         AliWarning(" stopped at not a minimum !");
+        break;
+     }
+
+     Double_t dd=dm;
+     for (Int_t div=1 ; ; div*=2) {
+        Evaluate(p1,t1+dt1,r1,g1,gg1);
+        Evaluate(p2,t2+dt2,r2,g2,gg2);
+        dx=r2[0]-r1[0]; dy=r2[1]-r1[1]; dz=r2[2]-r1[2];
+        dd=dx*dx/dx2 + dy*dy/dy2 + dz*dz/dz2;
+       if (dd<dm) break;
+        dt1*=0.5; dt2*=0.5;
+        if (div>512) {
+           AliWarning(" overshoot !"); break;
+        }   
+     }
+     dm=dd;
+
+     t1+=dt1;
+     t2+=dt2;
+
+  }
+
+  if (max<=0) AliWarning(" too many iterations !");
+
+  Double_t cs=TMath::Cos(GetAlpha());
+  Double_t sn=TMath::Sin(GetAlpha());
+  xthis=r1[0]*cs + r1[1]*sn;
+
+  cs=TMath::Cos(p->GetAlpha());
+  sn=TMath::Sin(p->GetAlpha());
+  xp=r2[0]*cs + r2[1]*sn;
+
+  return TMath::Sqrt(dm*TMath::Sqrt(dy2*dz2));
+}
+Double_t AliExternalTrackParam::
+PropagateToDCA(AliExternalTrackParam *p, Double_t b) {
+  //--------------------------------------------------------------
+  // Propagates this track and the argument track to the position of the
+  // distance of closest approach.
+  // Returns the (weighed !) distance of closest approach.
+  //--------------------------------------------------------------
+  Double_t xthis,xp;
+  Double_t dca=GetDCA(p,b,xthis,xp);
+
+  if (!PropagateTo(xthis,b)) {
+    //AliWarning(" propagation failed !");
+    return 1e+33;
+  }
+
+  if (!p->PropagateTo(xp,b)) {
+    //AliWarning(" propagation failed !";
+    return 1e+33;
+  }
+
+  return dca;
+}
+
+
+
 Bool_t Local2GlobalMomentum(Double_t p[3],Double_t alpha) {
   //----------------------------------------------------------------
   // This function performs local->global transformation of the
index 0549dbe6d67bfecb3a1181e12ca65bef34fb5146..0c7bf036f659c12c40e17f0c3379cff0c622006c 100644 (file)
@@ -39,12 +39,15 @@ class AliExternalTrackParam: public TObject {
 
   const Double_t *GetParameter() const {return fP;}
   const Double_t *GetCovariance() const {return fC;}
-  virtual Double_t GetX() const {return fX;}
-  virtual Double_t GetAlpha() const {return fAlpha;}
+  Double_t GetSigmaY2() const {return fC[0];}
+  Double_t GetSigmaZ2() const {return fC[2];}
+  Double_t GetX() const {return fX;}
+  Double_t GetAlpha() const {return fAlpha;}
   Double_t GetSign() const {return (fP[4]>0) ? 1 : -1;}
   Double_t GetP() const;
-  Double_t GetD(Double_t b, Double_t xv, Double_t yv) const; 
+  Double_t GetD(Double_t xv, Double_t yv, Double_t b) const; 
   Double_t GetLinearD(Double_t xv, Double_t yv) const; 
+  Bool_t CorrectForMaterial(Double_t d, Double_t x0, Double_t mass);
   Double_t GetPredictedChi2(Double_t p[2],Double_t cov[3]) const;
   Bool_t Update(Double_t p[2],Double_t cov[3]);
   Bool_t Rotate(Double_t alpha);
@@ -55,6 +58,11 @@ class AliExternalTrackParam: public TObject {
     return kFALSE;
   }
 
+  void GetHelixParameters(Double_t h[6], Double_t b) const;
+  Double_t GetDCA(const AliExternalTrackParam *p, Double_t b,
+    Double_t &xthis,Double_t &xp) const;
+  Double_t PropagateToDCA(AliExternalTrackParam *p, Double_t b);
+
   Bool_t GetPxPyPz(Double_t *p) const;
   Bool_t GetXYZ(Double_t *p) const;
   Bool_t GetCovarianceXYZPxPyPz(Double_t cv[21]) const;
index 90e6730416d51f4beb3571b68c38b08bff1c4a39..82663da53e2299b4cda985c724b69d3813ca1a91 100644 (file)
@@ -188,151 +188,6 @@ void AliKalmanTrack::External2Helix(Double_t helix[6]) const {
   helix[4]=helix[4]/GetLocalConvConst();  // C
 }
 
-static void Evaluate(const Double_t *h, Double_t t,
-                     Double_t r[3],  //radius vector
-                     Double_t g[3],  //first defivatives
-                     Double_t gg[3]) //second derivatives
-{
-  //--------------------------------------------------------------------
-  // Calculate position of a point on a track and some derivatives
-  //--------------------------------------------------------------------
-  Double_t phase=h[4]*t+h[2];
-  Double_t sn=TMath::Sin(phase), cs=TMath::Cos(phase);
-
-  r[0] = h[5] + (sn - h[6])/h[4];
-  r[1] = h[0] - (cs - h[7])/h[4];  
-  r[2] = h[1] + h[3]*t;
-
-  g[0] = cs; g[1]=sn; g[2]=h[3];
-  
-  gg[0]=-h[4]*sn; gg[1]=h[4]*cs; gg[2]=0.;
-}
-
-Double_t AliKalmanTrack::
-GetDCA(const AliKalmanTrack *p, Double_t &xthis, Double_t &xp) const {
-  //------------------------------------------------------------
-  // Returns the (weighed !) distance of closest approach between 
-  // this track and the track passed as the argument.
-  // Other returned values:
-  //   xthis, xt - coordinates of tracks' reference planes at the DCA 
-  //-----------------------------------------------------------
-  Double_t dy2=GetSigmaY2() + p->GetSigmaY2();
-  Double_t dz2=GetSigmaZ2() + p->GetSigmaZ2();
-  Double_t dx2=dy2; 
-
-  //dx2=dy2=dz2=1.;
-
-  Double_t p1[8]; External2Helix(p1);
-  p1[6]=TMath::Sin(p1[2]); p1[7]=TMath::Cos(p1[2]);
-  Double_t p2[8]; p->External2Helix(p2);
-  p2[6]=TMath::Sin(p2[2]); p2[7]=TMath::Cos(p2[2]);
-
-
-  Double_t r1[3],g1[3],gg1[3]; Double_t t1=0.;
-  Evaluate(p1,t1,r1,g1,gg1);
-  Double_t r2[3],g2[3],gg2[3]; Double_t t2=0.;
-  Evaluate(p2,t2,r2,g2,gg2);
-
-  Double_t dx=r2[0]-r1[0], dy=r2[1]-r1[1], dz=r2[2]-r1[2];
-  Double_t dm=dx*dx/dx2 + dy*dy/dy2 + dz*dz/dz2;
-
-  Int_t max=27;
-  while (max--) {
-     Double_t gt1=-(dx*g1[0]/dx2 + dy*g1[1]/dy2 + dz*g1[2]/dz2);
-     Double_t gt2=+(dx*g2[0]/dx2 + dy*g2[1]/dy2 + dz*g2[2]/dz2);
-     Double_t h11=(g1[0]*g1[0] - dx*gg1[0])/dx2 + 
-                  (g1[1]*g1[1] - dy*gg1[1])/dy2 +
-                  (g1[2]*g1[2] - dz*gg1[2])/dz2;
-     Double_t h22=(g2[0]*g2[0] + dx*gg2[0])/dx2 + 
-                  (g2[1]*g2[1] + dy*gg2[1])/dy2 +
-                  (g2[2]*g2[2] + dz*gg2[2])/dz2;
-     Double_t h12=-(g1[0]*g2[0]/dx2 + g1[1]*g2[1]/dy2 + g1[2]*g2[2]/dz2);
-
-     Double_t det=h11*h22-h12*h12;
-
-     Double_t dt1,dt2;
-     if (TMath::Abs(det)<1.e-33) {
-        //(quasi)singular Hessian
-        dt1=-gt1; dt2=-gt2;
-     } else {
-        dt1=-(gt1*h22 - gt2*h12)/det; 
-        dt2=-(h11*gt2 - h12*gt1)/det;
-     }
-
-     if ((dt1*gt1+dt2*gt2)>0) {dt1=-dt1; dt2=-dt2;}
-
-     //check delta(phase1) ?
-     //check delta(phase2) ?
-
-     if (TMath::Abs(dt1)/(TMath::Abs(t1)+1.e-3) < 1.e-4)
-     if (TMath::Abs(dt2)/(TMath::Abs(t2)+1.e-3) < 1.e-4) {
-        if ((gt1*gt1+gt2*gt2) > 1.e-4/dy2/dy2) 
-         AliWarning(" stopped at not a stationary point !");
-        Double_t lmb=h11+h22; lmb=lmb-TMath::Sqrt(lmb*lmb-4*det);
-        if (lmb < 0.) 
-         AliWarning(" stopped at not a minimum !");
-        break;
-     }
-
-     Double_t dd=dm;
-     for (Int_t div=1 ; ; div*=2) {
-        Evaluate(p1,t1+dt1,r1,g1,gg1);
-        Evaluate(p2,t2+dt2,r2,g2,gg2);
-        dx=r2[0]-r1[0]; dy=r2[1]-r1[1]; dz=r2[2]-r1[2];
-        dd=dx*dx/dx2 + dy*dy/dy2 + dz*dz/dz2;
-       if (dd<dm) break;
-        dt1*=0.5; dt2*=0.5;
-        if (div>512) {
-           AliWarning(" overshoot !"); break;
-        }   
-     }
-     dm=dd;
-
-     t1+=dt1;
-     t2+=dt2;
-
-  }
-
-  if (max<=0) AliWarning(" too many iterations !");
-
-  Double_t cs=TMath::Cos(GetAlpha());
-  Double_t sn=TMath::Sin(GetAlpha());
-  xthis=r1[0]*cs + r1[1]*sn;
-
-  cs=TMath::Cos(p->GetAlpha());
-  sn=TMath::Sin(p->GetAlpha());
-  xp=r2[0]*cs + r2[1]*sn;
-
-  return TMath::Sqrt(dm*TMath::Sqrt(dy2*dz2));
-}
-
-Double_t 
-AliKalmanTrack::PropagateToDCA(AliKalmanTrack *p, Double_t d, Double_t x0) {
-  //--------------------------------------------------------------
-  // Propagates this track and the argument track to the position of the
-  // distance of closest approach. 
-  // Returns the (weighed !) distance of closest approach.
-  //--------------------------------------------------------------
-  Double_t xthis,xp;
-  Double_t dca=GetDCA(p,xthis,xp);
-
-  if (!PropagateTo(xthis,d,x0)) {
-    //AliWarning(" propagation failed !");
-    return 1e+33;
-  }  
-
-  if (!p->PropagateTo(xp,d,x0)) {
-    //AliWarning(" propagation failed !";
-    return 1e+33;
-  }  
-
-  return dca;
-}
-
-
-
-
-
 Double_t AliKalmanTrack::MeanMaterialBudget(Double_t *start, Double_t *end, Double_t *mparam)
 {
   //
index 79fa2c119c663b6e2ccce131f4bc38542534fa8b..4596e3097aa494179d6dbbbd4ce633f6c9bf1daf 100644 (file)
@@ -42,9 +42,6 @@ public:
     return 0.;
   }
 
-  virtual Double_t GetDCA(const AliKalmanTrack *p,Double_t &xthis,Double_t &xp) const; 
-  virtual 
-  Double_t PropagateToDCA(AliKalmanTrack *p, Double_t d=0., Double_t x0=0.); 
   virtual Double_t GetAlpha() const {
     AliWarning("Method must be overloaded !\n");
     return 0.;
@@ -65,10 +62,9 @@ public:
 
   virtual Double_t GetPredictedChi2(const AliCluster *) const = 0;
   virtual Int_t PropagateTo(Double_t/*xr*/,Double_t/*x0*/,Double_t/*rho*/) = 0;
-  //virtual Int_t PropagateToVertex(Double_t /*d*/=0., Double_t /*x0*/=0.) = 0; 
+
   virtual Int_t Update(const AliCluster*, Double_t /*chi2*/, UInt_t) = 0;
 
-  //static Double_t GetConvConst();
   static Double_t MeanMaterialBudget(Double_t *start, Double_t *end, Double_t *mparam);
  
   // Time integration (S.Radomski@gsi.de)
diff --git a/STEER/AliV0vertexer.cxx b/STEER/AliV0vertexer.cxx
new file mode 100644 (file)
index 0000000..f7d4740
--- /dev/null
@@ -0,0 +1,141 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+
+//-------------------------------------------------------------------------
+//               Implementation of the V0 vertexer class
+//                  reads tracks writes out V0 vertices
+//                      fills the ESD with the V0s       
+//     Origin: Iouri Belikov, IReS, Strasbourg, Jouri.Belikov@cern.ch
+//-------------------------------------------------------------------------
+#include <TObjArray.h>
+#include <TTree.h>
+
+#include "AliESD.h"
+#include "AliESDv0.h"
+#include "AliESDtrack.h"
+#include "AliV0vertexer.h"
+
+ClassImp(AliV0vertexer)
+
+Int_t AliV0vertexer::Tracks2V0vertices(AliESD *event) {
+  //--------------------------------------------------------------------
+  //This function reconstructs V0 vertices
+  //--------------------------------------------------------------------
+
+   Int_t nentr=event->GetNumberOfTracks();
+   Double_t b=event->GetMagneticField();
+
+   TArrayI neg(nentr/2);
+   TArrayI pos(nentr/2);
+
+   Int_t nneg=0, npos=0, nvtx=0;
+
+   Int_t i;
+   for (i=0; i<nentr; i++) {
+     AliESDtrack *esdTrack=event->GetTrack(i);
+     UInt_t status=esdTrack->GetStatus();
+     UInt_t flags=AliESDtrack::kITSin|AliESDtrack::kTPCin|
+                  AliESDtrack::kTPCpid|AliESDtrack::kESDpid;
+
+     if ((status&AliESDtrack::kITSrefit)==0)
+        if (flags!=status) continue;
+
+     Double_t d=esdTrack->GetD(fX,fY,b);
+     if (TMath::Abs(d)<fDPmin) continue;
+     if (TMath::Abs(d)>fRmax) continue;
+
+     if (esdTrack->GetSign() < 0.) neg[nneg++]=i;
+     else pos[npos++]=i;
+   }   
+
+
+   for (i=0; i<nneg; i++) {
+      Int_t nidx=neg[i];
+      AliESDtrack *ntrk=event->GetTrack(nidx);
+
+      for (Int_t k=0; k<npos; k++) {
+         Int_t pidx=pos[k];
+        AliESDtrack *ptrk=event->GetTrack(pidx);
+
+         if (TMath::Abs(ntrk->GetD(fX,fY,b))<fDNmin)
+          if (TMath::Abs(ptrk->GetD(fX,fY,b))<fDNmin) continue;
+
+         Double_t xn, xp, dca=ntrk->GetDCA(ptrk,b,xn,xp);
+         if (dca > fDCAmax) continue;
+         if ((xn+xp) > 2*fRmax) continue;
+         if ((xn+xp) < 2*fRmin) continue;
+   
+         AliExternalTrackParam nt(*ntrk), pt(*ptrk);
+         Bool_t corrected=kFALSE;
+         if ((nt.GetX() > 3.) && (xn < 3.)) {
+          //correct for the beam pipe material
+           corrected=kTRUE;
+         }
+         if ((pt.GetX() > 3.) && (xp < 3.)) {
+          //correct for the beam pipe material
+           corrected=kTRUE;
+         }
+         if (corrected) {
+          dca=nt.GetDCA(&pt,b,xn,xp);
+           if (dca > fDCAmax) continue;
+           if ((xn+xp) > 2*fRmax) continue;
+           if ((xn+xp) < 2*fRmin) continue;
+        }
+
+         nt.PropagateTo(xn,b); pt.PropagateTo(xp,b);
+
+         AliESDv0 vertex(nt,nidx,pt,pidx);
+         if (vertex.GetChi2() > fChi2max) continue;
+        
+         /*  Think of something better here ! 
+         nt.PropagateToVertex(); if (TMath::Abs(nt.GetZ())<0.04) continue;
+         pt.PropagateToVertex(); if (TMath::Abs(pt.GetZ())<0.04) continue;
+        */
+
+         Double_t x,y,z; vertex.GetXYZ(x,y,z); 
+         Double_t px,py,pz; vertex.GetPxPyPz(px,py,pz);
+         Double_t p2=px*px+py*py+pz*pz;
+         Double_t cost=((x-fX)*px + (y-fY)*py + (z-fZ)*pz)/
+               TMath::Sqrt(p2*((x-fX)*(x-fX) + (y-fY)*(y-fY) + (z-fZ)*(z-fZ)));
+
+        //if (cost < (5*fCPAmax-0.9-TMath::Sqrt(r2)*(fCPAmax-1))/4.1) continue;
+         if (cost < fCPAmax) continue;
+        vertex.SetDcaDaughters(dca);
+         //vertex.ChangeMassHypothesis(); //default is Lambda0 
+
+         event->AddV0(&vertex);
+
+         nvtx++;
+      }
+   }
+
+   Info("Tracks2V0vertices","Number of reconstructed V0 vertices: %d",nvtx);
+
+   return 0;
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
+
similarity index 94%
rename from ITS/AliV0vertexer.h
rename to STEER/AliV0vertexer.h
index 6ddf03aa14faa3963d6fe3d80e67ca45e91b4362..c33420ef33a3a879fef9dc67b4d0feac84e45864 100644 (file)
@@ -12,7 +12,6 @@
 #include "TObject.h"
 
 class TTree;
-class AliITStrackV2;
 class AliESD;
 
 //_____________________________________________________________________________
@@ -25,9 +24,6 @@ public:
 
   Int_t Tracks2V0vertices(AliESD *event);
 
-  Int_t Tracks2V0vertices(TTree *in, TTree *out);
-  Double_t PropagateToDCA(AliITStrackV2 *nt, AliITStrackV2 *pt) const;
-
   void GetCuts(Double_t cuts[7]) const;
   void GetVertex(Double_t *vtx) const { vtx[0]=fX; vtx[1]=fY; vtx[2]=fZ; }
 
index 5a8e9cdaa8a167bd085b0eb3c613f02c180ce478..9feabbf47fb4d193e9ca9697020f7669a95f0043 100644 (file)
@@ -71,6 +71,8 @@
 #pragma link C++ class  AliReconstruction+;
 #pragma link C++ class  AliVertexGenFile+;
 #pragma link C++ class  AliVertexer+;
+#pragma link C++ class  AliV0vertexer+;
+#pragma link C++ class  AliCascadeVertexer+;
 
 #pragma link C++ class AliCDBPath;
 #pragma link C++ class AliCDBRunRange;
index e0157555d559229131b84a56bb839c1a72a3d3ed..d03fabab4724b2f2ac0e816928ebf6b89ab73545 100644 (file)
@@ -16,7 +16,7 @@ AliMergeCombi.cxx AliMagFMaps.cxx AliFieldMap.cxx \
 AliGausCorr.cxx AliTrackReference.cxx \
 AliTrackMap.cxx AliTrackMapper.cxx AliCollisionGeometry.cxx \
 AliMemoryWatcher.cxx \
-AliVertexer.cxx \
+AliVertexer.cxx  AliV0vertexer.cxx AliCascadeVertexer.cxx\
 AliMC.cxx AliSimulation.cxx AliReconstruction.cxx AliVertexGenFile.cxx \
 AliCDBEntry.cxx AliCDBId.cxx AliCDBMetaData.cxx \
 AliCDBPath.cxx AliCDBRunRange.cxx AliCDBManager.cxx\