PWG2/KINK/AliResonanceKinkLikeSign.cxx

   1 /**************************************************************************
   2  * Author: Paraskevi Ganoti, University of Athens (pganoti@phys.uoa.gr)   *
   3  * Contributors are mentioned in the code where appropriate.              *
   4  *                                                                        *
   5  * Permission to use, copy, modify and distribute this software and its   *
   6  * documentation strictly for non-commercial purposes is hereby granted   *
   7  * without fee, provided that the above copyright notice appears in all   *
   8  * copies and that both the copyright notice and this permission notice   *
   9  * appear in the supporting documentation. The authors make no claims     *
  10  * about the suitability of this software for any purpose. It is          *
  11  * provided "as is" without express or implied warranty.                  *
  12  **************************************************************************/
  13
  14 //----------------------------------------------------------------------------------------------------------------
  15 //                        class AliResonanceKinkLikeSign
  16 //        Example of an analysis task for producing a like-sign background for resonances having at least one
  17 //        kaon-kink in their decay products.
  18 //        Background is calculated from a positive kaon kink and a negative track.
  19 //-----------------------------------------------------------------------------------------------------------------
  20
  21 #include "TChain.h"
  22 #include "TTree.h"
  23 #include "TVector3.h"
  24 #include "TF1.h"
  25 #include "TH1D.h"
  26 #include "TH2D.h"
  27 #include <TDatabasePDG.h>
  28 #include <TParticlePDG.h>
  29 #include "TLorentzVector.h"
  30 #include "AliAnalysisTaskSE.h"
  31 #include "AliAnalysisManager.h"
  32
  33 #include "AliESDInputHandler.h"
  34
  35 #include "AliResonanceKinkLikeSign.h"
  36 #include "AliESDkink.h"
  37
  38 ClassImp(AliResonanceKinkLikeSign)
  39
  40 //________________________________________________________________________
  41 AliResonanceKinkLikeSign::AliResonanceKinkLikeSign()
  42   : AliAnalysisTaskSE(), fESD(0), fListOfHistos(0), f1(0), f2(0), fPosKaonLikeSign(0), fLikeSignInvmassPt(0)
  43 {
  44   // Constructor
  45 }
  46
  47 //________________________________________________________________________
  48 AliResonanceKinkLikeSign::AliResonanceKinkLikeSign(const char *name)
  49   : AliAnalysisTaskSE(name), fESD(0), fListOfHistos(0), f1(0), f2(0), fPosKaonLikeSign(0), fLikeSignInvmassPt(0)
  50
  51 {
  52   // Constructor
  53
  54   // Define input and output slots here
  55   // Input slot #0 works with a TChain
  56   DefineInput(0, TChain::Class());
  57   DefineOutput(1, TList::Class());
  58 }
  59
  60 //________________________________________________________________________
  61 void AliResonanceKinkLikeSign::ConnectInputData(Option_t *)
  62 {
  63   // Connect ESD or AOD here
  64   // Called once
  65
  66   TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
  67   if (!tree) {
  68     Printf("ERROR: Could not read chain from input slot 0");
  69   } else {
  70
  71     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
  72
  73     if (!esdH) {
  74       Printf("ERROR: Could not get ESDInputHandler");
  75     } else
  76       fESD = esdH->GetEvent();
  77   }
  78 }
  79
  80 //________________________________________________________________________
  81 void AliResonanceKinkLikeSign::CreateOutputObjects()
  82 {
  83   // Create histograms
  84   // Called once
  85
  86    f1=new TF1("f1","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
  87    f1->SetParameter(0,0.493677);
  88    f1->SetParameter(1,0.9127037);
  89    f1->SetParameter(2,TMath::Pi());
  90
  91    f2=new TF1("f2","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
  92    f2->SetParameter(0,0.13957018);
  93    f2->SetParameter(1,0.2731374);
  94    f2->SetParameter(2,TMath::Pi());
  95
  96    //OpenFile(1);  //uncomment for proof analysis
  97
  98     // for K*0(892)
  99    fPosKaonLikeSign=new TH1D("fPosKaonLikeSign"," ", 60,0.6,1.2);
 100    fLikeSignInvmassPt=new TH2D("fLikeSignInvmassPt"," ", 60,0.6,1.2, 100,0.0,10.0);
 101
 102    // for phi(1020)
 103  //  fPosKaonLikeSign=new TH1D("fPosKaonLikeSign"," ", 70,0.99,1.088);
 104
 105
 106    fListOfHistos=new TList();
 107    fListOfHistos->Add(fPosKaonLikeSign);
 108    fListOfHistos->Add(fLikeSignInvmassPt);
 109
 110 }
 111
 112 //________________________________________________________________________
 113 void AliResonanceKinkLikeSign::Exec(Option_t *)
 114 {
 115   // Main loop
 116   // Called for each event
 117
 118    Double_t daughter1Mass=TDatabasePDG::Instance()->GetParticle(kKPlus)->Mass();
 119    Double_t daughter2Mass=TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass();
 120
 121   if (!fESD) {
 122     Printf("ERROR: fESD not available");
 123     return;
 124   }
 125
 126   const AliESDVertex* vertex = GetEventVertex(fESD);
 127   if(!vertex) return;
 128
 129   Double_t ptrackpos[3], ptrackneg[3];
 130
 131   TLorentzVector p4pos;
 132   TLorentzVector p4neg;
 133   TLorentzVector p4comb;
 134
 135   for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
 136     AliESDtrack* trackpos = fESD->GetTrack(iTracks);
 137     if (!trackpos) {
 138       if (fDebug > 0) Printf("ERROR: Could not receive track %d", iTracks);
 139       continue;
 140     }
 141
 142     if (trackpos->GetSign() < 0) continue;
 143
 144     trackpos->GetPxPyPz(ptrackpos);
 145
 146     Float_t nSigmaToVertex = GetSigmaToVertex(trackpos);
 147
 148     Float_t bpos[2];
 149     Float_t bCovpos[3];
 150     trackpos->GetImpactParameters(bpos,bCovpos);
 151
 152     if (bCovpos[0]<=0 || bCovpos[2]<=0) {
 153      if (fDebug > 0) Printf("Estimated b resolution lower or equal zero!");
 154      bCovpos[0]=0; bCovpos[2]=0;
 155     }
 156
 157     Float_t dcaToVertexXYpos = bpos[0];
 158     Float_t dcaToVertexZpos = bpos[1];
 159
 160     if(nSigmaToVertex>=4) continue;
 161     if((dcaToVertexXYpos>3.0)||(dcaToVertexZpos>3.0)) continue;
 162
 163     TVector3 posTrackMom(ptrackpos[0],ptrackpos[1],ptrackpos[2]);
 164
 165     if(posTrackMom.Perp()<=0.25) continue;
 166
 167     //Uncomment the following block if the Like Sign is made of K+ kink + positive track
 168
 169     Int_t IndexKink=trackpos->GetKinkIndex(0);
 170     Int_t kaonKinkFlag=0;
 171     if(IndexKink<0){
 172
 173         AliESDkink *kink=fESD->GetKink(TMath::Abs(IndexKink)-1);
 174         const TVector3 motherMfromKink(kink->GetMotherP());
 175         const TVector3 daughterMKink(kink->GetDaughterP());
 176         Float_t qt=kink->GetQt();
 177
 178         Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
 179         Double_t maxDecAngpimu=f2->Eval(motherMfromKink.Mag(),0.,0.,0.);
 180
 181         Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
 182
 183         Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
 184         Float_t p1XM= motherMfromKink.Px();
 185         Float_t p1YM= motherMfromKink.Py();
 186         Float_t p1ZM= motherMfromKink.Pz();
 187         Float_t p2XM= daughterMKink.Px();
 188         Float_t p2YM= daughterMKink.Py();
 189         Float_t p2ZM= daughterMKink.Pz();
 190         Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
 191         Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
 192
 193         if((kinkAngle>maxDecAngpimu)&&(qt>0.05)&&(qt<0.25)&&((kink->GetR()>110.)&&(kink->GetR()<230.))&&(TMath::Abs(posTrackMom.Eta())<1.1)&&(invariantMassKmu<0.6)) {
 194
 195           if(posTrackMom.Mag()<=1.1) {
 196            kaonKinkFlag=1;
 197           }
 198           else
 199           if (kinkAngle<maxDecAngKmu) {
 200            kaonKinkFlag=1;
 201           }
 202         }
 203
 204     }  //End Kink Information
 205
 206     if(kaonKinkFlag==0) continue;
 207     if(kaonKinkFlag==1) p4pos.SetVectM(posTrackMom,daughter1Mass);
 208
 209       // Comment the following statements till the "for" if the Like Sign of K+ kink + positive track is needed
 210
 211 //       UInt_t status=trackpos->GetStatus();
 212 //       if((status&AliESDtrack::kTPCrefit)==0) continue;
 213 //       if(trackpos->GetTPCclusters(0)<50) continue;
 214 //       if((trackpos->GetTPCchi2()/trackpos->GetTPCclusters(0))>3.5) continue;
 215 //       Double_t extCovPos[15];
 216 //       trackpos->GetExternalCovariance(extCovPos);
 217 //       if(extCovPos[0]>2) continue;
 218 //       if(extCovPos[2]>2) continue;
 219 //       if(extCovPos[5]>0.5) continue;
 220 //       if(extCovPos[9]>0.5) continue;
 221 //       if(extCovPos[14]>2) continue;
 222 //
 223 //        p4pos.SetVectM(posTrackMom,daughter1Mass);
 224
 225       for (Int_t j=0; j<fESD->GetNumberOfTracks(); j++) {
 226         if(iTracks==j) continue;
 227         AliESDtrack* trackneg=fESD->GetTrack(j);
 228         if (trackneg->GetSign() < 0) continue;
 229
 230         trackneg->GetPxPyPz(ptrackneg);
 231         Float_t negSigmaToVertex = GetSigmaToVertex(trackneg);
 232
 233         Float_t bneg[2];
 234         Float_t bCovneg[3];
 235         trackneg->GetImpactParameters(bneg,bCovneg);
 236         if (bCovneg[0]<=0 || bCovneg[2]<=0) {
 237           if (fDebug > 0) Printf("Estimated b resolution lower or equal zero!");
 238           bCovneg[0]=0; bCovneg[2]=0;
 239         }
 240
 241         Float_t dcaToVertexXYneg = bneg[0];
 242         Float_t dcaToVertexZneg = bneg[1];
 243
 244         if(negSigmaToVertex>=4) continue;
 245         if((dcaToVertexXYneg>3.0)||(dcaToVertexZneg>3.0)) continue;
 246
 247         TVector3 negTrackMom(ptrackneg[0],ptrackneg[1],ptrackneg[2]);
 248
 249         if(negTrackMom.Perp()<=0.25) continue;
 250
 251         UInt_t statusneg=trackneg->GetStatus();
 252
 253         if((statusneg&AliESDtrack::kTPCrefit)==0) continue;
 254
 255         if(trackneg->GetTPCclusters(0)<50) continue;
 256         if((trackneg->GetTPCchi2()/trackneg->GetTPCclusters(0))>3.5) continue;
 257         Double_t extCovneg[15];
 258         trackneg->GetExternalCovariance(extCovneg);
 259         if(extCovneg[0]>2) continue;
 260         if(extCovneg[2]>2) continue;
 261         if(extCovneg[5]>0.5) continue;
 262         if(extCovneg[9]>0.5) continue;
 263         if(extCovneg[14]>2) continue;
 264
 265         p4neg.SetVectM(negTrackMom, daughter2Mass);
 266
 267         p4comb=p4pos;
 268         p4comb+=p4neg;
 269         fPosKaonLikeSign->Fill(p4comb.M());
 270         fLikeSignInvmassPt->Fill(p4comb.M(), p4comb.Vect().Pt());
 271
 272       } //inner track loop
 273
 274   } //outer track loop
 275
 276     PostData(1, fListOfHistos);
 277 }
 278
 279 //________________________________________________________________________
 280 void AliResonanceKinkLikeSign::Terminate(Option_t *)
 281 {
 282
 283 }
 284
 285 //____________________________________________________________________//
 286
 287 Float_t AliResonanceKinkLikeSign::GetSigmaToVertex(AliESDtrack* esdTrack) const {
 288   // Calculates the number of sigma to the vertex.
 289
 290   Float_t b[2];
 291   Float_t bRes[2];
 292   Float_t bCov[3];
 293
 294     esdTrack->GetImpactParameters(b,bCov);
 295
 296   if (bCov[0]<=0 || bCov[2]<=0) {
 297     //AliDebug(1, "Estimated b resolution lower or equal zero!");
 298     bCov[0]=0; bCov[2]=0;
 299   }
 300   bRes[0] = TMath::Sqrt(bCov[0]);
 301   bRes[1] = TMath::Sqrt(bCov[2]);
 302
 303   if (bRes[0] == 0 || bRes[1] ==0) return -1;
 304
 305   Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
 306
 307   if (TMath::Exp(-d * d / 2) < 1e-10) return 1000;
 308
 309   d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
 310
 311   return d;
 312 }
 313
 314 //____________________________________________________________________//
 315
 316 const AliESDVertex* AliResonanceKinkLikeSign::GetEventVertex(const AliESDEvent* esd) const
 317
 318 {
 319   // Get the vertex
 320
 321   const AliESDVertex* vertex = esd->GetPrimaryVertex();
 322
 323   if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
 324   else
 325   {
 326      vertex = esd->GetPrimaryVertexSPD();
 327      if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
 328      else
 329      return 0;
 330   }
 331 }
 332