From 4e539496ce2217c3700982049832b29b9ba911db Mon Sep 17 00:00:00 2001 From: hristov Date: Fri, 20 Apr 2007 08:50:20 +0000 Subject: [PATCH] D0 analysis moved to PWG3 (Andrea) --- ANALYSIS/AliD0toKpi.cxx | 535 -------------------------- ANALYSIS/AliD0toKpi.h | 195 ---------- ANALYSIS/AliD0toKpiAnalysis.cxx | 657 -------------------------------- ANALYSIS/AliD0toKpiAnalysis.h | 90 ----- ANALYSIS/AliD0toKpiPlots.C | 234 ------------ ANALYSIS/AliD0toKpiTest.C | 89 ----- ANALYSIS/AnalysisOldLinkDef.h | 2 - ANALYSIS/libAnalysisOld.pkg | 2 +- 8 files changed, 1 insertion(+), 1803 deletions(-) delete mode 100644 ANALYSIS/AliD0toKpi.cxx delete mode 100644 ANALYSIS/AliD0toKpi.h delete mode 100644 ANALYSIS/AliD0toKpiAnalysis.cxx delete mode 100644 ANALYSIS/AliD0toKpiAnalysis.h delete mode 100644 ANALYSIS/AliD0toKpiPlots.C delete mode 100644 ANALYSIS/AliD0toKpiTest.C diff --git a/ANALYSIS/AliD0toKpi.cxx b/ANALYSIS/AliD0toKpi.cxx deleted file mode 100644 index 31a4c94d156..00000000000 --- a/ANALYSIS/AliD0toKpi.cxx +++ /dev/null @@ -1,535 +0,0 @@ -/************************************************************************** - * 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. * - **************************************************************************/ - -//---------------------------------------------------------------------------- -// Implementation of the D0toKpi class -// for pp and PbPb interactions -// Note: the two decay tracks are labelled: 0 (positive track) -// 1 (negative track) -// Origin: A. Dainese andrea.dainese@lnl.infn.it -//---------------------------------------------------------------------------- - -#include -#include -#include -#include -#include - -#include "AliD0toKpi.h" - -ClassImp(AliD0toKpi) - -//---------------------------------------------------------------------------- -AliD0toKpi::AliD0toKpi() { - // Default constructor - - fSignal = kFALSE; - - fEvent = 0; - - fTrkNum[0] = 0; - fTrkNum[1] = 0; - - fV1x = 0.; - fV1y = 0.; - fV1z = 0.; - fV2x = 0.; - fV2y = 0.; - fV2z = 0.; - fDCA = 0.; - - fPx[0] = 0.; - fPy[0] = 0.; - fPz[0] = 0.; - fPx[1] = 0.; - fPy[1] = 0.; - fPz[1] = 0.; - - fd0[0] = 0.; - fd0[1] = 0.; - - fPdg[0] = 0; - fPdg[1] = 0; - fMum[0] = 0; - fMum[1] = 0; - - fTagPi[0] = 0.; - fTagPi[1] = 0.; - fTagKa[0] = 0.; - fTagKa[1] = 0.; - fTagNid[0] = 0.; - fTagNid[1] = 0.; - - fWgtAD0=fWgtAD0bar=fWgtBD0=fWgtBD0bar=fWgtCD0=fWgtCD0bar=fWgtDD0=fWgtDD0bar=0; - -} -//---------------------------------------------------------------------------- -AliD0toKpi::AliD0toKpi(Int_t ev,Int_t trkNum[2], - Double_t v1[3],Double_t v2[3], - Double_t dca, - Double_t mom[6],Double_t d0[2]) { - // Constructor - - fSignal = kFALSE; - - fEvent = ev; - fTrkNum[0] = trkNum[0]; - fTrkNum[1] = trkNum[1]; - - fV1x = v1[0]; - fV1y = v1[1]; - fV1z = v1[2]; - fV2x = v2[0]; - fV2y = v2[1]; - fV2z = v2[2]; - fDCA = dca; - - fPx[0] = mom[0]; - fPy[0] = mom[1]; - fPz[0] = mom[2]; - fPx[1] = mom[3]; - fPy[1] = mom[4]; - fPz[1] = mom[5]; - - fd0[0] = d0[0]; - fd0[1] = d0[1]; - - fPdg[0] = 0; - fPdg[1] = 0; - fMum[0] = 0; - fMum[1] = 0; - - fTagPi[0] = 0.; - fTagPi[1] = 0.; - fTagKa[0] = 0.; - fTagKa[1] = 0.; - fTagNid[0] = 0.; - fTagNid[1] = 0.; - - fWgtAD0=fWgtAD0bar=fWgtBD0=fWgtBD0bar=fWgtCD0=fWgtCD0bar=fWgtDD0=fWgtDD0bar=0; -} -//---------------------------------------------------------------------------- -AliD0toKpi::~AliD0toKpi() {} -//____________________________________________________________________________ -AliD0toKpi::AliD0toKpi( const AliD0toKpi& d0toKpi):TObject(d0toKpi) { - // dummy copy constructor -} -//---------------------------------------------------------------------------- -void AliD0toKpi::ApplyPID(TString pidScheme) { - // Applies particle identification - - if((!pidScheme.CompareTo("TOFparamPbPb") || !pidScheme.CompareTo("TOFparamPP")) && fPdg[0]==0) { - printf("AliD0toKpi::ApplyPID :\n Warning: TOF parameterized PID can be used only for simulation!\n"); - return; - } - - if(!pidScheme.CompareTo("TOFparamPbPb")) { - // tagging of the positive track - if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13 - || TMath::Abs(fPdg[0])==11) { // pion,muon,electron - fTagPi[0] = LinearInterpolation(PChild(0),kPiBinsPbPb,kPiBinWidthPbPb,kPiTagPiPbPb); - fTagNid[0] = 1.-fTagPi[0]; - fTagKa[0] = 0.; - fTagPr[0] = 0.; - } - if(TMath::Abs(fPdg[0])==321) { // kaon - fTagKa[0] = LinearInterpolation(PChild(0),kKBinsPbPb,kKBinWidthPbPb,kKTagKPbPb); - fTagNid[0] = LinearInterpolation(PChild(0),kKBinsPbPb,kKBinWidthPbPb,kKTagNidPbPb); - if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0]; - fTagPi[0] = 1.-fTagNid[0]-fTagKa[0]; - fTagPr[0] = 0.; - } - if(TMath::Abs(fPdg[0])==2212) { // proton - fTagPr[0] = LinearInterpolation(PChild(0),kPBinsPbPb,kPBinWidthPbPb,kPTagPPbPb); - fTagNid[0] = LinearInterpolation(PChild(0),kPBinsPbPb,kPBinWidthPbPb,kPTagNidPbPb); - if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0]; - fTagPi[0] = 1.-fTagNid[0]-fTagPr[0]; - fTagKa[0] = 0.; - } - // tagging of the negative track - if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13 - || TMath::Abs(fPdg[1])==11) { // pion,muon,electron - fTagPi[1] = LinearInterpolation(PChild(1),kPiBinsPbPb,kPiBinWidthPbPb,kPiTagPiPbPb); - fTagNid[1] = 1.-fTagPi[1]; - fTagKa[1] = 0.; - fTagPr[1] = 0.; - } - if(TMath::Abs(fPdg[1])==321) { // kaon - fTagKa[1] = LinearInterpolation(PChild(1),kKBinsPbPb,kKBinWidthPbPb,kKTagKPbPb); - fTagNid[1] = LinearInterpolation(PChild(1),kKBinsPbPb,kKBinWidthPbPb,kKTagNidPbPb); - if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1]; - fTagPi[1] = 1.-fTagNid[1]-fTagKa[1]; - fTagPr[1] = 0.; - } - if(TMath::Abs(fPdg[1])==2212) { // proton - fTagPr[1] = LinearInterpolation(PChild(1),kPBinsPbPb,kPBinWidthPbPb,kPTagPPbPb); - fTagNid[1] = LinearInterpolation(PChild(1),kPBinsPbPb,kPBinWidthPbPb,kPTagNidPbPb); - if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1]; - fTagPi[1] = 1.-fTagNid[1]-fTagPr[1]; - fTagKa[1] = 0.; - } - } - - - if(!pidScheme.CompareTo("TOFparamPP")) { - // tagging of the positive track - if(TMath::Abs(fPdg[0])==211 || TMath::Abs(fPdg[0])==13 - || TMath::Abs(fPdg[0])==11) { // pion,muon,electron - fTagPi[0] = LinearInterpolation(PChild(0),kPiBinsPP,kPiBinWidthPP,kPiTagPiPP); - fTagNid[0] = 1.-fTagPi[0]; - fTagKa[0] = 0.; - fTagPr[0] = 0.; - } - if(TMath::Abs(fPdg[0])==321) { // kaon - fTagKa[0] = LinearInterpolation(PChild(0),kKBinsPP,kKBinWidthPP,kKTagKPP); - fTagNid[0] = LinearInterpolation(PChild(0),kKBinsPP,kKBinWidthPP,kKTagNidPP); - if((fTagNid[0]+fTagKa[0])>1.) fTagNid[0] = 1.-fTagKa[0]; - fTagPi[0] = 1.-fTagNid[0]-fTagKa[0]; - fTagPr[0] = 0.; - } - if(TMath::Abs(fPdg[0])==2212) { // proton - fTagPr[0] = LinearInterpolation(PChild(0),kPBinsPP,kPBinWidthPP,kPTagPPP); - fTagNid[0] = LinearInterpolation(PChild(0),kPBinsPP,kPBinWidthPP,kPTagNidPP); - if((fTagNid[0]+fTagPr[0])>1.) fTagNid[0] = 1.-fTagPr[0]; - fTagPi[0] = 1.-fTagNid[0]-fTagPr[0]; - fTagKa[0] = 0.; - } - // tagging of the negative track - if(TMath::Abs(fPdg[1])==211 || TMath::Abs(fPdg[1])==13 - || TMath::Abs(fPdg[1])==11) { // pion,muon,electron - fTagPi[1] = LinearInterpolation(PChild(1),kPiBinsPP,kPiBinWidthPP,kPiTagPiPP); - fTagNid[1] = 1.-fTagPi[1]; - fTagKa[1] = 0.; - fTagPr[1] = 0.; - } - if(TMath::Abs(fPdg[1])==321) { // kaon - fTagKa[1] = LinearInterpolation(PChild(1),kKBinsPP,kKBinWidthPP,kKTagKPP); - fTagNid[1] = LinearInterpolation(PChild(1),kKBinsPP,kKBinWidthPP,kKTagNidPP); - if((fTagNid[1]+fTagKa[1])>1.) fTagNid[1] = 1.-fTagKa[1]; - fTagPi[1] = 1.-fTagNid[1]-fTagKa[1]; - fTagPr[1] = 0.; - } - if(TMath::Abs(fPdg[1])==2212) { // proton - fTagPr[1] = LinearInterpolation(PChild(1),kPBinsPP,kPBinWidthPP,kPTagPPP); - fTagNid[1] = LinearInterpolation(PChild(1),kPBinsPP,kPBinWidthPP,kPTagNidPP); - if((fTagNid[1]+fTagPr[1])>1.) fTagNid[1] = 1.-fTagPr[1]; - fTagPi[1] = 1.-fTagNid[1]-fTagPr[1]; - fTagKa[1] = 0.; - } - } - - return; -} -//---------------------------------------------------------------------------- -Double_t AliD0toKpi::ChildrenRelAngle() const { - // relative angle between K and pi - - TVector3 mom0(fPx[0],fPy[0],fPz[0]); - TVector3 mom1(fPx[1],fPy[1],fPz[1]); - - Double_t angle = mom0.Angle(mom1); - - return angle; -} -//---------------------------------------------------------------------------- -void AliD0toKpi::ComputeWgts() { - // calculate the weights for PID - - - // assignement of the weights from PID - fWgtAD0 = fTagKa[1]*(fTagPi[0]+fTagNid[0]); - fWgtAD0bar = fTagKa[0]*(fTagPi[1]+fTagNid[1]); - fWgtBD0 = fTagPi[0]*fTagNid[1]; - fWgtBD0bar = fTagPi[1]*fTagNid[0]; - fWgtCD0 = fTagNid[0]*fTagNid[1]; - fWgtCD0bar = fTagNid[0]*fTagNid[1]; - fWgtDD0 = 1.-fWgtAD0-fWgtBD0-fWgtCD0; - fWgtDD0bar = 1.-fWgtAD0bar-fWgtBD0bar-fWgtCD0bar; - - /* - cerr< 1.) { cerr<<"AliD0toKpi::CosThetaStar: > 1 "< 1.) { cerr<<"AliD0toKpi::CosThetaStar: > 1 "<0. ? absDD : -absDD); -} -//---------------------------------------------------------------------------- -void AliD0toKpi::InvMass(Double_t &mD0,Double_t &mD0bar) const { - // invariant mass as D0 and as D0bar - - Double_t energy[2]; - - // D0 -> K- Pi+ - energy[1] = TMath::Sqrt(kMK*kMK+PChild(1)*PChild(1)); - energy[0] = TMath::Sqrt(kMPi*kMPi+PChild(0)*PChild(0)); - - mD0 = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P()); - - - // D0bar -> K+ Pi- - energy[0] = TMath::Sqrt(kMK*kMK+PChild(0)*PChild(0)); - energy[1] = TMath::Sqrt(kMPi*kMPi+PChild(1)*PChild(1)); - - mD0bar = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-P()*P()); - - return; - -} -//---------------------------------------------------------------------------- -Double_t AliD0toKpi::Ql(Int_t child) const { - // longitudinal momentum of decay tracks w.r.t. to D0 momentum - - Double_t qL; - TVector3 mom(fPx[child],fPy[child],fPz[child]); - TVector3 momD(Px(),Py(),Pz()); - - qL = mom.Dot(momD)/momD.Mag(); - - return qL ; -} -//---------------------------------------------------------------------------- -Double_t AliD0toKpi::Qt() const { - // transverse momentum of decay tracks w.r.t. to D0 momentum - - TVector3 mom0(fPx[0],fPy[0],fPz[0]); - TVector3 momD(Px(),Py(),Pz()); - - return mom0.Perp(momD); -} -//---------------------------------------------------------------------------- -Bool_t AliD0toKpi::Select(const Double_t* cuts,Int_t& okD0,Int_t& okD0bar) - const { -// -// This function compares the D0 with a set of cuts: -// -// cuts[0] = inv. mass half width [GeV] -// cuts[1] = dca [micron] -// cuts[2] = cosThetaStar -// cuts[3] = pTK [GeV/c] -// cuts[4] = pTPi [GeV/c] -// cuts[5] = d0K [micron] upper limit! -// cuts[6] = d0Pi [micron] upper limit! -// cuts[7] = d0d0 [micron^2] -// cuts[8] = cosThetaPoint -// -// If the D0/D0bar doesn't pass the cuts it sets the weights to 0 -// If neither D0 nor D0bar pass the cuts return kFALSE -// - Double_t mD0,mD0bar,ctsD0,ctsD0bar; - okD0=1; okD0bar=1; - - if(PtChild(1) < cuts[3] || PtChild(0) < cuts[4]) okD0 = 0; - if(PtChild(0) < cuts[3] || PtChild(1) < cuts[4]) okD0bar = 0; - if(!okD0 && !okD0bar) return kFALSE; - - if(TMath::Abs(Getd0Child(1)) > cuts[5] || - TMath::Abs(Getd0Child(0)) > cuts[6]) okD0 = 0; - if(TMath::Abs(Getd0Child(0)) > cuts[6] || - TMath::Abs(Getd0Child(1)) > cuts[5]) okD0bar = 0; - if(!okD0 && !okD0bar) return kFALSE; - - if(GetDCA() > cuts[1]) { okD0 = okD0bar = 0; return kFALSE; } - - InvMass(mD0,mD0bar); - if(TMath::Abs(mD0-kMD0) > cuts[0]) okD0 = 0; - if(TMath::Abs(mD0bar-kMD0) > cuts[0]) okD0bar = 0; - if(!okD0 && !okD0bar) return kFALSE; - - CosThetaStar(ctsD0,ctsD0bar); - if(TMath::Abs(ctsD0) > cuts[2]) okD0 = 0; - if(TMath::Abs(ctsD0bar) > cuts[2]) okD0bar = 0; - if(!okD0 && !okD0bar) return kFALSE; - - if(ProdImpParams() > cuts[7]) { okD0 = okD0bar = 0; return kFALSE; } - - if(CosPointing() < cuts[8]) { okD0 = okD0bar = 0; return kFALSE; } - - return kTRUE; -} -//----------------------------------------------------------------------------- -void AliD0toKpi::SetPIDresponse(Double_t resp0[5],Double_t resp1[5]) { - // Set combined PID detector response probabilities - - fPIDrespEl[0] = resp0[0]; - fPIDrespEl[1] = resp1[0]; - fPIDrespMu[0] = resp0[1]; - fPIDrespMu[1] = resp1[1]; - fPIDrespPi[0] = resp0[2]; - fPIDrespPi[1] = resp1[2]; - fPIDrespKa[0] = resp0[3]; - fPIDrespKa[1] = resp1[3]; - fPIDrespPr[0] = resp0[4]; - fPIDrespPr[1] = resp1[4]; - - return; -} -//---------------------------------------------------------------------------- -Double_t AliD0toKpi::LinearInterpolation(Double_t p,Int_t nBins,Double_t Bin, - const Double_t *values) const { - // a linear interpolation method - - Double_t value=0; - Double_t slope; - - if(p<0.5*Bin) { - value = values[0]; - } else if(p>=(nBins-0.5)*Bin) { - slope = (2*values[nBins-1]-values[nBins-2]-values[nBins-3])/Bin/2; - value = values[nBins-2]+slope*(p-Bin*(nBins-1.5)); - } else { - for(Int_t i=0; i1.) value=1.; - - return value; -} -//---------------------------------------------------------------------------- - - - - - - - - - diff --git a/ANALYSIS/AliD0toKpi.h b/ANALYSIS/AliD0toKpi.h deleted file mode 100644 index 716dbb9d4dd..00000000000 --- a/ANALYSIS/AliD0toKpi.h +++ /dev/null @@ -1,195 +0,0 @@ -#ifndef AliD0toKpi_H -#define AliD0toKpi_H -/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * See cxx source for full Copyright notice */ - -//------------------------------------------------------------------------- -// Class AliD0toKpi -// Reconstructed D0 -> K^- pi^+ class -// -// Note: the two decay tracks are labelled: 0 (positive track) -// 1 (negative track) -// -// Origin: A. Dainese andrea.dainese@lnl.infn.it -//------------------------------------------------------------------------- - -#include -#include - -//---------------------------------------------------------------------------- -// Some constants (masses + parameterized TOF PID) -// -// particle masses -const Double_t kMD0 = 1.8645; // D0 mass -const Double_t kMK = 0.49368; // K+ mass -const Double_t kMPi = 0.13957; // pi+ mass - -// --- TOF tagging probabilities --- -// central HIJING -// B = 0.4 T -// tracking errors in TPC included -// With TRD -// -// *** Pb-Pb dNch/dy=6000 *** -// -// PIONS -const Int_t kPiBinsPbPb = 10; -const Double_t kPiBinWidthPbPb = 0.250; -const Double_t kPiTagPiPbPb[kPiBinsPbPb] = {0.211421,0.652184,0.624421,0.614727,0.610777,0.628015,0.631520,0.630324,0.637551,0.575235}; -const Double_t kPiTagNidPbPb[kPiBinsPbPb] = {0.788579,0.347816,0.375579,0.385273,0.389223,0.371985,0.368480,0.369676,0.362449,0.424765}; -// KAONS -const Int_t kKBinsPbPb = 10; -const Double_t kKBinWidthPbPb = 0.250; -const Double_t kKTagKPbPb[kKBinsPbPb] = {0.000000,0.101255,0.397662,0.467586,0.517008,0.555023,0.584185,0.519029,0.464117,0.247308}; -const Double_t kKTagPiPbPb[kKBinsPbPb] = {0.102049,0.289930,0.101930,0.057771,0.040286,0.028567,0.053108,0.094369,0.066302,0.247308}; -const Double_t kKTagNidPbPb[kKBinsPbPb] = {0.897951,0.608815,0.500408,0.474643,0.442705,0.416410,0.362707,0.386603,0.469580,0.505383}; -// PROTONS -const Int_t kPBinsPbPb = 9; -const Double_t kPBinWidthPbPb = 0.500; -const Double_t kPTagPPbPb[kPBinsPbPb] = {0.017940,0.350681,0.535286,0.583264,0.562935,0.560524,0.545992,0.598060,0.351245}; -const Double_t kPTagPiPbPb[kPBinsPbPb] = {0.195955,0.094949,0.039962,0.026039,0.007556,0.016986,0.030333,0.000000,0.000000}; -const Double_t kPTagNidPbPb[kPBinsPbPb] = {0.786105,0.554370,0.424751,0.390697,0.429508,0.422491,0.423675,0.401940,0.648755}; -// -// pp PYTHIA -// -// *** cuts for pp *** -// -// PIONS -const Int_t kPiBinsPP = 10; -const Double_t kPiBinWidthPP = 0.250; -const Double_t kPiTagPiPP[kPiBinsPP] = {0.194528,0.447097,0.603364,0.646413,0.647125,0.669157,0.688139,0.682564,0.689910,0.665710}; -const Double_t kPiTagNidPP[kPiBinsPP] = {0.805472,0.552903,0.396636,0.353587,0.352875,0.330843,0.311861,0.317436,0.310090,0.334290}; -// KAONS -const Int_t kKBinsPP = 10; -const Double_t kKBinWidthPP = 0.250; -const Double_t kKTagKPP[kKBinsPP] = {0.000000,0.173393,0.439690,0.519423,0.587025,0.605372,0.586021,0.650139,0.444444,0.299363}; -const Double_t kKTagPiPP[kKBinsPP] = {0.000000,0.001495,0.000000,-0.000000,-0.000000,0.000000,0.032258,0.060572,0.101449,0.242038}; -const Double_t kKTagNidPP[kKBinsPP] = {1.000000,0.825112,0.560310,0.480577,0.412975,0.394628,0.381720,0.289289,0.454106,0.458599}; -// PROTONS -const Int_t kPBinsPP = 9; -const Double_t kPBinWidthPP = 0.500; -const Double_t kPTagPPP[kPBinsPP] = {0.029404,0.438640,0.613710,0.665152,0.634961,0.657711,0.703704,0.685714,0.235294}; -const Double_t kPTagPiPP[kPBinsPP] = {0.000000,0.000000,0.000000,-0.000000,0.000000,0.000000,-0.000000,0.014286,-0.000000}; -const Double_t kPTagNidPP[kPBinsPP] = {0.970596,0.561360,0.386290,0.334848,0.365039,0.342289,0.296296,0.300000,0.764706}; - - - - -//----------------------------------------------------------------------------- -class AliD0toKpi : public TObject { - public: - // - AliD0toKpi(); - AliD0toKpi(Int_t ev,Int_t trkNum[2], - Double_t v1[3],Double_t v2[3],Double_t dca, - Double_t mom[6],Double_t d0[2]); - virtual ~AliD0toKpi(); - AliD0toKpi(const AliD0toKpi& d0toKpi); - - Double_t Alpha() const { return (Ql(0)-Ql(1))/(Ql(0)+Ql(1)); } - void ApplyPID(TString pidScheme="TOFparamPbPb"); - Double_t ChildrenRelAngle() const; - void ComputeWgts(); - void CorrectWgt4BR(Double_t factor); - Double_t CosPointing() const; - Double_t CosPointingXY() const; - void CosThetaStar(Double_t &ctsD0,Double_t &ctsD0bar) const; - Double_t Ct() const {return Length()*kMD0/P();} - Double_t Energy() const { return TMath::Sqrt(P()*P()+kMD0*kMD0); } - Double_t Eta() const; - Double_t EtaChild(Int_t child) const; - Int_t EventNo() const {return TMath::Abs(fEvent);} - Double_t GetDCA() const { return 10000.*fDCA; } - Int_t GetTrkNum(Int_t child) const { return fTrkNum[child]; } - Double_t Getd0Child(Int_t child) const { return fd0[child]; } - Int_t GetPdgChild(Int_t child) const { return fPdg[child]; } - Int_t GetPdgMum(Int_t child) const {return fMum[child]; } - void GetWgts(Double_t &WgtD0,Double_t &WgtD0bar,TString sample) const; - void GetPrimaryVtx(Double_t vtx[3]) const - { vtx[0]=fV1x; vtx[1]=fV1y; vtx[2]=fV1z; } - void GetSecondaryVtx(Double_t vtx[3]) const - { vtx[0]=fV2x; vtx[1]=fV2y; vtx[2]=fV2z; } - - Double_t ImpPar() const; - void InvMass(Double_t &mD0,Double_t &mD0bar) const; - Bool_t IsSignal() const { if(fSignal) return kTRUE; return kFALSE; } - Double_t Length() const - { return TMath::Sqrt((fV1x-fV2x)*(fV1x-fV2x) - +(fV1y-fV2y)*(fV1y-fV2y)+(fV1z-fV2z)*(fV1z-fV2z)); } - Double_t P() const { return TMath::Sqrt(Pt()*Pt()+Pz()*Pz()); } - Double_t PChild(Int_t child) const { return TMath::Sqrt(fPx[child]*fPx[child]+fPy[child]*fPy[child]+fPz[child]*fPz[child]); } - Double_t ProdImpParams() const { return fd0[0]*fd0[1]; } - Double_t Pt() const { return TMath::Sqrt(Px()*Px()+Py()*Py()); } - Double_t PtChild(Int_t child) const { return TMath::Sqrt(fPx[child]*fPx[child]+fPy[child]*fPy[child]); } - Double_t Px() const { return (fPx[0]+fPx[1]); } - Double_t Py() const { return (fPy[0]+fPy[1]); } - Double_t Pz() const { return (fPz[0]+fPz[1]); } - Double_t Ql(Int_t child) const; - Double_t Qt() const; - Double_t Rapidity() const { return 0.5*TMath::Log((Energy()+Pz())/(Energy()-Pz()+1.e-13)); } - Bool_t Select(const Double_t* cuts,Int_t& okD0,Int_t& okD0bar) const; - void SetPrimaryVtx(Double_t vtx[3]) - { fV1x=vtx[0]; fV1y=vtx[1]; fV1z=vtx[2]; } - void SetSignal() { fSignal = kTRUE; } - void SetTOFmasses(Double_t mass[2]) - { fTOFmass[0]=mass[0]; fTOFmass[1]=mass[1]; } - void SetPIDresponse(Double_t resp0[5],Double_t resp1[5]); - void SetPdgCodes(Int_t pdg[2]) {fPdg[0]=pdg[0];fPdg[1]=pdg[1]; } - void SetMumPdgCodes(Int_t mum[2]) {fMum[0]=mum[0];fMum[1]=mum[1]; } - Double_t LinearInterpolation(Double_t p,Int_t nBins,Double_t Bin, - const Double_t *values) const; - // - private: - // - Bool_t fSignal; // TRUE if signal, FALSE if background (for simulation) - Int_t fEvent; // number of the event this D0 comes from - // -1 if the D0 comes from ev. mixing - - Int_t fTrkNum[2]; // numbers of the two decay tracks - - Double_t fV1x; // X-position of the primary vertex of the event - Double_t fV1y; // Y-position of the primary vertex of the event - Double_t fV1z; // Z-position of the primary vertex of the event - Double_t fV2x; // X-position of the reconstructed secondary vertex - Double_t fV2y; // Y-position of the reconstructed secondary vertex - Double_t fV2z; // Z-position of the reconstructed secondary vertex - Double_t fDCA; // DCA of the two tracks - - Double_t fPx[2]; // X,Y,Z - Double_t fPy[2]; // momenta of the two tracks - Double_t fPz[2]; // at the reconstructed vertex - - Double_t fd0[2]; // impact parameters in the bending plane - - Int_t fPdg[2]; // PDG codes of the two tracks (for sim.) - Int_t fMum[2]; // PDG codes of the mothers (for sim.) - - Double_t fTagPi[2]; // probability to be tagged as pion - Double_t fTagKa[2]; // probability to be tagged as kaon - Double_t fTagPr[2]; // probability to be tagged as proton - Double_t fTagNid[2]; // probability to be tagged as "non-identified" - - Double_t fPIDrespEl[2]; // det. response to be electron - Double_t fPIDrespMu[2]; // det. response to be muon - Double_t fPIDrespPi[2]; // det. response to be pion - Double_t fPIDrespKa[2]; // det. response to be kaon - Double_t fPIDrespPr[2]; // det. response to be proton - Double_t fTOFmass[2]; // mass estimated by the TOF (-1000. if track not reached TOF) - - Double_t fWgtAD0,fWgtAD0bar; // weights for the 3 samples - Double_t fWgtBD0,fWgtBD0bar; // weights for the 3 samples - Double_t fWgtCD0,fWgtCD0bar; // A: (K,Pi)+(K,?) B: (?,Pi) C: (?,?) - Double_t fWgtDD0,fWgtDD0bar; // D: all other pairs - - ClassDef(AliD0toKpi,1) // Reconstructed D0 candidate class -}; - -#endif - - - - - - - - diff --git a/ANALYSIS/AliD0toKpiAnalysis.cxx b/ANALYSIS/AliD0toKpiAnalysis.cxx deleted file mode 100644 index e9294d2dd6c..00000000000 --- a/ANALYSIS/AliD0toKpiAnalysis.cxx +++ /dev/null @@ -1,657 +0,0 @@ -/************************************************************************** - * 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. * - **************************************************************************/ - -//---------------------------------------------------------------------------- -// Implementation of the D0toKpi reconstruction and analysis class -// Note: the two decay tracks are labelled: 0 (positive track) -// 1 (negative track) -// An example of usage can be found in the macro AliD0toKpiTest.C -// Origin: A. Dainese andrea.dainese@lnl.infn.it -//---------------------------------------------------------------------------- -#include -#include -#include -#include -#include -#include -#include "AliESD.h" -#include "AliMC.h" -#include "AliRun.h" -#include "AliRunLoader.h" -#include "AliVertexerTracks.h" -#include "AliESDVertex.h" -#include "AliESDv0.h" -#include "AliD0toKpi.h" -#include "AliD0toKpiAnalysis.h" -#include "AliLog.h" - -typedef struct { - Int_t lab; - Int_t pdg; - Int_t mumlab; - Int_t mumpdg; - Int_t mumprongs; - Float_t Vx,Vy,Vz; - Float_t Px,Py,Pz; -} REFTRACK; - -ClassImp(AliD0toKpiAnalysis) - -//---------------------------------------------------------------------------- -AliD0toKpiAnalysis::AliD0toKpiAnalysis() { - // Default constructor - - SetPtCut(); - Setd0Cut(); - SetMassCut(); - SetD0Cuts(); - SetVertex1(); - SetPID(); - fVertexOnTheFly = kFALSE; - fSim = kFALSE; - fOnlySignal = kFALSE; -} -//---------------------------------------------------------------------------- -AliD0toKpiAnalysis::~AliD0toKpiAnalysis() {} -//---------------------------------------------------------------------------- -void AliD0toKpiAnalysis::ApplySelection(const Char_t *inName,const Char_t *outName) const { - // select candidates that pass fD0Cuts and write them to a new file - - TFile *inFile = TFile::Open(inName); - - TTree *treeD0in=(TTree*)inFile->Get("TreeD0"); - AliD0toKpiAnalysis *inAnalysis = (AliD0toKpiAnalysis*)inFile->Get("D0toKpiAnalysis"); - printf("+++\n+++ I N P U T S T A T U S:\n+++\n"); - inAnalysis->PrintStatus(); - - - AliD0toKpi *d = 0; - treeD0in->SetBranchAddress("D0toKpi",&d); - Int_t entries = (Int_t)treeD0in->GetEntries(); - - printf("+++\n+++ Number of D0 in input tree: %d\n+++\n",entries); - - TTree *treeD0out = new TTree("TreeD0","Tree with selected D0 candidates"); - treeD0out->Branch("D0toKpi","AliD0toKpi",&d,200000,0); - - - Int_t okD0=0,okD0bar=0; - Int_t nSel = 0; - - for(Int_t i=0; iGetEvent(i); - - if(fSim && fOnlySignal && !d->IsSignal()) continue; - - // check if candidate passes selection (as D0 or D0bar) - if(d->Select(fD0Cuts,okD0,okD0bar)) { - nSel++; - treeD0out->Fill(); - } - - } - - AliD0toKpiAnalysis *outAnalysis = (AliD0toKpiAnalysis*)inAnalysis->Clone("D0toKpiAnalysis"); - outAnalysis->SetD0Cuts(fD0Cuts); - printf("------------------------------------------\n"); - printf("+++\n+++ O U T P U T S T A T U S:\n+++\n"); - outAnalysis->PrintStatus(); - - printf("+++\n+++ Number of D0 in output tree: %d\n+++\n",nSel); - - TFile* outFile = new TFile(outName,"recreate"); - treeD0out->Write(); - outAnalysis->Write(); - outFile->Close(); - - return; -} -//---------------------------------------------------------------------------- -Double_t AliD0toKpiAnalysis::CalculateTOFmass(Double_t mom,Double_t length, - Double_t time) const { - // calculated the mass from momentum, track length from vertex to TOF - // and time measured by the TOF - if(length==0.) return -1000.; - Double_t a = time*time/length/length; - if(a > 1.) { - a = TMath::Sqrt(a-1.); - } else { - a = -TMath::Sqrt(1.-a); - } - - return mom*a; -} -//---------------------------------------------------------------------------- -void AliD0toKpiAnalysis::FindCandidates(Int_t evFirst,Int_t evLast, - const Char_t *outName) { - // Find D0 candidates and calculate parameters - - - TString esdName="AliESDs.root"; - if(gSystem->AccessPathName(esdName.Data(),kFileExists)) { - printf("AliD0toKpiAnalysis::FindCandidatesESD(): No ESDs file found!\n"); - return; - } - - TString outName1=outName; - TString outName2="nTotEvents.dat"; - - Int_t nTotEv=0,nD0rec=0,nD0rec1ev=0; - Double_t dca; - Double_t v2[3],mom[6],d0[2]; - Int_t iTrkP,iTrkN,trkEntries; - Int_t nTrksP=0,nTrksN=0; - Int_t trkNum[2]; - Double_t tofmass[2]; - Int_t okD0=0,okD0bar=0; - AliESDtrack *postrack = 0; - AliESDtrack *negtrack = 0; - - // create the AliVertexerTracks object - // (it will be used only if fVertexOnTheFly=kTrue) - AliVertexerTracks *vertexer1 = new AliVertexerTracks; - if(fVertexOnTheFly) { - // open the mean vertex - TFile *invtx = new TFile("AliESDVertexMean.root"); - AliESDVertex *initVertex = (AliESDVertex*)invtx->Get("vtxmean"); - invtx->Close(); - vertexer1->SetVtxStart(initVertex); - delete invtx; - } - Int_t skipped[2]; - Bool_t goodVtx1; - - // create tree for reconstructed D0s - AliD0toKpi *ioD0toKpi=0; - TTree *treeD0 = new TTree("TreeD0","Tree with D0 candidates"); - treeD0->Branch("D0toKpi","AliD0toKpi",&ioD0toKpi,200000,0); - - // open file with tracks - TFile *esdFile = TFile::Open(esdName.Data()); - AliESD* event = new AliESD; - TTree* tree = (TTree*) esdFile->Get("esdTree"); - if(!tree) { - Error("FindCandidatesESD", "no ESD tree found"); - return; - } - tree->SetBranchAddress("ESD",&event); - - // loop on events in file - for(Int_t iEvent = evFirst; iEvent < tree->GetEntries(); iEvent++) { - if(iEvent > evLast) break; - tree->GetEvent(iEvent); - Int_t ev = (Int_t)event->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number. - printf("--- Finding D0 -> Kpi in event %d\n",ev); - // count the total number of events - nTotEv++; - - trkEntries = (Int_t)event->GetNumberOfTracks(); - printf(" Number of tracks: %d\n",trkEntries); - if(trkEntries<2) continue; - - // retrieve primary vertex from file - if(!event->GetPrimaryVertex()) { - printf(" No vertex\n"); - continue; - } - event->GetPrimaryVertex()->GetXYZ(fV1); - - // call function which applies sigle-track selection and - // separetes positives and negatives - TObjArray trksP(trkEntries/2); - Int_t *trkEntryP = new Int_t[trkEntries]; - TObjArray trksN(trkEntries/2); - Int_t *trkEntryN = new Int_t[trkEntries]; - TTree *trkTree = new TTree(); - SelectTracks(event,trksP,trkEntryP,nTrksP, - trksN,trkEntryN,nTrksN); - - printf(" pos. tracks: %d neg .tracks: %d\n",nTrksP,nTrksN); - - - nD0rec1ev = 0; - - // LOOP ON POSITIVE TRACKS - for(iTrkP=0; iTrkPGetMagneticField(); - AliESDtrack nt(*negtrack), pt(*postrack); - dca = nt.PropagateToDCA(&pt,b); - - // define the AliESDv0 object - AliESDv0 vertex2(nt,trkNum[0],pt,trkNum[1]); - - // get position of the secondary vertex - vertex2.GetXYZ(v2[0],v2[1],v2[2]); - vertex2.GetPPxPyPz(mom[0],mom[1],mom[2]); - vertex2.GetNPxPyPz(mom[3],mom[4],mom[5]); - // impact parameters of the tracks w.r.t. the primary vertex - d0[0] = 10000.*pt.GetD(fV1[0],fV1[1],b); - d0[1] = -10000.*nt.GetD(fV1[0],fV1[1],b); - goodVtx1 = kTRUE; - - // no vertexing if DeltaMass > fMassCut - if(fVertexOnTheFly) { - goodVtx1 = kFALSE; - if(SelectInvMass(mom)) { - // primary vertex from *other* tracks in the event - skipped[0] = trkEntryP[iTrkP]; - skipped[1] = trkEntryN[iTrkN]; - vertexer1->SetSkipTracks(2,skipped); - AliESDVertex *vertex1onfly = - (AliESDVertex*)vertexer1->FindPrimaryVertex(event); - if(vertex1onfly->GetNContributors()>0) goodVtx1 = kTRUE; - vertex1onfly->GetXYZ(fV1); - //vertex1onfly->PrintStatus(); - delete vertex1onfly; - } - } - - - // create the object AliD0toKpi - AliD0toKpi theD0(ev,trkNum,fV1,v2,dca,mom,d0); - // select D0s - if(goodVtx1 && theD0.Select(fD0Cuts,okD0,okD0bar)) { - // get PID info from ESD - AliESDtrack *t0 = (AliESDtrack*)event->GetTrack(trkNum[0]); - Double_t esdpid0[5]; - t0->GetESDpid(esdpid0); - if(t0->GetStatus()&AliESDtrack::kTOFpid) { - tofmass[0] = CalculateTOFmass(t0->GetP(), - t0->GetIntegratedLength(), - t0->GetTOFsignal()); - } else { - tofmass[0] = -1000.; - } - AliESDtrack *t1 = (AliESDtrack*)event->GetTrack(trkNum[1]); - Double_t esdpid1[5]; - t1->GetESDpid(esdpid1); - if(t1->GetStatus()&AliESDtrack::kTOFpid) { - tofmass[1] = CalculateTOFmass(t1->GetP(), - t1->GetIntegratedLength(), - t1->GetTOFsignal()); - } else { - tofmass[1] = -1000.; - } - - theD0.SetPIDresponse(esdpid0,esdpid1); - theD0.SetTOFmasses(tofmass); - - // fill the tree - ioD0toKpi=&theD0; - treeD0->Fill(); - - nD0rec++; nD0rec1ev++; - ioD0toKpi=0; - } - - negtrack = 0; - } // loop on negative tracks - postrack = 0; - } // loop on positive tracks - - delete [] trkEntryP; - delete [] trkEntryN; - delete trkTree; - - printf(" Number of D0 candidates: %d\n",nD0rec1ev); - } // loop on events in file - - - printf("\n+++\n+++ Total number of events: %d\n+++\n",nTotEv); - printf("\n+++\n+++ Total number of D0 candidates: %d\n+++\n",nD0rec); - - delete vertexer1; - - esdFile->Close(); - - // create a copy of this class to be written to output file - AliD0toKpiAnalysis *copy = (AliD0toKpiAnalysis*)this->Clone("D0toKpiAnalysis"); - - // add PDG codes to decay tracks in found candidates (in simulation mode) - // and store tree in the output file - if(!fSim) { - TFile *outroot = new TFile(outName1.Data(),"recreate"); - treeD0->Write(); - copy->Write(); - outroot->Close(); - delete outroot; - } else { - printf(" Now adding information from simulation (PDG codes) ...\n"); - TTree *treeD0sim = new TTree("TreeD0","Tree with D0 candidates"); - SimulationInfo(treeD0,treeD0sim); - delete treeD0; - TFile *outroot = new TFile(outName1.Data(),"recreate"); - treeD0sim->Write(); - copy->Write(); - outroot->Close(); - delete outroot; - } - - // write to a file the total number of events - FILE *outdat = fopen(outName2.Data(),"w"); - fprintf(outdat,"%d\n",nTotEv); - fclose(outdat); - - return; -} -//----------------------------------------------------------------------------- -void AliD0toKpiAnalysis::PrintStatus() const { - // Print parameters being used - - printf("Preselections:\n"); - printf(" fPtCut = %f GeV\n",fPtCut); - printf(" fd0Cut = %f micron\n",fd0Cut); - printf(" fMassCut = %f GeV\n",fMassCut); - if(fVertexOnTheFly) printf("Primary vertex on the fly\n"); - if(fSim) { - printf("Simulation mode\n"); - if(fOnlySignal) printf(" Only signal goes to file\n"); - } - printf("Cuts on candidates:\n"); - printf(" |M-MD0| [GeV] < %f\n",fD0Cuts[0]); - printf(" dca [micron] < %f\n",fD0Cuts[1]); - printf(" cosThetaStar < %f\n",fD0Cuts[2]); - printf(" pTK [GeV] > %f\n",fD0Cuts[3]); - printf(" pTpi [GeV] > %f\n",fD0Cuts[4]); - printf(" |d0K| [micron] < %f\n",fD0Cuts[5]); - printf(" |d0pi| [micron] < %f\n",fD0Cuts[6]); - printf(" d0d0 [micron^2] < %f\n",fD0Cuts[7]); - printf(" cosThetaPoint > %f\n",fD0Cuts[8]); - - return; -} -//----------------------------------------------------------------------------- -Bool_t AliD0toKpiAnalysis::SelectInvMass(const Double_t p[6]) const { - // Apply preselection in the invariant mass of the pair - - Double_t mD0 = 1.8645; - Double_t mPi = 0.13957; - Double_t mKa = 0.49368; - - Double_t energy[2]; - Double_t mom2[2],momTot2; - - mom2[0] = p[0]*p[0] + p[1]*p[1] + p[2]*p[2]; - mom2[1] = p[3]*p[3] + p[4]*p[4] + p[5]*p[5]; - - momTot2 = (p[0]+p[3])*(p[0]+p[3])+ - (p[1]+p[4])*(p[1]+p[4])+ - (p[2]+p[5])*(p[2]+p[5]); - - // D0 -> K- pi+ - energy[1] = TMath::Sqrt(mKa*mKa+mom2[1]); - energy[0] = TMath::Sqrt(mPi*mPi+mom2[0]); - - Double_t minvD0 = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2); - - // D0bar -> K+ pi- - energy[0] = TMath::Sqrt(mKa*mKa+mom2[0]); - energy[1] = TMath::Sqrt(mPi*mPi+mom2[1]); - - Double_t minvD0bar = TMath::Sqrt((energy[0]+energy[1])*(energy[0]+energy[1])-momTot2); - - if(TMath::Abs(minvD0-mD0) < fMassCut) return kTRUE; - if(TMath::Abs(minvD0bar-mD0) < fMassCut) return kTRUE; - return kFALSE; -} -//----------------------------------------------------------------------------- -void AliD0toKpiAnalysis::SelectTracks(AliESD *event, - TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP, - TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const { - // Create two TObjArrays with positive and negative tracks and - // apply single-track preselection - - nTrksP=0,nTrksN=0; - - Int_t entr = event->GetNumberOfTracks(); - - // transfer ITS tracks from ESD to arrays and to a tree - for(Int_t i=0; iGetTrack(i); - UInt_t status = esdtrack->GetStatus(); - - if(!(status&AliESDtrack::kITSin)) continue; - - // single track selection - if(!SingleTrkCuts(*esdtrack,event->GetMagneticField())) continue; - - if(esdtrack->GetSign()<0) { // negative track - trksN.AddLast(esdtrack); - trkEntryN[nTrksN] = i; - nTrksN++; - } else { // positive track - trksP.AddLast(esdtrack); - trkEntryP[nTrksP] = i; - nTrksP++; - } - - } // loop on ESD tracks - - return; -} -//----------------------------------------------------------------------------- -void AliD0toKpiAnalysis::SetD0Cuts(Double_t cut0,Double_t cut1, - Double_t cut2,Double_t cut3,Double_t cut4, - Double_t cut5,Double_t cut6, - Double_t cut7,Double_t cut8) { - // Set the cuts for D0 selection - fD0Cuts[0] = cut0; - fD0Cuts[1] = cut1; - fD0Cuts[2] = cut2; - fD0Cuts[3] = cut3; - fD0Cuts[4] = cut4; - fD0Cuts[5] = cut5; - fD0Cuts[6] = cut6; - fD0Cuts[7] = cut7; - fD0Cuts[8] = cut8; - - return; -} -//----------------------------------------------------------------------------- -void AliD0toKpiAnalysis::SetD0Cuts(const Double_t cuts[9]) { - // Set the cuts for D0 selection - - for(Int_t i=0; i<9; i++) fD0Cuts[i] = cuts[i]; - - return; -} -//----------------------------------------------------------------------------- -Bool_t -AliD0toKpiAnalysis::SingleTrkCuts(const AliESDtrack& trk, Double_t b) const { - // Check if track passes some kinematical cuts - // Magnetic field "b" (kG) - - if(TMath::Abs(1./trk.GetParameter()[4]) < fPtCut) - return kFALSE; - if(TMath::Abs(10000.*trk.GetD(fV1[0],fV1[1],b)) < fd0Cut) - return kFALSE; - - return kTRUE; -} -//---------------------------------------------------------------------------- -void AliD0toKpiAnalysis::MakeTracksRefFile(AliRun *gAlice, - Int_t evFirst,Int_t evLast) const { - // Create a file with simulation info for the reconstructed tracks - - TFile *outFile = TFile::Open("D0TracksRefFile.root","recreate"); - TFile *esdFile = TFile::Open("AliESDs.root"); - - AliMC *mc = gAlice->GetMCApp(); - - Int_t label; - TParticle *part; - TParticle *mumpart; - REFTRACK reftrk; - - AliESD* event = new AliESD; - TTree* tree = (TTree*) esdFile->Get("esdTree"); - tree->SetBranchAddress("ESD",&event); - // loop on events in file - for(Int_t iEvent=evFirst; iEventGetEntries(); iEvent++) { - if(iEvent>evLast) break; - tree->GetEvent(iEvent); - Int_t ev = (Int_t)event->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number. - - gAlice->GetEvent(ev); - - Int_t nentr=(Int_t)event->GetNumberOfTracks(); - - // Tree for true track parameters - char ttname[100]; - sprintf(ttname,"Tree_Ref_%d",ev); - TTree *reftree = new TTree(ttname,"Tree with true track params"); - reftree->Branch("rectracks",&reftrk,"lab/I:pdg:mumlab:mumpdg:Vx/F:Vy:Vz:Px:Py:Pz"); - - for(Int_t i=0; iGetTrack(i); - label = TMath::Abs(esdtrack->GetLabel()); - - part = (TParticle*)mc->Particle(label); - reftrk.lab = label; - reftrk.pdg = part->GetPdgCode(); - reftrk.mumlab = part->GetFirstMother(); - if(part->GetFirstMother()>=0) { - mumpart = (TParticle*)gAlice->GetMCApp()->Particle(part->GetFirstMother()); - reftrk.mumpdg = mumpart->GetPdgCode(); - reftrk.mumprongs = mumpart->GetNDaughters(); - } else { - reftrk.mumpdg=-1; - reftrk.mumprongs=-1; - } - reftrk.Vx = part->Vx(); - reftrk.Vy = part->Vy(); - reftrk.Vz = part->Vz(); - reftrk.Px = part->Px(); - reftrk.Py = part->Py(); - reftrk.Pz = part->Pz(); - - reftree->Fill(); - - } // loop on tracks - - outFile->cd(); - reftree->Write(); - - delete reftree; - } // loop on events - - esdFile->Close(); - outFile->Close(); - - return; -} -//----------------------------------------------------------------------------- -void AliD0toKpiAnalysis::SimulationInfo(TTree *treeD0in,TTree *treeD0out) const { - // add pdg codes to candidate decay tracks (for sim) - - TString refFileName("D0TracksRefFile.root"); - if(fSim && gSystem->AccessPathName(refFileName.Data(),kFileExists)) { - printf("AliD0toKpiAnalysis::SimulationInfo: no reference file found!\n"); - return; - } - TFile *refFile = TFile::Open(refFileName.Data()); - - Char_t refTreeName[100]; - Int_t event; - Int_t pdg[2],mumpdg[2],mumlab[2]; - REFTRACK reftrk; - - // read-in reference tree for event 0 (the only event for Pb-Pb) - sprintf(refTreeName,"Tree_Ref_%d",0); - TTree *refTree0 = (TTree*)refFile->Get(refTreeName); - refTree0->SetBranchAddress("rectracks",&reftrk); - - AliD0toKpi *theD0 = 0; - treeD0in->SetBranchAddress("D0toKpi",&theD0); - treeD0out->Branch("D0toKpi","AliD0toKpi",&theD0,200000,0); - - Int_t entries = (Int_t)treeD0in->GetEntries(); - - for(Int_t i=0; iGetEvent(i); - event = theD0->EventNo(); - - if(event==0) { // always true for Pb-Pb (avoid to read-in tree every time) - refTree0->GetEvent(theD0->GetTrkNum(0)); - pdg[0] = reftrk.pdg; - mumpdg[0] = reftrk.mumpdg; - mumlab[0] = reftrk.mumlab; - refTree0->GetEvent(theD0->GetTrkNum(1)); - pdg[1] = reftrk.pdg; - mumpdg[1] = reftrk.mumpdg; - mumlab[1] = reftrk.mumlab; - } else { - sprintf(refTreeName,"Tree_Ref_%d",event); - TTree *refTree = (TTree*)refFile->Get(refTreeName); - refTree->SetBranchAddress("rectracks",&reftrk); - refTree->GetEvent(theD0->GetTrkNum(0)); - pdg[0] = reftrk.pdg; - mumpdg[0] = reftrk.mumpdg; - mumlab[0] = reftrk.mumlab; - refTree->GetEvent(theD0->GetTrkNum(1)); - pdg[1] = reftrk.pdg; - mumpdg[1] = reftrk.mumpdg; - mumlab[1] = reftrk.mumlab; - delete refTree; - } - - theD0->SetPdgCodes(pdg); - theD0->SetMumPdgCodes(mumpdg); - - if(TMath::Abs(mumpdg[0])==421 && - TMath::Abs(mumpdg[1])==421 && - mumlab[0]==mumlab[1] && - reftrk.mumprongs==2 && - ((TMath::Abs(pdg[0])==211 && TMath::Abs(pdg[1])==321) || - (TMath::Abs(pdg[0])==321 && TMath::Abs(pdg[1])==211)) - ) theD0->SetSignal(); - - if(!fOnlySignal || theD0->IsSignal()) treeD0out->Fill(); - - } - - delete refTree0; - - refFile->Close(); - - return; -} - - - - - - diff --git a/ANALYSIS/AliD0toKpiAnalysis.h b/ANALYSIS/AliD0toKpiAnalysis.h deleted file mode 100644 index bc6ebcff56b..00000000000 --- a/ANALYSIS/AliD0toKpiAnalysis.h +++ /dev/null @@ -1,90 +0,0 @@ -#ifndef AliD0toKpiAnalysis_H -#define AliD0toKpiAnalysis_H -/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * See cxx source for full Copyright notice */ - -//------------------------------------------------------------------------- -// Class AliD0toKpiAnalysis -// Reconstruction and analysis D0 -> K^- pi^+ -// -// Origin: A. Dainese andrea.dainese@lnl.infn.it -//------------------------------------------------------------------------- - -#include -#include -#include "AliESD.h" -#include "AliRun.h" - -//----------------------------------------------------------------------------- -class AliD0toKpiAnalysis : public TNamed { - public: - // - AliD0toKpiAnalysis(); - virtual ~AliD0toKpiAnalysis(); - - void ApplySelection(const Char_t *inName="AliD0toKpi.root", - const Char_t *outName="AliD0toKpi_sele.root") const; - void FindCandidates(Int_t evFirst=0,Int_t evLast=0, - const Char_t *outName="AliD0toKpi.root"); - void MakeTracksRefFile(AliRun *gAlice,Int_t evFirst=0,Int_t evLast=0) const; - void PrintStatus() const; - void SetVertexOnTheFly() { fVertexOnTheFly=kTRUE; } - void SetSimulation() { fSim=kTRUE; } - void SetOnlySignal() { fOnlySignal=kTRUE; } - void SetPtCut(Double_t pt=0.) { fPtCut=pt; } - void Setd0Cut(Double_t d0=0.) { fd0Cut=d0; } - void SetMassCut(Double_t deltaM=1000.) { fMassCut=deltaM; } - void SetD0Cuts(Double_t cut0=1000.,Double_t cut1=100000., - Double_t cut2=1.1,Double_t cut3=0.,Double_t cut4=0., - Double_t cut5=100000.,Double_t cut6=100000., - Double_t cut7=100000000.,Double_t cut8=-1.1); - void SetD0Cuts(const Double_t cuts[9]); - void SetPID(const Char_t * pid="TOFparam_PbPb") { fPID=pid; } - // - private: - // - Bool_t fVertexOnTheFly; // flag for primary vertex reco on the fly - Bool_t fSim; // flag for the analysis of simulated events - Bool_t fOnlySignal; // write to file only signal candidates (for sim) - TString fPID; // PID scheme - - Double_t fV1[3]; // primary vertex position (in cm) - Double_t fPtCut; // minimum track pt (in GeV/c) - Double_t fd0Cut; // minimum track |rphi impact parameter| (in micron) - Double_t fMassCut; // maximum of |InvMass-MD0| (in GeV) - Double_t fD0Cuts[9]; // cuts on D0 candidates (see SetD0Cuts()) - // (to be passed to function AliD0toKpi::Select()) - // 0 = inv. mass half width [GeV] - // 1 = dca [micron] - // 2 = cosThetaStar - // 3 = pTK [GeV/c] - // 4 = pTPi [GeV/c] - // 5 = d0K [micron] upper limit! - // 6 = d0Pi [micron] upper limit! - // 7 = d0d0 [micron^2] - // 8 = cosThetaPoint - - // - Double_t CalculateTOFmass(Double_t mom,Double_t length,Double_t time) const; - Bool_t SelectInvMass(const Double_t p[6]) const; - void SelectTracks(AliESD *event, - TObjArray &trksP,Int_t *trkEntryP,Int_t &nTrksP, - TObjArray &trksN,Int_t *trkEntryN,Int_t &nTrksN) const; - void SetVertex1(Double_t x=0.,Double_t y=0.,Double_t z=0.) - { fV1[0]=x;fV1[1]=y;fV1[2]=z; } - void SimulationInfo(TTree *treeD0in,TTree *treeD0out) const; - Bool_t SingleTrkCuts(const AliESDtrack& trk, Double_t b) const; - // - ClassDef(AliD0toKpiAnalysis,3) // Reconstruction of D0 candidates class -}; - - -#endif - - - - - - - - diff --git a/ANALYSIS/AliD0toKpiPlots.C b/ANALYSIS/AliD0toKpiPlots.C deleted file mode 100644 index 7f68ca28484..00000000000 --- a/ANALYSIS/AliD0toKpiPlots.C +++ /dev/null @@ -1,234 +0,0 @@ -void AliD0toKpiPlots(const Char_t *inName="AliD0toKpi.root", - const Char_t *outName="D0histograms.root") { - //-------------------------------------------------------------------------- - // This macro histograms many variables of D0->Kpi candidates - // - // Andrea Dainese, andrea.dainese@lnl.infn.it - //-------------------------------------------------------------------------- - - gSystem->Load("libANALYSIS.so"); - - // set of cuts - Double_t D0Cuts[9] = {0.1, // mass [GeV] - 1000000., // dca [micron] - 1.1, // cosThetaStar - 0., // pT K [GeV/c] - 0., // pT Pi [GeV/c] - 100000., // d0K upper [micron] - 100000., // d0Pi upper [micron] - 10000000000., // d0d0 [micron^2] - -1.1}; // cosThetaPointing - - // number of events (for normalization) - Bool_t normalize = kFALSE; - Double_t events = 1.; - - - - // define histograms - TH1F *hptK = new TH1F("hptK","\"K\" p_{t} distribution",50,0,10); - hptK->SetXTitle("p_{t} [GeV]"); - - TH1F *hptPi = new TH1F("hptPi","\"#pi\" p_{t} distribution",50,0,10); - hptPi->SetXTitle("p_{t} [GeV]"); - - TH1F *hDCA = new TH1F("hDCA","DCA",50,0,1000); - hDCA->SetXTitle("dca [#mu m]"); - - TH1F *hptD0 = new TH1F("hptD0","D^{0} p_{t} distribution",40,0,40); - hptD0->SetXTitle("p_{t} [GeV]"); - - TH1F *hyD0 = new TH1F("hyD0","D^{0} rapidity distribution",50,-2,2); - hyD0->SetXTitle("y"); - - TH1F *hCPtaD0 = new TH1F("hCPtaD0","cosine of pointing angle distribution",100,-1,1); - hCPtaD0->SetXTitle("cos #theta_{point}"); - - TH1F *hCPtaXY = new TH1F("hCPtaXY","cosine of pointing angle in (x,y) plane",100,-1,1); - hCPtaXY->SetXTitle("cos #theta_{point}"); - - TH1F *hCts = new TH1F("hCts","cosine of decay angle",50,-1.2,1.2); - hCts->SetXTitle("cos #theta^{*}"); - - TH2F *hCtsVsPtK = new TH2F("hCtsVsPtK","cosine of decay angle VS \"K\" p_{t}",50,0,5,50,-1,1); - hCtsVsPtK->SetYTitle("cos #theta^{*}"); - hCtsVsPtK->SetXTitle("p_{t} [GeV]"); - - TH1F *hd0d0 = new TH1F("hd0d0","Product of the impact parameters",100,-100000,100000); - hd0d0->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]"); - - TH1F *hd0K = new TH1F("hd0K","Impact parameter of \"K\"",100,-5000,5000); - hd0K->SetXTitle("d_{0}^{K} [#mu m]"); - - TH1F *hd0Pi = new TH1F("hd0Pi","Impact parameter of \"#pi\"",100,-5000,5000); - hd0Pi->SetXTitle("d_{0}^{#pi} [#mu m]"); - - TH2F *hCPtaVsd0d0 = new TH2F("hCPtaVsd0d0","cos #theta_{point} vs d_{0}^{K} #times d_{0}^{#pi}",100,-100000,100000,100,-1,1); - hCPtaVsd0d0->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]"); - hCPtaVsd0d0->SetYTitle("cos #theta_{point}"); - - TH2F *hCPtaVsd0d0zoom = new TH2F("hCPtaVsd0d0zoom","cos #theta_{point} vs d_{0}^{K} #times d_{0}^{#pi}",100,-100000,0,100,.9,1); - hCPtaVsd0d0zoom->SetXTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]"); - hCPtaVsd0d0zoom->SetYTitle("cos #theta_{point}"); - - TH2F *hd0d0VSptD0 = new TH2F("hd0d0VSptD0","d_{0}^{K} #times d_{0}^{#pi} VS D^{0} p_{t}",50,0,25,100,-120000,120000); - hd0d0VSptD0->SetYTitle("d_{0}^{K} #times d_{0}^{#pi} [#mu m^{2}]"); - hd0d0VSptD0->SetXTitle("D^{0} p_{t} [GeV]"); - - TH1F *hMass = new TH1F("hMass","Invariant mass distribution",50,1.765,1.965); - hMass->SetXTitle("M[K,#pi] [GeV]"); - - TH2F *hArm = new TH2F("hArm","Armenteros plot",50,-2,2,50,0,1); - hArm->SetXTitle("#alpha"); - hArm->SetYTitle("q_{t}"); - - // open input file and get tree - TFile *inFile = TFile::Open(inName); - - TTree *treeD0 = (TTree*)inFile->Get("TreeD0"); - AliD0toKpi *D = 0; - treeD0->SetBranchAddress("D0toKpi",&D); - Int_t entries = (Int_t)treeD0->GetEntries(); - - printf("+++\n+++ Number of D0 in tree: %d\n+++\n",entries); - - Double_t MD0,MD0bar,ctsD0,ctsD0bar,ctsPiD0,ctsPiD0bar; - Double_t WgtD0,WgtD0bar; - Double_t sampleABC=0.; - Int_t okD0=0,okD0bar=0; - Int_t nSel = 0; - Int_t ptbin; - - // loop on D0 - for(Int_t i=0; iGetEvent(i); - //--- select the PID strategy & compute weights - // D->ApplyPID("TOFparam_PbPb"); - // D->ComputeWgts(); - // get weights for the three samples A+B+C - // D->GetWgts(WgtD0,WgtD0bar,"ABC"); - WgtD0 = 1.; WgtD0bar = 1.; - - // normalize to 1 event - if(normalize) { WgtD0 /= events; WgtD0bar /= events; } - - // check if candidate passes selection (as D0 or D0bar) - D->Select(D0Cuts,okD0,okD0bar); - - // set weights to 0 if the candidate doesn't pass selection - if(!okD0) WgtD0=0.; - if(!okD0bar) WgtD0bar=0.; - if(okD0 || okD0bar) nSel++; - - // count selected candidates - sampleABC += WgtD0 + WgtD0bar; - - // inv mass and cosThetaStar - D->InvMass(MD0,MD0bar); - D->CosThetaStar(ctsD0,ctsD0bar); - - // fill histograms - hptK->Fill(D->PtChild(1),WgtD0); - hptK->Fill(D->PtChild(0),WgtD0bar); - hptPi->Fill(D->PtChild(0),WgtD0); - hptPi->Fill(D->PtChild(1),WgtD0bar); - hd0K->Fill(D->Getd0Child(1),WgtD0); - hd0K->Fill(D->Getd0Child(0),WgtD0bar); - hd0Pi->Fill(D->Getd0Child(0),WgtD0); - hd0Pi->Fill(D->Getd0Child(1),WgtD0bar); - hMass->Fill(MD0,WgtD0); - hMass->Fill(MD0bar,WgtD0bar); - hCts->Fill(ctsD0,WgtD0); - hCts->Fill(ctsD0bar,WgtD0bar); - hCtsVsPtK->Fill(D->PtChild(1),ctsD0,WgtD0); - hCtsVsPtK->Fill(D->PtChild(0),ctsD0bar,WgtD0bar); - hDCA->Fill(D->GetDCA(),WgtD0+WgtD0bar); - hptD0->Fill(D->Pt(),WgtD0+WgtD0bar); - hyD0->Fill(D->Rapidity(),WgtD0+WgtD0bar); - hd0d0->Fill(D->ProdImpParams(),WgtD0+WgtD0bar); - hCPtaD0->Fill(D->CosPointing(),WgtD0+WgtD0bar); - hCPtaXY->Fill(D->CosPointingXY(),WgtD0+WgtD0bar); - hCPtaVsd0d0->Fill(D->ProdImpParams(),D->CosPointing(),WgtD0+WgtD0bar); - hd0d0VSptD0->Fill(D->Pt(),D->ProdImpParams(),WgtD0+WgtD0bar); - hCPtaVsd0d0zoom->Fill(D->ProdImpParams(),D->CosPointing(),WgtD0+WgtD0bar); - hArm->Fill(D->Alpha(),D->Qt(),WgtD0+WgtD0bar); - - - } // end loop on D0 candidates - - inFile->Close(); - - printf("\n\n --- Total number of candidates passing selection: %d\n\n --- Sum of weights sample A+B+C: %f\n\n",nSel,sampleABC); - - // draw histograms - TCanvas *c1 = new TCanvas("c1","pt K & pi",0,0,700,700); - c1->SetLogy(); - hptK->Draw(); - hptPi->Draw("same"); - - TCanvas *c2 = new TCanvas("c2","pt D0",0,0,700,700); - c2->SetLogy(); - hptD0->Draw(); - - TCanvas *c3 = new TCanvas("c3","rapidity D0",0,0,700,700); - hyD0->Draw(); - - TCanvas *c4 = new TCanvas("c4","pointing angle",0,0,700,700); - hCPtaD0->Draw(); - - TCanvas *c5 = new TCanvas("c5","d0 x d0",0,0,700,700); - c5->SetLogy(); - hd0d0->Draw(); - - TCanvas *c6 = new TCanvas("c6","pointing angle VS d0d0",0,0,700,700); - c6->SetLogz(); - hCPtaVsd0d0->Draw("box"); - - TCanvas *c7 = new TCanvas("c7","mass",0,0,700,700); - hMass->Draw(); - - TCanvas *c8 = new TCanvas("c8","armenteros",0,0,700,700); - hArm->Draw("box"); - - TCanvas *c9 = new TCanvas("c9","decay angle",0,0,700,700); - hCts->Draw(); - - TCanvas *c10 = new TCanvas("c10","dca",0,0,700,700); - c10->SetLogy(); - hDCA->Draw(); - - TCanvas *c11 = new TCanvas("c11","d0 K & pi",0,0,700,700); - c11->SetLogy(); - hd0K->Draw(); - hd0Pi->Draw("same"); - - // write all histograms to file - TFile *outFile = new TFile(outName,"recreate"); - hMass->Write(); - hDCA->Write(); - hCts->Write(); - hCtsVsPtK->Write(); - hArm->Write(); - hCPtaVsd0d0->Write(); - hd0d0VSptD0->Write(); - hCPtaVsd0d0zoom->Write(); - hd0d0->Write(); - hCPtaD0->Write(); - hCPtaXY->Write(); - hptK->Write(); - hptPi->Write(); - hptD0->Write(); - hyD0->Write(); - hd0K->Write(); - hd0Pi->Write(); - outFile->Close(); - - return; -} - - - - diff --git a/ANALYSIS/AliD0toKpiTest.C b/ANALYSIS/AliD0toKpiTest.C deleted file mode 100644 index 83a93f3c261..00000000000 --- a/ANALYSIS/AliD0toKpiTest.C +++ /dev/null @@ -1,89 +0,0 @@ -//-------------------------------------------------------------------------- -// Test macro for reconstruction and analysis of D0->Kpi -// -// Andrea Dainese, andrea.dainese@lnl.infn.it -//-------------------------------------------------------------------------- - -void AliD0toKpiReco() { - - gSystem->Load("libANALYSIS.so"); - - //============== R E C O N S T R U C T I O N ============================== - - Int_t evFirst = 0; - Int_t evLast = 1000000; - - // Get field from galice.root - if (gAlice) { - delete gAlice->GetRunLoader(); - delete gAlice; - gAlice=0; - } - AliRunLoader *rl = AliRunLoader::Open("galice.root"); - if (rl == 0x0) { - cerr<<"Can not open session"<LoadgAlice(); - if (retval) { - cerr<<"LoadgAlice returned error"<GetAliRun(); - AliMagF *fiel = (AliMagF*)gAlice->Field(); - // Set the conversion constant between curvature and Pt - AliTracker::SetFieldMap(fiel,kTRUE); - - AliD0toKpiAnalysis *analysis = new AliD0toKpiAnalysis(); - // set simulation to take info on PDG codes from kine tree - analysis->SetSimulation(); - rl->LoadKinematics(); - analysis->MakeTracksRefFile(gAlice,evFirst,evLast); - //--- set this is you want only signal candidates in output file - //analysis->SetOnlySignal(); - //--- set this if you want to compute primary vertex D0 by D0 using - // other tracks in the event (for pp, broad interaction region); - // it is time-consuming procedure, so it can be done after a - // preselection on invariant mass - //analysis->SetVertexOnTheFly(); - //analysis->SetMassCut(0.1); // GeV - //--- set single-track preselections - analysis->SetPtCut(0.); // GeV - analysis->Setd0Cut(0.); // micron - //--- set cuts on D0 candidates to be written to file - // (see AliD0toKpiAnalysis.h for a description and for the defaults) - //analysis->SetD0Cuts(0.1,1000.,1.1,0.,0.,10000.,10000.,0.,.5); - analysis->SetD0Cuts(); - - //--- check the settings - analysis->PrintStatus(); - - analysis->FindCandidates(evFirst,evLast,"AliD0toKpi.root"); - delete analysis; - - return; -} -//========================================================================== -void AliD0toKpiSele() { - - gSystem->Load("libANALYSIS.so"); - - //======================== S E L E C T I O N ============================ - - AliD0toKpiAnalysis *analysis = new AliD0toKpiAnalysis(); - analysis->SetSimulation(); - analysis->SetOnlySignal(); - analysis->SetD0Cuts(0.1,1000.,1.1,0.,0.,10000.,10000.,0.,.5); - analysis->ApplySelection("AliD0toKpi.root","AliD0toKpi_sele.root"); - delete analysis; - - return; -} -//========================================================================== - - - - - - diff --git a/ANALYSIS/AnalysisOldLinkDef.h b/ANALYSIS/AnalysisOldLinkDef.h index 16f15f5285c..9dee55c886b 100644 --- a/ANALYSIS/AnalysisOldLinkDef.h +++ b/ANALYSIS/AnalysisOldLinkDef.h @@ -6,8 +6,6 @@ #pragma link C++ class TGliteXmlEventlist+; -#pragma link C++ class AliD0toKpi+; -#pragma link C++ class AliD0toKpiAnalysis+; #pragma link C++ class AliRunAnalysis+; #pragma link C++ class AliAnalysis+; diff --git a/ANALYSIS/libAnalysisOld.pkg b/ANALYSIS/libAnalysisOld.pkg index e76eef4bfb4..90b662b1235 100644 --- a/ANALYSIS/libAnalysisOld.pkg +++ b/ANALYSIS/libAnalysisOld.pkg @@ -9,7 +9,7 @@ SRCS= TGliteXmlEventlist.cxx\ AliReader.cxx AliReaderAOD.cxx AliReaderKineTree.cxx \ AliReaderESD.cxx AliReaderESDTree.cxx \ AliTrackPoints.cxx AliClusterMap.cxx \ - AliD0toKpi.cxx AliD0toKpiAnalysis.cxx AliFlowAnalysis.cxx \ + AliFlowAnalysis.cxx \ AliMuonAnalysis.cxx AliAODv0.cxx AliAODxi.cxx AliAODevent.cxx \ AliAnalysisEventCuts.cxx AliAnalysisTrackCuts.cxx -- 2.31.1