#include "AliESDEvent.h"
#include "AliMCEvent.h"
-#include "AliESDInputHandler.h"
+#include "AliESDtrackCuts.h"
+#include "AliAnalysisFilter.h"
#include "AliInputEventHandler.h"
#include "AliVVertex.h"
-#include "AliAnalysisFilter.h"
#include "AliPID.h"
#include "AliPIDCombined.h"
#include "AliPIDResponse.h"
, fInputFromOtherTask(kFALSE)
, fDoPID(kTRUE)
, fDoEfficiency(kTRUE)
- , fDoPtResolution(kTRUE)
- , fDoDeDxCheck(kTRUE)
+ , fDoPtResolution(kFALSE)
+ , fDoDeDxCheck(kFALSE)
+ , fDoBinZeroStudy(kFALSE)
, fStoreCentralityPercentile(kFALSE)
, fStoreAdditionalJetInformation(kFALSE)
, fTakeIntoAccountMuons(kFALSE)
, fTOFmode(1)
, fEtaAbsCutLow(0.0)
, fEtaAbsCutUp(0.9)
+ , fPileUpRejectionType(AliAnalysisTaskPIDV0base::kPileUpRejectionOff)
, fDoAnySystematicStudiesOnTheExpectedSignal(kFALSE)
, fSystematicScalingSplinesThreshold(50.)
, fSystematicScalingSplinesBelowThreshold(1.0)
, fhEventsProcessed(0x0)
, fhEventsTriggerSel(0x0)
, fhEventsTriggerSelVtxCut(0x0)
- , fhSkippedTracksForSignalGeneration(0x0)
+ , fhEventsProcessedNoPileUpRejection(0x0)
+ , fChargedGenPrimariesTriggerSel(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCut(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCutZ(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCutZPileUpRej(0x0)
, fhMCgeneratedYieldsPrimaries(0x0)
, fh2FFJetPtRec(0x0)
, fh2FFJetPtGen(0x0)
, fh1Xsec(0x0)
, fh1Trials(0x0)
+ , fh1EvtsPtHardCut(0x0)
, fContainerEff(0x0)
+ , fQASharedCls(0x0)
, fDeDxCheck(0x0)
, fOutputContainer(0x0)
, fQAContainer(0x0)
, fInputFromOtherTask(kFALSE)
, fDoPID(kTRUE)
, fDoEfficiency(kTRUE)
- , fDoPtResolution(kTRUE)
- , fDoDeDxCheck(kTRUE)
+ , fDoPtResolution(kFALSE)
+ , fDoDeDxCheck(kFALSE)
+ , fDoBinZeroStudy(kFALSE)
, fStoreCentralityPercentile(kFALSE)
, fStoreAdditionalJetInformation(kFALSE)
, fTakeIntoAccountMuons(kFALSE)
, fTOFmode(1)
, fEtaAbsCutLow(0.0)
, fEtaAbsCutUp(0.9)
+ , fPileUpRejectionType(AliAnalysisTaskPIDV0base::kPileUpRejectionOff)
, fDoAnySystematicStudiesOnTheExpectedSignal(kFALSE)
, fSystematicScalingSplinesThreshold(50.)
, fSystematicScalingSplinesBelowThreshold(1.0)
, fhEventsProcessed(0x0)
, fhEventsTriggerSel(0x0)
, fhEventsTriggerSelVtxCut(0x0)
- , fhSkippedTracksForSignalGeneration(0x0)
+ , fhEventsProcessedNoPileUpRejection(0x0)
+ , fChargedGenPrimariesTriggerSel(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCut(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCutZ(0x0)
+ , fChargedGenPrimariesTriggerSelVtxCutZPileUpRej(0x0)
, fhMCgeneratedYieldsPrimaries(0x0)
, fh2FFJetPtRec(0x0)
, fh2FFJetPtGen(0x0)
, fh1Xsec(0x0)
, fh1Trials(0x0)
+ , fh1EvtsPtHardCut(0x0)
, fContainerEff(0x0)
+ , fQASharedCls(0x0)
, fDeDxCheck(0x0)
, fOutputContainer(0x0)
, fQAContainer(0x0)
fhMaxEtaVariation->SetBinContent(binY, maxAbs);
}
- printf("AliAnalysisTaskPID: Calculated max eta variation.\n");
+ printf("AliAnalysisTaskPID: Calculated max eta variation from map \"%s\".\n", hEta->GetTitle());
return kTRUE;
}
if(fDebug > 1)
printf("File: %s, Line: %d: UserCreateOutputObjects\n", (char*)__FILE__, __LINE__);
- SetUpPIDcombined();
-
- // Input handler
- AliAnalysisManager* man = AliAnalysisManager::GetAnalysisManager();
- AliInputEventHandler* inputHandler = dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
+ // Setup basic things, like PIDResponse
+ AliAnalysisTaskPIDV0base::UserCreateOutputObjects();
- if (!inputHandler)
- AliFatal("Input handler needed");
- else {
- // PID response object
- fPIDResponse = inputHandler->GetPIDResponse();
- if (!fPIDResponse)
- AliFatal("PIDResponse object was not created");
- }
+ if (!fPIDResponse)
+ AliFatal("PIDResponse object was not created");
- if(fDebug > 2)
- printf("File: %s, Line: %d: UserCreateOutputObjects -> Retrieved PIDresponse object\n", (char*)__FILE__, __LINE__);
+ SetUpPIDcombined();
OpenFile(1);
11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
- const Int_t nCentBins = 12;
- //-1 for pp; 90-100 has huge electro-magnetic impurities
- Double_t binsCent[nCentBins+1] = {-1, 0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100 };
+ const Bool_t useITSTPCtrackletsCentEstimatorWithNewBinning = fCentralityEstimator.CompareTo("ITSTPCtracklets", TString::kIgnoreCase) == 0
+ && fStoreCentralityPercentile;
+
+ const Int_t nCentBinsGeneral = 12;
+ const Int_t nCentBinsNewITSTPCtracklets = 17;
+
+ const Int_t nCentBins = useITSTPCtrackletsCentEstimatorWithNewBinning ? nCentBinsNewITSTPCtracklets : nCentBinsGeneral;
+
+ Double_t binsCent[nCentBins+1];
+ for (Int_t i = 0; i < nCentBins + 1; i++)
+ binsCent[i] = -1;
+
+ //-1 for pp (unless explicitely requested); 90-100 has huge electro-magnetic impurities
+ Double_t binsCentV0[nCentBinsGeneral+1] = {-1, 0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100 };
+
+ // These centrality estimators deal with integers! This implies that the ranges are always [lowlim, uplim - 1]
+ Double_t binsCentITSTPCTracklets[nCentBinsNewITSTPCtracklets+1] = { -9999, 0, 1, 4, 7, 10, 15, 20, 25, 30, 40, 50, 60, 70, 80, 90, 100, 9999 };
+ Double_t binsCentITSTPCTrackletsOldPreliminary[nCentBinsGeneral+1] = { 0, 7, 13, 20, 29, 40, 50, 60, 72, 83, 95, 105, 115 };
+
+ // Special centrality binning for pp
+ Double_t binsCentpp[nCentBinsGeneral+1] = { 0, 0.01, 0.1, 1, 5, 10, 15, 20, 30, 40, 50, 70, 100};
+
+ if (fCentralityEstimator.CompareTo("ITSTPCtrackletsOldPreliminaryBinning", TString::kIgnoreCase) == 0 && fStoreCentralityPercentile) {
+ // Special binning for this centrality estimator; but keep number of bins!
+ for (Int_t i = 0; i < nCentBinsGeneral+1; i++)
+ binsCent[i] = binsCentITSTPCTrackletsOldPreliminary[i];
+ }
+ else if (fCentralityEstimator.CompareTo("ITSTPCtracklets", TString::kIgnoreCase) == 0 && fStoreCentralityPercentile) {
+ // Special binning for this centrality estimator and different number of bins!
+ for (Int_t i = 0; i < nCentBinsNewITSTPCtracklets+1; i++)
+ binsCent[i] = binsCentITSTPCTracklets[i];
+ }
+ else if (fCentralityEstimator.Contains("ppMult", TString::kIgnoreCase) && fStoreCentralityPercentile) {
+ // Special binning for this pp centrality estimator; but keep number of bins!
+ for (Int_t i = 0; i < nCentBinsGeneral+1; i++)
+ binsCent[i] = binsCentpp[i];
+ }
+ else {
+ // Take default binning for VZERO
+ for (Int_t i = 0; i < nCentBinsGeneral+1; i++)
+ binsCent[i] = binsCentV0[i];
+ }
const Int_t nJetPtBins = 11;
Double_t binsJetPt[nJetPtBins+1] = {0, 2, 5, 10, 15, 20, 30, 40, 60, 80, 120, 200};
fhGenPr = new THnSparseD("hGenPr", "", nGenBins, genBins, genXmin, genXmax);
SetUpGenHist(fhGenPr, binsPt, deltaPrimeBins, binsCent, binsJetPt);
fOutputContainer->Add(fhGenPr);
-
- fhSkippedTracksForSignalGeneration = new TH2D("fhSkippedTracksForSignalGeneration",
- "Number of tracks skipped for the signal generation;p_{T}^{gen} (GeV/c);TPC signal N",
- nPtBins, binsPt, 161, -0.5, 160.5);
- fhSkippedTracksForSignalGeneration->Sumw2();
- fOutputContainer->Add(fhSkippedTracksForSignalGeneration);
}
fhEventsProcessed = new TH1D("fhEventsProcessed",
- "Number of events passing trigger selection, vtx and zvtx cuts;Centrality percentile",
+ "Number of events passing trigger selection, vtx and zvtx cuts and pile-up rejection;Centrality Percentile",
nCentBins, binsCent);
fhEventsProcessed->Sumw2();
fOutputContainer->Add(fhEventsProcessed);
fhEventsTriggerSelVtxCut = new TH1D("fhEventsTriggerSelVtxCut",
- "Number of events passing trigger selection and vtx cut;Centrality percentile",
+ "Number of events passing trigger selection and vtx cut;Centrality Percentile",
nCentBins, binsCent);
fhEventsTriggerSelVtxCut->Sumw2();
fOutputContainer->Add(fhEventsTriggerSelVtxCut);
fhEventsTriggerSel = new TH1D("fhEventsTriggerSel",
- "Number of events passing trigger selection;Centrality percentile",
+ "Number of events passing trigger selection;Centrality Percentile",
nCentBins, binsCent);
fOutputContainer->Add(fhEventsTriggerSel);
fhEventsTriggerSel->Sumw2();
+ fhEventsProcessedNoPileUpRejection = new TH1D("fhEventsProcessedNoPileUpRejection",
+ "Number of events passing trigger selection, vtx and zvtx cuts;Centrality Percentile",
+ nCentBins, binsCent);
+ fOutputContainer->Add(fhEventsProcessedNoPileUpRejection);
+ fhEventsProcessedNoPileUpRejection->Sumw2();
+
+
// Generated yields within acceptance
const Int_t nBinsGenYields = fStoreAdditionalJetInformation ? kGenYieldNumAxes : kGenYieldNumAxes - 3;
Int_t genYieldsBins[kGenYieldNumAxes] = { nMCPIDbins, nPtBins, nCentBins, nJetPtBins, nZBins, nXiBins,
fh1Trials->Sumw2();
fh1Trials->GetXaxis()->SetBinLabel(1, "#sum{ntrials}");
+ fh1EvtsPtHardCut = new TH1F("fh1EvtsPtHardCut", "#events before and after MC #it{p}_{T,hard} cut;;Events",2,0,2);
+ fh1EvtsPtHardCut->Sumw2();
+ fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(1, "All");
+ fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(2, "#it{p}_{T,hard}");
+
fOutputContainer->Add(fh1Xsec);
fOutputContainer->Add(fh1Trials);
+ fOutputContainer->Add(fh1EvtsPtHardCut);
if (fDoDeDxCheck || fDoPtResolution) {
OpenFile(3);
SetUpPtResHist(fPtResolution[i], pTbinsRes, binsJetPt, binsCent);
fQAContainer->Add(fPtResolution[i]);
}
+
+
+ // Besides the pT resolution, also perform check on shared clusters
+ const Int_t nBinsQASharedCls = kQASharedClsNumAxes;
+ Int_t qaSharedClsBins[kQASharedClsNumAxes] = { nJetPtBins, nPtBinsRes, 160, 160 };
+ Double_t qaSharedClsXmin[kQASharedClsNumAxes] = { binsJetPt[0], pTbinsRes[0], 0, -1 };
+ Double_t qaSharedClsXmax[kQASharedClsNumAxes] = { binsJetPt[nJetPtBins], pTbinsRes[nPtBinsRes], 160, 159 };
+
+ fQASharedCls = new THnSparseD("fQASharedCls", "QA shared clusters", nBinsQASharedCls, qaSharedClsBins, qaSharedClsXmin, qaSharedClsXmax);
+
+ SetUpSharedClsHist(fQASharedCls, pTbinsRes, binsJetPt);
+ fQAContainer->Add(fQASharedCls);
}
fQAContainer->Add(fDeDxCheck);
}
+ if (fDoBinZeroStudy) {
+ const Double_t etaLow = -0.9;
+ const Double_t etaUp = 0.9;
+ const Int_t nEtaBins = 18;
+
+ const Int_t nBinsBinZeroStudy = kBinZeroStudyNumAxes;
+ Int_t binZeroStudyBins[nBinsBinZeroStudy] = { nCentBins, nPtBins, nEtaBins };
+ Double_t binZeroStudyXmin[nBinsBinZeroStudy] = { binsCent[0], binsPt[0], etaLow };
+ Double_t binZeroStudyXmax[nBinsBinZeroStudy] = { binsCent[nCentBins], binsPt[nPtBins], etaUp };
+
+ fChargedGenPrimariesTriggerSel = new THnSparseD("fChargedGenPrimariesTriggerSel", "Trigger sel.", nBinsBinZeroStudy, binZeroStudyBins,
+ binZeroStudyXmin, binZeroStudyXmax);
+ SetUpBinZeroStudyHist(fChargedGenPrimariesTriggerSel, binsCent, binsPt);
+ fOutputContainer->Add(fChargedGenPrimariesTriggerSel);
+
+ fChargedGenPrimariesTriggerSelVtxCut = new THnSparseD("fChargedGenPrimariesTriggerSelVtxCut", "Vertex cut", nBinsBinZeroStudy,
+ binZeroStudyBins, binZeroStudyXmin, binZeroStudyXmax);
+ SetUpBinZeroStudyHist(fChargedGenPrimariesTriggerSelVtxCut, binsCent, binsPt);
+ fOutputContainer->Add(fChargedGenPrimariesTriggerSelVtxCut);
+
+ fChargedGenPrimariesTriggerSelVtxCutZ = new THnSparseD("fChargedGenPrimariesTriggerSelVtxCutZ", "Vertex #it{z} cut", nBinsBinZeroStudy,
+ binZeroStudyBins, binZeroStudyXmin, binZeroStudyXmax);
+ SetUpBinZeroStudyHist(fChargedGenPrimariesTriggerSelVtxCutZ, binsCent, binsPt);
+ fOutputContainer->Add(fChargedGenPrimariesTriggerSelVtxCutZ);
+
+ fChargedGenPrimariesTriggerSelVtxCutZPileUpRej = new THnSparseD("fChargedGenPrimariesTriggerSelVtxCutZPileUpRej", "Vertex #it{z} cut",
+ nBinsBinZeroStudy, binZeroStudyBins, binZeroStudyXmin, binZeroStudyXmax);
+ SetUpBinZeroStudyHist(fChargedGenPrimariesTriggerSelVtxCutZPileUpRej, binsCent, binsPt);
+ fOutputContainer->Add(fChargedGenPrimariesTriggerSelVtxCutZPileUpRej);
+ }
+
if(fDebug > 2)
printf("File: %s, Line: %d: UserCreateOutputObjects -> Posting output data\n", (char*)__FILE__, __LINE__);
PostData(1, fOutputContainer);
- PostData(2, fContainerEff);
- PostData(3, fQAContainer);
+ if (fDoEfficiency)
+ PostData(2, fContainerEff);
+ if (fDoDeDxCheck || fDoPtResolution)
+ PostData(3, fQAContainer);
if(fDebug > 2)
printf("File: %s, Line: %d: UserCreateOutputObjects -> Done\n", (char*)__FILE__, __LINE__);
if(fDebug > 1)
printf("File: %s, Line: %d: UserExec\n", (char*)__FILE__, __LINE__);
- Int_t run = InputEvent()->GetRunNumber();
-
- if (run != fRun){
- // If systematics on eta is investigated, need to calculate the maxEtaVariationMap
- if ((TMath::Abs(fSystematicScalingEtaCorrectionLowMomenta - 1.0) > fgkEpsilon) ||
- (TMath::Abs(fSystematicScalingEtaCorrectionHighMomenta - 1.0) > fgkEpsilon)) {
- if (!CalculateMaxEtaVariationMapFromPIDResponse())
- AliFatal("Systematics on eta correction requested, but failed to calculate max eta varation map!");
- }
- }
-
- fRun = run;
-
// No processing of event, if input is fed in directly from another task
if (fInputFromOtherTask)
return;
return;
}
+ ConfigureTaskForCurrentEvent(fEvent);
+
fMC = dynamic_cast<AliMCEvent*>(MCEvent());
if (!fPIDResponse || !fPIDcombined)
return;
Double_t centralityPercentile = -1;
- if (fStoreCentralityPercentile)
- centralityPercentile = fEvent->GetCentrality()->GetCentralityPercentile(fCentralityEstimator.Data());
-
- IncrementEventCounter(centralityPercentile, kTriggerSel);
+ if (fStoreCentralityPercentile) {
+ if (fCentralityEstimator.Contains("ITSTPCtracklets", TString::kIgnoreCase)) {
+ // Special pp centrality estimator
+ AliESDEvent* esdEvent = dynamic_cast<AliESDEvent*>(fEvent);
+ if (!esdEvent) {
+ AliError("Not esd event -> Cannot use tracklet multiplicity estimator!");
+ centralityPercentile = -1;
+ }
+ else
+ centralityPercentile = AliESDtrackCuts::GetReferenceMultiplicity(esdEvent, AliESDtrackCuts::kTrackletsITSTPC, fEtaAbsCutUp);
+ }
+ else if (fCentralityEstimator.Contains("ppMult", TString::kIgnoreCase)) {
+ // Another special pp centrality estimator
+ centralityPercentile = fAnaUtils->GetMultiplicityPercentile(fEvent, GetPPCentralityEstimator().Data());
+ }
+ else {
+ // Ordinary centrality estimator
+ centralityPercentile = fEvent->GetCentrality()->GetCentralityPercentile(fCentralityEstimator.Data());
+ }
+ }
// Check if vertex is ok, but don't apply cut on z position
- if (!GetVertexIsOk(fEvent, kFALSE))
- return;
+ const Bool_t passedVertexSelection = GetVertexIsOk(fEvent, kFALSE);
+ // Now check again, but also require z position to be in desired range
+ const Bool_t passedVertexZSelection = GetVertexIsOk(fEvent, kTRUE);
+ // Check pile-up
+ const Bool_t isPileUp = GetIsPileUp(fEvent, fPileUpRejectionType);
- fESD = dynamic_cast<AliESDEvent*>(fEvent);
- const AliVVertex* primaryVertex = fESD ? fESD->GetPrimaryVertexTracks() : fEvent->GetPrimaryVertex();
- if (!primaryVertex)
- return;
- if(primaryVertex->GetNContributors() <= 0)
+
+ if (fDoBinZeroStudy && fMC) {
+ for (Int_t iPart = 0; iPart < fMC->GetNumberOfTracks(); iPart++) {
+ AliMCParticle *mcPart = dynamic_cast<AliMCParticle*>(fMC->GetTrack(iPart));
+
+ if (!mcPart)
+ continue;
+
+ if (!fMC->IsPhysicalPrimary(iPart))
+ continue;
+
+ const Double_t etaGen = mcPart->Eta();
+ const Double_t ptGen = mcPart->Pt();
+
+ Double_t values[kBinZeroStudyNumAxes] = { 0. };
+ values[kBinZeroStudyCentrality] = centralityPercentile;
+ values[kBinZeroStudyGenPt] = ptGen;
+ values[kBinZeroStudyGenEta] = etaGen;
+
+ fChargedGenPrimariesTriggerSel->Fill(values);
+ if (passedVertexSelection) {
+ fChargedGenPrimariesTriggerSelVtxCut->Fill(values);
+ if (passedVertexZSelection) {
+ fChargedGenPrimariesTriggerSelVtxCutZ->Fill(values);
+ if (!isPileUp) {
+ fChargedGenPrimariesTriggerSelVtxCutZPileUpRej->Fill(values);
+ }
+ }
+ }
+ }
+ }
+
+
+
+ // Event counters for trigger selection, vertex cuts and pile-up rejection
+ IncrementEventCounter(centralityPercentile, kTriggerSel);
+
+ if (!passedVertexSelection)
return;
IncrementEventCounter(centralityPercentile, kTriggerSelAndVtxCut);
- // Now check again, but also require z position to be in desired range
- if (!GetVertexIsOk(fEvent, kTRUE))
+ if (!passedVertexZSelection)
+ return;
+
+ IncrementEventCounter(centralityPercentile, kTriggerSelAndVtxCutAndZvtxCutNoPileUpRejection);
+ // ATTENTION: Is this the right place for the pile-up rejection? Important to have still the proper bin-0 correction,
+ // which is done solely with sel and selVtx, since the zvtx selection does ~not change the spectra. The question is whether the pile-up
+ // rejection changes the spectra. If not, then it is perfectly fine to put it here and keep the usual histo for the normalisation to number
+ // of events. But if it does change the spectra, this must somehow be corrected for.
+ // NOTE: multiplicity >= 0 usually implies a properly reconstructed vertex. Hence, the bin-0 correction cannot be done in multiplicity bins.
+ // Furthermore, there seems to be no MC simulation with pile-up rejection, so the bin-0 correction cannot be extracted with it. Pile-up
+ // rejection has only a minor impact, so maybe there is no need to dig further.
+ if (isPileUp)
return;
IncrementEventCounter(centralityPercentile, kTriggerSelAndVtxCutAndZvtxCut);
// Only process tracks inside the desired eta window
if (!IsInAcceptedEtaRange(TMath::Abs(track->Eta()))) continue;
- if (fDoPID || fDoDeDxCheck)
+ if (fDoPID || fDoDeDxCheck || fDoPtResolution)
ProcessTrack(track, pdg, centralityPercentile, -1); // No jet information in this case -> Set jet pT to -1
if (fDoPtResolution) {
}
+//_____________________________________________________________________________
+void AliAnalysisTaskPID::ConfigureTaskForCurrentEvent(AliVEvent* event)
+{
+ // Configure the task for the current event. In particular, this is needed if the run number changes
+
+ if (!event) {
+ AliError("Could not set up task: no event!");
+ return;
+ }
+
+ Int_t run = event->GetRunNumber();
+
+ if (run != fRun){
+ // If systematics on eta is investigated, need to calculate the maxEtaVariationMap
+ if ((TMath::Abs(fSystematicScalingEtaCorrectionLowMomenta - 1.0) > fgkEpsilon) ||
+ (TMath::Abs(fSystematicScalingEtaCorrectionHighMomenta - 1.0) > fgkEpsilon)) {
+ if (!CalculateMaxEtaVariationMapFromPIDResponse())
+ AliFatal("Systematics on eta correction requested, but failed to calculate max eta varation map!");
+ }
+ }
+
+ fRun = run;
+}
+
+
//_____________________________________________________________________________
Int_t AliAnalysisTaskPID::PDGtoMCID(Int_t pdg)
{
}
if (absMotherPDG == 3122) { // Lambda
+ //if (absMotherPDG == 3122 || absMotherPDG == 3112 || absMotherPDG == 3222) { // Lambda / Sigma- / Sigma+
if (0.00 <= motherGenPt && motherGenPt < 0.20) fac = 0.645162;
else if(0.20 <= motherGenPt && motherGenPt < 0.40) fac = 0.627431;
else if(0.40 <= motherGenPt && motherGenPt < 0.60) fac = 0.457136;
// Get the (locally defined) particle type judged by TOF
if (!fPIDResponse) {
- Printf("ERROR: fEvent not available -> Cannot determine TOF type!");
+ Printf("ERROR: fPIDResponse not available -> Cannot determine TOF type!");
return kNoTOFinfo;
}
-
+
+ /*TODO still needs some further thinking how to implement it....
+ // TOF PID status (kTOFout, kTIME) automatically checked by return value of ComputeTPCProbability;
+ // also, probability array will be set there (no need to initialise here)
+ Double_t p[AliPID::kSPECIES];
+ const AliPIDResponse::EDetPidStatus tofStatus = fPIDResponse->ComputeTPCProbability(track, AliPID::kSPECIES, p);
+ if (tofStatus != AliPIDResponse::kDetPidOk)
+ return kNoTOFinfo;
+
+ // Do not consider muons
+ p[AliPID::kMuon] = 0.;
+
+ // Probabilities are not normalised yet
+ Double_t sum = 0.;
+ for (Int_t i = 0; i < AliPID::kSPECIES; i++)
+ sum += p[i];
+
+ if (sum <= 0.)
+ return kNoTOFinfo;
+
+ for (Int_t i = 0; i < AliPID::kSPECIES; i++)
+ p[i] /= sum;
+
+ Double_t probThreshold = -999.;
+
+ // If there is only one distribution, the threshold corresponds to...
+ if (tofMode == 0) {
+ probThreshold = ;
+ }
+ else if (tofMode == 1) { // default
+ probThreshold = 0.9973; // a 3-sigma inclusion cut
+ }
+ else if (tofMode == 2) {
+ inclusion = 3.;
+ exclusion = 3.5;
+ }
+ else {
+ Printf("ERROR: Bad TOF mode: %d!", tofMode);
+ return kNoTOFinfo;
+ }
+
+ */
+
+ ///* OLD: cut with n-sigmas (ATTENTION: Does not take into account properly the TOF tail!)
// Check kTOFout, kTIME, mismatch
const AliPIDResponse::EDetPidStatus tofStatus = fPIDResponse->CheckPIDStatus(AliPIDResponse::kTOF, track);
if (tofStatus != AliPIDResponse::kDetPidOk)
return kTOFkaon;
if (TMath::Abs(nsigma[kTOFpion]) > exclusion && TMath::Abs(nsigma[kTOFkaon]) > exclusion && TMath::Abs(nsigma[kTOFproton]) < inclusion)
return kTOFproton;
+ //*/
// There are no TOF electrons selected because the purity is rather bad, even for small momenta
// (also a small mismatch probability significantly affects electrons because their fraction
printf("\n");
+ printf("Pile up rejection type: %d\n", (Int_t)fPileUpRejectionType);
+
+ printf("\n");
+
printf("Use MC-ID for signal generation: %d\n", GetUseMCidForGeneration());
printf("Use ITS: %d\n", GetUseITS());
printf("Use TOF: %d\n", GetUseTOF());
printf("Do Efficiency: %d\n", fDoEfficiency);
printf("Do PtResolution: %d\n", fDoPtResolution);
printf("Do dEdxCheck: %d\n", fDoDeDxCheck);
+ printf("Do binZeroStudy: %d\n", fDoBinZeroStudy);
printf("\n");
if(fDebug > 1)
printf("File: %s, Line: %d: ProcessTrack\n", (char*)__FILE__, __LINE__);
- if (!fDoPID && !fDoDeDxCheck)
+ if (!fDoPID && !fDoDeDxCheck && !fDoPtResolution)
return kFALSE;
if(fDebug > 2)
Double_t dEdxTPC = tuneOnDataTPC ? fPIDResponse->GetTPCsignalTunedOnData(track) : track->GetTPCsignal();
if (dEdxTPC <= 0) {
- Printf("Skipping track with strange dEdx value: dEdx %f, pTPC %f, eta %f, ncl %d\n", track->GetTPCsignal(), pTPC,
- track->Eta(), track->GetTPCsignalN());
+ if (fDebug > 1)
+ Printf("Skipping track with strange dEdx value: dEdx %f, pTPC %f, eta %f, ncl %d\n", track->GetTPCsignal(), pTPC,
+ track->Eta(), track->GetTPCsignalN());
return kFALSE;
}
if ((pTPC >= 0.3 && (nSigmaPr > 10 && nSigmaEl > 10)) ||
(pTPC < 0.3 && (nSigmaPr > 15 && nSigmaEl > 15))) {
- Printf("Skipping track which seems to be a light nuclei: dEdx %f, pTPC %f, pT %f, eta %f, ncl %d, nSigmaPr %f, nSigmaEl %f\n",
- track->GetTPCsignal(), pTPC, pT, track->Eta(), track->GetTPCsignalN(), nSigmaPr, nSigmaEl);
+ if (fDebug > 1)
+ Printf("Skipping track which seems to be a light nucleus: dEdx %f, pTPC %f, pT %f, eta %f, ncl %d, nSigmaPr %f, nSigmaEl %f\n",
+ track->GetTPCsignal(), pTPC, pT, track->Eta(), track->GetTPCsignalN(), nSigmaPr, nSigmaEl);
return kFALSE;
}
printf("File: %s, Line: %d: ProcessTrack -> Generate Responses for dEdx check done\n", (char*)__FILE__, __LINE__);
}
+ if (fDoPtResolution) {
+ // Check shared clusters, which is done together with the pT resolution
+ Double_t qaEntry[kQASharedClsNumAxes];
+ qaEntry[kQASharedClsJetPt] = jetPt;
+ qaEntry[kQASharedClsPt] = pT;
+ qaEntry[kDeDxCheckP] = pTPC;
+ qaEntry[kQASharedClsNumSharedCls] = track->GetTPCSharedMapPtr()->CountBits();
+
+ Int_t iRowInd = -1;
+ // iRowInd == -1 for "all rows w/o multiple counting"
+ qaEntry[kQASharedClsPadRow] = iRowInd;
+ fQASharedCls->Fill(qaEntry);
+
+ // Fill hist for every pad row with shared cluster
+ for (iRowInd = 0; iRowInd < 159; iRowInd++) {
+ if (track->GetTPCSharedMapPtr()->TestBitNumber(iRowInd)) {
+ qaEntry[kQASharedClsPadRow] = iRowInd;
+ fQASharedCls->Fill(qaEntry);
+ }
+ }
+ }
+
if (!fDoPID)
return kTRUE;
}
if (errCode != kWarning) {
- fhSkippedTracksForSignalGeneration->Fill(track->GetTPCmomentum(), track->GetTPCsignalN());
return kFALSE;// Skip generated response in case of error
}
}
hist->GetAxis(kGenDeltaPrimeSpecies)->SetTitle("TPC #Delta'_{species} (arb. unit)");
- hist->GetAxis(kGenCentrality)->SetTitle(Form("Centrality Percentile (%s)", fCentralityEstimator.Data()));
+ hist->GetAxis(kGenCentrality)->SetTitle(Form("Centrality Percentile (%s)", GetPPCentralityEstimator().Data()));
if (fStoreAdditionalJetInformation) {
hist->GetAxis(kGenJetPt)->SetTitle("p_{T}^{jet} (GeV/c)");
// Set axes titles
hist->GetAxis(kGenYieldMCID)->SetTitle("MC PID");
hist->GetAxis(kGenYieldPt)->SetTitle("p_{T}^{gen} (GeV/c)");
- hist->GetAxis(kGenYieldCentrality)->SetTitle(Form("Centrality Percentile (%s)", fCentralityEstimator.Data()));
+ hist->GetAxis(kGenYieldCentrality)->SetTitle(Form("Centrality Percentile (%s)", GetPPCentralityEstimator().Data()));
if (fStoreAdditionalJetInformation) {
hist->GetAxis(kGenYieldJetPt)->SetTitle("p_{T}^{jet, gen} (GeV/c)");
hist->GetAxis(kDataDeltaPrimeSpecies)->SetTitle("TPC #Delta'_{species} (arb. unit)");
- hist->GetAxis(kDataCentrality)->SetTitle(Form("Centrality Percentile (%s)", fCentralityEstimator.Data()));
+ hist->GetAxis(kDataCentrality)->SetTitle(Form("Centrality Percentile (%s)", GetPPCentralityEstimator().Data()));
if (fStoreAdditionalJetInformation) {
hist->GetAxis(kDataJetPt)->SetTitle("p_{T}^{jet} (GeV/c)");
hist->GetAxis(kPtResRecPt)->SetTitle("p_{T}^{rec} (GeV/c)");
hist->GetAxis(kPtResCharge)->SetTitle("Charge (e_{0})");
- hist->GetAxis(kPtResCentrality)->SetTitle(Form("Centrality Percentile (%s)", fCentralityEstimator.Data()));
+ hist->GetAxis(kPtResCentrality)->SetTitle(Form("Centrality Percentile (%s)", GetPPCentralityEstimator().Data()));
}
+
+//________________________________________________________________________
+void AliAnalysisTaskPID::SetUpSharedClsHist(THnSparse* hist, Double_t* binsPt, Double_t* binsJetPt) const
+{
+ // Sets bin limits for axes which are not standard binned and the axes titles.
+
+ hist->SetBinEdges(kQASharedClsJetPt, binsJetPt);
+ hist->SetBinEdges(kQASharedClsPt, binsPt);
+
+ // Set axes titles
+ hist->GetAxis(kQASharedClsJetPt)->SetTitle("#it{p}_{T}^{jet} (GeV/#it{c})");
+ hist->GetAxis(kQASharedClsPt)->SetTitle("#it{p}_{T} (GeV/#it{c})");
+ hist->GetAxis(kQASharedClsNumSharedCls)->SetTitle("#it{N}_{shared}^{cls}");
+ hist->GetAxis(kQASharedClsPadRow)->SetTitle("Pad row");
+
+}
+
+
//________________________________________________________________________
void AliAnalysisTaskPID::SetUpDeDxCheckHist(THnSparse* hist, const Double_t* binsPt, const Double_t* binsJetPt, const Double_t* binsEtaAbs) const
{
hist->GetAxis(kDeDxCheckEtaAbs)->SetTitle("|#eta|");
hist->GetAxis(kDeDxCheckP)->SetTitle("p_{TPC} (GeV/c)");
hist->GetAxis(kDeDxCheckDeDx)->SetTitle("TPC dE/dx (arb. unit)");
+}
+
+
+//________________________________________________________________________
+void AliAnalysisTaskPID::SetUpBinZeroStudyHist(THnSparse* hist, const Double_t* binsCent, const Double_t* binsPt) const
+{
+ // Sets bin limits for axes which are not standard binned and the axes titles.
+ hist->SetBinEdges(kBinZeroStudyCentrality, binsCent);
+ hist->SetBinEdges(kBinZeroStudyGenPt, binsPt);
+
+ // Set axes titles
+ hist->GetAxis(kBinZeroStudyCentrality)->SetTitle(Form("Centrality Percentile (%s)", GetPPCentralityEstimator().Data()));
+ hist->GetAxis(kBinZeroStudyGenEta)->SetTitle("#it{#eta}^{gen}");
+ hist->GetAxis(kBinZeroStudyGenPt)->SetTitle("#it{p}_{T}^{gen} (GeV/#it{c})");
}
\ No newline at end of file
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff