[u/mrichter/AliRoot.git] / STEER / AliExternalTrackParam.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

///////////////////////////////////////////////////////////////////////////////
//                                                                           //
// track parameters in "external" format                                     //
//                                                                           //
// The track parameters are:                                                //
// - local y coordinate                                                      //
// - local z coordinate                                                      //
// - sin of azimutal angle                                                   //
// - tan of dip angle                                                        //
// - charge/pt                                                               //
// The parametrisation is given at the local x coordinate fX and the         //
// azimuthal angle fAlpha.                                                   //
//                                                                           //
// The external parametrisation can be used to exchange track parameters     //
// between different detectors.                                              //
//                                                                           //
///////////////////////////////////////////////////////////////////////////////

#include <TMath.h>
#include <TVector3.h>

#include "AliExternalTrackParam.h"
#include "AliKalmanTrack.h"
#include "AliTrackReference.h"
#include "AliLog.h"

ClassImp(AliExternalTrackParam)

const AliMagF *AliExternalTrackParam::fgkFieldMap=0;
Double_t AliExternalTrackParam::fgConvConst=0.;


//_____________________________________________________________________________
AliExternalTrackParam::AliExternalTrackParam() :
  fMass(-1),
  fX(0),
  fAlpha(0),
  fLocalConvConst(0)
{
  //
  // default constructor
  //
  for (Int_t i = 0; i < 5; i++) fParam[i] = 0;
  for (Int_t i = 0; i < 15; i++) fCovar[i] = 0;
}

//_____________________________________________________________________________
AliExternalTrackParam::AliExternalTrackParam(Double_t x, Double_t alpha, 
					     const Double_t param[5], 
					     const Double_t covar[15]) :
  fMass(-1),
  fX(x),
  fAlpha(alpha),
  fLocalConvConst(0)
{
  //
  // create external track parameters from given arguments
  //
  for (Int_t i = 0; i < 5; i++) fParam[i] = param[i];
  for (Int_t i = 0; i < 15; i++) fCovar[i] = covar[i];
}

//_____________________________________________________________________________
AliExternalTrackParam::AliExternalTrackParam(const AliKalmanTrack& track) :
  fMass(track.GetMass()),
  fX(0),
  fAlpha(track.GetAlpha()),
  fLocalConvConst(0)
{
  //
  //
  track.GetExternalParameters(fX,fParam);
  track.GetExternalCovariance(fCovar);
}


AliExternalTrackParam * AliExternalTrackParam::MakeTrack(const AliTrackReference *ref, Double_t mass)
{
  //
  // Make dummy track from the track reference 
  // negative mass means opposite charge 
  //
  Double_t xx[5];
  Double_t cc[15];
  for (Int_t i=0;i<15;i++) cc[i]=0;
  Double_t x = ref->X(), y = ref->Y(), z = ref->Z();
  Double_t alpha = TMath::ATan2(y,x);
  Double_t xr = TMath::Sqrt(x*x+y*y);
  xx[0] = 0;
  xx[1] = z;
  xx[3] = ref->Pz()/ref->Pt();
  Float_t b[3];
  Float_t xyz[3]={x,y,z};
  Float_t convConst = 0;
  (AliKalmanTrack::GetFieldMap())->Field(xyz,b);
  convConst=1000/0.299792458/(1e-13 - b[2]);
  xx[4] = 1./(convConst*ref->Pt()); // curvature rpresentation
  if (mass<0) xx[4]*=-1.;  // negative mass - negative direction
  Double_t alphap = TMath::ATan2(ref->Py(),ref->Px())-alpha;
  if (alphap> TMath::Pi()) alphap-=TMath::Pi();
  if (alphap<-TMath::Pi()) alphap+=TMath::Pi();
  xx[2] = TMath::Sin(alphap);
  xx[4]*=convConst;   // 1/pt representation 
  //  AliExternalTrackParam * track = new  AliExternalTrackParam(xx,cc,xr,alpha);
  AliExternalTrackParam * track = new  AliExternalTrackParam(xr,alpha,xx,cc);
  track->fMass = TMath::Abs(mass);
  //track->StartTimeIntegral();  
  track->SaveLocalConvConst(); 
  return track;
}


//_____________________________________________________________________________
const Double_t* AliExternalTrackParam::GetParameter() const
{
// get a pointer to the array of track parameters

  return fParam;
}

//_____________________________________________________________________________
const Double_t* AliExternalTrackParam::GetCovariance() const
{
// get a pointer to the array of the track parameter covariance matrix

  return fCovar;
}

//_____________________________________________________________________________
AliExternalTrackParam* AliExternalTrackParam::CreateExternalParam() const
{
// copy this instance

  return new AliExternalTrackParam(fX, fAlpha, fParam, fCovar);
}

//_____________________________________________________________________________
void AliExternalTrackParam::ResetCovariance(Double_t factor,
					    Bool_t clearOffDiagonal)
{
// reset the covariance matrix ("forget" track history)

  Int_t k = 0;
  for (Int_t i = 0; i < 5; i++) {
    for (Int_t j = 0; j < i; j++) {  // off diagonal elements
      if (clearOffDiagonal) {
	fCovar[k++] = 0;
      } else {
	fCovar[k++] *= factor;
      }
    }
    fCovar[k++] *= factor;     // diagonal elements
  }
}


//_____________________________________________________________________________
Bool_t AliExternalTrackParam::PropagateTo(Double_t xk, Double_t x0, Double_t rho)
{
  //
  // Propagate the track parameters to the given x coordinate assuming vacuum.
  // If length is not NULL, the change of track length is added to it.
  //
  
  Double_t lcc=GetLocalConvConst();  
  Double_t cur = fParam[4]/lcc;
  Double_t x1=fX, x2=xk, dx=x2-x1;
  Double_t f1=fParam[2], f2=f1 + cur*dx;
  if (TMath::Abs(f2) >= 0.98) {
    // MI change  - don't propagate highly inclined tracks
    //              covariance matrix distorted
    return kFALSE;
  }

  // old position [SR, GSI, 17.02.2003]
  Double_t oldX = fX, oldY = fParam[0], oldZ = fParam[1];
  Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);  
  fParam[0] += dx*(f1+f2)/(r1+r2);
  fParam[1] += dx*(f1+f2)/(f1*r2 + f2*r1)*fParam[3];
  fParam[2] += dx*cur;
  // transform error matrix to the curvature
  fCovar[10]/=lcc;
  fCovar[11]/=lcc;
  fCovar[12]/=lcc;
  fCovar[13]/=lcc;
  fCovar[14]/=lcc*lcc;

  //f = F - 1
  
  Double_t f02=    dx/(r1*r1*r1);
  Double_t f04=0.5*dx*dx/(r1*r1*r1);
  Double_t f12=    dx*fParam[3]*f1/(r1*r1*r1);
  Double_t f14=0.5*dx*dx*fParam[3]*f1/(r1*r1*r1);
  Double_t f13=    dx/r1;
  Double_t f24=    dx; 
  
  //b = C*ft
  Double_t b00=f02*fCovar[3] + f04*fCovar[10], b01=f12*fCovar[3] + f14*fCovar[10] + f13*fCovar[6];
  Double_t b02=f24*fCovar[10];
  Double_t b10=f02*fCovar[4] + f04*fCovar[11], b11=f12*fCovar[4] + f14*fCovar[11] + f13*fCovar[7];
  Double_t b12=f24*fCovar[11];
  Double_t b20=f02*fCovar[5] + f04*fCovar[12], b21=f12*fCovar[5] + f14*fCovar[12] + f13*fCovar[8];
  Double_t b22=f24*fCovar[12];
  Double_t b40=f02*fCovar[12] + f04*fCovar[14], b41=f12*fCovar[12] + f14*fCovar[14] + f13*fCovar[13];
  Double_t b42=f24*fCovar[14];
  Double_t b30=f02*fCovar[8] + f04*fCovar[13], b31=f12*fCovar[8] + f14*fCovar[13] + f13*fCovar[9];
  Double_t b32=f24*fCovar[13];
  
  //a = f*b = f*C*ft
  Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
  Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
  Double_t a22=f24*b42;

  //F*C*Ft = C + (b + bt + a)
  fCovar[0] += b00 + b00 + a00;
  fCovar[1] += b10 + b01 + a01; 
  fCovar[2] += b11 + b11 + a11;
  fCovar[3] += b20 + b02 + a02;
  fCovar[4] += b21 + b12 + a12;
  fCovar[5] += b22 + b22 + a22;
  fCovar[6] += b30;
  fCovar[7] += b31; 
  fCovar[8] += b32;
  fCovar[10] += b40;
  fCovar[11] += b41;
  fCovar[12] += b42;

  fX=x2;

  //Change of the magnetic field *************
  SaveLocalConvConst();
  // transform back error matrix from curvature to the 1/pt
  fCovar[10]*=lcc;
  fCovar[11]*=lcc;
  fCovar[12]*=lcc;
  fCovar[13]*=lcc;
  fCovar[14]*=lcc*lcc;

  Double_t dist = TMath::Sqrt((fX-oldX)*(fX-oldX)+(fParam[0]-oldY)*(fParam[0]-oldY)+
			      (fParam[1]-oldZ)*(fParam[1]-oldZ));
  if (!CorrectForMaterial(dist,x0,rho)) return 0;

  // Integrated Time [SR, GSI, 17.02.2003]
 //  if (IsStartedTimeIntegral() && fX>oldX) {
//     Double_t l2 = (fX-oldX)*(fX-oldX)+(fParam[0]-oldY)*(fParam[0]-oldY)+
//                   (fParam[1]-oldZ)*(fParam[1]-oldZ);
//     AddTimeStep(TMath::Sqrt(l2));
//   }
  //

  return kTRUE;
}

Bool_t     AliExternalTrackParam::PropagateToDCA(Double_t xd, Double_t yd,  Double_t x0, Double_t rho){
  //
  // Propagate the track parameters to the nearest point of given xv yv coordinate
  //
  Double_t a=fAlpha;
  Double_t cs=TMath::Cos(a),sn=TMath::Sin(a);

  Double_t xv= xd*cs + yd*sn;
  Double_t yv=-xd*sn + yd*cs;   // vertex position in local frame
  //  
  Double_t c=fParam[4]/GetLocalConvConst(), snp=fParam[2];
  //
  Double_t x=fX, y=fParam[1];
  Double_t tgfv=-(c*(x-xv)-snp)/(c*(y-yv) + TMath::Sqrt(1.-snp*snp));
  Double_t fv=TMath::ATan(tgfv);
  cs=TMath::Cos(fv); sn=TMath::Sin(fv);
  x = xv*cs + yv*sn;
  yv=-xv*sn + yv*cs; xv=x;
  RotateTo(fv+a);
  return PropagateTo(xv,x0,rho);  
}


//_____________________________________________________________________________
Bool_t AliExternalTrackParam::RotateTo(Double_t alp)
{
  // Rotate the reference axis for the parametrisation to the given angle.
  //
  Double_t  x=fX;
  Double_t p0=fParam[0];
  //
  if      (alp < -TMath::Pi()) alp += 2*TMath::Pi();
  else if (alp >= TMath::Pi()) alp -= 2*TMath::Pi();
  Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
  Double_t sf=fParam[2], cf=TMath::Sqrt(1.- fParam[2]*fParam[2]);
  // **** rotation **********************
  
  fAlpha = alp;
  fX =  x*ca + p0*sa;
  fParam[0]= -x*sa + p0*ca;
  fParam[2]=  sf*ca - cf*sa;
  Double_t rr=(ca+sf/cf*sa);  
  //
  fCovar[0] *= (ca*ca);
  fCovar[1] *= ca; 
  fCovar[3] *= ca*rr;
  fCovar[6] *= ca;
  fCovar[10] *= ca;
  fCovar[4] *= rr;
  fCovar[5] *= rr*rr;
  fCovar[7] *= rr;
  fCovar[11] *= rr;
  return kTRUE;
}

//_____________________________________________________________________________
Bool_t AliExternalTrackParam::CorrectForMaterial(Double_t d, Double_t x0, Double_t rho)
{
  //
  // Take into account material effects assuming:
  // x0  - mean rad length
  // rho - mean density

  //
  // multiple scattering
  //
  if (fMass<=0) {
    AliError("Non-positive mass");
    return kFALSE;
  }
  Double_t p2=(1.+ fParam[3]*fParam[3])/(fParam[4]*fParam[4]);
  Double_t beta2=p2/(p2 + fMass*fMass);
  Double_t theta2=14.1*14.1/(beta2*p2*1e6)*d/x0*rho;
  //
  fCovar[5] += theta2*(1.- fParam[2]*fParam[2])*(1. + fParam[3]*fParam[3]);
  fCovar[9] += theta2*(1. + fParam[3]*fParam[3])*(1. + fParam[3]*fParam[3]);
  fCovar[13] += theta2*fParam[3]*fParam[4]*(1. + fParam[3]*fParam[3]);
  fCovar[14] += theta2*fParam[3]*fParam[4]*fParam[3]*fParam[4];
  
  Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2+1e-10)) - beta2)*d*rho;  
  fParam[4] *=(1.- TMath::Sqrt(p2+fMass*fMass)/p2*dE);
  //
  Double_t sigmade = 0.02*TMath::Sqrt(TMath::Abs(dE));   // energy loss fluctuation 
  Double_t sigmac2 = sigmade*sigmade*fParam[4]*fParam[4]*(p2+fMass*fMass)/(p2*p2);
  fCovar[14] += sigmac2;
  //
  //
  

  return kTRUE;
}

//_____________________________________________________________________________
Bool_t AliExternalTrackParam::GetProlongationAt(Double_t xk, 
						Double_t& y, 
						Double_t& z) const
{
  //
  // Get the local y and z coordinates at the given x value
  //
  Double_t lcc=GetLocalConvConst();  
  Double_t cur = fParam[4]/lcc;
  Double_t x1=fX, x2=xk, dx=x2-x1; 
  Double_t f1=fParam[2], f2=f1 + cur*dx;
  //  
  if (TMath::Abs(f2) >= 0.98) {
    // MI change  - don't propagate highly inclined tracks
    //              covariance matrix distorted
    return kFALSE;
  }
  Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
  y = fParam[0] + dx*(f1+f2)/(r1+r2);
  z = fParam[1] + dx*(f1+f2)/(f1*r2 + f2*r1)*fParam[3];
  return kTRUE;
}

//_____________________________________________________________________________
Double_t AliExternalTrackParam::GetXAtVertex(Double_t /*x*/, 
					     Double_t /*y*/) const
{
// Get the x coordinate at the given vertex (x,y)
//
// NOT IMPLEMENTED for this class

  return 0;
}


// //_____________________________________________________________________________
// Double_t AliExternalTrackParam::GetPredictedChi2(const AliCluster* /*cluster*/)
// {
// // calculate the chi2 contribution of the given cluster
// //
// // NOT IMPLEMENTED for this class

//   return -1;
// }

// //_____________________________________________________________________________
// Bool_t AliExternalTrackParam::Update(const AliCluster* /*cluster*/)
// {
// // update the track parameters using the position and error 
// // of the given cluster
// //
// // NOT IMPLEMENTED for this class

//   return kFALSE;
// }


//_____________________________________________________________________________
Double_t AliExternalTrackParam::SigmaPhi() const
{
// get the error of the azimuthal angle

  return TMath::Sqrt(TMath::Abs(fCovar[5] / (1. - fParam[2]*fParam[2])));
}

//_____________________________________________________________________________
Double_t AliExternalTrackParam::SigmaTheta() const
{
// get the error of the polar angle

  return TMath::Sqrt(TMath::Abs(fCovar[9])) / (1. + fParam[3]*fParam[3]);
}

//_____________________________________________________________________________
Double_t AliExternalTrackParam::SigmaPt() const
{
// get the error of the transversal component of the momentum

  return TMath::Sqrt(fCovar[14]) / TMath::Abs(fParam[4]);
}

//_____________________________________________________________________________
TVector3 AliExternalTrackParam::Momentum() const
{
// get the momentum vector

  Double_t phi = TMath::ASin(fParam[2]) + fAlpha;
  Double_t pt = 1. / TMath::Abs(fParam[4]);
  return TVector3(pt * TMath::Cos(phi), 
		  pt * TMath::Sin(phi), 
		  pt * fParam[3]);
}

//_____________________________________________________________________________
TVector3 AliExternalTrackParam::Position() const
{
// get the current spatial position in global coordinates

  return TVector3(fX * TMath::Cos(fAlpha) - fParam[0] * TMath::Sin(fAlpha),
		  fX * TMath::Sin(fAlpha) + fParam[0] * TMath::Cos(fAlpha),
		  fParam[1]);
}


//_____________________________________________________________________________
void AliExternalTrackParam::Print(Option_t* /*option*/) const
{
// print the parameters and the covariance matrix

  printf("AliExternalTrackParam: x = %-12g  alpha = %-12g\n", fX, fAlpha);
  printf("  parameters: %12g %12g %12g %12g %12g\n",
	 fParam[0], fParam[1], fParam[2], fParam[3], fParam[4]);
  printf("  covariance: %12g\n", fCovar[0]);
  printf("              %12g %12g\n", fCovar[1], fCovar[2]);
  printf("              %12g %12g %12g\n", fCovar[3], fCovar[4], fCovar[5]);
  printf("              %12g %12g %12g %12g\n", 
	 fCovar[6], fCovar[7], fCovar[8], fCovar[9]);
  printf("              %12g %12g %12g %12g %12g\n", 
	 fCovar[10], fCovar[11], fCovar[12], fCovar[13], fCovar[14]);
}
Commit	Line	Data
51ad6848	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/* $Id$ */
	17
	18	///////////////////////////////////////////////////////////////////////////////
	19	// //
	20	// track parameters in "external" format //
	21	// //
	22	// The track parameters are: //
	23	// - local y coordinate //
	24	// - local z coordinate //
	25	// - sin of azimutal angle //
	26	// - tan of dip angle //
	27	// - charge/pt //
	28	// The parametrisation is given at the local x coordinate fX and the //
	29	// azimuthal angle fAlpha. //
	30	// //
	31	// The external parametrisation can be used to exchange track parameters //
	32	// between different detectors. //
	33	// //
	34	///////////////////////////////////////////////////////////////////////////////
	35
a0209ae4	36	#include <TMath.h>
a0209ae4	37	#include <TVector3.h>
51ad6848	38
	39	#include "AliExternalTrackParam.h"
	40	#include "AliKalmanTrack.h"
a5e407e9	41	#include "AliTrackReference.h"
0243e825	42	#include "AliLog.h"
51ad6848	43
	44	ClassImp(AliExternalTrackParam)
	45
a5e407e9	46	const AliMagF *AliExternalTrackParam::fgkFieldMap=0;
	47	Double_t AliExternalTrackParam::fgConvConst=0.;
	48
51ad6848	49
51ad6848	50	//_____________________________________________________________________________
90e48c0c	51	AliExternalTrackParam::AliExternalTrackParam() :
0243e825	52	fMass(-1),
90e48c0c	53	fX(0),
a5e407e9	54	fAlpha(0),
a5e407e9	55	fLocalConvConst(0)
51ad6848	56	{
90e48c0c	57	//
	58	// default constructor
	59	//
51ad6848	60	for (Int_t i = 0; i < 5; i++) fParam[i] = 0;
	61	for (Int_t i = 0; i < 15; i++) fCovar[i] = 0;
	62	}
	63
	64	//_____________________________________________________________________________
	65	AliExternalTrackParam::AliExternalTrackParam(Double_t x, Double_t alpha,
	66	const Double_t param[5],
90e48c0c	67	const Double_t covar[15]) :
0243e825	68	fMass(-1),
90e48c0c	69	fX(x),
a5e407e9	70	fAlpha(alpha),
a5e407e9	71	fLocalConvConst(0)
51ad6848	72	{
90e48c0c	73	//
	74	// create external track parameters from given arguments
	75	//
51ad6848	76	for (Int_t i = 0; i < 5; i++) fParam[i] = param[i];
	77	for (Int_t i = 0; i < 15; i++) fCovar[i] = covar[i];
	78	}
	79
90e48c0c	80	//_____________________________________________________________________________
90e48c0c	81	AliExternalTrackParam::AliExternalTrackParam(const AliKalmanTrack& track) :
0243e825	82	fMass(track.GetMass()),
90e48c0c	83	fX(0),
0243e825	84	fAlpha(track.GetAlpha()),
0243e825	85	fLocalConvConst(0)
51ad6848	86	{
	87	//
	88	//
	89	track.GetExternalParameters(fX,fParam);
	90	track.GetExternalCovariance(fCovar);
51ad6848	91	}
51ad6848	92
a5e407e9	93
	94	AliExternalTrackParam * AliExternalTrackParam::MakeTrack(const AliTrackReference *ref, Double_t mass)
	95	{
	96	//
	97	// Make dummy track from the track reference
	98	// negative mass means opposite charge
	99	//
	100	Double_t xx[5];
	101	Double_t cc[15];
	102	for (Int_t i=0;i<15;i++) cc[i]=0;
	103	Double_t x = ref->X(), y = ref->Y(), z = ref->Z();
	104	Double_t alpha = TMath::ATan2(y,x);
	105	Double_t xr = TMath::Sqrt(xx+yy);
	106	xx[0] = 0;
	107	xx[1] = z;
	108	xx[3] = ref->Pz()/ref->Pt();
	109	Float_t b[3];
	110	Float_t xyz[3]={x,y,z};
	111	Float_t convConst = 0;
	112	(AliKalmanTrack::GetFieldMap())->Field(xyz,b);
	113	convConst=1000/0.299792458/(1e-13 - b[2]);
	114	xx[4] = 1./(convConst*ref->Pt()); // curvature rpresentation
	115	if (mass<0) xx[4]*=-1.; // negative mass - negative direction
	116	Double_t alphap = TMath::ATan2(ref->Py(),ref->Px())-alpha;
	117	if (alphap> TMath::Pi()) alphap-=TMath::Pi();
	118	if (alphap<-TMath::Pi()) alphap+=TMath::Pi();
	119	xx[2] = TMath::Sin(alphap);
	120	xx[4]*=convConst; // 1/pt representation
	121	// AliExternalTrackParam * track = new AliExternalTrackParam(xx,cc,xr,alpha);
	122	AliExternalTrackParam * track = new AliExternalTrackParam(xr,alpha,xx,cc);
	123	track->fMass = TMath::Abs(mass);
	124	//track->StartTimeIntegral();
	125	track->SaveLocalConvConst();
	126	return track;
	127	}
	128
	129
	130
51ad6848	131	//_____________________________________________________________________________
	132	const Double_t* AliExternalTrackParam::GetParameter() const
	133	{
	134	// get a pointer to the array of track parameters
	135
	136	return fParam;
	137	}
	138
	139	//_____________________________________________________________________________
	140	const Double_t* AliExternalTrackParam::GetCovariance() const
	141	{
	142	// get a pointer to the array of the track parameter covariance matrix
	143
	144	return fCovar;
	145	}
	146
	147	//_____________________________________________________________________________
	148	AliExternalTrackParam* AliExternalTrackParam::CreateExternalParam() const
	149	{
	150	// copy this instance
	151
	152	return new AliExternalTrackParam(fX, fAlpha, fParam, fCovar);
	153	}
	154
	155	//_____________________________________________________________________________
	156	void AliExternalTrackParam::ResetCovariance(Double_t factor,
	157	Bool_t clearOffDiagonal)
	158	{
	159	// reset the covariance matrix ("forget" track history)
	160
	161	Int_t k = 0;
	162	for (Int_t i = 0; i < 5; i++) {
	163	for (Int_t j = 0; j < i; j++) { // off diagonal elements
	164	if (clearOffDiagonal) {
	165	fCovar[k++] = 0;
	166	} else {
	167	fCovar[k++] *= factor;
	168	}
	169	}
	170	fCovar[k++] *= factor; // diagonal elements
	171	}
	172	}
	173
	174
	175	//_____________________________________________________________________________
a5e407e9	176	Bool_t AliExternalTrackParam::PropagateTo(Double_t xk, Double_t x0, Double_t rho)
51ad6848	177	{
a5e407e9	178	//
	179	// Propagate the track parameters to the given x coordinate assuming vacuum.
	180	// If length is not NULL, the change of track length is added to it.
	181	//
	182
	183	Double_t lcc=GetLocalConvConst();
	184	Double_t cur = fParam[4]/lcc;
	185	Double_t x1=fX, x2=xk, dx=x2-x1;
	186	Double_t f1=fParam[2], f2=f1 + cur*dx;
	187	if (TMath::Abs(f2) >= 0.98) {
	188	// MI change - don't propagate highly inclined tracks
	189	// covariance matrix distorted
	190	return kFALSE;
	191	}
51ad6848	192
a5e407e9	193	// old position [SR, GSI, 17.02.2003]
	194	Double_t oldX = fX, oldY = fParam[0], oldZ = fParam[1];
	195	Double_t r1=sqrt(1.- f1f1), r2=sqrt(1.- f2f2);
	196	fParam[0] += dx*(f1+f2)/(r1+r2);
	197	fParam[1] += dx(f1+f2)/(f1r2 + f2r1)fParam[3];
	198	fParam[2] += dx*cur;
	199	// transform error matrix to the curvature
	200	fCovar[10]/=lcc;
	201	fCovar[11]/=lcc;
	202	fCovar[12]/=lcc;
	203	fCovar[13]/=lcc;
	204	fCovar[14]/=lcc*lcc;
	205
	206	//f = F - 1
	207
	208	Double_t f02= dx/(r1r1r1);
	209	Double_t f04=0.5dxdx/(r1r1r1);
	210	Double_t f12= dxfParam[3]f1/(r1r1r1);
	211	Double_t f14=0.5dxdxfParam[3]f1/(r1r1r1);
	212	Double_t f13= dx/r1;
	213	Double_t f24= dx;
	214
	215	//b = C*ft
	216	Double_t b00=f02fCovar[3] + f04fCovar[10], b01=f12fCovar[3] + f14fCovar[10] + f13*fCovar[6];
	217	Double_t b02=f24*fCovar[10];
	218	Double_t b10=f02fCovar[4] + f04fCovar[11], b11=f12fCovar[4] + f14fCovar[11] + f13*fCovar[7];
	219	Double_t b12=f24*fCovar[11];
	220	Double_t b20=f02fCovar[5] + f04fCovar[12], b21=f12fCovar[5] + f14fCovar[12] + f13*fCovar[8];
	221	Double_t b22=f24*fCovar[12];
	222	Double_t b40=f02fCovar[12] + f04fCovar[14], b41=f12fCovar[12] + f14fCovar[14] + f13*fCovar[13];
	223	Double_t b42=f24*fCovar[14];
	224	Double_t b30=f02fCovar[8] + f04fCovar[13], b31=f12fCovar[8] + f14fCovar[13] + f13*fCovar[9];
	225	Double_t b32=f24*fCovar[13];
	226
	227	//a = fb = fC*ft
	228	Double_t a00=f02b20+f04b40,a01=f02b21+f04b41,a02=f02b22+f04b42;
	229	Double_t a11=f12b21+f14b41+f13b31,a12=f12b22+f14b42+f13b32;
	230	Double_t a22=f24*b42;
	231
	232	//FCFt = C + (b + bt + a)
	233	fCovar[0] += b00 + b00 + a00;
	234	fCovar[1] += b10 + b01 + a01;
	235	fCovar[2] += b11 + b11 + a11;
	236	fCovar[3] += b20 + b02 + a02;
	237	fCovar[4] += b21 + b12 + a12;
	238	fCovar[5] += b22 + b22 + a22;
	239	fCovar[6] += b30;
	240	fCovar[7] += b31;
	241	fCovar[8] += b32;
	242	fCovar[10] += b40;
	243	fCovar[11] += b41;
	244	fCovar[12] += b42;
	245
	246	fX=x2;
	247
	248	//Change of the magnetic field *************
	249	SaveLocalConvConst();
	250	// transform back error matrix from curvature to the 1/pt
	251	fCovar[10]*=lcc;
	252	fCovar[11]*=lcc;
	253	fCovar[12]*=lcc;
	254	fCovar[13]*=lcc;
	255	fCovar[14]=lcclcc;
	256
257	Double_t dist = TMath::Sqrt((fX-oldX)(fX-oldX)+(fParam[0]-oldY)(fParam[0]-oldY)+
258	(fParam[1]-oldZ)*(fParam[1]-oldZ));
259	if (!CorrectForMaterial(dist,x0,rho)) return 0;
260
261	// Integrated Time [SR, GSI, 17.02.2003]
262	// if (IsStartedTimeIntegral() && fX>oldX) {
263	// Double_t l2 = (fX-oldX)(fX-oldX)+(fParam[0]-oldY)(fParam[0]-oldY)+
264	// (fParam[1]-oldZ)*(fParam[1]-oldZ);
265	// AddTimeStep(TMath::Sqrt(l2));
266	// }
267	//
268
269	return kTRUE;
270	}
271
272	Bool_t AliExternalTrackParam::PropagateToDCA(Double_t xd, Double_t yd, Double_t x0, Double_t rho){
273	//
274	// Propagate the track parameters to the nearest point of given xv yv coordinate
275	//
276	Double_t a=fAlpha;
277	Double_t cs=TMath::Cos(a),sn=TMath::Sin(a);
278
279	Double_t xv= xdcs + ydsn;
280	Double_t yv=-xdsn + ydcs; // vertex position in local frame
281	//
282	Double_t c=fParam[4]/GetLocalConvConst(), snp=fParam[2];
283	//
284	Double_t x=fX, y=fParam[1];
285	Double_t tgfv=-(c(x-xv)-snp)/(c(y-yv) + TMath::Sqrt(1.-snp*snp));
286	Double_t fv=TMath::ATan(tgfv);
287	cs=TMath::Cos(fv); sn=TMath::Sin(fv);
288	x = xvcs + yvsn;
289	yv=-xvsn + yvcs; xv=x;
290	RotateTo(fv+a);
291	return PropagateTo(xv,x0,rho);
51ad6848	292	}
51ad6848	293
a5e407e9	294
51ad6848	295	//_____________________________________________________________________________
a5e407e9	296	Bool_t AliExternalTrackParam::RotateTo(Double_t alp)
51ad6848	297	{
a5e407e9	298	// Rotate the reference axis for the parametrisation to the given angle.
	299	//
	300	Double_t x=fX;
	301	Double_t p0=fParam[0];
	302	//
	303	if (alp < -TMath::Pi()) alp += 2*TMath::Pi();
	304	else if (alp >= TMath::Pi()) alp -= 2*TMath::Pi();
	305	Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
	306	Double_t sf=fParam[2], cf=TMath::Sqrt(1.- fParam[2]*fParam[2]);
	307	// ** rotation ********************
	308
	309	fAlpha = alp;
	310	fX = xca + p0sa;
	311	fParam[0]= -xsa + p0ca;
	312	fParam[2]= sfca - cfsa;
	313	Double_t rr=(ca+sf/cf*sa);
	314	//
	315	fCovar[0] = (caca);
	316	fCovar[1] *= ca;
	317	fCovar[3] = carr;
	318	fCovar[6] *= ca;
	319	fCovar[10] *= ca;
	320	fCovar[4] *= rr;
	321	fCovar[5] = rrrr;
	322	fCovar[7] *= rr;
	323	fCovar[11] *= rr;
	324	return kTRUE;
51ad6848	325	}
	326
	327	//_____________________________________________________________________________
a5e407e9	328	Bool_t AliExternalTrackParam::CorrectForMaterial(Double_t d, Double_t x0, Double_t rho)
51ad6848	329	{
a5e407e9	330	//
	331	// Take into account material effects assuming:
	332	// x0 - mean rad length
	333	// rho - mean density
	334
	335	//
	336	// multiple scattering
	337	//
0243e825	338	if (fMass<=0) {
	339	AliError("Non-positive mass");
	340	return kFALSE;
	341	}
a5e407e9	342	Double_t p2=(1.+ fParam[3]fParam[3])/(fParam[4]fParam[4]);
	343	Double_t beta2=p2/(p2 + fMass*fMass);
	344	Double_t theta2=14.114.1/(beta2p21e6)d/x0*rho;
	345	//
	346	fCovar[5] += theta2(1.- fParam[2]fParam[2])(1. + fParam[3]fParam[3]);
	347	fCovar[9] += theta2(1. + fParam[3]fParam[3])(1. + fParam[3]fParam[3]);
	348	fCovar[13] += theta2fParam[3]fParam[4](1. + fParam[3]fParam[3]);
	349	fCovar[14] += theta2fParam[3]fParam[4]fParam[3]fParam[4];
	350
	351	Double_t dE=0.153e-3/beta2(log(5940beta2/(1-beta2+1e-10)) - beta2)drho;
	352	fParam[4] =(1.- TMath::Sqrt(p2+fMassfMass)/p2*dE);
	353	//
	354	Double_t sigmade = 0.02*TMath::Sqrt(TMath::Abs(dE)); // energy loss fluctuation
	355	Double_t sigmac2 = sigmadesigmadefParam[4]fParam[4](p2+fMassfMass)/(p2p2);
	356	fCovar[14] += sigmac2;
	357	//
	358	//
	359
51ad6848	360
a5e407e9	361
a5e407e9	362	return kTRUE;
51ad6848	363	}
	364
	365	//_____________________________________________________________________________
a5e407e9	366	Bool_t AliExternalTrackParam::GetProlongationAt(Double_t xk,
	367	Double_t& y,
	368	Double_t& z) const
51ad6848	369	{
a5e407e9	370	//
	371	// Get the local y and z coordinates at the given x value
	372	//
	373	Double_t lcc=GetLocalConvConst();
	374	Double_t cur = fParam[4]/lcc;
	375	Double_t x1=fX, x2=xk, dx=x2-x1;
	376	Double_t f1=fParam[2], f2=f1 + cur*dx;
	377	//
	378	if (TMath::Abs(f2) >= 0.98) {
	379	// MI change - don't propagate highly inclined tracks
	380	// covariance matrix distorted
	381	return kFALSE;
	382	}
	383	Double_t r1=sqrt(1.- f1f1), r2=sqrt(1.- f2f2);
	384	y = fParam[0] + dx*(f1+f2)/(r1+r2);
	385	z = fParam[1] + dx(f1+f2)/(f1r2 + f2r1)fParam[3];
	386	return kTRUE;
51ad6848	387	}
	388
	389	//_____________________________________________________________________________
	390	Double_t AliExternalTrackParam::GetXAtVertex(Double_t /x/,
	391	Double_t /y/) const
	392	{
	393	// Get the x coordinate at the given vertex (x,y)
	394	//
	395	// NOT IMPLEMENTED for this class
	396
	397	return 0;
	398	}
	399
	400
a5e407e9	401	// //_____________________________________________________________________________
	402	// Double_t AliExternalTrackParam::GetPredictedChi2(const AliCluster* /cluster/)
	403	// {
	404	// // calculate the chi2 contribution of the given cluster
	405	// //
	406	// // NOT IMPLEMENTED for this class
51ad6848	407
a5e407e9	408	// return -1;
a5e407e9	409	// }
51ad6848	410
a5e407e9	411	// //_____________________________________________________________________________
	412	// Bool_t AliExternalTrackParam::Update(const AliCluster* /cluster/)
	413	// {
	414	// // update the track parameters using the position and error
	415	// // of the given cluster
	416	// //
	417	// // NOT IMPLEMENTED for this class
51ad6848	418
a5e407e9	419	// return kFALSE;
a5e407e9	420	// }
51ad6848	421
	422
	423	//_____________________________________________________________________________
	424	Double_t AliExternalTrackParam::SigmaPhi() const
	425	{
	426	// get the error of the azimuthal angle
	427
	428	return TMath::Sqrt(TMath::Abs(fCovar[5] / (1. - fParam[2]*fParam[2])));
	429	}
	430
	431	//_____________________________________________________________________________
	432	Double_t AliExternalTrackParam::SigmaTheta() const
	433	{
	434	// get the error of the polar angle
	435
	436	return TMath::Sqrt(TMath::Abs(fCovar[9])) / (1. + fParam[3]*fParam[3]);
	437	}
	438
	439	//_____________________________________________________________________________
	440	Double_t AliExternalTrackParam::SigmaPt() const
	441	{
	442	// get the error of the transversal component of the momentum
	443
	444	return TMath::Sqrt(fCovar[14]) / TMath::Abs(fParam[4]);
	445	}
	446
	447	//_____________________________________________________________________________
	448	TVector3 AliExternalTrackParam::Momentum() const
	449	{
	450	// get the momentum vector
	451
	452	Double_t phi = TMath::ASin(fParam[2]) + fAlpha;
	453	Double_t pt = 1. / TMath::Abs(fParam[4]);
	454	return TVector3(pt * TMath::Cos(phi),
	455	pt * TMath::Sin(phi),
	456	pt * fParam[3]);
	457	}
	458
	459	//_____________________________________________________________________________
	460	TVector3 AliExternalTrackParam::Position() const
	461	{
	462	// get the current spatial position in global coordinates
	463
	464	return TVector3(fX * TMath::Cos(fAlpha) - fParam[0] * TMath::Sin(fAlpha),
	465	fX * TMath::Sin(fAlpha) + fParam[0] * TMath::Cos(fAlpha),
	466	fParam[1]);
	467	}
	468
	469
	470	//_____________________________________________________________________________
	471	void AliExternalTrackParam::Print(Option_t* /option/) const
	472	{
	473	// print the parameters and the covariance matrix
	474
	475	printf("AliExternalTrackParam: x = %-12g alpha = %-12g\n", fX, fAlpha);
	476	printf(" parameters: %12g %12g %12g %12g %12g\n",
	477	fParam[0], fParam[1], fParam[2], fParam[3], fParam[4]);
	478	printf(" covariance: %12g\n", fCovar[0]);
	479	printf(" %12g %12g\n", fCovar[1], fCovar[2]);
	480	printf(" %12g %12g %12g\n", fCovar[3], fCovar[4], fCovar[5]);
	481	printf(" %12g %12g %12g %12g\n",
	482	fCovar[6], fCovar[7], fCovar[8], fCovar[9]);
	483	printf(" %12g %12g %12g %12g %12g\n",
	484	fCovar[10], fCovar[11], fCovar[12], fCovar[13], fCovar[14]);
485	}