PWG1/AliESDRecV0Info.cxx

   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
  17 ///////////////////////////////////////////////////////////////////////////////
  18 //                                                                           //
  19 //                                                                           //
  20 //  Comparison class for V0 information                                      //
  21 //  responsible:
  22 //  marian.ivanov@cern.ch                                                    //
  23 //
  24 //
  25
  26
  27
  28
  29
  30 #include <stdio.h>
  31 #include <string.h>
  32 //ROOT includes
  33 #include "Rtypes.h"
  34 //
  35 //ALIROOT includes
  36 //
  37 #include "AliESDtrack.h"
  38 #include "AliTPCParam.h"
  39 #include "AliTrackReference.h"
  40 #include "AliTPCParamSR.h"
  41 #include "AliESD.h"
  42 #include "AliESDfriend.h"
  43 #include "AliESDtrack.h"
  44 #include "AliTPCseed.h"
  45 #include "AliITStrackMI.h"
  46 #include "AliTRDtrack.h"
  47 #include "AliHelix.h"
  48 #include "AliESDVertex.h"
  49 #include "AliExternalTrackParam.h"
  50 #include "AliESDkink.h"
  51 #include "AliESDv0.h"
  52 #include "AliV0.h"
  53 //
  54 #include "AliTreeDraw.h"
  55 #include "AliMCInfo.h"
  56 #include "AliGenKinkInfo.h"
  57 #include "AliGenV0Info.h"
  58
  59
  60 #include "AliESDRecV0Info.h"
  61
  62
  63
  64 ClassImp(AliESDRecV0Info)
  65
  66
  67 AliESDRecV0Info:: AliESDRecV0Info():
  68   TObject(),
  69   fT1(),      //track1
  70   fT2(),      //track2
  71   fDist1(0),    //info about closest distance according closest MC - linear DCA
  72   fDist2(0),    //info about closest distance parabolic DCA
  73   fInvMass(0),  //reconstructed invariant mass -
  74   //
  75   fDistMinR(0), // distance at minimal radius
  76   fRr(0),       // rec position of the vertex
  77   fPointAngleFi(0), //point angle fi
  78   fPointAngleTh(0), //point angle theta
  79   fPointAngle(0),   //point angle full
  80   fV0Status(0),       // status of the kink
  81   fV0tpc(0),           // Vo information from reconsturction according TPC
  82   fV0its(0),           // Vo information from reconsturction according ITS
  83   fV0rec(0),           // V0 information form the reconstruction
  84   fMultiple(0),     // how man times V0 was recostructed
  85   fV0Multiple(0),   // how man times was V0 reconstucted
  86   fRecStatus(0)    // status form the reconstuction
  87 {
  88   //
  89   // default constructor
  90   //
  91 }
  92
  93
  94 void  AliESDRecV0Info::Update(Float_t vertex[3])
  95 {
  96
  97   if ( (fT1.fStatus[1]>0)&& (fT2.fStatus[1]>0)){
  98     Float_t distance1,distance2;
  99     Double_t xx[3],pp[3];
 100     //
 101     Double_t xd[3],pd[3],signd;
 102     Double_t xm[3],pm[3],signm;
 103     //
 104     //
 105     if (fT1.fITSOn&&fT2.fITSOn){
 106       for (Int_t i=0;i<3;i++){
 107         xd[i] = fT2.fITSinR1[i];
 108         pd[i] = fT2.fITSinP1[i];
 109         xm[i] = fT1.fITSinR1[i];
 110         pm[i] = fT1.fITSinP1[i];
 111       }
 112     }
 113     else{
 114
 115       for (Int_t i=0;i<3;i++){
 116         xd[i] = fT2.fTPCinR1[i];
 117         pd[i] = fT2.fTPCinP1[i];
 118         xm[i] = fT1.fTPCinR1[i];
 119         pm[i] = fT1.fTPCinP1[i];
 120       }
 121     }
 122     //
 123     //
 124     signd =  fT2.fSign<0 ? -1:1;
 125     signm =  fT1.fSign<0 ? -1:1;
 126
 127     AliHelix dhelix1(xd,pd,signd);
 128     dhelix1.GetMomentum(0,pp,0);
 129     dhelix1.Evaluate(0,xx);
 130     //
 131     //  Double_t x2[3],p2[3];
 132     //
 133     AliHelix mhelix(xm,pm,signm);
 134     //
 135     //find intersection linear
 136     //
 137     Double_t phase[2][2],radius[2];
 138     Int_t  points = dhelix1.GetRPHIintersections(mhelix, phase, radius,200);
 139     Double_t delta1=10000,delta2=10000;
 140
 141     if (points==1){
 142       fRs[0] = TMath::Sqrt(radius[0]);
 143       fRs[1] = TMath::Sqrt(radius[0]);
 144     }
 145     if (points==2){
 146       fRs[0] =TMath::Min(TMath::Sqrt(radius[0]),TMath::Sqrt(radius[1]));
 147       fRs[1] =TMath::Max(TMath::Sqrt(radius[0]),TMath::Sqrt(radius[1]));
 148     }
 149
 150     if (points>0){
 151       dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
 152       dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
 153       dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
 154     }
 155     if (points==2){
 156       dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
 157       dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
 158       dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
 159     }
 160     if (points==1){
 161       fRs[0] = TMath::Sqrt(radius[0]);
 162       fRs[1] = TMath::Sqrt(radius[0]);
 163       fDistMinR = delta1;
 164     }
 165     if (points==2){
 166       if (radius[0]<radius[1]){
 167         fRs[0] = TMath::Sqrt(radius[0]);
 168         fRs[1] = TMath::Sqrt(radius[1]);
 169         fDistMinR = delta1;
 170       }
 171       else{
 172         fRs[0] = TMath::Sqrt(radius[1]);
 173         fRs[1] = TMath::Sqrt(radius[0]);
 174         fDistMinR = delta2;
 175       }
 176     }
 177     //
 178     //
 179     distance1 = TMath::Min(delta1,delta2);
 180     //
 181     //find intersection parabolic
 182     //
 183     points = dhelix1.GetRPHIintersections(mhelix, phase, radius);
 184     delta1=10000,delta2=10000;
 185
 186     if (points>0){
 187       dhelix1.ParabolicDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
 188     }
 189     if (points==2){
 190       dhelix1.ParabolicDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
 191     }
 192
 193     distance2 = TMath::Min(delta1,delta2);
 194     if (distance2>100) fDist2 =100;
 195     return;
 196     if (delta1<delta2){
 197       //get V0 info
 198       dhelix1.Evaluate(phase[0][0],fXr);
 199       dhelix1.GetMomentum(phase[0][0],fPdr);
 200       mhelix.GetMomentum(phase[0][1],fPm);
 201       dhelix1.GetAngle(phase[0][0],mhelix,phase[0][1],fAngle);
 202       fRr = TMath::Sqrt(radius[0]);
 203     }
 204     else{
 205       dhelix1.Evaluate(phase[1][0],fXr);
 206       dhelix1.GetMomentum(phase[1][0], fPdr);
 207       mhelix.GetMomentum(phase[1][1], fPm);
 208       dhelix1.GetAngle(phase[1][0],mhelix,phase[1][1],fAngle);
 209       fRr = TMath::Sqrt(radius[1]);
 210     }
 211     fDist1 = TMath::Sqrt(distance1);
 212     fDist2 = TMath::Sqrt(distance2);
 213
 214     if (fDist2<10.5){
 215       Double_t x,alpha,param[5],cov[15];
 216       //
 217       fT1.GetESDtrack()->GetInnerExternalParameters(alpha,x,param);
 218       fT1.GetESDtrack()->GetInnerExternalCovariance(cov);
 219       AliExternalTrackParam paramm(x,alpha,param,cov);
 220       //
 221       fT2.GetESDtrack()->GetInnerExternalParameters(alpha,x,param);
 222       fT2.GetESDtrack()->GetInnerExternalCovariance(cov);
 223       AliExternalTrackParam paramd(x,alpha,param,cov);
 224     }
 225     //
 226     //
 227
 228     Float_t v[3] = {fXr[0]-vertex[0],fXr[1]-vertex[1],fXr[2]-vertex[2]};
 229     Float_t p[3] = {fPdr[0]+fPm[0], fPdr[1]+fPm[1],fPdr[2]+fPm[2]};
 230
 231     Float_t vnorm2 = v[0]*v[0]+v[1]*v[1];
 232     Float_t vnorm3 = TMath::Sqrt(v[2]*v[2]+vnorm2);
 233     vnorm2 = TMath::Sqrt(vnorm2);
 234     Float_t pnorm2 = p[0]*p[0]+p[1]*p[1];
 235     Float_t pnorm3 = TMath::Sqrt(p[2]*p[2]+pnorm2);
 236     pnorm2 = TMath::Sqrt(pnorm2);
 237
 238     fPointAngleFi = (v[0]*p[0]+v[1]*p[1])/(vnorm2*pnorm2);
 239     fPointAngleTh = (v[2]*p[2]+vnorm2*pnorm2)/(vnorm3*pnorm3);
 240     fPointAngle   = (v[0]*p[0]+v[1]*p[1]+v[2]*p[2])/(vnorm3*pnorm3);
 241   }
 242 }
 243