#include <TMatrixDSym.h>
#include "AliExternalTrackParam.h"
#include "AliESDVertex.h"
+#include "TPolyMarker3D.h"
#include "AliLog.h"
ClassImp(AliExternalTrackParam)
for (Int_t i = 0; i < 15; i++) fC[i] = 0;
}
+//_____________________________________________________________________________
+void AliExternalTrackParam::AddCovariance(const Double_t c[15]) {
+ //
+ // Add "something" to the track covarince matrix.
+ // May be needed to account for unknown mis-calibration/mis-alignment
+ //
+ fC[0] +=c[0];
+ fC[1] +=c[1]; fC[2] +=c[2];
+ fC[3] +=c[3]; fC[4] +=c[4]; fC[5] +=c[5];
+ fC[6] +=c[6]; fC[7] +=c[7]; fC[8] +=c[8]; fC[9] +=c[9];
+ fC[10]+=c[10]; fC[11]+=c[11]; fC[12]+=c[12]; fC[13]+=c[13]; fC[14]+=c[14];
+}
+
+
Double_t AliExternalTrackParam::GetP() const {
//---------------------------------------------------------------------
// This function returns the track momentum
//Multiple scattering******************
if (xOverX0 != 0) {
Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(xOverX0);
+ if(theta2>TMath::Pi()*TMath::Pi()) return kFALSE;
//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);
Double_t e=TMath::Sqrt(p2 + mass*mass);
if ( TMath::Abs(dE) > 0.3*e ) return kFALSE; //30% energy loss is too much!
fP4*=(1.- e/p2*dE);
+ if (TMath::Abs(fP4)>100.) return kFALSE; // Do not track below 10 MeV/c
+
// Approximate energy loss fluctuation (M.Ivanov)
const Double_t knst=0.07; // To be tuned.
//Multiple scattering******************
if (d!=0) {
Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(d);
+ if(theta2>TMath::Pi()*TMath::Pi()) return kFALSE;
//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);
// the density effect taken into account at beta*gamma > 3.5
// (the approximation is reasonable only for solid materials)
//------------------------------------------------------------------
+ if (beta2 >= 1) return kVeryBig;
+
if (beta2/(1-beta2)>3.5*3.5)
return 0.153e-3/beta2*(log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2);
return kTRUE;
}
+Double_t *AliExternalTrackParam::GetResiduals(
+Double_t *p,Double_t *cov,Bool_t updated) const {
+ //------------------------------------------------------------------
+ // Returns the track residuals with the space point "p" having
+ // the covariance matrix "cov".
+ // If "updated" is kTRUE, the track parameters expected to be updated,
+ // otherwise they must be predicted.
+ //------------------------------------------------------------------
+ static Double_t res[2];
+
+ Double_t r00=cov[0], r01=cov[1], r11=cov[2];
+ if (updated) {
+ r00-=fC[0]; r01-=fC[1]; r11-=fC[2];
+ } else {
+ r00+=fC[0]; r01+=fC[1]; r11+=fC[2];
+ }
+ Double_t det=r00*r11 - r01*r01;
+
+ if (TMath::Abs(det) < kAlmost0) return 0;
+
+ Double_t tmp=r00; r00=r11/det; r11=tmp/det;
+
+ if (r00 < 0.) return 0;
+ if (r11 < 0.) return 0;
+
+ Double_t dy = fP[0] - p[0];
+ Double_t dz = fP[1] - p[1];
+
+ res[0]=dy*TMath::Sqrt(r00);
+ res[1]=dz*TMath::Sqrt(r11);
+
+ return res;
+}
+
Bool_t AliExternalTrackParam::Update(Double_t p[2], Double_t cov[3]) {
//------------------------------------------------------------------
// Update the track parameters with the space point "p" having
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 !");
+ AliDebug(1," 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 !");
+ AliDebug(1," stopped at not a minimum !");
break;
}
if (dd<dm) break;
dt1*=0.5; dt2*=0.5;
if (div>512) {
- AliWarning(" overshoot !"); break;
+ AliDebug(1," overshoot !"); break;
}
}
dm=dd;
}
- if (max<=0) AliWarning(" too many iterations !");
+ if (max<=0) AliDebug(1," too many iterations !");
Double_t cs=TMath::Cos(GetAlpha());
Double_t sn=TMath::Sin(GetAlpha());
}
-
-
-Bool_t AliExternalTrackParam::PropagateToDCA(const AliESDVertex *vtx, Double_t b, Double_t maxd){
+Bool_t AliExternalTrackParam::PropagateToDCA(const AliESDVertex *vtx,
+Double_t b, Double_t maxd, Double_t dz[2], Double_t covar[3]) {
//
- // Try to relate this track to the vertex "vtx",
+ // Propagate this track to the DCA to vertex "vtx",
// if the (rough) transverse impact parameter is not bigger then "maxd".
// Magnetic field is "b" (kG).
//
Double_t sn=TMath::Sin(alpha), cs=TMath::Cos(alpha);
Double_t x=GetX(), y=GetParameter()[0], snp=GetParameter()[2];
Double_t xv= vtx->GetXv()*cs + vtx->GetYv()*sn;
- Double_t yv=-vtx->GetXv()*sn + vtx->GetYv()*cs;
+ Double_t yv=-vtx->GetXv()*sn + vtx->GetYv()*cs, zv=vtx->GetZv();
x-=xv; y-=yv;
//Estimate the impact parameter neglecting the track curvature
if (d > maxd) return kFALSE;
//Propagate to the DCA
- Double_t crv=0.299792458e-3*b*GetParameter()[4];
+ Double_t crv=kB2C*b*GetParameter()[4];
+ if (TMath::Abs(b) < kAlmost0Field) crv=0.;
+
Double_t tgfv=-(crv*x - snp)/(crv*y + TMath::Sqrt(1.-snp*snp));
sn=tgfv/TMath::Sqrt(1.+ tgfv*tgfv); cs=TMath::Sqrt(1.- sn*sn);
+ if (TMath::Abs(tgfv)>0.) cs = sn/tgfv;
+ else cs=1.;
x = xv*cs + yv*sn;
yv=-xv*sn + yv*cs; xv=x;
if (!Propagate(alpha+TMath::ASin(sn),xv,b)) return kFALSE;
- return kTRUE;
-}
-
-
+ if (dz==0) return kTRUE;
+ dz[0] = GetParameter()[0] - yv;
+ dz[1] = GetParameter()[1] - zv;
+
+ if (covar==0) return kTRUE;
+ Double_t cov[6]; vtx->GetCovMatrix(cov);
-Bool_t Local2GlobalMomentum(Double_t p[3],Double_t alpha) {
- //----------------------------------------------------------------
- // This function performs local->global transformation of the
- // track momentum.
- // When called, the arguments are:
- // p[0] = 1/pt of the track;
- // p[1] = sine of local azim. angle of the track momentum;
- // p[2] = tangent of the track momentum dip angle;
- // alpha - rotation angle.
- // The result is returned as:
- // p[0] = px
- // p[1] = py
- // p[2] = pz
- // Results for (nearly) straight tracks are meaningless !
- //----------------------------------------------------------------
- if (TMath::Abs(p[0])<=kAlmost0) return kFALSE;
- if (TMath::Abs(p[1])> kAlmost1) return kFALSE;
-
- Double_t pt=1./TMath::Abs(p[0]);
- Double_t cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
- Double_t r=TMath::Sqrt(1 - p[1]*p[1]);
- p[0]=pt*(r*cs - p[1]*sn); p[1]=pt*(p[1]*cs + r*sn); p[2]=pt*p[2];
+ //***** Improvements by A.Dainese
+ alpha=GetAlpha(); sn=TMath::Sin(alpha); cs=TMath::Cos(alpha);
+ Double_t s2ylocvtx = cov[0]*sn*sn + cov[2]*cs*cs - 2.*cov[1]*cs*sn;
+ covar[0] = GetCovariance()[0] + s2ylocvtx; // neglecting correlations
+ covar[1] = GetCovariance()[1]; // between (x,y) and z
+ covar[2] = GetCovariance()[2] + cov[5]; // in vertex's covariance matrix
+ //*****
return kTRUE;
}
-Bool_t Local2GlobalPosition(Double_t r[3],Double_t alpha) {
- //----------------------------------------------------------------
- // This function performs local->global transformation of the
- // track position.
- // When called, the arguments are:
- // r[0] = local x
- // r[1] = local y
- // r[2] = local z
- // alpha - rotation angle.
- // The result is returned as:
- // r[0] = global x
- // r[1] = global y
- // r[2] = global z
- //----------------------------------------------------------------
- Double_t cs=TMath::Cos(alpha), sn=TMath::Sin(alpha), x=r[0];
- r[0]=x*cs - r[1]*sn; r[1]=x*sn + r[1]*cs;
-
- return kTRUE;
-}
void AliExternalTrackParam::GetDirection(Double_t d[3]) const {
//----------------------------------------------------------------
d[2]=fP[3]/norm;
}
-Bool_t AliExternalTrackParam::GetPxPyPz(Double_t *p) const {
+Bool_t AliExternalTrackParam::GetPxPyPz(Double_t p[3]) const {
//---------------------------------------------------------------------
// This function returns the global track momentum components
// Results for (nearly) straight tracks are meaningless !
return p[2];
}
+Double_t AliExternalTrackParam::Xv() const {
+ //---------------------------------------------------------------------
+ // Returns x-component of first track point
+ //---------------------------------------------------------------------
+
+ Double_t r[3]={0.,0.,0.};
+ GetXYZ(r);
+
+ return r[0];
+}
+
+Double_t AliExternalTrackParam::Yv() const {
+ //---------------------------------------------------------------------
+ // Returns y-component of first track point
+ //---------------------------------------------------------------------
+
+ Double_t r[3]={0.,0.,0.};
+ GetXYZ(r);
+
+ return r[1];
+}
+
+Double_t AliExternalTrackParam::Zv() const {
+ //---------------------------------------------------------------------
+ // Returns z-component of first track point
+ //---------------------------------------------------------------------
+
+ Double_t r[3]={0.,0.,0.};
+ GetXYZ(r);
+
+ return r[2];
+}
+
Double_t AliExternalTrackParam::Theta() const {
// return theta angle of momentum
- return TMath::ATan2(Pt(), Pz());
+ return 0.5*TMath::Pi() - TMath::ATan(fP[3]);
}
Double_t AliExternalTrackParam::Phi() const {
return kTRUE;
}
+
+
+//
+// Draw functionality.
+// Origin: Marian Ivanov, Marian.Ivanov@cern.ch
+//
+
+void AliExternalTrackParam::DrawTrack(Float_t magf, Float_t minR, Float_t maxR, Float_t stepR){
+ //
+ // Draw track line
+ //
+ if (minR>maxR) return ;
+ if (stepR<=0) return ;
+ Int_t npoints = TMath::Nint((maxR-minR)/stepR)+1;
+ if (npoints<1) return;
+ TPolyMarker3D *polymarker = new TPolyMarker3D(npoints);
+ FillPolymarker(polymarker, magf,minR,maxR,stepR);
+ polymarker->Draw();
+}
+
+//
+void AliExternalTrackParam::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
+ //
+ // Fill points in the polymarker
+ //
+ Int_t counter=0;
+ for (Double_t r=minR; r<maxR; r+=stepR){
+ Double_t point[3];
+ GetXYZAt(r,magF,point);
+ pol->SetPoint(counter,point[0],point[1], point[2]);
+ printf("xyz\t%f\t%f\t%f\n",point[0], point[1],point[2]);
+ counter++;
+ }
+}