* provided "as is" without express or implied warranty. *
**************************************************************************/
-/* $Id:$ */
-
#include "TChain.h"
#include "TList.h"
#include "TFile.h"
#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"
#include "AliPIDCombined.h"
+
ClassImp(AliAnalysisTaskDptDptCorrelations)
AliAnalysisTaskDptDptCorrelations::AliAnalysisTaskDptDptCorrelations()
_nClusterMin ( 80),
_trackFilterBit (0),
fNSigmaCut (3.),
-//_min_eta_1(0),
-//_max_eta_1 (0),
-//_min_eta_2 (0),
-//_max_eta_2 (0),
-
+_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),
_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
-//_nBins_eta_1(18),
-//_min_eta_1(-0.9), _max_eta_1(0.9), _width_eta_1(0.1),
- _nBins_eta_1(0),
- _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
+_nBins_eta_1(0), _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
_nBins_etaPhi_1(0),
_nBins_etaPhiPt_1(0),
_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
-//_nBins_eta_2(18),
-//_min_eta_2(-0.9), _max_eta_2(0.9), _width_eta_2(0.1),
- _nBins_eta_2(0),
- _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
-
-
+_nBins_eta_2(0), _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
_nBins_etaPhi_2(0),
_nBins_etaPhiPt_2(0),
_n1_1_vsEtaVsPhi ( 0),
_s1pt_1_vsEtaVsPhi ( 0),
_n1_1_vsZVsEtaVsPhiVsPt ( 0),
-_n1_1_vsM ( 0), // w/ weight
+_n1_1_vsM ( 0),
_s1pt_1_vsM ( 0),
-_n1Nw_1_vsM ( 0), // w/o weight
+_n1Nw_1_vsM ( 0),
_s1ptNw_1_vsM ( 0),
_dedxVsP_1 ( 0),
_corrDedxVsP_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),
_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
- //_nBins_eta_1(18),
- //_min_eta_1(-0.9), _max_eta_1(0.9), _width_eta_1(0.1),
- _nBins_eta_1(0),
- _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
+_nBins_eta_1(0), _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
_nBins_etaPhi_1(0),
_nBins_etaPhiPt_1(0),
_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
- //_nBins_eta_2(18),
- //_min_eta_2(-0.9), _max_eta_2(0.9), _width_eta_2(0.1),
- _nBins_eta_2(0),
- _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
-
+_nBins_eta_2(0), _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
_nBins_etaPhi_2(0),
_nBins_etaPhiPt_2(0),
_n1_1_vsEtaVsPhi ( 0),
_s1pt_1_vsEtaVsPhi ( 0),
_n1_1_vsZVsEtaVsPhiVsPt ( 0),
-_n1_1_vsM ( 0), // w/ weight
+_n1_1_vsM ( 0),
_s1pt_1_vsM ( 0),
-_n1Nw_1_vsM ( 0), // w/o weight
+_n1Nw_1_vsM ( 0),
_s1ptNw_1_vsM ( 0),
_dedxVsP_1 ( 0),
_corrDedxVsP_1 ( 0),
AliAnalysisTaskDptDptCorrelations::~AliAnalysisTaskDptDptCorrelations()
{
- /*
- delete _id_1;
- delete _charge_1;
- delete _iEtaPhi_1;
- delete _iPt_1;
- delete _pt_1;
- delete _px_1;
- delete _py_1;
- delete _pz_1;
- delete _correction_1;
- delete _dedx_1;
- delete __n1_1_vsPt;
- delete __n1_1_vsEtaPhi;
- delete __s1pt_1_vsEtaPhi;
- delete __n1_1_vsZEtaPhiPt;
- if (_correctionWeight_1) delete _correctionWeight_1;
-
- if (!_sameFilter)
- {
- delete _id_2;
- delete _charge_2;
- delete _iEtaPhi_2;
- delete _iPt_2;
- delete _pt_2;
- delete _px_2;
- delete _py_2;
- delete _pz_2;
- delete _correction_2;
- delete _dedx_2;
- delete __n1_2_vsPt;
- delete __n1_2_vsEtaPhi;
- delete __s1pt_2_vsEtaPhi;
- delete __n1_2_vsZEtaPhiPt;
- if (_correctionWeight_2) delete _correctionWeight_2;
- }
-
- if (!_singlesOnly)
- {
- delete __n2_12_vsPtPt;
- delete __n2_12_vsEtaPhi;
- delete __s2ptpt_12_vsEtaPhi;
- delete __s2PtN_12_vsEtaPhi;
- delete __s2NPt_12_vsEtaPhi;
- }
- */
+
}
void AliAnalysisTaskDptDptCorrelations::UserCreateOutputObjects()
{
- //cout<< "AliAnalysisTaskDptDptCorrelations::CreateOutputObjects() Starting " << endl;
OpenFile(0);
_outputHistoList = new TList();
_outputHistoList->SetOwner();
- //if (_useWeights) DefineInput(2, TList::Class());
-
- //Setup the parameters of the histograms
_nBins_M0 = 500; _min_M0 = 0.; _max_M0 = 5000.; _width_M0 = (_max_M0-_min_M0)/_nBins_M0;
_nBins_M1 = 500; _min_M1 = 0.; _max_M1 = 5000.; _width_M1 = (_max_M1-_min_M1)/_nBins_M1;
_nBins_M2 = 500; _min_M2 = 0.; _max_M2 = 5000.; _width_M2 = (_max_M2-_min_M2)/_nBins_M2;
_nBins_etaPhiPt_2 = _nBins_etaPhi_2 * _nBins_pt_2;
_nBins_zEtaPhiPt_2 = _nBins_vertexZ * _nBins_etaPhiPt_2;
_nBins_etaPhi_12 = _nBins_etaPhi_1 * _nBins_etaPhi_2;
-
- //setup the work arrays
-
+
_id_1 = new int[arraySize];
_charge_1 = new int[arraySize];
- //_iPhi_1 = new int[arraySize];
- //_iEta_1 = new int[arraySize];
_iEtaPhi_1 = new int[arraySize];
_iPt_1 = new int[arraySize];
_pt_1 = new float[arraySize];
_px_1 = new float[arraySize];
_py_1 = new float[arraySize];
_pz_1 = new float[arraySize];
- //_phi_1 = new float[arraySize];
- //_eta_1 = new float[arraySize];
_correction_1 = new float[arraySize];
- //_dedx_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.);
__n1_1_vsZEtaPhiPt = getFloatArray(_nBins_zEtaPhiPt_1, 0.);
- //cout << "==========================================================================================" << endl;
- //cout << "=============== Booking for particle 1 done." << endl;
- //cout << "_requestedCharge_1: " << _requestedCharge_1 << endl;
- //cout << "_requestedCharge_2: " << _requestedCharge_2 << endl;
-
+
if (_requestedCharge_2!=_requestedCharge_1)
{
- //cout << " creating arrays for particle 2 with size: " << arraySize << endl;
- _sameFilter = 0;
+ _sameFilter = 0;
//particle 2
_id_2 = new int[arraySize];
_charge_2 = new int[arraySize];
- //_iPhi_2 = new int[arraySize];
- //_iEta_2 = new int[arraySize];
_iEtaPhi_2 = new int[arraySize];
_iPt_2 = new int[arraySize];
_pt_2 = new float[arraySize];
_px_2 = new float[arraySize];
_py_2 = new float[arraySize];
_pz_2 = new float[arraySize];
- //_phi_2 = new float[arraySize];
- //_eta_2 = new float[arraySize];
_correction_2 = new float[arraySize];
- //_dedx_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 = "eventAccounting";
- // bin index : what it is...
- // 0 : number of event submitted
- // 1 : number accepted by centrality cut
- // 2 : number accepted by centrality cut and z cut
- // 3 : total number of particles that satisfied filter 1
- // 4 : total number of particles that satisfied filter 2
- _eventAccounting = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
+ _eventAccounting = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
name = "m0"; _m0 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m0, _title_counts);
name = "m1"; _m1 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m1, _title_counts);
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,_nBins_vertexZ, _min_vertexZ, _max_vertexZ, "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 = "Phi"; _phidis = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","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 = "dedxVsP_1"; _dedxVsP_1 = createHisto2D(name,name,1000,-10.,10.,1000,0.,1000.,"p (GeV/c)", "dedx", "counts");
- name = "dedxVsP_2"; _dedxVsP_2 = createHisto2D(name,name,1000,-10.,10.,1000,0.,1000.,"p (GeV/c)", "dedx", "counts");
- name = "corrDedxVsP_1"; _corrDedxVsP_1 = createHisto2D(name,name,1000,-10.,10.,1000,0.,500,"p (GeV/c)", "dedx", "counts");
- name = "corrDedxVsP_2"; _corrDedxVsP_2 = createHisto2D(name,name,1000,-10.,10.,1000,0.,500,"p (GeV/c)", "dedx", "counts");
- */
+
name = "Nclus1"; _Ncluster1 = createHisto1F(name,name, 200, 0, 200, "Ncluster1","counts");
- name = "Nclus2"; _Ncluster2 = createHisto1F(name,name, 200, 0, 200, "Ncluster2","counts");
+ //name = "Nclus2"; _Ncluster2 = createHisto1F(name,name, 200, 0, 200, "Ncluster2","counts");
if (_singlesOnly)
{
name = n1Name+part_1_Name+vsPt; _n1_1_vsPt = createHisto1F(name,name, _nBins_pt_1, _min_pt_1, _max_pt_1, _title_pt_1, _title_AvgN_1);
name = n1Name+part_1_Name+vsZ+vsEtaPhi+vsPt; _n1_1_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_1, 0., double(_nBins_etaPhi_1), _nBins_pt_1, _min_pt_1, _max_pt_1, "zVertex", _title_etaPhi_1, _title_pt_1);
- //name = "dedxVsP_1"; _dedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
- //name = "corrDedxVsP_1"; _corrDedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
- //name = "betaVsP_1"; _betaVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
name = n1Name+part_2_Name+vsPt; _n1_2_vsPt = createHisto1F(name,name, _nBins_pt_2, _min_pt_2, _max_pt_2, _title_pt_2, _title_AvgN_2);
name = n1Name+part_2_Name+vsZ+vsEtaPhi+vsPt; _n1_2_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_2, 0., double(_nBins_etaPhi_2), _nBins_pt_2, _min_pt_2, _max_pt_2, "zVertex", _title_etaPhi_2, _title_pt_2);
- //name = "dedxVsP_2"; _dedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
- //name = "corrDedxVsP_2"; _corrDedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
- //name = "betaVsP_2"; _betaVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
}
else
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;
- //float p, dedx, dedx_1, dedx_2;
- //double nSigma;
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();
if (!manager) {
- //cout<<"ERROR: Analysis manager not found."<<endl;
return;
}
- //coneect to the inputHandler------------
AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler());
if (!inputHandler) {
- //cout<<"ERROR: Input handler not found."<<endl;
return;
}
-
+
fAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());
+ //AliAODEvent* fAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());
+
if (!fAODEvent)
{
- //cout<< "ERROR 01: AOD not found " <<endl;
return;
}
fPIDResponse =inputHandler->GetPIDResponse();
AliFatal("This Task needs the PID response attached to the inputHandler");
return;
}
-
+
// count all events looked at here
_eventCount++;
-
- // A. Adare - Fix from above to avoid coverity complaint about null dereference
+
if (_eventAccounting)
{
_eventAccounting->Fill(0);// count all calls to this function
}
else
{
-
- return;
+
+ return;
}
-
- //fPIDResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetPIDResponse();
-//fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
-
-
_eventAccounting->Fill(1);// count all calls to this function with a valid pointer
//reset single particle counters
k1 = k2 = 0;
float vertexX = -999;
float vertexY = -999;
float vertexZ = -999;
- float vertexXY = -999;
- //float dcaX = -999;
- //float dcaY = -999;
- float dcaZ = -999;
- float dcaXY = -999;
+ //float vertexXY = -999;
+ //float dcaZ = -999;
+ //float dcaXY = -999;
float centrality = -999;
- //Double_t nSigma =-999;
-
+
if(fAODEvent)
{
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 5" << endl;
-
- //Centrality
- AliCentrality* centralityObject = fAODEvent->GetHeader()->GetCentralityP();
- if (centralityObject)
- {
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 6" << endl;
+ //Centrality
+ AliCentrality* centralityObject = ((AliVAODHeader*)fAODEvent->GetHeader())->GetCentralityP();
+ if (centralityObject)
+ {
+ //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 6" << endl;
+
+ v0Centr = centralityObject->GetCentralityPercentile("V0M");
+ v0ACentr = centralityObject->GetCentralityPercentile("V0A");
+ trkCentr = centralityObject->GetCentralityPercentile("TRK");
+ spdCentr = centralityObject->GetCentralityPercentile("CL1");
+
+ }
- v0Centr = centralityObject->GetCentralityPercentile("V0M");
- v0ACentr = centralityObject->GetCentralityPercentile("V0A");
- trkCentr = centralityObject->GetCentralityPercentile("TRK");
- spdCentr = centralityObject->GetCentralityPercentile("CL1");
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 7" << endl;
+ _nTracks =fAODEvent->GetNumberOfTracks();//NEW Test
- }
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 8" << endl;
-
- // _nTracks = fAODEvent->GetNTracks(); //OLD
- _nTracks =fAODEvent->GetNumberOfTracks();//NEW Test
-
+ _mult3 = _nTracks;
+ _mult4 = v0Centr;
+ _mult4a = v0ACentr;
+ _mult5 = trkCentr;
+ _mult6 = spdCentr;
+ _field = fAODEvent->GetMagneticField();
-
- //_mult0 = 0;
- //_mult1 = 0;
- //_mult2 = 0;
- _mult3 = _nTracks;
- _mult4 = v0Centr;
- _mult4a = v0ACentr;
- _mult5 = trkCentr;
- _mult6 = spdCentr;
- _field = fAODEvent->GetMagneticField();
-
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - _field:" << _field << endl;
-
- //_centralityMethod
- switch (_centralityMethod)
- {
- case 0: centrality = _mult0; break;
- case 1: centrality = _mult1; break;
- case 2: centrality = _mult2; break;
- case 3: centrality = _mult3; break;
- case 4: centrality = _mult4; break;
- case 5: centrality = _mult5; break;
- case 6: centrality = _mult6; break;
- case 7: centrality = _mult4a; break;
- }
-
-
- if ( centrality < _centralityMin ||
- centrality > _centralityMax ||
- fabs(v0Centr-trkCentr)>5.0)
- {
-
- return;
- }
- _eventAccounting->Fill(2);// count all events with right centrality
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 12" << endl;
-
- // filter on z and xy vertex
- vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertexSPD();
- if (!vertex || vertex->GetNContributors()<1)
- {
- 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);
- }
- if (!vertex ||
- vertexZ < _vertexZMin ||
- vertexZ > _vertexZMax ||
- vertexXY < _vertexXYMin ||
- vertexXY > _vertexXYMax)
- return;
-
- iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
- iVertexP1 = iVertex*_nBins_etaPhiPt_1;
- iVertexP2 = iVertex*_nBins_etaPhiPt_2;
- if (iVertex<0 || iVertex>=_nBins_vertexZ)
- {
- AliError("AliAnalysisTaskDptDptCorrelations::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
- return;
- }
- _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
- //======================
-
- TExMap *trackMap = new TExMap();//Mapping matrix----
-
- //1st loop track
- for(Int_t i = 0; i < _nTracks; i++)
- {
- AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(i));
- if(!aodTrack) {
- AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
- continue;
+ //_centralityMethod
+ switch (_centralityMethod)
+ {
+ case 0: centrality = _mult0; break;
+ case 1: centrality = _mult1; break;
+ case 2: centrality = _mult2; break;
+ case 3: centrality = _mult3; break;
+ case 4: centrality = _mult4; break;
+ case 5: centrality = _mult5; break;
+ case 6: centrality = _mult6; break;
+ case 7: centrality = _mult4a; break;
}
- //Double_t tpcSignalAll = aodTrack->GetTPCsignal();
- //_dedxVsP_1->Fill(aodTrack->GetTPCmomentum(),tpcSignalAll);
- Int_t gID = aodTrack->GetID();
- if (aodTrack->TestFilterBit(1)) trackMap->Add(gID, i);//Global tracks
- }
-
-
- //--------------------------
-
- AliAODTrack* newAodTrack;
-
- //2nd loop track starts here
- for (int iTrack=0; iTrack< _nTracks; iTrack++)
- {
- //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 16 iTrack:" << iTrack << endl;
- AliAODTrack * t = (AliAODTrack *) fAODEvent->GetTrack(iTrack);
- if (!t)
- {
- AliError(Form("AliAnalysisTaskDptDptCorrelations::Exec(Option_t *option) No track ofr iTrack=%d", iTrack));
- continue;
- }
- bitOK = t->TestFilterBit(_trackFilterBit);
- if (!bitOK) continue;
-
- Int_t gID = t->GetID();
- newAodTrack = gID >= 0 ?t : fAODEvent->GetTrack(trackMap->GetValue(-1-gID));
-
-
- q = t->Charge();
- charge = int(q);
- phi = t->Phi();
- //p = t->P();
- pt = t->Pt();
- px = t->Px();
- py = t->Py();
- pz = t->Pz();
- eta = t->Eta();
- //dedx = t->GetTPCsignal();
- nClus = t->GetTPCNcls();
- Double_t nclus2 = t->GetTPCClusterInfo(2,1);
+ if ( centrality < _centralityMin || centrality > _centralityMax )
+ {
+ return;
+ }
+ _eventAccounting->Fill(2);// count all events with right centrality
+ // filter on z and xy vertex
+ vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertex();
+ // Double_t V[2];
+ //vertex->GetXYZ(V);
+
+ if(vertex)
+ {
+ 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
+ }
+ }
+ }
+ _vertexZ->Fill(vertexZ);
-
- //Float_t nsigmaTPCPID = -999.;
- //Float_t nsigmaTOFPID = -999.;
+ iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
+ iVertexP1 = iVertex*_nBins_etaPhiPt_1;
+ iVertexP2 = iVertex*_nBins_etaPhiPt_2;
+ if (iVertex<0 || iVertex>=_nBins_vertexZ)
+ {
+ AliError("AliAnalysisTaskDptDptCorrelations::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
+ return;
+ }
+ _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
+ //======================
- 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));
+ //*********************************************************
+ 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));
+ if(!aodTrack) {
+ AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
+ continue;
+ }
+ Int_t gID = aodTrack->GetID();
+ if (aodTrack->TestFilterBit(1)) trackMap->Add(gID, i);//Global tracks
+ }
+
+ AliAODTrack* newAodTrack;
- //Fill QA before the PID
- //_dedxVsP_1 -> Fill(newAodTrack->P()*newAodTrack->Charge(),newAodTrack->GetTPCsignal());
- //_corrDedxVsP_1 -> Fill(newAodTrack->P()*newAodTrack->Charge(),nsigmaelectron);
- //end of QA-before pid
-
- //-----------------------------------------
- //nsigma cut to reject electron (NEW)
- if(nsigmaelectron < fNSigmaCut
- && nsigmapion > fNSigmaCut
- && nsigmakaon > fNSigmaCut
- && nsigmaproton > fNSigmaCut ) continue;
- //------------------------------------------
-
- //Fill QA after the PID
- // _dedxVsP_2 -> Fill(newAodTrack->P()*newAodTrack->Charge(),newAodTrack->GetTPCsignal());
- //_corrDedxVsP_2 -> Fill(newAodTrack->P()*newAodTrack->Charge(),nsigmaelectron);
- //-----------------
-
- //Particle 1
- if (_requestedCharge_1 == charge &&
- pt >= _min_pt_1 &&
- pt < _max_pt_1 &&
- eta >= _min_eta_1 &&
- eta < _max_eta_1)
- {
-
+ //Track Loop starts here
+ for (int iTrack=0; iTrack< _nTracks; iTrack++)
+ {
+ AliAODTrack* t = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(iTrack));
+ if (!t) {
+ AliError(Form("Could not receive track %d", iTrack));
+ continue;
+ }
+
+ bitOK = t->TestFilterBit(_trackFilterBit);
+ 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();
+ pt = t->Pt();
+ px = t->Px();
+ py = t->Py();
+ pz = t->Pz();
+ eta = t->Eta();
+ dedx = t->GetTPCsignal();
+ //dcaXY = t->DCA();
+ //dcaZ = t->ZAtDCA();
+ nClus = t->GetTPCNcls();
+
+ if ( nClus<_nClusterMin ) continue;
+
+ _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));
+
+ //nsigma cut to reject electron
+
+ if(nsigmaelectron < fNSigmaCut
+ && nsigmapion > fNSigmaCut
+ && nsigmakaon > fNSigmaCut
+ && nsigmaproton > fNSigmaCut ) continue;
- iPhi = int( phi/_width_phi_1);
-
- if (iPhi<0 || iPhi>=_nBins_phi_1 )
- {
- AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
- return;
- }
-
- iEta = int((eta-_min_eta_1)/_width_eta_1);
- if (iEta<0 || iEta>=_nBins_eta_1)
- {
- AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
- continue;
- }
- iPt = int((pt -_min_pt_1 )/_width_pt_1 );
- if (iPt<0 || iPt >=_nBins_pt_1)
- {
- AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
- continue;
- }
- iEtaPhi = iEta*_nBins_phi_1+iPhi;
- iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
-
- if (_correctionWeight_1)
- corr = _correctionWeight_1[iZEtaPhiPt];
- else
- corr = 1;
- if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
- {
- AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
- continue;
- }
-
- if (_singlesOnly)
- {
-
- __n1_1_vsPt[iPt] += corr; //cout << "step 15" << endl;
- __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
-
- }
- else
- {
- corrPt = corr*pt;
- _id_1[k1] = iTrack;
- _charge_1[k1] = charge;
- _iEtaPhi_1[k1] = iEtaPhi;
- _iPt_1[k1] = iPt;
- _pt_1[k1] = pt;
- _px_1[k1] = px;
- _py_1[k1] = py;
- _pz_1[k1] = pz;
- _correction_1[k1] = corr;
- __n1_1 += corr;
- __n1_1_vsEtaPhi[iEtaPhi] += corr;
- __s1pt_1 += corrPt;
- __s1pt_1_vsEtaPhi[iEtaPhi] += corrPt;
- __n1Nw_1 += 1;
- __s1ptNw_1 += pt;
- ++k1;
- if (k1>=arraySize)
- {
- AliError(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
- return;
- }
- }
- }
-
- if (!_sameFilter &&
- _requestedCharge_2 == charge &&
- pt >= _min_pt_2 &&
- pt < _max_pt_2 &&
- eta >= _min_eta_2 &&
- eta < _max_eta_2)
- {
-
- dcaXY = t->DCA(); //new change Prabhat
- dcaZ = t->ZAtDCA(); //new change Prabhat
+ 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;
- //Check dca cuts for systematics
- //Comment this portion when not required
- /*if (dcaZ < _dcaZMin ||
- dcaZ > _dcaZMax ||
- dcaXY < _dcaXYMin ||
- dcaXY > _dcaXYMax ) 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;
+ //---------------------------
+
//==== QA ===========================
- _dcaz->Fill(dcaZ);
- _dcaxy->Fill(dcaXY);
+ //_dcaz->Fill(DCAZ);
+ //_dcaxy->Fill(DCAXY);
_etadis->Fill(eta);
- _phidis->Fill(phi);
- _Ncluster1->Fill(nClus);
- _Ncluster2->Fill(nclus2);
-
+ //_phidis->Fill(phi);
//===================================
-
- iPhi = int( phi/_width_phi_2);
+ //*************************************************
+
+ //Particle 1
+ if (_requestedCharge_1 == charge && dedx >= _dedxMin && dedx < _dedxMax)
+ {
+ iPhi = int( phi/_width_phi_1);
+
+ if (iPhi<0 || iPhi>=_nBins_phi_1 )
+ {
+ AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
+ return;
+ }
+
+ iEta = int((eta-_min_eta_1)/_width_eta_1);
+ if (iEta<0 || iEta>=_nBins_eta_1)
+ {
+ AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
+ continue;
+ }
+ iPt = int((pt -_min_pt_1 )/_width_pt_1 );
+ if (iPt<0 || iPt >=_nBins_pt_1)
+ {
+ AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
+ continue;
+ }
+ iEtaPhi = iEta*_nBins_phi_1+iPhi;
+ iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
+
+ if (_correctionWeight_1)
+ corr = _correctionWeight_1[iZEtaPhiPt];
+ else
+ corr = 1;
+ if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
+ {
+ AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
+ continue;
+ }
+
+
+ if (_singlesOnly)
+ {
+
+ __n1_1_vsPt[iPt] += corr; //cout << "step 15" << endl;
+ __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
+
+ }
+ else
+ {
+ corrPt = corr*pt;
+ _id_1[k1] = iTrack;
+ _charge_1[k1] = charge;
+ _iEtaPhi_1[k1] = iEtaPhi;
+ _iPt_1[k1] = iPt;
+ _pt_1[k1] = pt;
+ _px_1[k1] = px;
+ _py_1[k1] = py;
+ _pz_1[k1] = pz;
+ _correction_1[k1] = corr;
+ __n1_1 += corr;
+ __n1_1_vsEtaPhi[iEtaPhi] += corr;
+ __s1pt_1 += corrPt;
+ __s1pt_1_vsEtaPhi[iEtaPhi] += corrPt;
+ __n1Nw_1 += 1;
+ __s1ptNw_1 += pt;
+ ++k1;
+ if (k1>=arraySize)
+ {
+ AliError(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
+ return;
+ }
+ }
+ }
- if (iPhi<0 || iPhi>=_nBins_phi_2 )
- {
- AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
- return;
- }
-
- iEta = int((eta-_min_eta_2)/_width_eta_2);
- if (iEta<0 || iEta>=_nBins_eta_2)
- {
- AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
- continue;
- }
- iPt = int((pt -_min_pt_2 )/_width_pt_2 );
- if (iPt<0 || iPt >=_nBins_pt_2)
- {
- AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
- continue;
- }
-
- iEtaPhi = iEta*_nBins_phi_2+iPhi;
- iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
- if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
- {
- AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
- continue;
- }
-
-
- if (_correctionWeight_2)
- corr = _correctionWeight_2[iZEtaPhiPt];
- else
- corr = 1;
- //dpt = pt - (charge>0) ? _avgPt_vsEtaPhi_2p[iEtaPhi] : _avgPt_vsEtaPhi_2m[iEtaPhi];
-
- if (_singlesOnly)
- {
- //_dedxVsP_2->Fill(p*q,dedx);
- __n1_2_vsPt[iPt] += corr; //cout << "step 15" << endl;
- __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
- }
- else
- {
- corrPt = corr*pt;
- _id_2[k2] = iTrack; //cout << "step 1" << endl;
- _charge_2[k2] = charge; //cout << "step 2" << endl;
- _iEtaPhi_2[k2] = iEtaPhi; //cout << "step 3" << endl;
- _iPt_2[k2] = iPt; //cout << "step 4" << endl;
- _pt_2[k2] = pt; //cout << "step 5" << endl;
- _px_2[k2] = px; //cout << "step 6" << endl;
- _py_2[k2] = py; //cout << "step 7" << endl;
- _pz_2[k2] = pz; //cout << "step 8" << endl;
- _correction_2[k2] = corr; //cout << "step 9" << endl;
- __n1_2 += corr; //cout << "step 10" << endl;
- __s1pt_2 += corrPt; //cout << "step 13" << endl;
- __n1Nw_2 += 1;
- __n1_2_vsEtaPhi[iEtaPhi] += corr; //cout << "step 11" << endl;
- __s1pt_2_vsEtaPhi[iEtaPhi] += corrPt; //cout << "step 14" << endl;
- __s1ptNw_2 += pt;
- ++k2;
- if (k2>=arraySize)
- {
- AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
- return;
- }
- }
-
- //cout << "done with track" << endl;
- } //iTrack
- } //aod
+ if (!_sameFilter && _requestedCharge_2 == charge &&
+ dedx >= _dedxMin && dedx < _dedxMax)
+
+ {
+
+ iPhi = int( phi/_width_phi_2);
+
+ if (iPhi<0 || iPhi>=_nBins_phi_2 )
+ {
+ AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
+ return;
+ }
+
+ iEta = int((eta-_min_eta_2)/_width_eta_2);
+ if (iEta<0 || iEta>=_nBins_eta_2)
+ {
+ AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
+ continue;
+ }
+ iPt = int((pt -_min_pt_2 )/_width_pt_2 );
+ if (iPt<0 || iPt >=_nBins_pt_2)
+ {
+ AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
+ continue;
+ }
+
+ iEtaPhi = iEta*_nBins_phi_2+iPhi;
+ iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
+ if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
+ {
+ AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
+ continue;
+ }
+
+
+ if (_correctionWeight_2)
+ corr = _correctionWeight_2[iZEtaPhiPt];
+ else
+ corr = 1;
+
+ if (_singlesOnly)
+ {
+ __n1_2_vsPt[iPt] += corr; //cout << "step 15" << endl;
+ __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
+ }
+ else
+ {
+ corrPt = corr*pt;
+ _id_2[k2] = iTrack; //cout << "step 1" << endl;
+ _charge_2[k2] = charge; //cout << "step 2" << endl;
+ _iEtaPhi_2[k2] = iEtaPhi; //cout << "step 3" << endl;
+ _iPt_2[k2] = iPt; //cout << "step 4" << endl;
+ _pt_2[k2] = pt; //cout << "step 5" << endl;
+ _px_2[k2] = px; //cout << "step 6" << endl;
+ _py_2[k2] = py; //cout << "step 7" << endl;
+ _pz_2[k2] = pz; //cout << "step 8" << endl;
+ _correction_2[k2] = corr; //cout << "step 9" << endl;
+ __n1_2 += corr; //cout << "step 10" << endl;
+ __s1pt_2 += corrPt; //cout << "step 13" << endl;
+ __n1Nw_2 += 1;
+ __n1_2_vsEtaPhi[iEtaPhi] += corr; //cout << "step 11" << endl;
+ __s1pt_2_vsEtaPhi[iEtaPhi] += corrPt; //cout << "step 14" << endl;
+ __s1ptNw_2 += pt;
+ ++k2;
+ if (k2>=arraySize)
+ {
+ AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
+ return;
+ }
+ }
+
+ //cout << "done with track" << endl;
+ } //iTrack
+ } //aod
}
-
+
//cout << "Filling histograms now" << endl;
_m0->Fill(_mult0);
_m4->Fill(_mult4);
_m5->Fill(_mult5);
_m6->Fill(_mult6);
- _vertexZ->Fill(vertexZ);
+ //_vertexZ->Fill(vertexZ);
if (_singlesOnly)
{
- // nothing to do here.
+ // nothing to do here.
}
else
{
- if (_sameFilter)
- {
+ if (_sameFilter)
+ {
_n1_1_vsM->Fill(centrality, __n1_1);
_s1pt_1_vsM->Fill(centrality, __s1pt_1);
_n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
__n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
if (_field>0)
{
- for (int i1=0; i1<k1; i1++)
- {
- ////cout << " i1:" << i1 << endl;
- id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
- q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
- iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
- 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;
- //1 and 2
- for (int i2=i1+1; i2<k1; i2++)
- {
- ////cout << " i2:" << i2 << endl;
- id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
- if (id_1!=id_2)
- {
- q_2 = _charge_1[i2]; ////cout << " q_1:" << q_1 << endl;
- iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
- 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;
- 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))
- {
- ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
- }
- else // swap particles
- {
- ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
- }
-
- __n2_12 += corr;
- __n2_12_vsEtaPhi[ij] += corr;
- ptpt = pt_1*pt_2;
- __s2ptpt_12 += corr*ptpt;
- __s2PtN_12 += corr*pt_1;
- __s2NPt_12 += corr*pt_2;
- __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
- __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
- __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
- __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
-
- __n2Nw_12 += 1;
- __s2ptptNw_12 += ptpt;
- __s2PtNNw_12 += pt_1;
- __s2NPtNw_12 += pt_2;
-
- }
- } //i2
- } //i1
+ for (int i1=0; i1<k1; i1++)
+ {
+ ////cout << " i1:" << i1 << endl;
+ id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
+ q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
+ iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
+ 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;
+ //1 and 2
+ for (int i2=i1+1; i2<k1; i2++)
+ {
+ ////cout << " i2:" << i2 << endl;
+ id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
+ if (id_1!=id_2)
+ {
+ q_2 = _charge_1[i2]; ////cout << " q_1:" << q_1 << endl;
+ iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
+ 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;
+ 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))
+ {
+ ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
+ }
+ else // swap particles
+ {
+ ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
+ }
+
+ __n2_12 += corr;
+ __n2_12_vsEtaPhi[ij] += corr;
+ ptpt = pt_1*pt_2;
+ __s2ptpt_12 += corr*ptpt;
+ __s2PtN_12 += corr*pt_1;
+ __s2NPt_12 += corr*pt_2;
+ __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
+ __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
+ __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
+ __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
+
+ __n2Nw_12 += 1;
+ __s2ptptNw_12 += ptpt;
+ __s2PtNNw_12 += pt_1;
+ __s2NPtNw_12 += pt_2;
+
+ }
+ } //i2
+ } //i1
}
else // field<0
{
- for (int i1=0; i1<k1; i1++)
- {
- ////cout << " i1:" << i1 << endl;
- id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
- q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
- iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
- 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;
- //1 and 2
- for (int i2=i1+1; i2<k1; i2++)
- {
- ////cout << " i2:" << i2 << endl;
- id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
- if (id_1!=id_2)
- {
- q_2 = _charge_1[i2]; ////cout << " q_2:" << q_2 << endl;
- iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
- 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;
- 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))
- {
- ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
- }
- else // swap particles
- {
- ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
- }
-
- __n2_12 += corr;
- __n2_12_vsEtaPhi[ij] += corr;
- ptpt = pt_1*pt_2;
- __s2ptpt_12 += corr*ptpt;
- __s2PtN_12 += corr*pt_1;
- __s2NPt_12 += corr*pt_2;
- __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
- __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
- __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
- __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
-
- __n2Nw_12 += 1;
- __s2ptptNw_12 += ptpt;
- __s2PtNNw_12 += pt_1;
- __s2NPtNw_12 += pt_2;
-
- }
- } //i2
- } //i1
+ for (int i1=0; i1<k1; i1++)
+ {
+ ////cout << " i1:" << i1 << endl;
+ id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
+ q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
+ iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
+ 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;
+ //1 and 2
+ for (int i2=i1+1; i2<k1; i2++)
+ {
+ ////cout << " i2:" << i2 << endl;
+ id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
+ if (id_1!=id_2)
+ {
+ q_2 = _charge_1[i2]; ////cout << " q_2:" << q_2 << endl;
+ iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
+ 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;
+ 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))
+ {
+ ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
+ }
+ else // swap particles
+ {
+ ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
+ }
+
+ __n2_12 += corr;
+ __n2_12_vsEtaPhi[ij] += corr;
+ ptpt = pt_1*pt_2;
+ __s2ptpt_12 += corr*ptpt;
+ __s2PtN_12 += corr*pt_1;
+ __s2NPt_12 += corr*pt_2;
+ __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
+ __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
+ __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
+ __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
+
+ __n2Nw_12 += 1;
+ __s2ptptNw_12 += ptpt;
+ __s2PtNNw_12 += pt_1;
+ __s2NPtNw_12 += pt_2;
+
+ }
+ } //i2
+ } //i1
}
- }
- else // filter 1 and 2 are different -- must do all particle pairs...
- {
- _n1_1_vsM->Fill(centrality, __n1_1);
- _s1pt_1_vsM->Fill(centrality, __s1pt_1);
- _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
- _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
- _n1_2_vsM->Fill(centrality, __n1_2);
- _s1pt_2_vsM->Fill(centrality, __s1pt_2);
- _n1Nw_2_vsM->Fill(centrality, __n1Nw_2);
- _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_2);
- // reset pair counters
- __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
- __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
- for (int i1=0; i1<k1; i1++)
- {
- ////cout << " i1:" << i1 << endl;
- id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
- q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
- iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
- 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;
- 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;
-
- //1 and 2
- for (int i2=0; i2<k2; i2++)
- {
- ////cout << " i2:" << i2 << endl;
- id_2 = _id_2[i2]; ////cout << " id_2:" << id_2 << endl;
- if (id_1!=id_2) // exclude auto correlation
- {
- q_2 = _charge_2[i2]; ////cout << " q_2:" << q_2 << endl;
- iEtaPhi_2 = _iEtaPhi_2[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
- iPt_2 = _iPt_2[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
- corr_2 = _correction_2[i2]; ////cout << " corr_2:" << corr_2 << endl;
- pt_2 = _pt_2[i2]; ////cout << " pt_2:" << pt_2 << endl;
- 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;
-
-
- if (_rejectPairConversion)
- {
- float e1Sq = massElecSq + pt_1*pt_1 + pz_1*pz_1;
- float e2Sq = massElecSq + pt_2*pt_2 + pz_2*pz_2;
- 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.5)
- //{
- // if (dedx_1>75. && dedx_2>75.)
- //{
- _invMassElec->Fill(mInv);
- //if (mInv<0.05) continue;
- //}
- //}
- }
- corr = corr_1*corr_2;
- ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
- __n2_12 += corr;
- __n2_12_vsEtaPhi[ij] += corr;
- ptpt = pt_1*pt_2;
- __s2ptpt_12 += corr*ptpt;
- __s2PtN_12 += corr*pt_1;
- __s2NPt_12 += corr*pt_2;
- __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
- __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
- __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
- __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
- __n2Nw_12 += 1;
- __s2ptptNw_12 += ptpt;
- __s2PtNNw_12 += pt_1;
- __s2NPtNw_12 += pt_2;
-
- }
- } //i2
- } //i1
- }
-
- _n2_12_vsM->Fill(centrality, __n2_12);
- _s2PtPt_12_vsM->Fill(centrality, __s2ptpt_12);
- _s2PtN_12_vsM->Fill(centrality, __s2NPt_12);
- _s2NPt_12_vsM->Fill(centrality, __s2PtN_12);
-
- _n2Nw_12_vsM->Fill(centrality, __n2Nw_12);
- _s2PtPtNw_12_vsM->Fill(centrality, __s2ptptNw_12);
- _s2PtNNw_12_vsM->Fill(centrality, __s2NPtNw_12);
- _s2NPtNw_12_vsM->Fill(centrality, __s2PtNNw_12);
-
+ }
+ else // filter 1 and 2 are different -- must do all particle pairs...
+ {
+ _n1_1_vsM->Fill(centrality, __n1_1);
+ _s1pt_1_vsM->Fill(centrality, __s1pt_1);
+ _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
+ _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
+ _n1_2_vsM->Fill(centrality, __n1_2);
+ _s1pt_2_vsM->Fill(centrality, __s1pt_2);
+ _n1Nw_2_vsM->Fill(centrality, __n1Nw_2);
+ _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_2);
+ // reset pair counters
+ __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
+ __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
+ for (int i1=0; i1<k1; i1++)
+ {
+ ////cout << " i1:" << i1 << endl;
+ id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
+ q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
+ iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
+ 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;
+ 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;
+
+ //1 and 2
+ for (int i2=0; i2<k2; i2++)
+ {
+ ////cout << " i2:" << i2 << endl;
+ id_2 = _id_2[i2]; ////cout << " id_2:" << id_2 << endl;
+ if (id_1!=id_2) // exclude auto correlation
+ {
+ q_2 = _charge_2[i2]; ////cout << " q_2:" << q_2 << endl;
+ iEtaPhi_2 = _iEtaPhi_2[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
+ iPt_2 = _iPt_2[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
+ corr_2 = _correction_2[i2]; ////cout << " corr_2:" << corr_2 << endl;
+ pt_2 = _pt_2[i2]; ////cout << " pt_2:" << pt_2 << endl;
+ 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;
+
+
+ if (_rejectPairConversion)
+ {
+ float e1Sq = massElecSq + pt_1*pt_1 + pz_1*pz_1;
+ float e2Sq = massElecSq + pt_2*pt_2 + pz_2*pz_2;
+ 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;
+ ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
+ __n2_12 += corr;
+ __n2_12_vsEtaPhi[ij] += corr;
+ ptpt = pt_1*pt_2;
+ __s2ptpt_12 += corr*ptpt;
+ __s2PtN_12 += corr*pt_1;
+ __s2NPt_12 += corr*pt_2;
+ __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
+ __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
+ __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
+ __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
+ __n2Nw_12 += 1;
+ __s2ptptNw_12 += ptpt;
+ __s2PtNNw_12 += pt_1;
+ __s2NPtNw_12 += pt_2;
+
+ }
+ } //i2
+ } //i1
+ }
+
+ _n2_12_vsM->Fill(centrality, __n2_12);
+ _s2PtPt_12_vsM->Fill(centrality, __s2ptpt_12);
+ _s2PtN_12_vsM->Fill(centrality, __s2NPt_12);
+ _s2NPt_12_vsM->Fill(centrality, __s2PtN_12);
+
+ _n2Nw_12_vsM->Fill(centrality, __n2Nw_12);
+ _s2PtPtNw_12_vsM->Fill(centrality, __s2ptptNw_12);
+ _s2PtNNw_12_vsM->Fill(centrality, __s2NPtNw_12);
+ _s2NPtNw_12_vsM->Fill(centrality, __s2PtNNw_12);
+
}