[u/mrichter/AliRoot.git] / PWGPP / TRD / info / AliTRDtrackInfo.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: AliTRDtrackInfo.cxx 27496 2008-07-22 08:35:45Z cblume $ */

////////////////////////////////////////////////////////////////////////////
//                                                                        //
//  Reconstruction QA                                                     //
//                                                                        //
//  Authors:                                                              //
//    Alex Bercuci <A.Bercuci@gsi.de>                                     //
//    Markus Fasel <M.Fasel@gsi.de>                                       //
//                                                                        //
////////////////////////////////////////////////////////////////////////////

#include "TDatabasePDG.h"
#include "TPDGCode.h"
#include "TVectorT.h"

#include "AliTrackReference.h"
#include "AliTrackPointArray.h"
#include "AliExternalTrackParam.h"

#include "AliTRDseedV1.h"
#include "AliTRDtrackV1.h"
#include "AliTRDgeometry.h"
#include "AliTRDtrackerV1.h"

#include "AliTRDtrackInfo.h"

ClassImp(AliTRDtrackInfo)
ClassImp(AliTRDtrackInfo::AliMCinfo)
ClassImp(AliTRDtrackInfo::AliESDinfo)
Double_t AliTRDtrackInfo::AliMCinfo::fgKalmanStep = 2.;
Bool_t AliTRDtrackInfo::AliMCinfo::fgKalmanUpdate = kTRUE;

//___________________________________________________
AliTRDtrackInfo::AliTRDtrackInfo():
  TObject()
  ,fNClusters(0)
  ,fTRDtrack(NULL)
  ,fMC(NULL)
  ,fESD()
{
  //
  // Default constructor
  //
}


//___________________________________________________
AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
  TObject((const TObject&)trdInfo)  
  ,fNClusters(trdInfo.fNClusters)
  ,fTRDtrack(NULL)
  ,fMC(NULL)
  ,fESD(trdInfo.fESD)
{
  //
  // copy Entries
  //

  if(trdInfo.fMC) fMC = new AliMCinfo(*trdInfo.fMC);
  SetTrack(trdInfo.fTRDtrack);
}

//___________________________________________________
AliTRDtrackInfo::AliMCinfo::AliMCinfo()
  :fLabel(0)
  ,fTRDlabel(0)
  ,fPDG(0)
  ,fNTrackRefs(0)
{
  // Set 0-Pointers
  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
}

//___________________________________________________
AliTRDtrackInfo::AliMCinfo::AliMCinfo(const AliMCinfo &mc)
  :fLabel(mc.fLabel)
  ,fTRDlabel(mc.fTRDlabel)
  ,fPDG(mc.fPDG)
  ,fNTrackRefs(mc.fNTrackRefs)
{
  //
  // Constructor
  //

  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
  for(Int_t ien = 0; ien < 12; ien++){
    if(mc.fTrackRefs[ien])
      fTrackRefs[ien] = new AliTrackReference(*(mc.fTrackRefs[ien]));
  }
}

//________________________________________________________
Int_t AliTRDtrackInfo::AliMCinfo::GetPID() const
{
// Translate pdg code to PID index

  switch(fPDG){
  case kElectron: 
  case kPositron: return AliPID::kElectron;  
  case kMuonPlus:
  case kMuonMinus: return AliPID::kMuon;  
  case kPiPlus: 
  case kPiMinus: return AliPID::kPion;  
  case kKPlus: 
  case kKMinus: return AliPID::kKaon;
  case kProton: 
  case kProtonBar: return AliPID::kProton;
  } 
  return -1;
}

//___________________________________________________
AliTRDtrackInfo::AliESDinfo::AliESDinfo()
  :fHasV0(0)
  ,fId(-1)
  ,fStatus(0)
  ,fKinkIndex(0)
  ,fTPCncls(0)
  ,fTPCdedx(0.)
  ,fTOFbc(0)
  ,fTRDpidQuality(0)
  ,fTRDnSlices(0)
  ,fPt(0.)
  ,fPhi(-999.)
  ,fEta(-999.)
  ,fTRDslices(NULL)
  ,fOP(NULL)
  ,fTPCout(NULL)
  ,fITSout(NULL)
  ,fTPArray(NULL)
{
  //
  // Constructor
  //

  memset(fTRDr, 0, AliPID::kSPECIES*sizeof(Double32_t));
  memset(fTRDv0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
}

//___________________________________________________
AliTRDtrackInfo::AliESDinfo::AliESDinfo(const AliESDinfo &esd)
  :fHasV0(esd.fHasV0)
  ,fId(esd.fId)
  ,fStatus(esd.fStatus)
  ,fKinkIndex(esd.fKinkIndex)
  ,fTPCncls(esd.fTPCncls)
  ,fTPCdedx(esd.fTPCdedx)
  ,fTOFbc(esd.fTOFbc)
  ,fTRDpidQuality(esd.fTRDpidQuality)
  ,fTRDnSlices(esd.fTRDnSlices)
  ,fPt(esd.fPt)
  ,fPhi(esd.fPhi)
  ,fEta(esd.fEta)
  ,fTRDslices(NULL)
  ,fOP(NULL)
  ,fTPCout(NULL)
  ,fITSout(NULL)
  ,fTPArray(NULL)
{
  //
  // Constructor
  //

  memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
  memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));

  if(fTRDnSlices){
    fTRDslices = new Double32_t[fTRDnSlices];
    memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
  }
  SetOuterParam(esd.fOP);
  SetTPCoutParam(esd.fTPCout);
  SetITSoutParam(esd.fITSout);
  SetTrackPointArray(esd.fTPArray);
}


//___________________________________________________
AliTRDtrackInfo::~AliTRDtrackInfo()
{
  //
  // Destructor
  //

  if(fMC) delete fMC;
  if(fTRDtrack) delete fTRDtrack;
}

//___________________________________________________
AliTRDtrackInfo::AliMCinfo::~AliMCinfo()
{
  //
  // Destructor
  //

  fNTrackRefs = 0;
  for(Int_t ien = 0; ien < 12; ien++){
    if(fTrackRefs[ien]) delete fTrackRefs[ien];
    fTrackRefs[ien] = NULL;
  }
}

//___________________________________________________
AliTRDtrackInfo::AliESDinfo::~AliESDinfo()
{
  //
  // Destructor
  //

  if(fTRDnSlices){
    delete [] fTRDslices;
    fTRDslices = NULL;
    fTRDnSlices = 0;
  }
  if(fOP) delete fOP; fOP = NULL;
  if(fTPCout) delete fTPCout; fTPCout = NULL;
  if(fITSout) delete fITSout; fITSout = NULL;
  if(fTPArray) delete fTPArray;
}

//___________________________________________________
void AliTRDtrackInfo::AliESDinfo::Delete(const Option_t *){
  //
  // Delete Pointer members 
  // 
  if(fTRDnSlices){
    delete [] fTRDslices;
    fTRDslices = NULL;
    fTRDnSlices = 0;
  }
  if(fOP) delete fOP; fOP = NULL;
  if(fTPCout) delete fTPCout; fTPCout = NULL;
  if(fITSout) delete fITSout; fITSout = NULL;
  if(fTPArray) delete fTPArray; fTPArray = NULL;
}


//___________________________________________________
AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
{
  //
  // = Operator
  //
  if(this == &trdInfo) return *this;

  fNClusters  = trdInfo.fNClusters;
  fESD        = trdInfo.fESD;

  if(trdInfo.fMC){
    if(!fMC) fMC = new AliMCinfo(*trdInfo.fMC);
    else{
      fMC->~AliMCinfo();
      new(fMC) AliMCinfo(*trdInfo.fMC);
    }
  } else {if(fMC) delete fMC; fMC = NULL;}

  SetTrack(trdInfo.fTRDtrack);

  return *this;
}

//___________________________________________________
AliTRDtrackInfo::AliMCinfo& AliTRDtrackInfo::AliMCinfo::operator=(const AliMCinfo &mc)
{
  //
  // Assignment operator
  //

  if(this == &mc) return *this;
  fLabel      = mc.fLabel;
  fTRDlabel   = mc.fTRDlabel;
  fPDG        = mc.fPDG;
  fNTrackRefs = mc.fNTrackRefs;

  AliTrackReference **itr = &fTrackRefs[0];
  AliTrackReference* const *jtr = &mc.fTrackRefs[0];
  for(Int_t ien = 0; ien < 2*AliTRDgeometry::kNlayer; ien++, itr++, jtr++){
    if((*jtr)){
      if(!(*itr)) (*itr) = new AliTrackReference(*(*jtr));
      else{
        (*itr)->~AliTrackReference();
        new(&(*itr)) AliTrackReference(*(*jtr));
      }
    } else {if((*itr)) delete (*itr); (*itr) = NULL;}
  }
  return *this;
}

//___________________________________________________
AliTRDtrackInfo::AliESDinfo& AliTRDtrackInfo::AliESDinfo::operator=(const AliESDinfo &esd)
{
  //
  // Assignment operator
  //

  if(this == &esd) return *this;
  fHasV0       = esd.fHasV0;
  fId          = esd.fId;
  fStatus      = esd.fStatus;
  fKinkIndex   = esd.fKinkIndex;
  fTPCncls     = esd.fTPCncls;
  fTPCdedx     = esd.fTPCdedx;
  fTOFbc       = esd.fTOFbc;
  fTRDpidQuality= esd.fTRDpidQuality;
  fTRDnSlices  = esd.fTRDnSlices;
  fPt          = esd.fPt;
  fPhi         = esd.fPhi;
  fEta         = esd.fEta;
  
  memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
  memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));

  if(fTRDnSlices){
    if(!fTRDslices) fTRDslices = new Double32_t[fTRDnSlices];
    memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
  }
  SetOuterParam(esd.fOP);
  SetTPCoutParam(esd.fTPCout);
  SetITSoutParam(esd.fITSout);
  SetTrackPointArray(esd.fTPArray);

  return *this;
}

//___________________________________________________
void AliTRDtrackInfo::Delete(const Option_t *)
{
  //
  // Delete
  //

  AliDebug(2, Form("track[%p] mc[%p]", (void*)fTRDtrack, (void*)fMC));
  fNClusters  = 0;
  if(fMC) delete fMC; fMC = NULL;
  fESD.Delete(NULL);
  if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
}

//___________________________________________________
void AliTRDtrackInfo::SetTrack(const AliTRDtrackV1 *track)
{
  //
  // Set the TRD track
  //

  if(track){
    if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*track);
    else{
      fTRDtrack->~AliTRDtrackV1();
      new(fTRDtrack) AliTRDtrackV1(*track);
    }
    if(track->IsOwner()) fTRDtrack->SetOwner();
  } else {
    if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
  }
}

//___________________________________________________
void AliTRDtrackInfo::AliESDinfo::SetTrackPointArray(const AliTrackPointArray *tps)
{
  //
  // Set the track point array for alignment task
  //


  if(tps){
    if(!fTPArray) fTPArray = new AliTrackPointArray(*tps);
    else{
      fTPArray->~AliTrackPointArray();
      new(fTPArray) AliTrackPointArray(*tps);
    }
  } else {
    if(fTPArray) delete fTPArray; fTPArray = NULL;
  }
}

//___________________________________________________
void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
{
  //
  // Add track reference
  //

  if(fMC->fNTrackRefs >= 2*AliTRDgeometry::kNlayer){ 
    SetCurved();
    return;
  }
  // Make a copy for the object in order to avoid ownership problems
  fMC->fTrackRefs[fMC->fNTrackRefs++] = new AliTrackReference(*tref);
}

//___________________________________________________
AliTrackReference* AliTRDtrackInfo::GetTrackRef(Int_t idx) const
{
//
// Returns a track reference
//
  if(!fMC) return NULL;
  return (idx>=0 && idx < 12) ? fMC->fTrackRefs[idx] : NULL;
}

//___________________________________________________
AliTrackReference* AliTRDtrackInfo::GetTrackRef(const AliTRDseedV1* const tracklet) const
{
//
// Returns a track reference
//
  if(!fMC) return NULL;
  Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
  AliTrackReference * const* jtr = &(fMC->fTrackRefs[0]);
  for(Int_t itr = 0; itr < fMC->fNTrackRefs; itr++, ++jtr){
    if(!(*jtr)) break;   
    if(TMath::Abs(tracklet->GetX0() - (*jtr)->LocalX()) < cw) return (*jtr);
  }
  return NULL;
}

//___________________________________________________
Int_t AliTRDtrackInfo::GetNumberOfClusters() const
{
  //
  // Returns the number of clusters
  //

  Int_t n = 0;
  if(!fTRDtrack) return 0;
  if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
  AliTRDseedV1 *tracklet = NULL;
  for(Int_t ip=0; ip<AliTRDgeometry::kNlayer; ip++){
    if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
    n+=tracklet->GetN();
  }
  return n;
}


//___________________________________________________
void  AliTRDtrackInfo::AliESDinfo::SetOuterParam(const AliExternalTrackParam *op)
{
  //
  // Set outer track parameters
  //

  if(op){
    if(fOP){
      fOP->~AliExternalTrackParam();
      // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
      new(fOP) AliExternalTrackParam(*op);
    } else fOP = new AliExternalTrackParam(*op);
  } else {
    if(fOP) delete fOP; fOP = NULL;
  }
}

//___________________________________________________
void  AliTRDtrackInfo::AliESDinfo::SetITSoutParam(const AliExternalTrackParam *op)
{
  //
  // Set TPCout track parameters
  //

  if(op){
    if(fITSout){
      fITSout->~AliExternalTrackParam();
      // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
      new(fITSout) AliExternalTrackParam(*op);
    } else fITSout = new AliExternalTrackParam(*op);
  } else {
    if(fITSout) delete fITSout; fITSout = NULL;
  }
}

//___________________________________________________
void  AliTRDtrackInfo::AliESDinfo::SetTPCoutParam(const AliExternalTrackParam *op)
{
  //
  // Set TPCout track parameters
  //

  if(op){
    if(fTPCout){
      fTPCout->~AliExternalTrackParam();
      // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
      new(fTPCout) AliExternalTrackParam(*op);
    } else fTPCout = new AliExternalTrackParam(*op);
  } else {
    if(fTPCout) delete fTPCout; fTPCout = NULL;
  }
}

//___________________________________________________
Int_t AliTRDtrackInfo::GetNTracklets() const
{
  //
  // Return the number of tracklets
  //

  if(!fTRDtrack) return 0;
  return fTRDtrack->GetNumberOfTracklets();
}

//___________________________________________________
void AliTRDtrackInfo::SetSlices(Int_t n, Double32_t *s)
{
  //
  // Set the slices
  //
  if(fESD.fTRDnSlices != n){
    fESD.fTRDnSlices = 0;
    delete [] fESD.fTRDslices;
    fESD.fTRDslices = NULL;
  }

  if(!fESD.fTRDnSlices){
    fESD.fTRDnSlices = n;
    fESD.fTRDslices = new Double32_t[fESD.fTRDnSlices];
  }

  memcpy(fESD.fTRDslices, s, n*sizeof(Double32_t));
}
 
//___________________________________________________
Bool_t AliTRDtrackInfo::AliMCinfo::GetDirections(Float_t &x0, Float_t &y0, Float_t &z0, Float_t &dydx, Float_t &dzdx, Float_t &pt, Float_t &eta, Float_t &phi, UChar_t &status) const
{
// Check for 2 track ref for the tracklet defined bythe radial position x0
// The "status" is a bit map and gives a more informative output in case of failure:
//   - 0 : everything is OK
//   - BIT(0) : 0 track refs found
//   - BIT(1) : 1 track reference found
//   - BIT(2) : dx <= 0 between track references
//   - BIT(3) : dx > 0 && dx < 3.7 - tangent tracks 

  status = 0;
  Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
  Int_t nFound = 0;
  AliTrackReference *tr[2] = {NULL, NULL};
  AliTrackReference * const* jtr = &fTrackRefs[0];
  for(Int_t itr = 0; itr < fNTrackRefs; itr++, ++jtr){
    if(!(*jtr)) break;
/*
    if(fDebugLevel>=5) printf("\t\tref[%2d] x[%6.3f]\n", itr, (*jtr)->LocalX());*/
    if(TMath::Abs(x0 - (*jtr)->LocalX()) > cw) continue;
    tr[nFound++] = (*jtr);
    if(nFound == 2) break;
  } 
  if(nFound < 2){ 
    //AliDebug(1, Form("Missing track ref x0[%6.3f] nref[%d]", x0, nFound));
    if(!nFound) SETBIT(status, 0);
    else SETBIT(status, 1);
    return kFALSE;
  }
  pt=tr[1]->Pt();
  if(pt < 1.e-3) return kFALSE;

  Double_t dx = tr[1]->LocalX() - tr[0]->LocalX();
  if(dx <= 0.){
    AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Track ref with wrong radial distances refX0[%6.3f] refX1[%6.3f]", tr[0]->LocalX(), tr[1]->LocalX()));
    SETBIT(status, 2);
    return kFALSE;
  }
  if(TMath::Abs(dx-AliTRDgeometry::CamHght()-AliTRDgeometry::CdrHght())>1.E-3) SETBIT(status, 3); 

  dydx = (tr[1]->LocalY() - tr[0]->LocalY()) / dx;
  dzdx = (tr[1]->Z() - tr[0]->Z()) / dx;
  //Float_t dx0 = tr[1]->LocalX() - x0;
  y0   =  tr[1]->LocalY()/* - dydx*dx0*/;
  z0   =  tr[1]->Z()/* - dzdx*dx0*/;
  x0   =  tr[1]->LocalX();
  eta  =  -TMath::Log(TMath::Tan(0.5 * tr[1]->Theta()));
  phi  =  TMath::ATan2(tr[1]->Y(), tr[1]->X());
  return kTRUE;
}

//___________________________________________________
Bool_t AliTRDtrackInfo::AliMCinfo::PropagateKalman(
      TVectorD *x, TVectorD *y, TVectorD *z,
      TVectorD *dx, TVectorD *dy, TVectorD *dz,
      TVectorD *pt, TVectorD *dpt, TVectorD *budget, TVectorD *c, Double_t mass) const
{
// Propagate Kalman from the first TRD track reference to 
// last one and save residuals in the y, z and pt.
// 
// This is to calibrate the dEdx and MS corrections

  if(!fNTrackRefs) return kFALSE;
  for(Int_t itr=kNTrackRefs; itr--;){
    (*x)[itr] = 0.;(*y)[itr] = 0.;(*z)[itr] = 0.;
    (*dx)[itr] = -1.; (*dy)[itr] = 100.; (*dz)[itr] = 100.; (*dpt)[itr] = 100.;
  }

  // Initialize TRD track to the first track reference
  AliTrackReference *tr(NULL);
  Int_t itr(0); while(!(tr = fTrackRefs[itr])) itr++;
  if(tr->Pt()<1.e-3) return kFALSE;

  AliTRDtrackV1 tt;
  Double_t xyz[3]={tr->X(),tr->Y(),tr->Z()};
  Double_t pxyz[3]={tr->Px(),tr->Py(),tr->Pz()};
  Double_t var[6] = {1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4};
  Double_t cov[21]={
    var[0],  0.,  0.,  0.,  0.,  0.,
         var[1],  0.,  0.,  0.,  0.,
              var[2],  0.,  0.,  0.,
                   var[3],  0.,  0.,
                        var[4],  0.,
                             var[5]
  };
  TDatabasePDG db;
  const TParticlePDG *pdg=db.GetParticle(fPDG);
  if(!pdg){
    AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::PropagateKalman()", Form("PDG entry missing for code %d. References for track %d", fPDG, fNTrackRefs));
    return kFALSE;
  }
  tt.Set(xyz, pxyz, cov, Short_t(pdg->Charge()));
  if(mass<0){ // mass 0 use PDG
    tt.SetMass(pdg->Mass());
  } else { // use rec value
    tt.SetMass(mass);
  }

//  Double_t bg(tr->P()/pdg->Mass());
//  printf("\n\nNEW track PDG[%d] bg[%f] x[%f]\n", fPDG, bg, TMath::Log(bg));
  Double_t x0(tr->LocalX());
  const Double_t *cc(NULL);
  for(Int_t ip=0; itr<fNTrackRefs; itr++){
    if(!(tr = fTrackRefs[itr])) continue;
//    printf("ip[%d] det[%d]\n", ip, tr->DetectorId());
    if(!AliTRDtrackerV1::PropagateToX(tt, tr->LocalX(), fgKalmanStep)) continue;

    //if(update) ...
    tt.GetXYZ(xyz);
    (*x)[ip]   = xyz[0];
    (*y)[ip]   = xyz[1];
    (*z)[ip]   = xyz[2];
    (*dx)[ip]  = tt.GetX() - x0;
    (*dy)[ip]  = tt.GetY() - tr->LocalY();
    (*dz)[ip]  = tt.GetZ() - tr->Z();
    (*pt)[ip]  = tr->Pt();
    (*dpt)[ip] = tt.Pt()- tr->Pt();
//    printf("pt : kalman[%e] MC[%e]\n", tt.Pt(), tr->Pt());
    (*budget)[ip] = tt.GetBudget(0);
    cc = tt.GetCovariance();
//    printf("dx[%5.2f] sy[%f]\n", (*dx)[ip], TMath::Sqrt(cc[0]));
    for(Int_t ic(0), jp(ip*15); ic<15; ic++, jp++) (*c)[jp]=cc[ic];
    ip++;
  }
  return kTRUE;
}
Commit	Line	Data
1ee39b3a	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: AliTRDtrackInfo.cxx 27496 2008-07-22 08:35:45Z cblume $ */
	17
	18	////////////////////////////////////////////////////////////////////////////
	19	// //
	20	// Reconstruction QA //
	21	// //
	22	// Authors: //
	23	// Alex Bercuci <A.Bercuci@gsi.de> //
	24	// Markus Fasel <M.Fasel@gsi.de> //
	25	// //
	26	////////////////////////////////////////////////////////////////////////////
	27
	28	#include "TDatabasePDG.h"
b9ddd472	29	#include "TPDGCode.h"
4226db3e	30	#include "TVectorT.h"
1ee39b3a	31
1ee39b3a	32	#include "AliTrackReference.h"
eb05d549	33	#include "AliTrackPointArray.h"
1ee39b3a	34	#include "AliExternalTrackParam.h"
	35
	36	#include "AliTRDseedV1.h"
	37	#include "AliTRDtrackV1.h"
	38	#include "AliTRDgeometry.h"
	39	#include "AliTRDtrackerV1.h"
	40
	41	#include "AliTRDtrackInfo.h"
	42
	43	ClassImp(AliTRDtrackInfo)
	44	ClassImp(AliTRDtrackInfo::AliMCinfo)
	45	ClassImp(AliTRDtrackInfo::AliESDinfo)
5066aa9a	46	Double_t AliTRDtrackInfo::AliMCinfo::fgKalmanStep = 2.;
3ceb45ae	47	Bool_t AliTRDtrackInfo::AliMCinfo::fgKalmanUpdate = kTRUE;
1ee39b3a	48
	49	//___________________________________________________
	50	AliTRDtrackInfo::AliTRDtrackInfo():
	51	TObject()
	52	,fNClusters(0)
f232df0d	53	,fTRDtrack(NULL)
f232df0d	54	,fMC(NULL)
1ee39b3a	55	,fESD()
	56	{
	57	//
	58	// Default constructor
	59	//
	60	}
	61
	62
	63	//___________________________________________________
	64	AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
	65	TObject((const TObject&)trdInfo)
	66	,fNClusters(trdInfo.fNClusters)
f232df0d	67	,fTRDtrack(NULL)
a340fe39	68	,fMC(NULL)
1ee39b3a	69	,fESD(trdInfo.fESD)
	70	{
	71	//
	72	// copy Entries
	73	//
	74
	75	if(trdInfo.fMC) fMC = new AliMCinfo(*trdInfo.fMC);
eb05d549	76	SetTrack(trdInfo.fTRDtrack);
1ee39b3a	77	}
	78
	79	//___________________________________________________
	80	AliTRDtrackInfo::AliMCinfo::AliMCinfo()
	81	:fLabel(0)
eb05d549	82	,fTRDlabel(0)
1ee39b3a	83	,fPDG(0)
	84	,fNTrackRefs(0)
	85	{
	86	// Set 0-Pointers
	87	memset(fTrackRefs, 0, sizeof(AliTrackReference ) 12);
	88	}
	89
	90	//___________________________________________________
	91	AliTRDtrackInfo::AliMCinfo::AliMCinfo(const AliMCinfo &mc)
	92	:fLabel(mc.fLabel)
eb05d549	93	,fTRDlabel(mc.fTRDlabel)
1ee39b3a	94	,fPDG(mc.fPDG)
	95	,fNTrackRefs(mc.fNTrackRefs)
	96	{
	97	//
	98	// Constructor
	99	//
	100
	101	memset(fTrackRefs, 0, sizeof(AliTrackReference ) 12);
	102	for(Int_t ien = 0; ien < 12; ien++){
	103	if(mc.fTrackRefs[ien])
	104	fTrackRefs[ien] = new AliTrackReference(*(mc.fTrackRefs[ien]));
	105	}
	106	}
	107
b9ddd472	108	//________________________________________________________
	109	Int_t AliTRDtrackInfo::AliMCinfo::GetPID() const
	110	{
61f6b45e	111	// Translate pdg code to PID index
61f6b45e	112
b9ddd472	113	switch(fPDG){
	114	case kElectron:
	115	case kPositron: return AliPID::kElectron;
	116	case kMuonPlus:
	117	case kMuonMinus: return AliPID::kMuon;
	118	case kPiPlus:
	119	case kPiMinus: return AliPID::kPion;
	120	case kKPlus:
	121	case kKMinus: return AliPID::kKaon;
	122	case kProton:
	123	case kProtonBar: return AliPID::kProton;
	124	}
	125	return -1;
	126	}
	127
1ee39b3a	128	//___________________________________________________
1ee39b3a	129	AliTRDtrackInfo::AliESDinfo::AliESDinfo()
b9ddd472	130	:fHasV0(0)
b9ddd472	131	,fId(-1)
1ee39b3a	132	,fStatus(0)
	133	,fKinkIndex(0)
	134	,fTPCncls(0)
2b33fd9c	135	,fTPCdedx(0.)
3ceb45ae	136	,fTOFbc(0)
1ee39b3a	137	,fTRDpidQuality(0)
1ee39b3a	138	,fTRDnSlices(0)
eb05d549	139	,fPt(0.)
	140	,fPhi(-999.)
	141	,fEta(-999.)
f232df0d	142	,fTRDslices(NULL)
f232df0d	143	,fOP(NULL)
35983729	144	,fTPCout(NULL)
b9058d72	145	,fITSout(NULL)
eb05d549	146	,fTPArray(NULL)
1ee39b3a	147	{
	148	//
	149	// Constructor
	150	//
	151
	152	memset(fTRDr, 0, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00	153	memset(fTRDv0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a	154	}
	155
	156	//___________________________________________________
	157	AliTRDtrackInfo::AliESDinfo::AliESDinfo(const AliESDinfo &esd)
b9ddd472	158	:fHasV0(esd.fHasV0)
b9ddd472	159	,fId(esd.fId)
1ee39b3a	160	,fStatus(esd.fStatus)
	161	,fKinkIndex(esd.fKinkIndex)
	162	,fTPCncls(esd.fTPCncls)
2b33fd9c	163	,fTPCdedx(esd.fTPCdedx)
3ceb45ae	164	,fTOFbc(esd.fTOFbc)
1ee39b3a	165	,fTRDpidQuality(esd.fTRDpidQuality)
1ee39b3a	166	,fTRDnSlices(esd.fTRDnSlices)
eb05d549	167	,fPt(esd.fPt)
	168	,fPhi(esd.fPhi)
	169	,fEta(esd.fEta)
f232df0d	170	,fTRDslices(NULL)
f232df0d	171	,fOP(NULL)
35983729	172	,fTPCout(NULL)
b9058d72	173	,fITSout(NULL)
eb05d549	174	,fTPArray(NULL)
1ee39b3a	175	{
	176	//
	177	// Constructor
	178	//
	179
	180	memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00	181	memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a	182
	183	if(fTRDnSlices){
	184	fTRDslices = new Double32_t[fTRDnSlices];
	185	memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
	186	}
eb05d549	187	SetOuterParam(esd.fOP);
	188	SetTPCoutParam(esd.fTPCout);
	189	SetITSoutParam(esd.fITSout);
	190	SetTrackPointArray(esd.fTPArray);
1ee39b3a	191	}
	192
	193
	194	//___________________________________________________
	195	AliTRDtrackInfo::~AliTRDtrackInfo()
	196	{
	197	//
	198	// Destructor
	199	//
	200
	201	if(fMC) delete fMC;
	202	if(fTRDtrack) delete fTRDtrack;
	203	}
	204
	205	//___________________________________________________
	206	AliTRDtrackInfo::AliMCinfo::~AliMCinfo()
	207	{
	208	//
	209	// Destructor
	210	//
	211
	212	fNTrackRefs = 0;
	213	for(Int_t ien = 0; ien < 12; ien++){
	214	if(fTrackRefs[ien]) delete fTrackRefs[ien];
f232df0d	215	fTrackRefs[ien] = NULL;
1ee39b3a	216	}
	217	}
	218
	219	//___________________________________________________
	220	AliTRDtrackInfo::AliESDinfo::~AliESDinfo()
	221	{
	222	//
	223	// Destructor
	224	//
	225
	226	if(fTRDnSlices){
	227	delete [] fTRDslices;
f232df0d	228	fTRDslices = NULL;
1ee39b3a	229	fTRDnSlices = 0;
1ee39b3a	230	}
f232df0d	231	if(fOP) delete fOP; fOP = NULL;
35983729	232	if(fTPCout) delete fTPCout; fTPCout = NULL;
b9058d72	233	if(fITSout) delete fITSout; fITSout = NULL;
eb05d549	234	if(fTPArray) delete fTPArray;
1ee39b3a	235	}
1ee39b3a	236
8778aca3	237	//___________________________________________________
	238	void AliTRDtrackInfo::AliESDinfo::Delete(const Option_t *){
	239	//
	240	// Delete Pointer members
	241	//
	242	if(fTRDnSlices){
	243	delete [] fTRDslices;
	244	fTRDslices = NULL;
	245	fTRDnSlices = 0;
	246	}
	247	if(fOP) delete fOP; fOP = NULL;
35983729	248	if(fTPCout) delete fTPCout; fTPCout = NULL;
b9058d72	249	if(fITSout) delete fITSout; fITSout = NULL;
eb05d549	250	if(fTPArray) delete fTPArray; fTPArray = NULL;
8778aca3	251	}
8778aca3	252
1ee39b3a	253
	254	//___________________________________________________
	255	AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
	256	{
	257	//
	258	// = Operator
	259	//
01ccc21a	260	if(this == &trdInfo) return *this;
1ee39b3a	261
1ee39b3a	262	fNClusters = trdInfo.fNClusters;
1881bb2b	263	fESD = trdInfo.fESD;
1ee39b3a	264
1ee39b3a	265	if(trdInfo.fMC){
2c33fb46	266	if(!fMC) fMC = new AliMCinfo(*trdInfo.fMC);
	267	else{
	268	fMC->~AliMCinfo();
1ee39b3a	269	new(fMC) AliMCinfo(*trdInfo.fMC);
2c33fb46	270	}
1881bb2b	271	} else {if(fMC) delete fMC; fMC = NULL;}
1ee39b3a	272
eb05d549	273	SetTrack(trdInfo.fTRDtrack);
1ee39b3a	274
	275	return *this;
	276	}
	277
	278	//___________________________________________________
	279	AliTRDtrackInfo::AliMCinfo& AliTRDtrackInfo::AliMCinfo::operator=(const AliMCinfo &mc)
	280	{
	281	//
	282	// Assignment operator
	283	//
	284
01ccc21a	285	if(this == &mc) return *this;
1ee39b3a	286	fLabel = mc.fLabel;
eb05d549	287	fTRDlabel = mc.fTRDlabel;
1ee39b3a	288	fPDG = mc.fPDG;
	289	fNTrackRefs = mc.fNTrackRefs;
	290
	291	AliTrackReference **itr = &fTrackRefs[0];
	292	AliTrackReference* const *jtr = &mc.fTrackRefs[0];
eb05d549	293	for(Int_t ien = 0; ien < 2*AliTRDgeometry::kNlayer; ien++, itr++, jtr++){
1ee39b3a	294	if((*jtr)){
1ee39b3a	295	if(!(itr)) (itr) = new AliTrackReference((jtr));
2c33fb46	296	else{
	297	(*itr)->~AliTrackReference();
	298	new(&(itr)) AliTrackReference((*jtr));
	299	}
1881bb2b	300	} else {if((itr)) delete (itr); (*itr) = NULL;}
1ee39b3a	301	}
	302	return *this;
	303	}
	304
	305	//___________________________________________________
	306	AliTRDtrackInfo::AliESDinfo& AliTRDtrackInfo::AliESDinfo::operator=(const AliESDinfo &esd)
	307	{
	308	//
	309	// Assignment operator
	310	//
	311
01ccc21a	312	if(this == &esd) return *this;
01ccc21a	313	fHasV0 = esd.fHasV0;
1ee39b3a	314	fId = esd.fId;
01ccc21a	315	fStatus = esd.fStatus;
	316	fKinkIndex = esd.fKinkIndex;
	317	fTPCncls = esd.fTPCncls;
2b33fd9c	318	fTPCdedx = esd.fTPCdedx;
01ccc21a	319	fTOFbc = esd.fTOFbc;
01ccc21a	320	fTRDpidQuality= esd.fTRDpidQuality;
1ee39b3a	321	fTRDnSlices = esd.fTRDnSlices;
eb05d549	322	fPt = esd.fPt;
	323	fPhi = esd.fPhi;
	324	fEta = esd.fEta;
178d284a	325
1ee39b3a	326	memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
01ccc21a	327	memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a	328
1ee39b3a	329	if(fTRDnSlices){
56838f9d	330	if(!fTRDslices) fTRDslices = new Double32_t[fTRDnSlices];
1ee39b3a	331	memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
1ee39b3a	332	}
eb05d549	333	SetOuterParam(esd.fOP);
	334	SetTPCoutParam(esd.fTPCout);
	335	SetITSoutParam(esd.fITSout);
	336	SetTrackPointArray(esd.fTPArray);
1ee39b3a	337
	338	return *this;
	339	}
	340
	341	//___________________________________________________
	342	void AliTRDtrackInfo::Delete(const Option_t *)
	343	{
	344	//
	345	// Delete
	346	//
	347
2c33fb46	348	AliDebug(2, Form("track[%p] mc[%p]", (void)fTRDtrack, (void)fMC));
1ee39b3a	349	fNClusters = 0;
f232df0d	350	if(fMC) delete fMC; fMC = NULL;
8778aca3	351	fESD.Delete(NULL);
f232df0d	352	if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
1ee39b3a	353	}
	354
	355	//___________________________________________________
	356	void AliTRDtrackInfo::SetTrack(const AliTRDtrackV1 *track)
	357	{
	358	//
	359	// Set the TRD track
	360	//
	361
eb05d549	362	if(track){
	363	if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*track);
	364	else{
	365	fTRDtrack->~AliTRDtrackV1();
	366	new(fTRDtrack) AliTRDtrackV1(*track);
	367	}
	368	if(track->IsOwner()) fTRDtrack->SetOwner();
	369	} else {
	370	if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
	371	}
	372	}
	373
	374	//___________________________________________________
	375	void AliTRDtrackInfo::AliESDinfo::SetTrackPointArray(const AliTrackPointArray *tps)
	376	{
	377	//
	378	// Set the track point array for alignment task
	379	//
	380
	381
	382	if(tps){
	383	if(!fTPArray) fTPArray = new AliTrackPointArray(*tps);
	384	else{
	385	fTPArray->~AliTrackPointArray();
	386	new(fTPArray) AliTrackPointArray(*tps);
	387	}
	388	} else {
	389	if(fTPArray) delete fTPArray; fTPArray = NULL;
2c33fb46	390	}
1ee39b3a	391	}
	392
	393	//___________________________________________________
	394	void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
	395	{
	396	//
	397	// Add track reference
	398	//
	399
eb05d549	400	if(fMC->fNTrackRefs >= 2*AliTRDgeometry::kNlayer){
1ee39b3a	401	SetCurved();
	402	return;
	403	}
	404	// Make a copy for the object in order to avoid ownership problems
	405	fMC->fTrackRefs[fMC->fNTrackRefs++] = new AliTrackReference(*tref);
	406	}
	407
	408	//___________________________________________________
	409	AliTrackReference* AliTRDtrackInfo::GetTrackRef(Int_t idx) const
	410	{
	411	//
	412	// Returns a track reference
	413	//
f232df0d	414	if(!fMC) return NULL;
f232df0d	415	return (idx>=0 && idx < 12) ? fMC->fTrackRefs[idx] : NULL;
1ee39b3a	416	}
	417
	418	//___________________________________________________
61f6b45e	419	AliTrackReference* AliTRDtrackInfo::GetTrackRef(const AliTRDseedV1* const tracklet) const
1ee39b3a	420	{
	421	//
	422	// Returns a track reference
	423	//
f232df0d	424	if(!fMC) return NULL;
1ee39b3a	425	Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
	426	AliTrackReference * const* jtr = &(fMC->fTrackRefs[0]);
	427	for(Int_t itr = 0; itr < fMC->fNTrackRefs; itr++, ++jtr){
	428	if(!(*jtr)) break;
	429	if(TMath::Abs(tracklet->GetX0() - (jtr)->LocalX()) < cw) return (jtr);
	430	}
f232df0d	431	return NULL;
1ee39b3a	432	}
	433
	434	//___________________________________________________
	435	Int_t AliTRDtrackInfo::GetNumberOfClusters() const
	436	{
	437	//
	438	// Returns the number of clusters
	439	//
	440
	441	Int_t n = 0;
	442	if(!fTRDtrack) return 0;
	443	if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
f232df0d	444	AliTRDseedV1 *tracklet = NULL;
eb05d549	445	for(Int_t ip=0; ip<AliTRDgeometry::kNlayer; ip++){
1ee39b3a	446	if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
	447	n+=tracklet->GetN();
	448	}
	449	return n;
	450	}
	451
	452
	453	//___________________________________________________
eb05d549	454	void AliTRDtrackInfo::AliESDinfo::SetOuterParam(const AliExternalTrackParam *op)
1ee39b3a	455	{
	456	//
	457	// Set outer track parameters
	458	//
	459
eb05d549	460	if(op){
	461	if(fOP){
	462	fOP->~AliExternalTrackParam();
	463	// RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
	464	new(fOP) AliExternalTrackParam(*op);
	465	} else fOP = new AliExternalTrackParam(*op);
	466	} else {
	467	if(fOP) delete fOP; fOP = NULL;
	468	}
1ee39b3a	469	}
1ee39b3a	470
b9058d72	471	//___________________________________________________
eb05d549	472	void AliTRDtrackInfo::AliESDinfo::SetITSoutParam(const AliExternalTrackParam *op)
b9058d72	473	{
	474	//
	475	// Set TPCout track parameters
	476	//
	477
eb05d549	478	if(op){
	479	if(fITSout){
	480	fITSout->~AliExternalTrackParam();
	481	// RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
	482	new(fITSout) AliExternalTrackParam(*op);
	483	} else fITSout = new AliExternalTrackParam(*op);
	484	} else {
	485	if(fITSout) delete fITSout; fITSout = NULL;
	486	}
b9058d72	487	}
b9058d72	488
35983729	489	//___________________________________________________
eb05d549	490	void AliTRDtrackInfo::AliESDinfo::SetTPCoutParam(const AliExternalTrackParam *op)
35983729	491	{
	492	//
	493	// Set TPCout track parameters
	494	//
	495
eb05d549	496	if(op){
	497	if(fTPCout){
	498	fTPCout->~AliExternalTrackParam();
	499	// RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
	500	new(fTPCout) AliExternalTrackParam(*op);
	501	} else fTPCout = new AliExternalTrackParam(*op);
	502	} else {
	503	if(fTPCout) delete fTPCout; fTPCout = NULL;
	504	}
35983729	505	}
35983729	506
1ee39b3a	507	//___________________________________________________
	508	Int_t AliTRDtrackInfo::GetNTracklets() const
	509	{
	510	//
	511	// Return the number of tracklets
	512	//
	513
e2927481	514	if(!fTRDtrack) return 0;
1ee39b3a	515	return fTRDtrack->GetNumberOfTracklets();
	516	}
	517
	518	//___________________________________________________
	519	void AliTRDtrackInfo::SetSlices(Int_t n, Double32_t *s)
	520	{
	521	//
	522	// Set the slices
	523	//
1ee39b3a	524	if(fESD.fTRDnSlices != n){
	525	fESD.fTRDnSlices = 0;
	526	delete [] fESD.fTRDslices;
f232df0d	527	fESD.fTRDslices = NULL;
1ee39b3a	528	}
	529
	530	if(!fESD.fTRDnSlices){
	531	fESD.fTRDnSlices = n;
f232df0d	532	fESD.fTRDslices = new Double32_t[fESD.fTRDnSlices];
1ee39b3a	533	}
	534
	535	memcpy(fESD.fTRDslices, s, n*sizeof(Double32_t));
	536	}
	537
	538	//___________________________________________________
f429b017	539	Bool_t AliTRDtrackInfo::AliMCinfo::GetDirections(Float_t &x0, Float_t &y0, Float_t &z0, Float_t &dydx, Float_t &dzdx, Float_t &pt, Float_t &eta, Float_t &phi, UChar_t &status) const
1ee39b3a	540	{
	541	// Check for 2 track ref for the tracklet defined bythe radial position x0
	542	// The "status" is a bit map and gives a more informative output in case of failure:
	543	// - 0 : everything is OK
	544	// - BIT(0) : 0 track refs found
	545	// - BIT(1) : 1 track reference found
	546	// - BIT(2) : dx <= 0 between track references
	547	// - BIT(3) : dx > 0 && dx < 3.7 - tangent tracks
	548
	549	status = 0;
	550	Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
	551	Int_t nFound = 0;
f232df0d	552	AliTrackReference *tr[2] = {NULL, NULL};
1ee39b3a	553	AliTrackReference * const* jtr = &fTrackRefs[0];
	554	for(Int_t itr = 0; itr < fNTrackRefs; itr++, ++jtr){
	555	if(!(*jtr)) break;
	556	/*
	557	if(fDebugLevel>=5) printf("\t\tref[%2d] x[%6.3f]\n", itr, (jtr)->LocalX());/
	558	if(TMath::Abs(x0 - (*jtr)->LocalX()) > cw) continue;
	559	tr[nFound++] = (*jtr);
	560	if(nFound == 2) break;
	561	}
	562	if(nFound < 2){
e02d9eb3	563	//AliDebug(1, Form("Missing track ref x0[%6.3f] nref[%d]", x0, nFound));
1ee39b3a	564	if(!nFound) SETBIT(status, 0);
	565	else SETBIT(status, 1);
	566	return kFALSE;
	567	}
61f6b45e	568	pt=tr[1]->Pt();
61f6b45e	569	if(pt < 1.e-3) return kFALSE;
1ee39b3a	570
	571	Double_t dx = tr[1]->LocalX() - tr[0]->LocalX();
	572	if(dx <= 0.){
	573	AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Track ref with wrong radial distances refX0[%6.3f] refX1[%6.3f]", tr[0]->LocalX(), tr[1]->LocalX()));
	574	SETBIT(status, 2);
	575	return kFALSE;
	576	}
	577	if(TMath::Abs(dx-AliTRDgeometry::CamHght()-AliTRDgeometry::CdrHght())>1.E-3) SETBIT(status, 3);
	578
	579	dydx = (tr[1]->LocalY() - tr[0]->LocalY()) / dx;
	580	dzdx = (tr[1]->Z() - tr[0]->Z()) / dx;
	581	//Float_t dx0 = tr[1]->LocalX() - x0;
	582	y0 = tr[1]->LocalY()/* - dydxdx0/;
	583	z0 = tr[1]->Z()/* - dzdxdx0/;
	584	x0 = tr[1]->LocalX();
3ceb45ae	585	eta = -TMath::Log(TMath::Tan(0.5 * tr[1]->Theta()));
01ccc21a	586	phi = TMath::ATan2(tr[1]->Y(), tr[1]->X());
1ee39b3a	587	return kTRUE;
	588	}
	589
	590	//___________________________________________________
3ceb45ae	591	Bool_t AliTRDtrackInfo::AliMCinfo::PropagateKalman(
	592	TVectorD x, TVectorD y, TVectorD *z,
	593	TVectorD dx, TVectorD dy, TVectorD *dz,
	594	TVectorD pt, TVectorD dpt, TVectorD budget, TVectorD c, Double_t mass) const
1ee39b3a	595	{
	596	// Propagate Kalman from the first TRD track reference to
	597	// last one and save residuals in the y, z and pt.
	598	//
	599	// This is to calibrate the dEdx and MS corrections
	600
3ceb45ae	601	if(!fNTrackRefs) return kFALSE;
1ee39b3a	602	for(Int_t itr=kNTrackRefs; itr--;){
3ceb45ae	603	(x)[itr] = 0.;(y)[itr] = 0.;(*z)[itr] = 0.;
4226db3e	604	(dx)[itr] = -1.; (dy)[itr] = 100.; (dz)[itr] = 100.; (dpt)[itr] = 100.;
1ee39b3a	605	}
1ee39b3a	606
1ee39b3a	607	// Initialize TRD track to the first track reference
4226db3e	608	AliTrackReference *tr(NULL);
4226db3e	609	Int_t itr(0); while(!(tr = fTrackRefs[itr])) itr++;
3ceb45ae	610	if(tr->Pt()<1.e-3) return kFALSE;
1ee39b3a	611
	612	AliTRDtrackV1 tt;
	613	Double_t xyz[3]={tr->X(),tr->Y(),tr->Z()};
	614	Double_t pxyz[3]={tr->Px(),tr->Py(),tr->Pz()};
	615	Double_t var[6] = {1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4};
	616	Double_t cov[21]={
	617	var[0], 0., 0., 0., 0., 0.,
	618	var[1], 0., 0., 0., 0.,
	619	var[2], 0., 0., 0.,
	620	var[3], 0., 0.,
	621	var[4], 0.,
	622	var[5]
	623	};
	624	TDatabasePDG db;
	625	const TParticlePDG *pdg=db.GetParticle(fPDG);
	626	if(!pdg){
	627	AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::PropagateKalman()", Form("PDG entry missing for code %d. References for track %d", fPDG, fNTrackRefs));
3ceb45ae	628	return kFALSE;
1ee39b3a	629	}
1ee39b3a	630	tt.Set(xyz, pxyz, cov, Short_t(pdg->Charge()));
3ceb45ae	631	if(mass<0){ // mass 0 use PDG
	632	tt.SetMass(pdg->Mass());
	633	} else { // use rec value
	634	tt.SetMass(mass);
	635	}
	636
	637	// Double_t bg(tr->P()/pdg->Mass());
	638	// printf("\n\nNEW track PDG[%d] bg[%f] x[%f]\n", fPDG, bg, TMath::Log(bg));
4226db3e	639	Double_t x0(tr->LocalX());
	640	const Double_t *cc(NULL);
	641	for(Int_t ip=0; itr<fNTrackRefs; itr++){
	642	if(!(tr = fTrackRefs[itr])) continue;
3ceb45ae	643	// printf("ip[%d] det[%d]\n", ip, tr->DetectorId());
5066aa9a	644	if(!AliTRDtrackerV1::PropagateToX(tt, tr->LocalX(), fgKalmanStep)) continue;
1ee39b3a	645
1ee39b3a	646	//if(update) ...
3ceb45ae	647	tt.GetXYZ(xyz);
	648	(*x)[ip] = xyz[0];
	649	(*y)[ip] = xyz[1];
	650	(*z)[ip] = xyz[2];
4226db3e	651	(*dx)[ip] = tt.GetX() - x0;
	652	(*dy)[ip] = tt.GetY() - tr->LocalY();
	653	(*dz)[ip] = tt.GetZ() - tr->Z();
5066aa9a	654	(*pt)[ip] = tr->Pt();
4226db3e	655	(*dpt)[ip] = tt.Pt()- tr->Pt();
3ceb45ae	656	// printf("pt : kalman[%e] MC[%e]\n", tt.Pt(), tr->Pt());
3ceb45ae	657	(*budget)[ip] = tt.GetBudget(0);
1ee39b3a	658	cc = tt.GetCovariance();
3ceb45ae	659	// printf("dx[%5.2f] sy[%f]\n", (*dx)[ip], TMath::Sqrt(cc[0]));
3ceb45ae	660	for(Int_t ic(0), jp(ip15); ic<15; ic++, jp++) (c)[jp]=cc[ic];
1ee39b3a	661	ip++;
1ee39b3a	662	}
3ceb45ae	663	return kTRUE;
1ee39b3a	664	}