#include "TH3D.h"
#include "THnSparse.h"
#include "TCanvas.h"
+#include "TRandom.h"
#include <TROOT.h>
#include <TChain.h>
#include <TH1D.h>
#include <TH2D.h>
#include <TH3D.h>
-#include <TRandom.h>
#include "AliAnalysisManager.h"
#include "AliAODHandler.h"
_nClusterMin ( 80),
_trackFilterBit (0),
fNSigmaCut (3.),
+_tpcnclus ( 50),
+_chi2ndf (5.),
_field ( 1.),
_nTracks ( 0 ),
_mult0 ( 0 ),
_mult4a ( 0 ),
_mult5 ( 0 ),
_mult6 ( 0 ),
-arraySize ( 2000),
+arraySize ( 2500),
_id_1(0),
_charge_1(0),
_iEtaPhi_1(0),
_nClusterMin ( 80),
_trackFilterBit ( 0),
fNSigmaCut ( 3.),
+_tpcnclus ( 50),
+_chi2ndf (5.),
_field ( 1.),
_nTracks ( 0 ),
_mult0 ( 0 ),
_mult4a ( 0 ),
_mult5 ( 0 ),
_mult6 ( 0 ),
-arraySize ( 2000),
+arraySize ( 2500),
_id_1(0),
_charge_1(0),
_iEtaPhi_1(0),
_py_1 = new float[arraySize];
_pz_1 = new float[arraySize];
_correction_1 = new float[arraySize];
-
+ _dedx_1 = new float[arraySize];
__n1_1_vsPt = getDoubleArray(_nBins_pt_1, 0.);
__n1_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
__s1pt_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
_py_2 = new float[arraySize];
_pz_2 = new float[arraySize];
_correction_2 = new float[arraySize];
-
+ _dedx_2 = new float[arraySize];
__n1_2_vsPt = getDoubleArray(_nBins_pt_2, 0.);
__n1_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
__s1pt_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
name = "m4"; _m4 = createHisto1D(name,name,_nBins_M4, _min_M4, _max_M4, _title_m4, _title_counts);
name = "m5"; _m5 = createHisto1D(name,name,_nBins_M5, _min_M5, _max_M5, _title_m5, _title_counts);
name = "m6"; _m6 = createHisto1D(name,name,_nBins_M6, _min_M6, _max_M6, _title_m6, _title_counts);
- name = "zV"; _vertexZ = createHisto1D(name,name,1200, -12, 12, "z-Vertex (cm)", _title_counts);
+ name = "zV"; _vertexZ = createHisto1D(name,name,100, -10, 10, "z-Vertex (cm)", _title_counts);
name = "Eta"; _etadis = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
name = "DCAz"; _dcaz = createHisto1F(name,name, 500, -5.0, 5.0, "dcaZ","counts");
name = "DCAxy"; _dcaxy = createHisto1F(name,name, 500, -5.0, 5.0, "dcaXY","counts");
- //name = "Nclus1"; _Ncluster1 = createHisto1F(name,name, 200, 0, 200, "Ncluster1","counts");
+ name = "Nclus1"; _Ncluster1 = createHisto1F(name,name, 200, 0, 200, "Ncluster1","counts");
//name = "Nclus2"; _Ncluster2 = createHisto1F(name,name, 200, 0, 200, "Ncluster2","counts");
if (_singlesOnly)
name = s2PtNNwName+pair_12_Name + vsM; _s2PtNNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
name = s2NPtNwName+pair_12_Name + vsM; _s2NPtNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
- name = "mInv"; _invMass = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
+ name = "mInv"; _invMass = createHisto1F(name,name, 50, 0.41, 0.55, "M_{inv}","counts");
name = "mInvElec"; _invMassElec = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
}
int k1,k2;
int iPhi, iEta, iEtaPhi, iPt, charge;
- float q, phi, pt, eta, corr, corrPt, px, py, pz;
+ float q, phi, pt, eta, corr, corrPt, px, py, pz, dedx;
int ij;
int id_1, q_1, iEtaPhi_1, iPt_1;
- float pt_1, px_1, py_1, pz_1, corr_1;
+ float pt_1, px_1, py_1, pz_1, corr_1, dedx_1;
int id_2, q_2, iEtaPhi_2, iPt_2;
- float pt_2, px_2, py_2, pz_2, corr_2;
+ float pt_2, px_2, py_2, pz_2, corr_2, dedx_2;
float ptpt;
int iVertex, iVertexP1, iVertexP2;
int iZEtaPhiPt;
- float massElecSq = 2.5e-7;
+ float massElecSq = 1.94797849000000016e-02;
//double b[2];
//double bCov[3];
const AliAODVertex* vertex;
+ int nClus;
bool bitOK;
AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
float vertexX = -999;
float vertexY = -999;
float vertexZ = -999;
- float vertexXY = -999;
- float dcaZ = -999;
- float dcaXY = -999;
+ //float vertexXY = -999;
+ //float dcaZ = -999;
+ //float dcaXY = -999;
float centrality = -999;
if(fAODEvent)
{
//Centrality
- AliCentrality* centralityObject = fAODEvent->GetHeader()->GetCentralityP();
+ AliCentrality* centralityObject = ((AliVAODHeader*)fAODEvent->GetHeader())->GetCentralityP();
if (centralityObject)
{
//cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 6" << endl;
_eventAccounting->Fill(2);// count all events with right centrality
// filter on z and xy vertex
- vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertexSPD();
- if (!vertex || vertex->GetNContributors()<1)
+ vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertex();
+ // Double_t V[2];
+ //vertex->GetXYZ(V);
+
+ if(vertex)
{
- vertexZ = -999;
- vertexXY = -999;
- AliInfo("AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - No valid vertex object or poor vertex");
- }
- else
- {
- vertexX = vertex->GetX();
- vertexY = vertex->GetY();
- vertexZ = vertex->GetZ();
- vertexXY = sqrt(vertexX*vertexX+vertexY*vertexY);
+ Double32_t fCov[6];
+ vertex->GetCovarianceMatrix(fCov);
+ if(vertex->GetNContributors() > 0)
+ {
+ if(fCov[5] != 0)
+ {
+ vertexX = vertex->GetX();
+ vertexY = vertex->GetY();
+ vertexZ = vertex->GetZ();
+
+ if(TMath::Abs(vertexZ) > 10)
+ {
+ return;
+ } // Z-Vertex Cut
+ }
+ }
}
- if (!vertex ||
- vertexZ < _vertexZMin ||
- vertexZ > _vertexZMax ||
- vertexXY < _vertexXYMin ||
- vertexXY > _vertexXYMax)
- return;
_vertexZ->Fill(vertexZ);
//======================
//*********************************************************
- TExMap *trackMap = new TExMap();//Mapping matrix----
- //1st loop track
+ TExMap *trackMap = new TExMap();//Mapping matrix----
+
+ //1st loop track for Global tracks
for(Int_t i = 0; i < _nTracks; i++)
{
AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(i));
}
Int_t gID = aodTrack->GetID();
if (aodTrack->TestFilterBit(1)) trackMap->Add(gID, i);//Global tracks
- }
+ }
AliAODTrack* newAodTrack;
}
bitOK = t->TestFilterBit(_trackFilterBit);
- if (!bitOK) continue; //128bit
-
- Int_t gID = t->GetID();
- newAodTrack = gID >= 0 ?t : fAODEvent->GetTrack(trackMap->GetValue(-1-gID));
+ if (!bitOK) continue; //128bit or 272bit
+ Int_t gID = t->GetID();
+ newAodTrack = gID >= 0 ?t : dynamic_cast<AliAODTrack*>(fAODEvent->GetTrack(trackMap->GetValue(-1-gID)));
+ if(!newAodTrack) AliFatal("Not a standard AOD?");
+
q = t->Charge();
charge = int(q);
phi = t->Phi();
py = t->Py();
pz = t->Pz();
eta = t->Eta();
- dcaXY = t->DCA();
- dcaZ = t->ZAtDCA();
+ dedx = t->GetTPCsignal();
+ //dcaXY = t->DCA();
+ //dcaZ = t->ZAtDCA();
+ nClus = t->GetTPCNcls();
+
+ if ( nClus<_nClusterMin ) continue;
- Double_t nsigmaelectron = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kElectron));
+ _Ncluster1->Fill(nClus);
+
+ /*
+ //cuts on more than 0 shared cluster (suggested by Michael)
+ if(t->GetTPCnclsS() > 0){
+ continue;
+ }*/
+
+ //for Global tracks
+ Double_t nsigmaelectron = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kElectron));
Double_t nsigmapion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kPion));
Double_t nsigmakaon = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kKaon));
Double_t nsigmaproton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kProton));
&& nsigmakaon > fNSigmaCut
&& nsigmaproton > fNSigmaCut ) continue;
- if(charge == 0) continue;
-
- /*if (newAodTrack->PropagateToDCA(vertex, _field, 100., b, bCov))
- {
- AliAODTrack* clone = (AliAODTrack*)newAodTrack->Clone("trk_clone");
- if (clone->PropagateToDCA(vertex, _field, 100., b, bCov))
- {
- dcaXY = b[0];
- dcaZ = b[1];
- }
- else{
- dcaXY = -999999;
- dcaZ = -999999;
- }
- delete clone;
- }*/
-
+ if(charge == 0) continue;
// Kinematics cuts used
- if( pt < _min_pt_1 || pt > _max_pt_1) continue;
- if( eta < _min_eta_1 || eta > _max_eta_1) continue;
+ if( pt < _min_pt_1 || pt > _max_pt_1) continue;
+ if( eta < _min_eta_1 || eta > _max_eta_1) continue;
+ /* Double_t pos[3];
+ newAodTrack->GetXYZ(pos);
+
+ Double_t DCAX = pos[0] - vertexX;
+ Double_t DCAY = pos[1] - vertexY;
+ Double_t DCAZ = pos[2] - vertexZ;
+
+ Double_t DCAXY = TMath::Sqrt((DCAX*DCAX) + (DCAY*DCAY));
+
+ if (DCAZ < _dcaZMin ||
+ DCAZ > _dcaZMax ||
+ DCAXY > _dcaXYMax ) continue;
+ */
+
+ //------- Eff. test---------- //just for checking
+ //Double_t yy = (1 - 0.7)/1.8;
+ //Double_t zz = (pt - 0.2);
+ //Double_t effValue = 0.7 + yy*zz;
+ //Double_t R = gRandom->Rndm();
+ //if(R > effValue) continue;
+ //---------------------------
- /*Double_t pos[3];
-
- // if not constrained track the DCA is stored as position (at the vertex)
- if(!t->GetXYZ(pos)){
- dcaXY = TMath::Sign(1.0,pos[0])*TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);
- dcaZ = pos[2];
- }*/
-
-
- if (dcaZ < _dcaZMin ||
- dcaZ > _dcaZMax ||
- dcaXY < _dcaXYMin ||
- dcaXY > _dcaXYMax
- )continue;
//==== QA ===========================
- _dcaz->Fill(dcaZ);
- _dcaxy->Fill(dcaXY);
+ //_dcaz->Fill(DCAZ);
+ //_dcaxy->Fill(DCAXY);
_etadis->Fill(eta);
- _phidis->Fill(phi);
+ //_phidis->Fill(phi);
//===================================
//*************************************************
//Particle 1
- if (_requestedCharge_1 == charge)
+ if (_requestedCharge_1 == charge && dedx >= _dedxMin && dedx < _dedxMax)
{
-
iPhi = int( phi/_width_phi_1);
if (iPhi<0 || iPhi>=_nBins_phi_1 )
}
}
- if (!_sameFilter && _requestedCharge_2 == charge)
+ if (!_sameFilter && _requestedCharge_2 == charge &&
+ dedx >= _dedxMin && dedx < _dedxMax)
+
{
iPhi = int( phi/_width_phi_2);
iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
- //dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
+ dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
//1 and 2
for (int i2=i1+1; i2<k1; i2++)
{
iPt_2 = _iPt_1[i2]; ////cout << " iPt_1:" << iPt_1 << endl;
corr_2 = _correction_1[i2]; ////cout << " corr_1:" << corr_1 << endl;
pt_2 = _pt_1[i2]; ////cout << " pt_1:" << pt_1 << endl;
- //dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
+ dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
corr = corr_1*corr_2;
if (q_2>q_1 || (q_1>0 && q_2>0 && pt_2<=pt_1) || (q_1<0 && q_2<0 && pt_2>=pt_1))
{
iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
- //dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
+ dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
//1 and 2
for (int i2=i1+1; i2<k1; i2++)
{
iPt_2 = _iPt_1[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
corr_2 = _correction_1[i2]; ////cout << " corr_2:" << corr_2 << endl;
pt_2 = _pt_1[i2]; ////cout << " pt_2:" << pt_2 << endl;
- //dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
+ dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
corr = corr_1*corr_2;
if ( q_2<q_1 || (q_1>0 && q_2>0 && pt_2>=pt_1) || (q_1<0 && q_2<0 && pt_2<=pt_1))
{
px_1 = _px_1[i1]; ////cout << " px_1:" << px_1 << endl;
py_1 = _py_1[i1]; ////cout << " py_1:" << py_1 << endl;
pz_1 = _pz_1[i1]; ////cout << " pz_1:" << pz_1 << endl;
- //dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
+ dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
//1 and 2
for (int i2=0; i2<k2; i2++)
px_2 = _px_2[i2]; ////cout << " px_2:" << px_2 << endl;
py_2 = _py_2[i2]; ////cout << " py_2:" << py_2 << endl;
pz_2 = _pz_2[i2]; ////cout << " pz_2:" << pz_2 << endl;
- //dedx_2 = _dedx_2[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
+ dedx_2 = _dedx_2[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
if (_rejectPairConversion)
float mInvSq = 2*(massElecSq + sqrt(e1Sq*e2Sq) - px_1*px_2 - py_1*py_2 - pz_1*pz_2 );
float mInv = sqrt(mInvSq);
_invMass->Fill(mInv);
+ if (mInv<0.51)
+ {
+ if (dedx_1>75. && dedx_2>75.)
+ {
+ //_invMassElec->Fill(mInv);
+ if (mInv<0.05) continue;
+ }
+ }
}
corr = corr_1*corr_2;