// calculate same- and mixed-event correlations, and fill THnSparse
// output. -A. Adare, Apr 2011
-#include "TCanvas.h"
-#include "TChain.h"
-#include "TFormula.h"
-#include "TH1.h"
-#include "TH2.h"
-#include "TAxis.h"
-#include "TProfile2D.h"
-#include "TROOT.h"
-#include "TTree.h"
+#include <TAxis.h>
+#include <TCanvas.h>
+#include <TChain.h>
+#include <TFormula.h>
+#include <TH1.h>
+#include <TH2.h>
+#include <TH3.h>
+#include <THn.h>
+#include <TProfile2D.h>
+#include <TROOT.h>
+#include <TTree.h>
+#include "AliAnalysisUtils.h"
#include "AliAODEvent.h"
#include "AliAODInputHandler.h"
#include "AliAnalysisManager.h"
#include "AliDhcTask.h"
#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
+#include "AliESDMuonTrack.h"
#include "AliESDtrackCuts.h"
#include "AliPool.h"
#include "AliVParticle.h"
ClassImp(AliDhcTask)
//________________________________________________________________________
-AliDhcTask::AliDhcTask(const char *name)
-: AliAnalysisTaskSE(name), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15),
- fTrackDepth(1000), fPoolSize(200), fTracksName(), fDoWeights(0),
- fESD(0), fAOD(0), fOutputList(0), fHEvt(0), fHTrk(0), fHPtAss(0), fHPtTrg(0), fHPtTrg_Evt(0),
- fHCent(0), fHZvtx(0), fNbins(0), fHSs(0), fHMs(0), fHPts(0), fIndex(0),
+AliDhcTask::AliDhcTask()
+: AliAnalysisTaskSE(), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15),
+ fTrackDepth(1000), fPoolSize(200), fTracksName(), fDoWeights(kFALSE), fFillMuons(kFALSE),
+ fRequireMuon(kFALSE), fReqPtLo(0.0), fReqPtHi(1000.0),
+ fPtTACrit(kTRUE), fAllTAHists(kFALSE), fMixInEtaT(kFALSE),
+ fUseMuonUtils(kFALSE), fMuonCutMask(0), fMuonTrackCuts(0x0),
+ fEtaTLo(-1.0), fEtaTHi(1.0), fEtaALo(-1.0), fEtaAHi(1.0), fOmitFirstEv(kTRUE),
+ fDoFillSame(kFALSE), fDoMassCut(kFALSE), fCheckVertex(kTRUE), fClassName(),
+ fCentMethod("V0M"), fNBdeta(20), fNBdphi(36), fTriggerMatch(kTRUE),
+ fBPtT(0x0), fBPtA(0x0), fBCent(0x0), fBZvtx(0x0),
+ fMixBCent(0x0), fMixBZvtx(0x0), fHEffT(0x0), fHEffA(0x0),
+ fESD(0x0), fAOD(0x0), fOutputList(0x0), fHEvt(0x0), fHEvtWTr(0x0), fHTrk(0x0), fHPoolReady(0x0),
+ fHPtAss(0x0), fHPtTrg(0x0), fHPtTrgEvt(0x0),
+ fHPtTrgNorm1S(0x0), fHPtTrgNorm1M(0x0), fHPtTrgNorm2S(0x0), fHPtTrgNorm2M(0x0),
+ fHCent(0x0), fHZvtx(0x0), fNbins(0), fHSs(0x0), fHMs(0x0), fHPts(0x0), fHSMass(0x0), fHMMass(0x0),
+ fHQATp(0x0), fHQAAp(0x0), fHQATpCorr(0x0), fHQAApCorr(0x0),
+ fHQATm(0x0), fHQAAm(0x0), fHQATmCorr(0x0), fHQAAmCorr(0x0),
+ fHPtCentT(0x0), fHPtCentA(0x0), fIndex(0x0),
fCentrality(99), fZVertex(99), fEsdTPCOnly(0), fPoolMgr(0),
- fCentMethod("V0M"), fNBdeta(20), fNBdphi(36),
+ fUtils(0x0)
+{
+ // Constructor
+}
+
+//________________________________________________________________________
+AliDhcTask::AliDhcTask(const char *name, Bool_t def)
+: AliAnalysisTaskSE(name), fVerbosity(0), fEtaMax(1), fZVtxMax(10), fPtMin(0.25), fPtMax(15),
+ fTrackDepth(1000), fPoolSize(200), fTracksName(), fDoWeights(kFALSE), fFillMuons(kFALSE),
+ fRequireMuon(kFALSE), fReqPtLo(0.0), fReqPtHi(1000.0),
+ fPtTACrit(kTRUE), fAllTAHists(kFALSE), fMixInEtaT(kFALSE),
+ fUseMuonUtils(kFALSE), fMuonCutMask(0), fMuonTrackCuts(0x0),
+ fEtaTLo(-1.0), fEtaTHi(1.0), fEtaALo(-1.0), fEtaAHi(1.0), fOmitFirstEv(kTRUE),
+ fDoFillSame(kFALSE), fDoMassCut(kFALSE), fCheckVertex(kTRUE), fClassName(),
+ fCentMethod("V0M"), fNBdeta(20), fNBdphi(36), fTriggerMatch(kTRUE),
fBPtT(0x0), fBPtA(0x0), fBCent(0x0), fBZvtx(0x0),
- fMixBCent(0x0), fMixBZvtx(0x0)
+ fMixBCent(0x0), fMixBZvtx(0x0), fHEffT(0x0), fHEffA(0x0),
+ fESD(0x0), fAOD(0x0), fOutputList(0x0), fHEvt(0x0), fHEvtWTr(0x0), fHTrk(0x0), fHPoolReady(0x0),
+ fHPtAss(0x0), fHPtTrg(0x0), fHPtTrgEvt(0x0),
+ fHPtTrgNorm1S(0x0), fHPtTrgNorm1M(0x0), fHPtTrgNorm2S(0x0), fHPtTrgNorm2M(0x0),
+ fHCent(0x0), fHZvtx(0x0), fNbins(0), fHSs(0x0), fHMs(0x0), fHPts(0x0), fHSMass(0x0), fHMMass(0x0),
+ fHQATp(0x0), fHQAAp(0x0), fHQATpCorr(0x0), fHQAApCorr(0x0),
+ fHQATm(0x0), fHQAAm(0x0), fHQATmCorr(0x0), fHQAAmCorr(0x0),
+ fHPtCentT(0x0), fHPtCentA(0x0), fIndex(0x0),
+ fCentrality(99), fZVertex(99), fEsdTPCOnly(0), fPoolMgr(0),
+ fUtils(0x0)
{
// Constructor
fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.,SPDVertex.,TPCVertex.,Tracks "
"AOD:header,tracks,vertices,";
- Double_t ptt[4] = {0.25, 1.0, 2.0, 15.0};
- fBPtT = new TAxis(3,ptt);
- Double_t pta[4] = {0.25, 1.0, 2.0, 15.0};
- fBPtA = new TAxis(3,pta);
- Double_t cent[2] = {-100.0, 100.0};
- fBCent = new TAxis(1,cent);
- Double_t zvtx[2] = {-10, 10};
- fBZvtx = new TAxis(1,zvtx);
- Double_t centmix[2] = {-100.0, 100.0};
- fMixBCent = new TAxis(1,centmix);
- Double_t zvtxmix[9] = {-10,-6,-4,-2,0,2,4,6,10};
- fMixBZvtx = new TAxis(8,zvtxmix);
+ if (def) {
+ Double_t ptt[4] = {0.25, 1.0, 2.0, 15.0};
+ fBPtT = new TAxis(3,ptt);
+ Double_t pta[4] = {0.25, 1.0, 2.0, 15.0};
+ fBPtA = new TAxis(3,pta);
+ Double_t cent[2] = {-100.0, 100.0};
+ fBCent = new TAxis(1,cent);
+ Double_t zvtx[2] = {-10, 10};
+ fBZvtx = new TAxis(1,zvtx);
+ Double_t centmix[2] = {-100.0, 100.0};
+ fMixBCent = new TAxis(1,centmix);
+ Double_t zvtxmix[9] = {-10,-6,-4,-2,0,2,4,6,10};
+ fMixBZvtx = new TAxis(8,zvtxmix);
+ }
}
//________________________________________________________________________
{
// Create histograms
// Called once (per slave on PROOF!)
+ PrintDhcSettings();
fOutputList = new TList();
fOutputList->SetOwner(1);
- fEsdTPCOnly = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
- //fEsdTPCOnly->SetMinNClustersTPC(70);
- fEsdTPCOnly->SetMinNCrossedRowsTPC(70);
- fEsdTPCOnly->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
+ fUtils = new AliAnalysisUtils();
+ if (fUseMuonUtils) {
+ fMuonTrackCuts = new AliMuonTrackCuts("StdMuonCuts","StdMuonCuts");
+ fMuonTrackCuts->SetCustomParamFromRun(197388,"muon_pass2"); // for LHC13b,c,d,e,f ?
+ fMuonTrackCuts->SetFilterMask(fMuonCutMask);
+ AliInfo(Form(" using muon track cuts with bit %u\n", fMuonCutMask));
+ }
BookHistos();
InitEventMixer();
PostData(1, fOutputList);
}
+//________________________________________________________________________
+void AliDhcTask::PrintDhcSettings()
+{
+ AliInfo(Form("Dhc Task %s settings",fName.Data()));
+ AliInfo(Form(" centrality estimator %s", fCentMethod.Data()));
+ AliInfo(Form(" using tracks named %s", fTracksName.Data()));
+ AliInfo(Form(" efficiency correct triggers? %d", fHEffT!=0));
+ if (fHEffT!=0) {
+ AliInfo(Form(" %d dimensions (t)", fHEffT->GetNdimensions()));
+ }
+ AliInfo(Form(" efficiency correct associates? %d", fHEffA!=0));
+ if (fHEffA!=0) {
+ AliInfo(Form(" %d dimensions (a)", fHEffA->GetNdimensions()));
+ }
+ AliInfo(Form(" fill muons? %d", fFillMuons));
+ AliInfo(Form(" require muon even if not filling them? %d", fRequireMuon));
+ if (fRequireMuon) AliInfo(Form(" require muon with %f < pt < %f", fReqPtLo, fReqPtHi));
+ AliInfo(Form(" use pTT > pTA criterion? %d", fPtTACrit));
+ AliInfo(Form(" create all pTT, pTA hists? %d", fAllTAHists));
+ AliInfo(Form(" Mix in eta_T bins instead of z_vertex? %d", fMixInEtaT));
+ AliInfo(Form(" trigger eta range %f .. %f", fEtaTLo, fEtaTHi));
+ AliInfo(Form(" associate eta range %f .. %f", fEtaALo, fEtaAHi));
+ AliInfo(Form(" fill same event in any case? %d", fDoFillSame));
+ AliInfo(Form(" do invariant mass cut? %d", fDoMassCut));
+ AliInfo(Form(" omit first event? %d", fOmitFirstEv));
+ AliInfo(Form(" check the vertex? %d", fCheckVertex));
+ AliInfo(Form(" use the muon PWG muon cuts? %d", fUseMuonUtils));
+}
+
//________________________________________________________________________
void AliDhcTask::BookHistos()
{
for (Int_t i=1; i<=nZvtx; i++) {
zvtx[i] = fBZvtx->GetBinUpEdge(i);
}
-
+
fHEvt = new TH2F("fHEvt", "Event-level variables; Zvtx; Cent", nZvtx, zvtx, nCent, cent);
fOutputList->Add(fHEvt);
- fHTrk = new TH2F("fHTrk", "Track-level variables",
- 100, 0, TMath::TwoPi(), 100, -fEtaMax, +fEtaMax);
+ fHEvtWTr = new TH2F("fHEvtWTr", "Events with tracks; Zvtx; Cent", nZvtx, zvtx, nCent, cent);
+ fOutputList->Add(fHEvtWTr);
+ fHTrk = new TH2F("fHTrk", "Track-level variables",
+ 100, 0, TMath::TwoPi(), 100, -fEtaMax, +fEtaMax);
fOutputList->Add(fHTrk);
+ // Centrality mixing axis
+ Int_t nCentBins=fMixBCent->GetNbins();
+ Double_t centBins[nCentBins+1];
+ centBins[0] = fMixBCent->GetBinLowEdge(1);
+ for (Int_t i=1; i<=nCentBins; i++) {
+ centBins[i] = fMixBCent->GetBinUpEdge(i);
+ }
+ // Z-vertex mixing axis
+ Int_t nZvtxBins=fMixBZvtx->GetNbins();
+ Double_t zvtxbin[nZvtxBins+1];
+ zvtxbin[0] = fMixBZvtx->GetBinLowEdge(1);
+ for (Int_t i=1; i<=nZvtxBins; i++) {
+ zvtxbin[i] = fMixBZvtx->GetBinUpEdge(i);
+ }
+ fHPoolReady = new TH2F("fHPoolReady","mixing started", nZvtxBins, zvtxbin, nCentBins, centBins);
+ fOutputList->Add(fHPoolReady);
+
fHPtAss = new TH1F("fHPtAss","PtAssoc;P_{T} (GeV/c) [GeV/c]",nPtAssc,pta);
fOutputList->Add(fHPtAss);
fHPtTrg = new TH1F("fHPtTrg","PtTrig;P_{T} (GeV/c) [GeV/c]",nPtTrig,ptt);
fOutputList->Add(fHPtTrg);
- fHPtTrg_Evt = (TH1*) fHPtTrg->Clone("fHPtTrg_Evt");
+ fHPtTrgEvt = new TH1F("fHPtTrgEvt","PtTrig;P_{T} (GeV/c) [GeV/c]",nPtTrig,ptt);
+ fHPtTrgNorm1S = new TH3F("fHPtTrgNorm1S","PtTrig;P_{T} (GeV/c) [GeV/c];centrality;z_{vtx}",nPtTrig,ptt,nCent,cent,nZvtx,zvtx);
+ fOutputList->Add(fHPtTrgNorm1S);
+ fHPtTrgNorm1M = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm1M");
+ fOutputList->Add(fHPtTrgNorm1M);
+ fHPtTrgNorm2S = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm2S");
+ fOutputList->Add(fHPtTrgNorm2S);
+ fHPtTrgNorm2M = (TH3*) fHPtTrgNorm1S->Clone("fHPtTrgNorm2M");
+ fOutputList->Add(fHPtTrgNorm2M);
+
fHCent = new TH1F("fHCent","Cent;bins",nCent,cent);
fOutputList->Add(fHCent);
fHZvtx = new TH1F("fHZvtx","Zvertex;bins",nZvtx,zvtx);
fOutputList->Add(fHZvtx);
-
- fNbins = nPtTrig*nPtAssc*nCent*nZvtx;
- fHSs = new TH2*[fNbins];
- fHMs = new TH2*[fNbins];
- fHPts = new TH1*[fNbins];
-
- fIndex = new TFormula("GlobIndex","(t-1)*[0]*[1]*[2]+(z-1)*[0]*[1]+(x-1)*[0]+(y-1)+0*[4]");
- fIndex->SetParameters(nPtTrig,nPtAssc,nZvtx,nCent);
- fIndex->SetParNames("NTrigBins","NAssocBins", "NZvertexBins", "NCentBins");
- //fOutputList->Add(fIndex);
-
+
+ fHQATp = new TH3F("fHQATp","QA trigger;p_{T} (GeV/c);#eta;#phi",
+ 100,0.0,10.0,
+ 40,fEtaTLo,fEtaTHi,
+ 36,0.0,TMath::TwoPi());
+ fOutputList->Add(fHQATp);
+ fHQAAp = new TH3F("fHQAAp","QA associated;p_{T} (GeV/c);#eta;#phi",
+ 100,0.0,10.0,
+ 40,fEtaALo,fEtaAHi,
+ 36,0.0,TMath::TwoPi());
+ fOutputList->Add(fHQAAp);
+ fHQATpCorr = (TH3 *) fHQATp->Clone("fHQATpCorr");
+ fOutputList->Add(fHQATpCorr);
+ fHQAApCorr = (TH3 *) fHQAAp->Clone("fHQAApCorr");
+ fOutputList->Add(fHQAApCorr);
+ fHQATm = (TH3 *) fHQATp->Clone("fHQATm");
+ fOutputList->Add(fHQATm);
+ fHQAAm = (TH3 *) fHQAAp->Clone("fHQAAm");
+ fOutputList->Add(fHQAAm);
+ fHQATmCorr = (TH3 *) fHQATm->Clone("fHQATmCorr");
+ fOutputList->Add(fHQATmCorr);
+ fHQAAmCorr = (TH3 *) fHQAAm->Clone("fHQAAmCorr");
+ fOutputList->Add(fHQAAmCorr);
+
+ fHPtCentT = new TH2F("fHPtCentT",Form("trigger particles;p_{T} (GeV/c);centrality (%s)",fCentMethod.Data()),
+ 100,0.0,10.0,
+ 100,cent[0],cent[nCent]);
+ fOutputList->Add(fHPtCentT);
+ fHPtCentA = new TH2F("fHPtCentA",Form("associated particles;p_{T} (GeV/c);centrality (%s)",fCentMethod.Data()),
+ 100,0.0,10.0,
+ 100,cent[0],cent[nCent]);
+ fOutputList->Add(fHPtCentA);
+
+ fNbins = nPtTrig*nPtAssc*nCent*nZvtx;
+ fHSs = new TH2*[fNbins];
+ fHMs = new TH2*[fNbins];
+ fHPts = new TH1*[fNbins];
+ fHSMass = new TH1*[fNbins];
+ fHMMass = new TH1*[fNbins];
+
+ fIndex = new TFormula("GlobIndex","(t-1)*[0]*[1]*[2]+(z-1)*[0]*[1]+(x-1)*[1]+(y-1)");
+ fIndex->SetParameters(nPtTrig,nPtAssc,nZvtx);
+ fIndex->SetParNames("NTrigBins","NAssocBins","NZvertexBins");
+ fOutputList->Add(fIndex);
+
+ Double_t dEtaMin = fEtaTLo - fEtaAHi;
+ Double_t dEtaMax = fEtaTHi - fEtaALo;
+
Int_t count = 0;
for (Int_t c=1; c<=nCent; ++c) {
for (Int_t z=1; z<=nZvtx; ++z) {
for (Int_t t=1; t<=nPtTrig; ++t) {
- for (Int_t a=1; a<=nPtAssc; ++a) {
- fHSs[count] = 0;
- fHMs[count] = 0;
- fHPts[count] = 0;
- if (a>t) {
- ++count;
- continue;
- }
+ for (Int_t a=1; a<=nPtAssc; ++a) {
+ fHSs[count] = 0;
+ fHMs[count] = 0;
+ fHPts[count] = 0;
+ fHSMass[count] = 0;
+ fHMMass[count] = 0;
+ if ((a>t)&&!fAllTAHists) {
+ ++count;
+ continue;
+ }
+ if (z==1 && t==1 && a==1) {
+ TString title(Form("cen=%d (%.1f to %.1f)",
+ c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c)));
+ fHSMass[count] = new TH1F(Form("hSMass%d",count), Form("Mass Same Event %s;m (GeV)",title.Data()), 10000, 0,10);
+ fOutputList->Add(fHSMass[count]);
+ fHMMass[count] = new TH1F(Form("hMMass%d",count), Form("Mass Mixed Event %s;m (GeV)",title.Data()), 10000, 0,10);
+ fOutputList->Add(fHMMass[count]);
+ }
if (t==1 && a==1) {
TString title(Form("cen=%d (%.1f to %.1f), zVtx=%d (%.1f to %.1f)",
- c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c),
+ c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c),
z, fBZvtx->GetBinLowEdge(z), fBZvtx->GetBinUpEdge(z)));
fHPts[count] = new TH1F(Form("hPt%d",count), Form("Ptdist %s;p_{T} (GeV/c)",title.Data()), 200,0,20);
fOutputList->Add(fHPts[count]);
}
- TString title(Form("cen=%d (%.1f to %.1f), zVtx=%d (%.1f to %.1f), trig=%d (%.1f to %.1f), assoc=%d (%.1f to %.1f)",
- c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c),
+ TString title(Form("cen=%d (%.1f to %.1f), zVtx=%d (%.1f to %.1f), trig=%d (%.1f to %.1f), assoc=%d (%.1f to %.1f)",
+ c, fBCent->GetBinLowEdge(c), fBCent->GetBinUpEdge(c),
z, fBZvtx->GetBinLowEdge(z), fBZvtx->GetBinUpEdge(z),
- t, fBPtT->GetBinLowEdge(t), fBPtT->GetBinUpEdge(t),
+ t, fBPtT->GetBinLowEdge(t), fBPtT->GetBinUpEdge(t),
a, fBPtA->GetBinLowEdge(a), fBPtA->GetBinUpEdge(a)));
- fHSs[count] = new TH2F(Form("hS%d",count), Form("Signal %s;#Delta #varphi;#Delta #eta",title.Data()),
- fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,-2*fEtaMax,2*fEtaMax);
- fHMs[count] = new TH2F(Form("hM%d",count), Form("Mixed %s;#Delta #varphi;#Delta #eta",title.Data()),
- fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,-2*fEtaMax,2*fEtaMax);
- fOutputList->Add(fHSs[count]);
- fOutputList->Add(fHMs[count]);
+ fHSs[count] = new TH2F(Form("hS%d",count), Form("Signal %s;#Delta #varphi;#Delta #eta",title.Data()),
+ fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,dEtaMin,dEtaMax);
+ fHMs[count] = new TH2F(Form("hM%d",count), Form("Mixed %s;#Delta #varphi;#Delta #eta",title.Data()),
+ fNBdphi,-0.5*TMath::Pi(),1.5*TMath::Pi(),fNBdeta,dEtaMin,dEtaMax);
+ fOutputList->Add(fHSs[count]);
+ fOutputList->Add(fHMs[count]);
Int_t tcount = (Int_t)fIndex->Eval(t,a,z,c);
- if (fVerbosity>5)
- cout << count << " " << tcount << ": " << title << endl;
+ if (fVerbosity>5)
+ cout << count << " " << tcount << ": " << title << endl;
if (count != tcount)
AliFatal(Form("Index mismatch: %d %d", count, tcount));
- ++count;
- }
+ ++count;
+ }
}
}
}
}
+//________________________________________________________________________
+void AliDhcTask::SetAnaMode(Int_t iAna)
+{
+ if (iAna==kHH) {
+ SetFillMuons(kFALSE);
+ fEtaTLo = -1.6;
+ fEtaTHi = +1.6;
+ fEtaALo = -1.6;
+ fEtaAHi = +1.6;
+ } else if (iAna==kMuH) {
+ SetFillMuons(kTRUE);
+ fEtaTLo = -5.0;
+ fEtaTHi = -2.0;
+ fEtaALo = -1.6;
+ fEtaAHi = +1.6;
+ } else if (iAna==kHMu) {
+ SetFillMuons(kTRUE);
+ fEtaTLo = -1.6;
+ fEtaTHi = +1.6;
+ fEtaALo = -5.0;
+ fEtaAHi = -2.0;
+ } else if (iAna==kMuMu) {
+ SetFillMuons(kTRUE);
+ fEtaTLo = -5.0;
+ fEtaTHi = -2.0;
+ fEtaALo = -5.0;
+ fEtaAHi = -2.0;
+ } else if (iAna==kPSide) {
+ SetFillMuons(kFALSE);
+ fEtaTLo = -1.6;
+ fEtaTHi = -0.465;
+ fEtaALo = -1.6;
+ fEtaAHi = +1.6;
+ } else if (iAna==kASide) {
+ SetFillMuons(kFALSE);
+ fEtaTLo = +0.465;
+ fEtaTHi = +1.6;
+ fEtaALo = -1.6;
+ fEtaAHi = +1.6;
+ } else {
+ AliInfo(Form("Unrecognized analysis option: %d", iAna));
+ }
+}
+
//________________________________________________________________________
void AliDhcTask::InitEventMixer()
{
{
// Main loop, called for each event.
- static int iEvent = -1; ++iEvent;
-
- if (fVerbosity>2) {
- if (iEvent % 10 == 0)
- cout << iEvent << endl;
- }
-
- Int_t dType = -1; // Will be set to kESD or kAOD.
- MiniEvent* sTracks = new MiniEvent(0); // deleted by pool mgr.
+ if (fVerbosity > 10)
+ AliInfo(Form("======= Dhc Task %s start next event",fName.Data()));
LoadBranches();
+ if (fOmitFirstEv) {
+ if (fUtils->IsFirstEventInChunk(InputEvent()))
+ return;
+ }
+
// Get event pointers, check for signs of life
+ Int_t dType = -1; // Will be set to kESD or kAOD.
fESD = dynamic_cast<AliESDEvent*>(InputEvent());
fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
- if (fESD)
+ if (fESD) {
dType = kESD;
- else if (fAOD)
+ } else if (fAOD) {
dType = kAOD;
- else {
+ } else {
AliError("Neither fESD nor fAOD available");
return;
}
+ // select specific trigger classes?
+ if (fClassName.Length()>0) {
+ TString strFiredClass;
+ if (fESD)
+ strFiredClass = fESD->GetFiredTriggerClasses();
+ else
+ strFiredClass = fAOD->GetFiredTriggerClasses();
+
+ if (fVerbosity > 10) {
+ AliInfo(Form("Trigger class selection: This event has classes %s", strFiredClass.Data()));
+ AliInfo(Form("Trigger class selection: selecting classes %s", fClassName.Data()));
+ }
+
+ TObjArray *arrClass = fClassName.Tokenize(",");
+ Int_t nClass = arrClass->GetEntries();
+
+ TString strOneClass;
+ Bool_t bAccEvent = kFALSE;
+ for (Int_t iClass=0; iClass<nClass; iClass++) {
+ strOneClass = arrClass->At(iClass)->GetName();
+ if (strFiredClass.Contains(strOneClass))
+ bAccEvent = kTRUE;
+ }
+
+ if (!bAccEvent)
+ return;
+
+ if (fVerbosity > 10)
+ AliInfo(Form("Passed trigger class selection, this event has classes %s", strFiredClass.Data()));
+ }
+
Bool_t mcgen = 0;
if (fTracksName.Contains("Gen"))
mcgen = 1;
if (mcgen) {
fZVertex = 0;
TList *list = InputEvent()->GetList();
- TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
+ TClonesArray *tcaTracks = 0;
+ if (list)
+ tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
if (tcaTracks)
fCentrality = tcaTracks->GetEntries();
} else {
if (dType == kESD) {
- if (!VertexOk(fESD)) {
+ const AliESDVertex* vertex = fESD->GetPrimaryVertex();
+ fZVertex = vertex->GetZ();
+ if (fCheckVertex && !VertexOk()) {
if (fVerbosity > 1)
AliInfo(Form("Event REJECTED (ESD vertex not OK). z = %.1f", fZVertex));
return;
}
- const AliESDVertex* vertex = fESD->GetPrimaryVertex();
- fZVertex = vertex->GetZ();
if(fESD->GetCentrality()) {
fCentrality =
fESD->GetCentrality()->GetCentralityPercentile(fCentMethod);
}
- }
- if (dType == kAOD) {
+ } else if (dType == kAOD) {
const AliAODVertex* vertex = fAOD->GetPrimaryVertex();
fZVertex = vertex->GetZ();
- if (!VertexOk(fAOD)) {
+ if (fCheckVertex && !VertexOk()) {
if (fVerbosity > 1)
Info("Exec()", "Event REJECTED (AOD vertex not OK). z = %.1f", fZVertex);
return;
}
- const AliCentrality *aodCent = fAOD->GetHeader()->GetCentralityP();
+ AliAODHeader * header = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
+ if(!header) AliFatal("Not a standard AOD");
+
+ const AliCentrality *aodCent = header->GetCentralityP();
if (aodCent) {
fCentrality = aodCent->GetCentralityPercentile(fCentMethod);
}
}
}
- // Fill Event histogram
+ // Fill event histogram
fHEvt->Fill(fZVertex, fCentrality);
if (fCentrality > fBCent->GetXmax() || fCentrality < fBCent->GetXmin()) {
- AliInfo(Form("Event REJECTED (centrality out of range). fCentrality = %.1f", fCentrality));
+ if (fVerbosity > 1)
+ AliInfo(Form("Event REJECTED (centrality out of range). fCentrality = %.1f", fCentrality));
return;
}
-
- // Get array of selected tracks
- if (dType == kESD) {
- GetESDTracks(sTracks);
+ if (fZVertex > fBZvtx->GetXmax() || fZVertex < fBZvtx->GetXmin()) {
+ if (fVerbosity > 1)
+ AliInfo(Form("Event REJECTED (z_vertex out of range). fZVertex = %.1f", fZVertex));
+ return;
}
- if (dType == kAOD) {
- GetAODTracks(sTracks);
+
+ // reject events without muon if required
+ if (fRequireMuon) {
+ Bool_t bHasMuon = kFALSE;
+ if (dType == kESD) {
+ bHasMuon = HasMuonESD();
+ } else if (dType == kAOD) {
+ bHasMuon = HasMuonAOD();
+ }
+ if (!bHasMuon) {
+ if (fVerbosity > 1)
+ AliInfo(Form("Event REJECTED (no muon). fCentrality = %.1f", fCentrality));
+ return;
+ }
}
// Get pool containing tracks from other events like this one
return;
}
+ if (fVerbosity > 10)
+ AliInfo("--- prepare to get tracks");
+
+ MiniEvent* sTracks = new MiniEvent(0); // deleted by pool mgr.
+ if (dType == kESD) {
+ GetESDTracks(sTracks);
+ } else if (dType == kAOD) {
+ GetAODTracks(sTracks);
+ }
+
+ if (fVerbosity > 10)
+ AliInfo(Form("--- got a track array with %lu tracks",sTracks->size()));
+
+ if (sTracks->size()==0) {
+ AliWarning(Form("Track array empty"));
+ delete sTracks;
+ return;
+ }
+
+ fHEvtWTr->Fill(fZVertex, fCentrality);
+
if (!pool->IsReady()) {
pool->UpdatePool(sTracks);
+ if (fDoFillSame) { // fill same event right away if requested
+ Correlate(*sTracks, *sTracks, kSameEvt);
+ }
+ PostData(1, fOutputList);
return;
}
if (pool->IsFirstReady()) {
+ fHPoolReady->Fill(fZVertex,fCentrality);
// recover events missed before the pool is ready
Int_t nEvs = pool->GetCurrentNEvents();
if (nEvs>1) {
for (Int_t i=0; i<nEvs; ++i) {
- MiniEvent* evI = pool->GetEvent(i);
- for (Int_t j=0; j<nEvs; ++j) {
- MiniEvent* evJ = pool->GetEvent(j);
- if (i==j) {
- Correlate(*evI, *evJ, kSameEvt);
- } else {
- Correlate(*evI, *evJ, kDiffEvt);
- }
- }
+ MiniEvent* evI = pool->GetEvent(i);
+ for (Int_t j=0; j<nEvs; ++j) {
+ MiniEvent* evJ = pool->GetEvent(j);
+ if (i==j) {
+ if (!fDoFillSame) { // do not fill the same events twice
+ Correlate(*evI, *evJ, kSameEvt);
+ }
+ } else {
+ Correlate(*evI, *evJ, kDiffEvt);
+ }
+ }
}
}
} else { /* standard case: same event, then mix*/
//________________________________________________________________________
void AliDhcTask::GetESDTracks(MiniEvent* miniEvt)
{
+
// Loop twice: 1. Count sel. tracks. 2. Fill vector.
-
if (fTracksName.IsNull()) {
const AliESDVertex *vtxSPD = fESD->GetPrimaryVertexSPD();
if (!vtxSPD)
TObjArray arr(nTrax);
arr.SetOwner(1);
+ if (!fEsdTPCOnly) {
+ fEsdTPCOnly = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+ //fEsdTPCOnly->SetMinNClustersTPC(70);
+ fEsdTPCOnly->SetMinNCrossedRowsTPC(70);
+ fEsdTPCOnly->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
+ }
+
for (Int_t i = 0; i < nTrax; ++i) {
AliESDtrack* esdtrack = fESD->GetTrack(i);
if (!esdtrack) {
arr.Add(newtrack);
nSelTrax++;
}
- return;
- }
- TList *list = InputEvent()->GetList();
- TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
- const Int_t ntracks = tcaTracks->GetEntries();
- if (miniEvt)
- miniEvt->reserve(ntracks);
- else {
- AliError("Ptr to miniEvt zero");
- return;
+ for(Int_t itrack = 0; itrack < nSelTrax; itrack++) {
+ AliVTrack *esdtrack = static_cast<AliESDtrack*>(arr.At(itrack));
+ if(!esdtrack) {
+ AliError(Form("ERROR: Could not retrieve esdtrack %d",itrack));
+ continue;
+ }
+ Double_t pt = esdtrack->Pt();
+ Double_t eta = esdtrack->Eta();
+ Double_t phi = esdtrack->Phi();
+ Int_t sign = esdtrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
+ }
+ } else {
+ TList *list = InputEvent()->GetList();
+ TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
+
+ if (!tcaTracks){
+ AliError("Ptr to tcaTracks zero");
+ return;
+ }
+
+ const Int_t ntracks = tcaTracks->GetEntries();
+ Int_t nGoodTracks = 0;
+ // count good tracks
+ for (Int_t itrack = 0; itrack < ntracks; itrack++) {
+ AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack));
+ if (!vtrack) {
+ AliError(Form("ERROR: Could not retrieve vtrack %d",itrack));
+ continue;
+ }
+ Double_t pt = vtrack->Pt();
+ Double_t eta = vtrack->Eta();
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ nGoodTracks++;
+ }
+ if (miniEvt)
+ miniEvt->reserve(nGoodTracks);
+ else {
+ AliError("Ptr to miniEvt zero");
+ return;
+ }
+ // fill good tracks into minievent
+ for (Int_t itrack = 0; itrack < ntracks; itrack++) {
+ AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack));
+ if (!vtrack) {
+ AliError(Form("ERROR: Could not retrieve vtrack %d",itrack));
+ continue;
+ }
+ Double_t pt = vtrack->Pt();
+ Double_t eta = vtrack->Eta();
+ Double_t phi = vtrack->Phi();
+ Int_t sign = vtrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
+ }
}
- for(Int_t itrack = 0; itrack < ntracks; itrack++) {
- AliVTrack *esdtrack = static_cast<AliESDtrack*>(tcaTracks->At(itrack));
- if(!esdtrack) {
- AliError(Form("ERROR: Could not retrieve esdtrack %d",itrack));
- continue;
+ // count and fill muons if required
+ if (fFillMuons) {
+ // count muons
+ Int_t nGoodMuons = 0;
+ for (Int_t iMu = 0; iMu<fESD->GetNumberOfMuonTracks(); iMu++) {
+ AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iMu);
+ if (muonTrack) {
+ if ( !IsGoodMUONtrack(*muonTrack) ) continue;
+ Double_t ptMu = muonTrack->Pt();
+ Double_t etaMu = muonTrack->Eta();
+ if (IsTrigger(etaMu,ptMu)||IsAssociated(etaMu,ptMu))
+ nGoodMuons++;
+ }
+ }
+ // fill them into the mini event
+ miniEvt->reserve(miniEvt->size()+nGoodMuons);
+ for (Int_t iMu = 0; iMu<fESD->GetNumberOfMuonTracks(); iMu++) {
+ AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iMu);
+ if (muonTrack) {
+ if ( !IsGoodMUONtrack(*muonTrack) ) continue;
+ Double_t ptMu = muonTrack->Pt();
+ Double_t etaMu = muonTrack->Eta();
+ Double_t phiMu = muonTrack->Phi();
+ Int_t signMu = muonTrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(etaMu,ptMu)||IsAssociated(etaMu,ptMu))
+ miniEvt->push_back(AliMiniTrack(ptMu, etaMu, phiMu, signMu));
+ }
}
- Double_t pt = esdtrack->Pt();
- Double_t eta = esdtrack->Eta();
- Double_t phi = esdtrack->Phi();
- Int_t sign = esdtrack->Charge() > 0 ? 1 : -1;
- miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
}
+
}
//________________________________________________________________________
void AliDhcTask::GetAODTracks(MiniEvent* miniEvt)
{
+
// Loop twice: 1. Count sel. tracks. 2. Fill vector.
+ if (fTracksName.IsNull()) {
+ Int_t nTrax = fAOD->GetNumberOfTracks();
+ Int_t nSelTrax = 0;
- Int_t nTrax = fAOD->GetNumberOfTracks();
- Int_t nSelTrax = 0;
-
- if (fVerbosity > 2)
- AliInfo(Form("%d tracks in event",nTrax));
+ if (fVerbosity > 2)
+ AliInfo(Form("%d tracks in event",nTrax));
- // Loop 1.
- for (Int_t i = 0; i < nTrax; ++i) {
- AliAODTrack* aodtrack = fAOD->GetTrack(i);
- if (!aodtrack) {
- AliError(Form("Couldn't get AOD track %d\n", i));
- continue;
+ // Loop 1.
+ for (Int_t i = 0; i < nTrax; ++i) {
+ AliAODTrack* aodtrack = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(i));
+ if(!aodtrack) AliFatal("Not a standard AOD");
+ if (!aodtrack) {
+ AliError(Form("Couldn't get AOD track %d\n", i));
+ continue;
+ }
+ // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
+ UInt_t tpcOnly = 1 << 7;
+ Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
+ if (!trkOK)
+ continue;
+ Double_t pt = aodtrack->Pt();
+ Double_t eta = aodtrack->Eta();
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ nSelTrax++;
}
- // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
- UInt_t tpcOnly = 1 << 7;
- Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
- if (!trkOK)
- continue;
- Double_t pt = aodtrack->Pt();
- Bool_t ptOK = pt >= fPtMin && pt < fPtMax;
- if (!ptOK)
- continue;
- Double_t eta = aodtrack->Eta();
- if (TMath::Abs(eta) > fEtaMax)
- continue;
- nSelTrax++;
- }
- if (miniEvt)
- miniEvt->reserve(nSelTrax);
- else {
- AliError("!miniEvt");
- return;
- }
-
- // Loop 2.
- for (Int_t i = 0; i < nTrax; ++i) {
- AliAODTrack* aodtrack = fAOD->GetTrack(i);
- if (!aodtrack) {
- AliError(Form("Couldn't get AOD track %d\n", i));
- continue;
+ if (miniEvt)
+ miniEvt->reserve(nSelTrax);
+ else {
+ AliError("!miniEvt");
+ return;
}
+
+ // Loop 2.
+ for (Int_t i = 0; i < nTrax; ++i) {
+ AliAODTrack* aodtrack = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(i));
+ if(!aodtrack) AliFatal("Not a standard AOD");
+ if (!aodtrack) {
+ AliError(Form("Couldn't get AOD track %d\n", i));
+ continue;
+ }
- // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
- UInt_t tpcOnly = 1 << 7;
- Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
- if (!trkOK)
- continue;
- Double_t pt = aodtrack->Pt();
- Bool_t ptOK = pt >= fPtMin && pt < fPtMax;
- if (!ptOK)
- continue;
- Double_t eta = aodtrack->Eta();
- if (TMath::Abs(eta) > fEtaMax)
- continue;
+ // See $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C
+ UInt_t tpcOnly = 1 << 7;
+ Bool_t trkOK = aodtrack->TestFilterBit(tpcOnly);
+ if (!trkOK)
+ continue;
+ Double_t pt = aodtrack->Pt();
+ Double_t eta = aodtrack->Eta();
+ Double_t phi = aodtrack->Phi();
+ Int_t sign = aodtrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
+ }
+ } else {
+ TList *list = InputEvent()->GetList();
+ TClonesArray *tcaTracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
+
+ if (!tcaTracks){
+ AliError("Ptr to tcaTracks zero");
+ return;
+ }
+
+ const Int_t ntracks = tcaTracks->GetEntries();
+ Int_t nGoodTracks = 0;
+ // count good tracks
+ for (Int_t itrack = 0; itrack < ntracks; itrack++) {
+ AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack));
+ if (!vtrack) {
+ AliError(Form("ERROR: Could not retrieve vtrack %d",itrack));
+ continue;
+ }
+ Double_t pt = vtrack->Pt();
+ Double_t eta = vtrack->Eta();
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ nGoodTracks++;
+ }
+ if (miniEvt)
+ miniEvt->reserve(nGoodTracks);
+ else {
+ AliError("Ptr to miniEvt zero");
+ return;
+ }
+ // fill good tracks into minievent
+ for (Int_t itrack = 0; itrack < ntracks; itrack++) {
+ AliVTrack *vtrack = static_cast<AliVTrack*>(tcaTracks->At(itrack));
+ if (!vtrack) {
+ AliError(Form("ERROR: Could not retrieve vtrack %d",itrack));
+ continue;
+ }
+ Double_t pt = vtrack->Pt();
+ Double_t eta = vtrack->Eta();
+ Double_t phi = vtrack->Phi();
+ Int_t sign = vtrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(eta,pt)||IsAssociated(eta,pt))
+ miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
+ }
+ }
- Double_t phi = aodtrack->Phi();
- Int_t sign = aodtrack->Charge() > 0 ? 1 : -1;
- miniEvt->push_back(AliMiniTrack(pt, eta, phi, sign));
+ // count and fill muons if required
+ if (fFillMuons) {
+ // count muons
+ Int_t nGoodMuons = 0;
+ for (Int_t iMu = 0; iMu<fAOD->GetNumberOfTracks(); iMu++) {
+ AliAODTrack* muonTrack = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iMu));
+ if(!muonTrack) AliFatal("Not a standard AOD");
+ if (muonTrack) {
+ if ( !IsGoodMUONtrack(*muonTrack) ) continue;
+ Double_t ptMu = muonTrack->Pt();
+ Double_t etaMu = muonTrack->Eta();
+ if (IsTrigger(etaMu,ptMu)||IsAssociated(etaMu,ptMu))
+ nGoodMuons++;
+ }
+ }
+ // fill them into the mini event
+ miniEvt->reserve(miniEvt->size()+nGoodMuons);
+ for (Int_t iMu = 0; iMu<fAOD->GetNumberOfTracks(); iMu++) {
+ AliAODTrack* muonTrack = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iMu));
+ if(!muonTrack) AliFatal("Not a standard AOD");
+ if (muonTrack) {
+ if ( !IsGoodMUONtrack(*muonTrack) ) continue;
+ Double_t ptMu = muonTrack->Pt();
+ Double_t etaMu = muonTrack->Eta();
+ Double_t phiMu = muonTrack->Phi();
+ Int_t signMu = muonTrack->Charge() > 0 ? 1 : -1;
+ if (IsTrigger(etaMu,ptMu)||IsAssociated(etaMu,ptMu))
+ miniEvt->push_back(AliMiniTrack(ptMu, etaMu, phiMu, signMu));
+ }
+ }
}
+
}
//________________________________________________________________________
TH2 **hist = fHMs;
if (pairing == kSameEvt) {
hist = fHSs;
- fHCent->AddBinContent(cbin);
- fHZvtx->AddBinContent(zbin);
+ fHCent->Fill(fCentrality);
+ fHZvtx->Fill(fZVertex);
}
Int_t nZvtx = fHZvtx->GetNbinsX();
Int_t nPtAssc = fHPtAss->GetNbinsX();
Int_t globIndex = (cbin-1)*nZvtx*nPtTrig*nPtAssc+(zbin-1)*nPtTrig*nPtAssc;
+ Int_t ptindex = (Int_t)fIndex->Eval(1,1,zbin,cbin);
+ Int_t mindex = (Int_t)fIndex->Eval(1,1,1,cbin);
- Int_t ptindex = (Int_t)fIndex->Eval(1,1,zbin,cbin);
- fHPtTrg_Evt->Reset();
- if (fDoWeights) { // Count trigger particles in this event
- for (Int_t i=0; i<iMax; ++i) {
- const AliMiniTrack &a(evt1.at(i));
- Float_t pta = a.Pt();
- fHPtTrg_Evt->Fill(pta);
+ fHPtTrgEvt->Reset();
+ for (Int_t i=0; i<iMax; ++i) {
+ const AliMiniTrack &a(evt1.at(i));
+ Float_t pta = a.Pt();
+ fHPtTrgEvt->Fill(pta);
+ if (pairing == kSameEvt) {
fHPts[ptindex]->Fill(pta);
}
}
- for (Int_t i=0; i<iMax; ++i) {
+ Bool_t bCountTrg; // only count the trigger if an associated particle is found
+ for (Int_t i=0; i<iMax; ++i) {
// Trigger particles
const AliMiniTrack &a(evt1.at(i));
Float_t pta = a.Pt();
- Int_t abin = fHPtTrg->FindBin(pta);
- if (fHPtTrg->IsBinOverflow(abin) ||
- fHPtTrg->IsBinUnderflow(abin))
+ Float_t etaa = a.Eta();
+ Float_t phia = a.Phi();
+ Int_t sgna = a.Sign();
+
+ // brief intermezzo in the trigger particle loop: do associated particle QA here in order to avoid double counting
+ if (pairing == kSameEvt) {
+ if (IsAssociated(etaa,pta)) {
+ Double_t aQAWght = 1.0;
+ if (fHEffA) {
+ const Int_t nEffDimA = fHEffA->GetNdimensions();
+ Int_t effBinA[nEffDimA];
+ effBinA[0] = fHEffA->GetAxis(0)->FindBin(etaa);
+ effBinA[1] = fHEffA->GetAxis(1)->FindBin(pta);
+ effBinA[2] = fHEffA->GetAxis(2)->FindBin(fCentrality);
+ effBinA[3] = fHEffA->GetAxis(3)->FindBin(fZVertex);
+ if (nEffDimA>4) {
+ effBinA[4] = fHEffA->GetAxis(4)->FindBin(phia);
+ if (nEffDimA>5) {
+ effBinA[5] = fHEffA->GetAxis(5)->FindBin(sgna);
+ }
+ }
+ aQAWght = fHEffA->GetBinContent(effBinA);
+ }
+ // fill every associated particle once unweighted, once weighted
+ if (sgna>0.0) {
+ fHQAAp->Fill(pta,etaa,phia);
+ fHQAApCorr->Fill(pta,etaa,phia,aQAWght);
+ } else {
+ fHQAAm->Fill(pta,etaa,phia);
+ fHQAAmCorr->Fill(pta,etaa,phia,aQAWght);
+ }
+ fHPtCentA->Fill(pta,fCentrality);
+ }
+ }
+
+ // back to triggers
+ if (!IsTrigger(etaa,pta))
continue;
+ Int_t abin = fHPtAss->FindBin(pta);
+
+ // efficiency weighting
+ Double_t effWtT = 1.0;
+ if (fHEffT) {
+ // trigger particle
+ const Int_t nEffDimT = fHEffT->GetNdimensions();
+ Int_t effBinT[nEffDimT];
+ effBinT[0] = fHEffT->GetAxis(0)->FindBin(etaa);
+ effBinT[1] = fHEffT->GetAxis(1)->FindBin(pta);
+ effBinT[2] = fHEffT->GetAxis(2)->FindBin(fCentrality);
+ effBinT[3] = fHEffT->GetAxis(3)->FindBin(fZVertex);
+ if (nEffDimT>4) {
+ effBinT[4] = fHEffT->GetAxis(4)->FindBin(phia);
+ if (nEffDimT>5) {
+ effBinT[5] = fHEffT->GetAxis(5)->FindBin(sgna);
+ }
+ }
+ effWtT = fHEffT->GetBinContent(effBinT);
+ }
+
if (pairing == kSameEvt) {
- fHTrk->Fill(a.Phi(),a.Eta());
- fHPtTrg->AddBinContent(abin);
+ fHTrk->Fill(phia,etaa);
+ if (sgna>0.0) {
+ fHQATp->Fill(pta,etaa,phia);
+ fHQATpCorr->Fill(pta,etaa,phia,effWtT);
+ } else {
+ fHQATm->Fill(pta,etaa,phia);
+ fHQATmCorr->Fill(pta,etaa,phia,effWtT);
+ }
+ fHPtCentT->Fill(pta,fCentrality);
+ fHPtTrg->Fill(pta);
+ fHPtTrgNorm1S->Fill(pta,fCentrality,fZVertex,effWtT);
+ } else {
+ fHPtTrgNorm1M->Fill(pta,fCentrality,fZVertex,effWtT);
}
+ bCountTrg = kFALSE;
for (Int_t j=0; j<jMax; ++j) {
// Associated particles
if (pairing == kSameEvt && i==j)
- continue;
+ continue;
const AliMiniTrack &b(evt2.at(j));
Float_t ptb = b.Pt();
- if (pta < ptb)
- continue;
+ Float_t etab = b.Eta();
+ Float_t phib = b.Phi();
+ Int_t sgnb = b.Sign();
+ if (fPtTACrit&&(pta < ptb)) {
+ continue;
+ }
+
+ if (!IsAssociated(etab,ptb))
+ continue;
Int_t bbin = fHPtAss->FindBin(ptb);
- if (fHPtAss->IsBinOverflow(bbin) ||
- fHPtAss->IsBinUnderflow(bbin))
- continue;
- Float_t dphi = DeltaPhi(a.Phi(), b.Phi());
- Float_t deta = a.Eta() - b.Eta();
+ Float_t dphi = DeltaPhi(phia, phib);
+ Float_t deta = etaa - etab;
+ // invariant mass under mu mu assumption
+ Float_t muMass2 = 0.1056583715*0.1056583715;
+ Float_t ea2 = muMass2 + pta*pta*TMath::CosH(etaa)*TMath::CosH(etaa);
+ Float_t eb2 = muMass2 + ptb*ptb*TMath::CosH(etab)*TMath::CosH(etab);
+ Float_t papb = pta*TMath::Cos(phia)*ptb*TMath::Cos(phib) +
+ pta*TMath::Sin(phia)*ptb*TMath::Sin(phib) +
+ pta*TMath::SinH(etaa)*ptb*TMath::SinH(etab);
+ Float_t mass = TMath::Sqrt(2.0*(muMass2 + TMath::Sqrt(ea2*eb2) - papb));
+ Int_t q2 = sgna + sgnb;
+ if ((q2==0) && fDoMassCut) {
+ if (mass>3.0 && mass<3.2)
+ continue;
+ if (mass>9.2&&mass<9.8)
+ continue;
+ }
Int_t index = globIndex+(abin-1)*nPtAssc+(bbin-1);
Double_t weight = 1.0;
- if (fDoWeights)
- weight = fHPtTrg_Evt->GetBinContent(abin);
+ // number of trigger weight event-by-event (CMS method)
+ if (fDoWeights) {
+ Double_t nTrgs = fHPtTrgEvt->GetBinContent(abin);
+ if (nTrgs==0.0) {
+ AliError(Form("Trying to work with number of triggers weight = %f",nTrgs));
+ return 0;
+ }
+ weight *= 1./nTrgs;
+ }
- if (weight==0.0) {
- AliError(Form("Trying to work with weight = %f",weight));
- } else {
- hist[index]->Fill(dphi,deta,1./weight);
+ // efficiency weighting
+ if (fHEffT) {
+ // trigger particle
+ weight *= effWtT;
}
+ if (fHEffA) {
+ // associated particle
+ const Int_t nEffDimA = fHEffA->GetNdimensions();
+ Int_t effBinA[nEffDimA];
+ effBinA[0] = fHEffA->GetAxis(0)->FindBin(etab);
+ effBinA[1] = fHEffA->GetAxis(1)->FindBin(ptb);
+ effBinA[2] = fHEffA->GetAxis(2)->FindBin(fCentrality);
+ effBinA[3] = fHEffA->GetAxis(3)->FindBin(fZVertex);
+ if (nEffDimA>4) {
+ effBinA[4] = fHEffA->GetAxis(4)->FindBin(phib);
+ if (nEffDimA>5) {
+ effBinA[5] = fHEffA->GetAxis(5)->FindBin(sgnb);
+ }
+ }
+ weight *= fHEffA->GetBinContent(effBinA);
+ }
+
+ if (hist[index]) { // check if this histogram exists, relevant in the fPtTACrit==kFALSE case
+ hist[index]->Fill(dphi,deta,weight);
+ bCountTrg = kTRUE;
+ if (pairing == kSameEvt) {
+ fHPtAss->Fill(ptb); // fill every associated particle every time it is used
+ if (q2==0)
+ fHSMass[mindex]->Fill(mass);
+ } else {
+ if (q2==0)
+ fHMMass[mindex]->Fill(mass);
+ }
+ }
+ }
+ if (bCountTrg) {
if (pairing == kSameEvt) {
- fHPtAss->AddBinContent(bbin);
+ fHPtTrgNorm2S->Fill(pta,fCentrality,fZVertex,effWtT);
+ } else {
+ fHPtTrgNorm2M->Fill(pta,fCentrality,fZVertex,effWtT);
}
}
}
- return 0;
+ return 1;
+}
+
+
+//________________________________________________________________________
+Bool_t AliDhcTask::HasMuonESD()
+{
+ // does this ESD event have a good muon?
+ for (Int_t iMu = 0; iMu<fESD->GetNumberOfMuonTracks(); iMu++) {
+ AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iMu);
+ if (muonTrack) {
+ if ( IsGoodMUONtrack(*muonTrack) ) {
+ Double_t ptMu = muonTrack->Pt();
+ if (ptMu>fReqPtLo && ptMu<fReqPtHi)
+ return kTRUE;
+ }
+ }
+ }
+ return kFALSE;
}
+//________________________________________________________________________
+Bool_t AliDhcTask::HasMuonAOD()
+{
+ // does this AOD event have a good muon?
+ for (Int_t iMu = 0; iMu<fAOD->GetNumberOfTracks(); iMu++) {
+ AliAODTrack* muonTrack = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iMu));
+ if(!muonTrack) AliFatal("Not a standard AOD");
+ if (muonTrack) {
+ if ( IsGoodMUONtrack(*muonTrack) ) {
+ Double_t ptMu = muonTrack->Pt();
+ if (ptMu>fReqPtLo && ptMu<fReqPtHi)
+ return kTRUE;
+ }
+ }
+ }
+ return kFALSE;
+}
//________________________________________________________________________
void AliDhcTask::Terminate(Option_t *)
{
// Draw result to the screen
// Called once at the end of the query
-
- delete fPoolMgr;
- fHCent->SetEntries(fHCent->Integral());
- fHZvtx->SetEntries(fHZvtx->Integral());
- fHPtTrg->SetEntries(fHPtTrg->Integral());
- fHPtAss->SetEntries(fHPtAss->Integral());
}
//________________________________________________________________________
-Bool_t AliDhcTask::VertexOk(TObject* obj) const
+Bool_t AliDhcTask::VertexOk() const
{
// Modified from AliAnalyseLeadingTrackUE::VertexSelection()
Int_t nContributors = 0;
Double_t zVertex = 999;
TString name("");
-
- if (obj->InheritsFrom("AliESDEvent")) {
- AliESDEvent* esdevt = (AliESDEvent*) obj;
- const AliESDVertex* vtx = esdevt->GetPrimaryVertex();
+
+ Int_t runno = InputEvent()->GetRunNumber();
+ if (runno>=176326 && runno<=197692) { // year 12 and 13
+ if (!fUtils->IsVertexSelected2013pA(InputEvent()))
+ return 0;
+ }
+
+ if (fESD) {
+ const AliESDVertex* vtx = fESD->GetPrimaryVertex();
if (!vtx)
return 0;
nContributors = vtx->GetNContributors();
zVertex = vtx->GetZ();
name = vtx->GetName();
- }
- else if (obj->InheritsFrom("AliAODEvent")) {
- AliAODEvent* aodevt = (AliAODEvent*) obj;
- if (aodevt->GetNumberOfVertices() < 1)
+ } else {
+ if (fAOD->GetNumberOfVertices() < 1)
return 0;
- const AliAODVertex* vtx = aodevt->GetPrimaryVertex();
+ const AliAODVertex* vtx = fAOD->GetPrimaryVertex();
nContributors = vtx->GetNContributors();
zVertex = vtx->GetZ();
name = vtx->GetName();
return kTRUE;
}
+
+//________________________________________________________________________
+Bool_t AliDhcTask::IsGoodMUONtrack(AliESDMuonTrack &track)
+{
+ // Applying track cuts for MUON tracks
+
+ if (fUseMuonUtils) { // muon cut using official class
+ return fMuonTrackCuts->IsSelected(&track);
+ } else { // manual muon cut
+ if (!track.ContainTrackerData())
+ return kFALSE;
+ Double_t thetaTrackAbsEnd = TMath::ATan(track.GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
+ if ((thetaTrackAbsEnd < 2.) || (thetaTrackAbsEnd > 10.))
+ return kFALSE;
+ Double_t eta = track.Eta();
+ if ((eta < -4.) || (eta > -2.5))
+ return kFALSE;
+ if (fTriggerMatch) {
+ if (!track.ContainTriggerData())
+ return kFALSE;
+ if (track.GetMatchTrigger() < 0.5)
+ return kFALSE;
+ }
+ return kTRUE;
+ }
+}
+
+//________________________________________________________________________
+Bool_t AliDhcTask::IsGoodMUONtrack(AliAODTrack &track)
+{
+ // Applying track cuts for MUON tracks
+
+ if (fUseMuonUtils) { // muon cut using official class
+ return fMuonTrackCuts->IsSelected(&track);
+ } else { // manual muon cut
+ if (!track.IsMuonTrack())
+ return kFALSE;
+ Double_t dThetaAbs = TMath::ATan(track.GetRAtAbsorberEnd()/505.)* TMath::RadToDeg();
+ if ((dThetaAbs<2.) || (dThetaAbs>10.))
+ return kFALSE;
+ Double_t dEta = track.Eta();
+ if ((dEta<-4.) || (dEta>-2.5))
+ return kFALSE;
+ if (fTriggerMatch) {
+ if (track.GetMatchTrigger()<0.5)
+ return kFALSE;
+ }
+ return kTRUE;
+ }
+}
+
+//________________________________________________________________________
+Bool_t AliDhcTask::IsTrigger(Double_t eta, Double_t pt)
+{
+ // is in trigger eta,pt range?
+ Int_t bin = fHPtTrg->FindBin(pt);
+ if (fHPtTrg->IsBinOverflow(bin) || fHPtTrg->IsBinUnderflow(bin))
+ return kFALSE;
+
+ if (eta<fEtaTLo || eta>fEtaTHi)
+ return kFALSE;
+
+ return kTRUE;
+}
+
+//________________________________________________________________________
+Bool_t AliDhcTask::IsAssociated(Double_t eta, Double_t pt)
+{
+ // is in associated eta,pt range?
+ Int_t bin = fHPtAss->FindBin(pt);
+ if (fHPtAss->IsBinOverflow(bin) || fHPtAss->IsBinUnderflow(bin))
+ return kFALSE;
+
+ if (eta<fEtaALo || eta>fEtaAHi)
+ return kFALSE;
+
+ return kTRUE;
+}
+
+//________________________________________________________________________
+AliDhcTask::~AliDhcTask()
+{
+ //Destructor
+ if (fPoolMgr)
+ delete fPoolMgr;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff