From b2cf720741bc21cc992e6b827fb2fe6d6b9bfac1 Mon Sep 17 00:00:00 2001 From: esicking Date: Fri, 7 Oct 2011 13:00:15 +0000 Subject: [PATCH] Added analysis task for effiency studies (original by Veronica Canoa Roman). Now, different track types can be tested (global, TPC, ITS). It can be used to study nuclei and anti-nuclei transport, too. --- .../scripts/efficiency/AddTaskEfficiency.C | 163 +++ .../efficiency/AliAnalysisTaskEfficiency.cxx | 1085 +++++++++++++++++ .../efficiency/AliAnalysisTaskEfficiency.h | 149 +++ .../scripts/efficiency/CreateAnalysisPlugin.C | 145 +++ .../scripts/efficiency/SetupAnalysis.C | 148 +++ test/vmctest/scripts/efficiency/doubleRatio.C | 310 +++++ test/vmctest/scripts/efficiency/macrosalida.C | 290 +++++ test/vmctest/scripts/efficiency/runALICE.C | 65 + 8 files changed, 2355 insertions(+) create mode 100644 test/vmctest/scripts/efficiency/AddTaskEfficiency.C create mode 100644 test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.cxx create mode 100644 test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.h create mode 100644 test/vmctest/scripts/efficiency/CreateAnalysisPlugin.C create mode 100644 test/vmctest/scripts/efficiency/SetupAnalysis.C create mode 100644 test/vmctest/scripts/efficiency/doubleRatio.C create mode 100644 test/vmctest/scripts/efficiency/macrosalida.C create mode 100644 test/vmctest/scripts/efficiency/runALICE.C diff --git a/test/vmctest/scripts/efficiency/AddTaskEfficiency.C b/test/vmctest/scripts/efficiency/AddTaskEfficiency.C new file mode 100644 index 00000000000..c20ade42ce9 --- /dev/null +++ b/test/vmctest/scripts/efficiency/AddTaskEfficiency.C @@ -0,0 +1,163 @@ +AliAnalysisTaskEfficiency * AddTaskEfficiency(Int_t runNumber) + +{ + // Creates a QA task exploiting simple symmetries phi, eta +/-, charge ... + + // Get the pointer to the existing analysis manager via the static access method. + //============================================================================== + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) { + ::Error("AddTaskEfficiency", "No analysis manager to connect to."); + return NULL; + } + + // Check the analysis type using the event handlers connected to the analysis manager. + //============================================================================== + if (!mgr->GetInputEventHandler()) { + ::Error("AddTaskEfficiency", "This task requires an input event handler"); + return NULL; + } + + TString inputDataType = mgr->GetInputEventHandler()->GetDataType(); + // can be "ESD" or "AOD" + + // Configure analysis + //=========================================================================== + + //TASKS + //Task for global tracks + AliAnalysisTaskEfficiency *task0 = + new AliAnalysisTaskEfficiency("AliAnalysisTaskEfficiency_Global"); + task0->SelectCollisionCandidates(); // default setting: kMB = min bias trigger + task0->SetTrackType(0);// global tracks + + // Task for ITS tracks + AliAnalysisTaskEfficiency *task1 = + new AliAnalysisTaskEfficiency("AliAnalysisTaskEfficiency_ITS"); + task1->SelectCollisionCandidates(); + task1->SetTrackType(1);// ITS tracks + + //Task for ITS_SA tracks + AliAnalysisTaskEfficiency *task1sa = + new AliAnalysisTaskEfficiency("AliAnalysisTaskEfficiency_ITS_SA"); + task1sa->SelectCollisionCandidates(); + task1sa->SetTrackType(1);// ITS tracks + + //Task for TPC tracks + AliAnalysisTaskEfficiency *task2 = + new AliAnalysisTaskEfficiency("AliAnalysisTaskEfficiency_TPC"); + task2->SelectCollisionCandidates(); + task2->SetTrackType(2);// TPC only tracks + + //CUTS + //cuts for global tracks + AliESDtrackCuts* esdTrackCutsL0 = new AliESDtrackCuts("AliESDtrackCuts0","Global"); + esdTrackCutsL0->SetMinNClustersTPC(70); + esdTrackCutsL0->SetRequireTPCRefit(kTRUE); + esdTrackCutsL0->SetRequireITSRefit(kTRUE); + esdTrackCutsL0->SetMaxDCAToVertexXY(3.); + esdTrackCutsL0->SetMaxDCAToVertexZ(3.); + esdTrackCutsL0->SetAcceptKinkDaughters(kFALSE); + + //cuts for ITS tracks + AliESDtrackCuts* esdTrackCutsL1 = new AliESDtrackCuts("AliESDtrackCuts1","ITS"); + esdTrackCutsL1->SetMaxDCAToVertexXY(3.); + esdTrackCutsL1->SetMaxDCAToVertexZ(3.); + esdTrackCutsL1->SetAcceptKinkDaughters(kFALSE); + esdTrackCutsL1->SetRequireITSRefit(kTRUE); + esdTrackCutsL1->SetRequireITSStandAlone(kTRUE); + + //cuts for ITS tracks SA + AliESDtrackCuts* esdTrackCutsL1sa = new AliESDtrackCuts("AliESDtrackCuts1sa","ITS_SA"); + esdTrackCutsL1sa->SetMaxDCAToVertexXY(3.); + esdTrackCutsL1sa->SetMaxDCAToVertexZ(3.); + esdTrackCutsL1sa->SetAcceptKinkDaughters(kFALSE); + esdTrackCutsL1sa->SetRequireITSRefit(kTRUE); + esdTrackCutsL1sa->SetRequireITSPureStandAlone(kTRUE); + + //cuts for TPC tracks + AliESDtrackCuts* esdTrackCutsL2 = new AliESDtrackCuts("AliESDtrackCuts2","TPC"); + esdTrackCutsL2->SetRequireTPCRefit(kFALSE); + esdTrackCutsL2->SetAcceptKinkDaughters(kFALSE); + //jacek's cuts: + esdTrackCutsL2->SetMinNClustersTPC(70); + // cut on max ncl=160 in Task + esdTrackCutsL2->SetMaxDCAToVertexXY(3.); + esdTrackCutsL2->SetMaxDCAToVertexZ(3.); + esdTrackCutsL2->SetMaxChi2PerClusterTPC(3.999); + //cut minChi=0 in task + //esdTrackCutsL2->SetPRange(0.15,16); // not needed for QA + //esdTrackCutsL2->SetEtaRange(-0.8, 0.7999); // not needed for QA + //cuts for ITS tracks + + //add cuts to tasks + task0->SetCuts(esdTrackCutsL0); + task1->SetCuts(esdTrackCutsL1); + task1sa->SetCuts(esdTrackCutsL1sa); + task2->SetCuts(esdTrackCutsL2); + + // add tasks to manager + mgr->AddTask(task0); + mgr->AddTask(task1); + mgr->AddTask(task1sa); + mgr->AddTask(task2); + + //output container for tasks + AliAnalysisDataContainer *cout0 = 0; + AliAnalysisDataContainer *cout1 = 0; + AliAnalysisDataContainer *cout1sa = 0; + AliAnalysisDataContainer *cout2 = 0; + + + if(runNumber>0){ + cout0 = mgr->CreateContainer("QAHists_Global",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("run%d.root",runNumber)); + cout1 = mgr->CreateContainer("QAHists_ITS",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("run%d.root",runNumber)); + cout1sa = mgr->CreateContainer("QAHists_ITS_SA",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("run%d.root",runNumber)); + cout2 = mgr->CreateContainer("QAHists_TPC",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("run%d.root",runNumber)); + } + + else{ + cout0 = mgr->CreateContainer("QAHists_Global",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("%s:QAHists",AliAnalysisManager:: + GetCommonFileName())); + cout1 = mgr->CreateContainer("QAHists_ITS",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("%s:QAHists",AliAnalysisManager:: + GetCommonFileName())); + cout1sa = mgr->CreateContainer("QAHists_ITS_SA",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("%s:QAHists",AliAnalysisManager:: + GetCommonFileName())); + cout2 = mgr->CreateContainer("QAHists_TPC",TList::Class(), + AliAnalysisManager::kOutputContainer, + Form("%s:QAHists",AliAnalysisManager:: + GetCommonFileName())); + } + + //connect input to task + mgr->ConnectInput (task0, 0, mgr->GetCommonInputContainer()); + mgr->ConnectInput (task1, 0, mgr->GetCommonInputContainer()); + mgr->ConnectInput (task1sa, 0, mgr->GetCommonInputContainer()); + mgr->ConnectInput (task2, 0, mgr->GetCommonInputContainer()); + + //connect output to task + mgr->ConnectOutput (task0, 1, cout0); + mgr->ConnectOutput (task1, 1, cout1); + mgr->ConnectOutput (task1sa, 1, cout1sa); + mgr->ConnectOutput (task2, 1, cout2); + + + return task0; + +} + + diff --git a/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.cxx b/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.cxx new file mode 100644 index 00000000000..037f0231620 --- /dev/null +++ b/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.cxx @@ -0,0 +1,1085 @@ +#include "TChain.h" +#include "TTree.h" +#include "TH1F.h" +#include "TH2F.h" +#include "TProfile.h" + +#include "TCanvas.h" +#include "TList.h" +#include "TParticle.h" +#include "TParticlePDG.h" +#include "TProfile.h" + +#include "AliAnalysisTask.h" +#include "AliAnalysisManager.h" + +#include "AliESDEvent.h" +#include "AliStack.h" +#include "AliESDVertex.h" +#include "AliESDInputHandler.h" +#include "AliESDtrackCuts.h" +#include "AliMultiplicity.h" + +#include "AliMCParticle.h" +#include "AliMCEvent.h" +#include "AliAnalysisTaskEfficiency.h" +#include "AliExternalTrackParam.h" +#include "AliTrackReference.h" +#include "AliStack.h" +#include "AliHeader.h" +#include "AliGenEventHeader.h" +#include "AliGenCocktailEventHeader.h" +#include "AliGenPythiaEventHeader.h" +#include "AliGenDPMjetEventHeader.h" + +// Analysis Task for basic QA on the ESD +// Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing, Veronica Canoa, AM, Eva Sicking + +ClassImp(AliAnalysisTaskEfficiency) + +//________________________________________________________________________ + AliAnalysisTaskEfficiency::AliAnalysisTaskEfficiency(const char *name) + : AliAnalysisTaskSE(name) + ,fFieldOn(kTRUE) + ,fHists(0) + ,fHistRECpt(0) + ,fHistMCpt(0) + ,fHistMCNRpt(0) + ,fHistFAKEpt(0) + ,fHistMULTpt(0) + ,fHistRECeta(0) + ,fHistMCeta(0) + ,fHistMCNReta(0) + ,fHistFAKEeta(0) + ,fHistMULTeta(0) + ,fHistRECphi(0) + ,fHistMCphi(0) + ,fHistMCNRphi(0) + ,fHistFAKEphi(0) + ,fHistMULTphi(0) + ,fHistRECHPTeta(0) + ,fHistMCHPTeta(0) + ,fHistMCNRHPTeta(0) + ,fHistFAKEHPTeta(0) + ,fHistRECHPTphi(0) + ,fHistMCHPTphi(0) + ,fHistMCNRHPTphi(0) + ,fHistFAKEHPTphi(0) + ,fHistRecMult(0) + ,fHistNCluster(0) + ,fh1VtxEff(0) + ,fh2VtxTpcSpd(0) + ,fh2EtaCorrelation(0) + ,fh2EtaCorrelationShift(0) + ,fh2PhiCorrelation(0) + ,fh2PhiCorrelationShift(0) + ,fh2PtCorrelation(0) + ,fh2PtCorrelationShift(0) + ,fHistDeltaphiprimaries(0) + ,fHistDeltaphisecondaries(0) + ,fHistDeltaphireject(0) + ,fHistTPCITSdeltax(0) + ,fHistTPCITSdeltay(0) + ,fHistTPCITSdeltaz(0) + ,fHistTPCITSdeltar(0) + ,fHistTPCTRDdeltax(0) + ,fHistTPCTRDdeltay(0) + ,fHistTPCTRDdeltaz(0) + ,fHistTPCTRDdeltar(0) + ,fHistTPCTRDdeltaphi(0) + ,fPtPullsPos(0) + ,fPtPullsNeg(0) + ,fPtShiftPos(0) + ,fPtShiftNeg(0) + ,fTrackType(0) + ,fCuts(0) + ,fVertexRvsZ(0) + ,fVertexRvsZC(0) + ,fEtaMultiFluc(0) + ,fEtaMulti(0) + ,fEtaMultiH(0) + ,fPhiMultiFluc(0) + ,fPhiMulti(0) + ,fPhiMultiH(0) +{ + // + // Constructor + // + + for(int i = 0;i< 2;++i){ + fh2MultSpdChips[i] = 0; + } + + for(int i = 0;i < kTotalAC;++i){ + fh2PhiPadRow[i] = fh2PhiLayer[i] = 0; + } + + for(int i = 0; i < kTotalVtx; ++i){ + fh2VertexCorrelation[i] = + fh2VertexCorrelationShift[i] = 0; + fh1VertexShift[i] = + fh1VertexShiftNorm[i] = 0; + } + + for(int i = 0;i< 8;++i){ + fHistRECptCharge[i]=0; + fHistMCptCharge[i]=0; + fHistMCNRptCharge[i]=0; + fHistFAKEptCharge[i]=0; + + fHistRECetaCharge[i]=0; + fHistMCetaCharge[i]=0; + fHistMCNRetaCharge[i]=0; + fHistFAKEetaCharge[i]=0; + + fHistRECphiCharge[i]=0; + fHistMCphiCharge[i]=0; + fHistMCNRphiCharge[i]=0; + fHistFAKEphiCharge[i]=0; + + } + + DefineOutput(1, TList::Class()); +} + + +//________________________________________________________________________ +void AliAnalysisTaskEfficiency::UserCreateOutputObjects() +{ + // Create histograms + + + fHists = new TList(); + fHistRECpt = new TH1F("fHistRECpt", " p_{T} distribution: all reconstructed", 100, 0., 20.); + fHistMCpt = new TH1F("fHistMCpt", " p_{T} distribution: all MC", 100, 0., 20.); + fHistMCNRpt = new TH1F("fHistMCNRpt", " p_{T} distribution: all not-reconstructed MC", 100, 0., 20.); + fHistFAKEpt = new TH1F("fHistFAKEpt", " p_{T} distribution: all Fake", 100, 0., 20.); + fHistMULTpt = new TH1F("fHistMULTpt", " p_{T} distribution: multiply rec.", 100, 0., 20.); + + fHistRECeta = new TH1F("fHistRECeta", " #eta distribution: all reconstructed", 100,-1.0,1.0); + fHistMCeta = new TH1F("fHistMCeta", " #eta distribution: all MC", 100,-1.0,1.0); + fHistMCNReta = new TH1F("fHistMCNReta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0); + fHistFAKEeta = new TH1F("fHistFAKEeta", " #eta distribution: all Fake", 100,-1.0,1.0); + fHistMULTeta = new TH1F("fHistMULTeta", " #eta distribution: multiply rec.", 100,-1.0,1.0); + + fHistRECphi = new TH1F("fHistRECphi", " #phi distribution: all reconstructed", 314, 0., 6.28); + fHistMCphi = new TH1F("fHistMCphi", " #phi distribution: all MC", 314, 0., 6.28); + fHistMCNRphi = new TH1F("fHistMCNRphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28); + fHistFAKEphi = new TH1F("fHistFAKEphi", " #phi distribution: all Fake", 314, 0., 6.28); + fHistMULTphi = new TH1F("fHistMULTphi", " #phi distribution: multipli rec.", 314, 0., 6.28); + + fHistRECHPTeta = new TH1F("fHistRECHPTeta", " #eta distribution: all reconstructed", 100,-1.0,1.0); + fHistMCHPTeta = new TH1F("fHistMCHPTeta", " #eta distribution: all MC", 100,-1.0,1.0); + fHistMCNRHPTeta = new TH1F("fHistMCNRHPTeta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0); + fHistFAKEHPTeta = new TH1F("fHistFAKEHPTeta", " #eta distribution: all Fake", 100,-1.0,1.0); + + fHistRECHPTphi = new TH1F("fHistRECHPTphi", " #phi distribution: all reconstructed", 314, 0., 6.28); + fHistMCHPTphi = new TH1F("fHistMCHPTphi", " #phi distribution: all MC", 314, 0., 6.28); + fHistMCNRHPTphi = new TH1F("fHistMCNRHPTphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28); + fHistFAKEHPTphi = new TH1F("fHistFAKEHPTphi", " #phi distribution: all Fake", 314, 0., 6.28); + + fHistRecMult = new TH1F("fHistRecMult", "Multiple reconstructed tracks", 50, 0., 50.); + fHistNCluster = new TH1F("fHistNCluster", "Number of clusters for suspicious tracks", 300, 0., 300.); + + // CKB + + fh1VtxEff = new TH1F("fh1VtxEff","VtxEff TPC",4,-0.5,3.5); + fh1VtxEff->GetXaxis()->SetBinLabel(1,"NO TPC VTX"); + fh1VtxEff->GetXaxis()->SetBinLabel(2,"TPC VTX"); + fh1VtxEff->GetXaxis()->SetBinLabel(3,"NO SPD VTX"); + fh1VtxEff->GetXaxis()->SetBinLabel(4,"SPD VTX"); + + TString labels[kTotalAC]; + labels[kNegA]="NegA"; + labels[kPosA]="PosA"; + labels[kNegC]="NegC"; + labels[kPosC]="PosC"; + + for(int i = 0;i< kTotalAC;++i){ + fh2PhiPadRow[i] = new TH2F(Form("fh2PhiPadRow_%d",i),Form("Padrow vs phi AC %s;phi;padrow",labels[i].Data()),360,0.,360.,159,-0.5,158.5); + fh2PhiLayer[i] = new TH2F(Form("fh2PhiLayer_%d",i),Form("layer vs phi AC %s;phi;layer",labels[i].Data()),360,0.,360.,6,-0.5,5.5); + } + for(int i = 0;i<2;++i){ + fh2MultSpdChips[i] = new TH2F(Form("fh2ChipsSpdMult_%d",i),"mult SPD vs chips;chips;Mult",201,-1.5,199.5,201,-1.5,199.5); + } + fh2VtxTpcSpd = new TH2F("fh2VtxTpcSpd","SPD vtx vs TPC;tpc-z;spd-z",200,-20.,20.,200,-20,20.); + + TString labelsVtx[kTotalVtx]; + labelsVtx[kTPC] = "TPC"; + labelsVtx[kSPD] = "SPD"; + for(int i = 0; i < kTotalVtx; ++i){ + fh2VertexCorrelation[i] = new TH2F(Form("fh2VertexCorrelation_%d",i),Form("VertexCorrelation %s;MC z-vtx;ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 120, -30, 30); + fh2VertexCorrelationShift[i] = new TH2F(Form("fh2VertexCorrelationShift_%d",i), Form("VertexCorrelationShift %s;MC z-vtx;MC z-vtx - ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 100, -1, 1); + fh1VertexShift[i] = new TH1F(Form("fh1VertexShift_%d",i), Form("VertexShift %s;(MC z-vtx - ESD z-vtx);Entries",labelsVtx[i].Data()), 201, -2, 2); + fh1VertexShiftNorm[i] = new TH1F(Form("fh1VertexShiftNorm_%d",i), Form("VertexShiftNorm %s;(MC z-vtx - ESD z-vtx) / #sigma_{ESD z-vtx};Entries",labelsVtx[i].Data()), 200, -100, 100); + } + + + fh2EtaCorrelation = new TH2F("fh2EtaCorrelation", "EtaCorrelation;MC #eta;ESD #eta", 120, -3, 3, 120, -3, 3); + fh2EtaCorrelationShift = new TH2F("fh2EtaCorrelationShift", "EtaCorrelationShift;ESD #eta;MC #eta - ESD #eta", 120, -3, 3, 100, -0.1, 0.1); + fh2PhiCorrelation = new TH2F("fh2PhiCorrelation", "PhiCorrelation;MC #phi;ESD #phi", 120, 0, 360, 120, 0, 360); + fh2PhiCorrelationShift = new TH2F("fh2PhiCorrelationShift", "PhiCorrelationShift;ESD #phi;MC #phi - ESD #phi", 120, 0, 360, 100, -10, 10); + fh2PtCorrelation = new TH2F("fh2PtCorrelation", "PtCorrelation;MC p_T;ESD p_T", 200, 0., 200., 200, 0, 200); + fh2PtCorrelationShift = new TH2F("fh2PtCorrelationShift", "PtCorrelationShift;ESD p_T;MC p_T - ESD #p_T", 120, 0, 10, 100, -2, 2); + + fHistDeltaphiprimaries =new TH1F("fHistDeltaphiprimaries", " #Delta#phi distribution: primaries", 314, -0.2,0.2); + fHistDeltaphisecondaries=new TH1F("fHistDeltaphisecondaries", "#Delta#phi distribution: secondaries", 314, -0.2,0.2); + fHistDeltaphireject =new TH1F("fHistDeltaphireject", " #Delta#phi distribution: reject trackled", 314, -0.2,0.2); + + fHistTPCITSdeltax =new TH1F("fHistTPCITSdeltax", "TPC-ITS matching:#Delta x", 100,-20,20); + fHistTPCITSdeltay =new TH1F("fHistTPCITSdeltay", "TPC-ITS matching:#Delta y", 100,-20,20); + fHistTPCITSdeltaz =new TH1F("fHistTPCITSdeltaz", "TPC-ITS matching:#Delta z", 100,-20,20); + fHistTPCITSdeltar =new TH1F("fHistTPCITSdeltar", "TPC-ITS matching:#Delta r", 100,0,20); + + fHistTPCTRDdeltax = new TH1F("fHistTPCTRDdeltax", "TPC-TRD matching:#Delta x", 400,-400, 400); + fHistTPCTRDdeltay = new TH1F("fHistTPCTRDdeltay", "TPC-TRD matching:#Delta y", 400,-400, 400); + fHistTPCTRDdeltaz = new TH1F("fHistTPCTRDdeltaz", "TPC-TRD matching:#Delta z", 400,-400, 400); + fHistTPCTRDdeltar = new TH1F("fHistTPCTRDdeltar", "TPC-TRD matching:#Delta r", 100,0,20); + fHistTPCTRDdeltaphi = new TH1F("fHistTPCTRDdeltarphi", "TPC-TRD matching:#Delta #phi", 128,-3.14 , 3.14); + // Pulls + fPtPullsPos = new TProfile("fPtPullsPos", "Pt Pulls for pos. primaries", 50 , 0., 10., -5., 5.,"S"); + fPtPullsNeg = new TProfile("fPtPullsNeg", "Pt Pulls for neg. primaries", 50 , 0., 10., -5., 5.,"S"); + fPtShiftPos = new TProfile("fPtShiftPos", "Pt Shift for pos. primaries", 50 , 0., 10., -0.2, 0.2,"S"); + fPtShiftNeg = new TProfile("fPtShiftNeg", "Pt Shift for neg. primaries", 50 , 0., 10., -0.2, 0.2,"S"); + + fEtaMultiFluc = new TProfile("fEtaMultiFluc", "eta multiplicity fluctuations", 10, -2., 2., 0., 600., "S"); + fEtaMulti = new TH1F("fEtaMulti", "eta multiplicity fluctuations", 10, -2., 2.); + fEtaMultiH = new TH1F("fEtaMultiH", "eta multiplicity fluctuations", 600, 0., 600.); + + fPhiMultiFluc = new TProfile("fPhiMultiFluc", "phi multiplicity fluctuations", 64, 0., 6.4, 0., 100., "S"); + fPhiMulti = new TH1F("fPhiMulti", "phi multiplicity fluctuations", 64, 0., 6.4); + fPhiMultiH = new TH1F("fPhiMultiH", "phi multiplicity fluctuations", 100, 0., 100.); + + // SPD Vertex + fVertexRvsZ = new TH2F("fVertexRvsZ", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.); + fVertexRvsZC = new TH2F("fVertexRvsZC", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.); + + TString charge[8]={"Deu", "AntiDeu", "Tri", "AntiTri", "He3", "AntiHe3", "He4", "AntiHe4"}; + for(Int_t i=0;i<8;i++){ + fHistRECptCharge[i] = new TH1F(Form("fHistRECptCharge%s",charge[i].Data()), + "p_{T} distribution: all reconstructed", + 100, 0., 20.); + fHistMCptCharge[i] = new TH1F(Form("fHistMCptCharge%s",charge[i].Data()), + "p_{T} distribution: all MC", + 100, 0., 20.); + fHistMCNRptCharge[i] = new TH1F(Form("fHistMCNRptCharge%s",charge[i].Data()), + "p_{T} distribution: all not-reconstructed MC", + 100, 0., 20.); + fHistFAKEptCharge[i] = new TH1F( Form("fHistFAKEptCharge%s",charge[i].Data()), + "p_{T} distribution: all Fake", + 100, 0., 20.); + + fHistRECetaCharge[i] = new TH1F(Form("fHistRECetaCharge%s",charge[i].Data()), + "p_{T} distribution: all reconstructed", + 100, 0., 20.); + fHistMCetaCharge[i] = new TH1F(Form("fHistMCetaCharge%s",charge[i].Data()), + "p_{T} distribution: all MC", + 100, 0., 20.); + fHistMCNRetaCharge[i] = new TH1F(Form("fHistMCNRetaCharge%s",charge[i].Data()), + "p_{T} distribution: all not-reconstructed MC", + 100, 0., 20.); + fHistFAKEetaCharge[i] = new TH1F( Form("fHistFAKEetaCharge%s",charge[i].Data()), + "p_{T} distribution: all Fake", + 100, 0., 20.); + + fHistRECphiCharge[i] = new TH1F(Form("fHistRECphiCharge%s",charge[i].Data()), + "p_{T} distribution: all reconstructed", + 100, 0., 20.); + fHistMCphiCharge[i] = new TH1F(Form("fHistMCphiCharge%s",charge[i].Data()), + "p_{T} distribution: all MC", + 100, 0., 20.); + fHistMCNRphiCharge[i] = new TH1F(Form("fHistMCNRphiCharge%s",charge[i].Data()), + "p_{T} distribution: all not-reconstructed MC", + 100, 0., 20.); + fHistFAKEphiCharge[i] = new TH1F( Form("fHistFAKEphiCharge%s",charge[i].Data()), + "p_{T} distribution: all Fake", + 100, 0., 20.); + + } + // BKC + + + fHists->SetOwner(); + + fHists->Add(fHistRECpt); + fHists->Add(fHistMCpt); + fHists->Add(fHistMCNRpt); + fHists->Add(fHistFAKEpt); + fHists->Add(fHistMULTpt); + + fHists->Add(fHistRECeta); + fHists->Add(fHistMCeta); + fHists->Add(fHistMCNReta); + fHists->Add(fHistFAKEeta); + fHists->Add(fHistMULTeta); + + fHists->Add(fHistRECphi); + fHists->Add(fHistMCphi); + fHists->Add(fHistMCNRphi); + fHists->Add(fHistFAKEphi); + fHists->Add(fHistMULTphi); + + fHists->Add(fHistRECHPTeta); + fHists->Add(fHistMCHPTeta); + fHists->Add(fHistMCNRHPTeta); + fHists->Add(fHistFAKEHPTeta); + + fHists->Add(fHistRECHPTphi); + fHists->Add(fHistMCHPTphi); + fHists->Add(fHistMCNRHPTphi); + fHists->Add(fHistFAKEHPTphi); + + // CKB + + fHists->Add(fh1VtxEff); + for(int i = 0;i < kTotalAC;++i){ + fHists->Add(fh2PhiPadRow[i]); + fHists->Add(fh2PhiLayer[i]); + } + for(int i = 0;i<2;++i){ + fHists->Add(fh2MultSpdChips[i]); + } + fHists->Add(fh2VtxTpcSpd); + + for(int i = 0; i < kTotalVtx; ++i){ + fHists->Add(fh2VertexCorrelation[i]); + fHists->Add(fh2VertexCorrelationShift[i]); + fHists->Add(fh1VertexShift[i]); + fHists->Add(fh1VertexShiftNorm[i]); + } + + + fHists->Add(fh2EtaCorrelation); + fHists->Add(fh2EtaCorrelationShift); + fHists->Add(fh2PhiCorrelation); + fHists->Add(fh2PhiCorrelationShift); + fHists->Add(fh2PtCorrelation); + fHists->Add(fh2PtCorrelationShift); + + fHists->Add(fHistDeltaphiprimaries); + fHists->Add(fHistDeltaphisecondaries); + fHists->Add(fHistDeltaphireject); + + fHists->Add(fHistTPCITSdeltax); + fHists->Add(fHistTPCITSdeltay); + fHists->Add(fHistTPCITSdeltaz); + fHists->Add(fHistTPCITSdeltar); + + fHists->Add(fHistTPCTRDdeltax); + fHists->Add(fHistTPCTRDdeltay); + fHists->Add(fHistTPCTRDdeltaz); + fHists->Add(fHistTPCTRDdeltar); + fHists->Add(fHistTPCTRDdeltaphi); + + fHists->Add(fHistRecMult); + fHists->Add(fHistNCluster); + + + fHists->Add(fPtPullsPos); + fHists->Add(fPtPullsNeg); + fHists->Add(fPtShiftPos); + fHists->Add(fPtShiftNeg); + + fHists->Add(fVertexRvsZ); + fHists->Add(fVertexRvsZC); + fHists->Add(fEtaMultiFluc); + fHists->Add(fEtaMulti); + fHists->Add(fEtaMultiH); + fHists->Add(fPhiMultiFluc); + fHists->Add(fPhiMulti); + fHists->Add(fPhiMultiH); + + for(Int_t i=0;i<8;i++){ + fHists->Add(fHistRECptCharge[i]); + fHists->Add(fHistMCptCharge[i]); + fHists->Add(fHistMCNRptCharge[i]); + fHists->Add(fHistFAKEptCharge[i]); + + fHists->Add(fHistRECetaCharge[i]); + fHists->Add(fHistMCetaCharge[i]); + fHists->Add(fHistMCNRetaCharge[i]); + fHists->Add(fHistFAKEetaCharge[i]); + + fHists->Add(fHistRECphiCharge[i]); + fHists->Add(fHistMCphiCharge[i]); + fHists->Add(fHistMCNRphiCharge[i]); + fHists->Add(fHistFAKEphiCharge[i]); + } + + for (Int_t i=0; iGetEntries(); ++i) { + TH1 *h1 = dynamic_cast(fHists->At(i)); + if (h1){ + h1->Sumw2(); + } + } + // BKC + + + // Post output data. + PostData(1, fHists); + +} + +//________________________________________________________________________ +void AliAnalysisTaskEfficiency::UserExec(Option_t *) +{ + // Event loop + // Analysis of MC true particles and reconstructed tracks + // Different track types (Global, TPC, ITS) can be selected + + if (!fInputEvent) { + Printf("ERROR: fESD not available"); + return; + } + static Int_t nfc = 0; + + //AliESDInputHandler* esdI = (AliESDInputHandler*) + //(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler(); + //AliESDEvent* hltEvent = esdI->GetHLTEvent(); + + AliStack* stack = MCEvent()->Stack(); + AliESDEvent* esdE = (AliESDEvent*) fInputEvent; + + // CKB + // Fetch the SPD vertex: + Double_t vtxSPD[3]; + const AliESDVertex* spdVertex = esdE->GetPrimaryVertexSPD(); + if (spdVertex) { + if(spdVertex->GetNContributors()<=0){ + spdVertex = 0; + } + else{ + spdVertex->GetXYZ(vtxSPD); + } + } + + // Fetch the TPC vertex + Double_t vtxTPC[3]; + const AliESDVertex* tpcVertex = esdE->GetPrimaryVertexTPC(); + if (tpcVertex) { + if(tpcVertex->GetNContributors()<=0){ + tpcVertex = 0; + } + else{ + tpcVertex->GetXYZ(vtxTPC); + } + } + // SPD Vertex + if (spdVertex) { + Double_t x,y,z; + x = spdVertex->GetX() + 0.07; + y = spdVertex->GetY() - 0.25; + z = spdVertex->GetZ(); + if (TMath::Abs(z) < 10.) { + fVertexRvsZ->Fill(x,y); + if (TMath::Sqrt(x*x + y*y) > 5. * 0.028) + fVertexRvsZC->Fill(x,y); + } + } + // BKC + //Printf("%s:%d %5d",(char*)__FILE__,__LINE__, esdE->GetNumberOfTracks()); + TArrayI labels(esdE->GetNumberOfTracks()); + Int_t igood = 0; + // Track loop to fill a pT spectrum + //printf("ESD track loop \n"); + + AliESDtrack *tpcP = 0x0; + + for (Int_t iTracks = 0; iTracks < esdE->GetNumberOfTracks(); iTracks++) { + + //prevent mem leak for TPConly track + if(fTrackType==2&&tpcP){ + delete tpcP; + tpcP = 0; + } + + AliVParticle *track = esdE->GetTrack(iTracks); + AliESDtrack *esdtrack = dynamic_cast(track); + esdtrack->PropagateToDCA(esdE->GetPrimaryVertex(), + esdE->GetMagneticField(), 10000.); + + if (!track) { + Printf("ERROR: Could not receive track %d", iTracks); + continue; + } + //__________ + // run Task for global tracks or ITS tracks or TPC tracks + const AliExternalTrackParam *tpcPin = 0x0; + //Double_t phiIn=0.; + Float_t phi = 0.; + Float_t eta = 0.; + Float_t pt = 0.; + Int_t indexAC = GetIndexAC(esdtrack); + + //select in the steering macro, which track type should be analysed. + // Four track types are selectable: + // 0. Global tracks + // 1. ITS tracks (SA or Pure SA) + // 2. TPC only tracks + // Track selection copied from PWG1/AliAnalysisTaskQAsym + + if(fTrackType==0){ + //Fill all histograms with global tracks + tpcP = esdtrack; + if (!tpcP) continue; + if (!fCuts->AcceptTrack(tpcP)) continue; + phi = tpcP->Phi(); + eta = tpcP->Eta(); + pt = tpcP->Pt(); + } + else if(fTrackType==1){ + //Fill all histograms with ITS tracks + tpcP = esdtrack; + if (!tpcP) continue; + if (!fCuts->AcceptTrack(tpcP)) continue; + phi = tpcP->Phi(); + eta = tpcP->Eta(); + pt = tpcP->Pt(); + } + else if(fTrackType==2){ + //Fill all histograms with TPC track information + tpcPin = esdtrack->GetInnerParam(); + if (!tpcPin) continue; + tpcP = AliESDtrackCuts:: + GetTPCOnlyTrack(dynamic_cast(esdE), + esdtrack->GetID()); + if (!tpcP) continue; + if (!fCuts->AcceptTrack(tpcP)) continue; + if(tpcP->GetNcls(1)>160)continue;//jacek's track cut + if(tpcP->GetConstrainedChi2TPC()<0)continue; // jacek's track cut + phi=tpcPin->Phi(); + eta=tpcPin->Eta(); + pt=tpcPin->Pt(); + } + else{ + Printf("ERROR: wrong track type \n"); + continue; + } + //___________ + // + + + labels.AddAt(-1, iTracks); + + + // Selected + if (TMath::Abs(eta) > 0.9) { + } else { + + //Int_t charge=track->Charge(); + // Within acceptance + Int_t ind = TMath::Abs(esdtrack->GetLabel()); + AliMCParticle* mcPtest = (AliMCParticle*) MCEvent()->GetTrack(ind); + if (stack->IsPhysicalPrimary(mcPtest->Label())){ + if(TMath::Abs(mcPtest->PdgCode())>=1000020030){//all helium (+-1000020030,+-1000020040) + pt*=2;// reconstruction takes charge=1 -> for particles with charge=2, pT,rec need to be corrected + } + + Int_t index=ConvertePDG(mcPtest->PdgCode()); + if(fDebug>1)Printf("PDG=%d, index=%d", mcPtest->PdgCode(), index); + //fill tracks comming from d,t,3He, 4He + if(index<8){ + fHistRECptCharge [index] ->Fill(pt); + fHistRECetaCharge[index]->Fill(eta); + fHistRECphiCharge[index]->Fill(phi); + } + + } + fHistRECpt ->Fill(pt); + fHistRECeta->Fill(eta); + fHistRECphi->Fill(phi); + + if (pt > 2.) { + fHistRECHPTeta->Fill(eta); + fHistRECHPTphi->Fill(phi); + } + + // Int_t ind = TMath::Abs(esdtrack->GetLabel()); + AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(ind); + if (!(stack->IsPhysicalPrimary(mcP->Label()))) { + fHistFAKEpt ->Fill(pt); + fHistFAKEeta->Fill(eta); + fHistFAKEphi->Fill(phi); + + if (pt > 2.) { + fHistFAKEHPTeta->Fill(eta); + fHistFAKEHPTphi->Fill(phi); + } + + } + + labels.AddAt(TMath::Abs(esdtrack->GetLabel()), iTracks); + igood++; + + + Float_t phiDeg = TMath::RadToDeg()*track->Phi(); + + //TPC-ITS matching + //TPC-TRD matching + + if (tpcP){ + Int_t labeltpcits = esdtrack->GetLabel(); + AliExternalTrackParam * tpcPCopy = new AliExternalTrackParam(*tpcP); + Double_t xk = 43.6; // cm + Double_t bz = 5.0; // kG + if(tpcPCopy->PropagateTo(xk,bz)){ + Double_t alpha=tpcPCopy->GetAlpha(); + // if(tpcPCopy->Rotate(0.)){ + Float_t xtpc=tpcPCopy->GetX(); + Float_t ytpc=tpcPCopy->GetY(); + Float_t ztpc=tpcPCopy->GetZ(); + Float_t xpr,ypr ; + xpr=xtpc*TMath::Cos(alpha)-ytpc*TMath::Sin(alpha); + ypr=xtpc*TMath::Sin(alpha)+ytpc*TMath::Cos(alpha); + + AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpcits)); + Int_t ntref = mcPart->GetNumberOfTrackReferences(); + + for (Int_t k = 0; k < ntref; k++) { + AliTrackReference* ref = mcPart->GetTrackReference(k); + Float_t xhits = ref->X(); + Float_t yhits = ref->Y(); + Float_t radio = TMath::Sqrt(xhits*xhits+yhits*yhits); + + if(ref->DetectorId() == 0 && (radio > 42.8 && radio < 43.7)) { + + Float_t xits = ref->X(); + Float_t yits = ref->Y(); + Float_t zits = ref->Z(); + + Float_t difx=(xits-xpr); + fHistTPCITSdeltax->Fill(difx); + Float_t dify=(yits-ypr); + fHistTPCITSdeltay->Fill(dify); + Float_t difz=(zits-ztpc); + fHistTPCITSdeltaz->Fill(difz); + Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz); + fHistTPCITSdeltar->Fill(deltar); + break; + } + // } + } // trackrefs + } + } //ITS-TPC maching + //TPC-TRD maching + + const AliExternalTrackParam *trd = esdtrack->GetOuterParam(); + + if (trd){Int_t labeltpctrd = track->GetLabel(); + + AliExternalTrackParam * trdCopy = new AliExternalTrackParam(*trd); + Double_t xktrd=296.0; + Double_t bztrd=5.0; + if(trdCopy->PropagateTo(xktrd,bztrd)){ + Float_t xtpc2=trdCopy->GetX(); + Float_t ytpc2=trdCopy->GetY(); + Float_t ztpc2=trdCopy->GetZ(); + Double_t alpha=trdCopy->GetAlpha(); + Float_t xpr,ypr ; + xpr=xtpc2*TMath::Cos(alpha)-ytpc2*TMath::Sin(alpha); + ypr=xtpc2*TMath::Sin(alpha)+ytpc2*TMath::Cos(alpha); + AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpctrd)); + Int_t ntreftrd = mcPart->GetNumberOfTrackReferences(); + + for (Int_t k = 0; k < ntreftrd; k++) { + AliTrackReference* ref = mcPart->GetTrackReference(k); + if(ref->DetectorId() == 3 && ref->Label() == labeltpctrd){ + + Float_t xtrd = ref->X(); + Float_t ytrd = ref->Y(); + Float_t ztrd = ref->Z(); + + Float_t difx=(xtrd-xpr); + fHistTPCTRDdeltax->Fill(difx); + Float_t dify=(ytrd-ypr); + fHistTPCTRDdeltay->Fill(dify); + Float_t difz=(ztrd-ztpc2); + fHistTPCTRDdeltaz->Fill(difz); + Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz); + fHistTPCTRDdeltar->Fill(deltar); + Float_t phi_tpc = TMath::ATan2(ypr, xpr); + Float_t phi_trd = TMath::ATan2(ytrd, xtrd); + Float_t dphi = phi_trd - phi_tpc; + if (dphi > TMath::Pi()) dphi -= 2. * TMath::Pi(); + if (dphi < - TMath::Pi()) dphi += 2. * TMath::Pi(); + + fHistTPCTRDdeltaphi->Fill(dphi); + break; + } + + } + } + + }//TRD-TPC maching + + // CKB + if(pt>2.&&fFieldOn){ + TBits bits(esdtrack->GetTPCClusterMap()); + for(unsigned int ip = 0;ip < bits.GetNbits();++ip){ + if(bits[ip]){ + fh2PhiPadRow[indexAC]->Fill(phiDeg,ip); + } + } + for(int i = 0;i < 6;++i){ + if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i); + } + } + else if(!fFieldOn){ // field not on all momenta are set to MPV + TBits bits(esdtrack->GetTPCClusterMap()); + for(unsigned int ip = 0;ip < bits.GetNbits();++ip){ + if(bits[ip]){ + fh2PhiPadRow[indexAC]->Fill(phiDeg,ip); + } + for(int i = 0;i < 6;++i){ + if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i); + } + } + + } + // Fill the correlation + if(MCEvent()){ + TParticle *part = MCEvent()->Stack()->Particle(TMath::Abs(esdtrack->GetLabel())); + if(part){ + Float_t mcEta = part->Eta(); + Float_t mcPhi = TMath::RadToDeg()*part->Phi(); + Float_t mcPt = part->Pt(); + // if (pt - mcPt > 20.) { + fh2EtaCorrelation->Fill(mcEta,eta); + fh2EtaCorrelationShift->Fill(eta,mcEta-eta); + fh2PhiCorrelation->Fill(mcPhi,phiDeg); + fh2PhiCorrelationShift->Fill(phiDeg,mcPhi-phiDeg); + fh2PtCorrelation->Fill(mcPt,pt); + fh2PtCorrelationShift->Fill(pt,mcPt-pt); + //} + Double_t sigmaPt = TMath::Sqrt(esdtrack->GetSigma1Pt2()); + if (track->Charge() > 0.) { + fPtPullsPos->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt); + fPtShiftPos->Fill(mcPt, (1./pt - 1./mcPt)); + } else { + fPtPullsNeg->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt); + fPtShiftNeg->Fill(mcPt, (1./pt - 1./mcPt)); + } + } + // BKC + } + } + + } //track loop + + //prevent mem leak for TPConly track + if(fTrackType==2&&tpcP){ + delete tpcP; + tpcP = 0; + } + + //Printf("%s:%d",(char*)__FILE__,__LINE__); + // CKB + // Vertex eff + if(tpcVertex){ + fh1VtxEff->Fill("TPC VTX",1); + } + else{ + fh1VtxEff->Fill("NO TPC VTX",1); + } + + if(spdVertex){ + fh1VtxEff->Fill("SPD VTX",1); + } + else{ + fh1VtxEff->Fill("NO SPD VTX",1); + } + + + // Vertex correlation SPD vs. TPC + if(tpcVertex&&spdVertex){ + fh2VtxTpcSpd->Fill(vtxTPC[2],vtxSPD[2]); + } + // Printf("%s:%d",(char*)__FILE__,__LINE__); + // Multiplicity checks in the SPD + const AliMultiplicity *spdMult = esdE->GetMultiplicity(); + // Multiplicity Analysis + Int_t nt = spdMult->GetNumberOfTracklets(); + nfc++; + if (nfc == 520) { + nfc = 0; + for (Int_t ib = 1; ib <= 10; ib++) { + Int_t ic = Int_t(fEtaMulti->GetBinContent(ib)); + Float_t eta = -1.8 + Float_t(ib-1) * 0.4; + fEtaMultiFluc->Fill(eta, Float_t(ic)); + if (ib == 5) fEtaMultiH->Fill(Float_t(ic)); + } + + for (Int_t ib = 1; ib <= 64; ib++) { + Int_t ic = Int_t(fPhiMulti->GetBinContent(ib)); + Float_t phi = 0.05 + Float_t(ib-1) * 0.1; + fPhiMultiFluc->Fill(phi, Float_t(ic)); + if (ib == 2) fPhiMultiH->Fill(Float_t(ic)); + } + + fEtaMulti->Reset(); + fPhiMulti->Reset(); + } + + for (Int_t j = 0; j < nt; j++) { + fEtaMulti->Fill(spdMult->GetEta(j)); + fPhiMulti->Fill(spdMult->GetPhi(j)); + } + + Short_t chips0 = spdMult->GetNumberOfFiredChips(0); + Short_t chips1 = spdMult->GetNumberOfFiredChips(1); + //Printf("%s:%d",(char*)__FILE__,__LINE__); + Int_t inputCount = 0; + + if(spdVertex){ + // get multiplicity from ITS tracklets + for (Int_t i=0; iGetNumberOfTracklets(); ++i){ + Float_t deltaPhi = spdMult->GetDeltaPhi(i); + // prevent values to be shifted by 2 Pi() + if (deltaPhi < -TMath::Pi()) + deltaPhi += TMath::Pi() * 2; + if (deltaPhi > TMath::Pi()) + deltaPhi -= TMath::Pi() * 2; + if (TMath::Abs(deltaPhi) > 1) + printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, spdMult->GetTheta(i), spdMult->GetPhi(i), deltaPhi); + //trackled Deltaphi + Int_t label1=spdMult->GetLabel(i,0); + Int_t label2=spdMult->GetLabel(i,1); + if ( label1==label2){ + if(stack->IsPhysicalPrimary(label1) == kTRUE) + fHistDeltaphiprimaries->Fill(deltaPhi); + else + fHistDeltaphisecondaries->Fill(deltaPhi); + } + else{ + fHistDeltaphireject->Fill(deltaPhi); + } + ++inputCount; + } + // Printf("%s:%d",(char*)__FILE__,__LINE__); + fh2MultSpdChips[0]->Fill(chips0,inputCount); + fh2MultSpdChips[1]->Fill(chips1,inputCount); + //Printf("%s:%d",(char*)__FILE__,__LINE__); + } + // BKC + + // MC analysis + Int_t nmc = MCEvent()->GetNumberOfTracks(); + AliGenEventHeader* header = MCEvent()->GenEventHeader(); + // some test + AliGenCocktailEventHeader* cH = dynamic_cast (header); + AliGenPythiaEventHeader* pH; + if (cH == 0) + { + header->Dump(); + pH = dynamic_cast(header); + } else { + TObject* entry = cH->GetHeaders()->FindObject("Pythia"); + pH = dynamic_cast(entry); + } + Int_t iproc = -1; + + if (pH) iproc = pH->ProcessType(); + + TArrayF mcV; + header->PrimaryVertex(mcV); + Float_t vzMC = mcV[2]; + + // CKB + // Fill the correlation with MC + if(spdVertex&&MCEvent()){ + Float_t diff = vzMC-vtxSPD[2]; + fh2VertexCorrelation[kSPD]->Fill(vzMC,vtxSPD[2]); + fh2VertexCorrelationShift[kSPD]->Fill(vzMC,diff); + fh1VertexShift[kSPD]->Fill(diff); + if(spdVertex->GetZRes()>0)fh1VertexShiftNorm[kSPD]->Fill(diff/spdVertex->GetZRes()); + } + if(tpcVertex&&MCEvent()){ + Float_t diff = vzMC-vtxTPC[2]; + fh2VertexCorrelation[kTPC]->Fill(vzMC,vtxTPC[2]); + fh2VertexCorrelationShift[kTPC]->Fill(vzMC,diff); + fh1VertexShift[kTPC]->Fill(diff); + if(tpcVertex->GetZRes()>0)fh1VertexShiftNorm[kTPC]->Fill(diff/tpcVertex->GetZRes()); + } + // BKC + + + // MC event loop + //printf("MC particle loop \n"); + for (Int_t i = 0; i < nmc; i++) { + AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(i); + //printf("MC particle loop %5d \n", i); + // Primaries only + if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue; + if (mcP->Particle()->GetPDG()->Charge() == 0) continue; + // Int_t charge=Int_t(mcP->Particle()->GetPDG()->Charge()); + Float_t phi = mcP->Phi(); + Float_t eta = mcP->Eta(); + Float_t pt = mcP->Pt(); + if (TMath::Abs(eta) > 0.9) continue; + Int_t ntr = mcP->GetNumberOfTrackReferences(); + Int_t nITS = 0; + Int_t nTPC = 0; + Int_t nFRA = 0; + Float_t x = 0.; + Float_t y = 0.; + + for (Int_t j = 0; j < ntr; j++) { + + AliTrackReference* ref = mcP->GetTrackReference(j); + + if (ref->DetectorId() == 0) nITS++; + if (ref->DetectorId() == 1) nTPC++; + if (ref->DetectorId() == 2) nFRA++; + if (nTPC == 1) { + x = ref->X(); + y = ref->Y(); + break; + } + } + + fHistMCpt ->Fill(pt); + fHistMCeta->Fill(eta); + fHistMCphi->Fill(phi); + Int_t index=ConvertePDG(mcP->PdgCode()); + if(index<8){ + fHistMCptCharge [index] ->Fill(pt); + fHistMCetaCharge[index]->Fill(eta); + fHistMCphiCharge[index]->Fill(phi); + } + + + if (pt > 2.) { + fHistMCHPTeta->Fill(eta); + fHistMCHPTphi->Fill(phi); + } + + Bool_t reco = kFALSE; + Int_t multR = 0; + Int_t jold = -1; + + for (Int_t j = 0; j < esdE->GetNumberOfTracks(); j++) { + if (i == labels.At(j)) { + reco = kTRUE; + multR++; + //AliESDtrack* test = esdE->GetTrack(j); + if (multR > 1) { + Int_t nclus = 0; + AliESDtrack* track = esdE->GetTrack(jold); + nclus = track->GetTPCclusters(0); + fHistNCluster->Fill(nclus); + + + track = esdE->GetTrack(j); + nclus = track->GetTPCclusters(0); + fHistNCluster->Fill(nclus); + + + } + jold = j; + } + } + + fHistRecMult->Fill(Float_t(multR), 1.); + if (multR == 0) { + fHistMCNRpt ->Fill(pt); + fHistMCNReta->Fill(eta); + fHistMCNRphi->Fill(phi); + if(index<8){ + fHistMCNRptCharge [index] ->Fill(pt); + fHistMCNRetaCharge[index]->Fill(eta); + fHistMCNRphiCharge[index]->Fill(phi); + } + + + if (pt > 2.) { + fHistMCNRHPTeta->Fill(eta); + fHistMCNRHPTphi->Fill(phi); + } + } + else if (multR > 1) + { + fHistMULTpt ->Fill(pt); + fHistMULTeta->Fill(eta); + fHistMULTphi->Fill(phi); + } + } // MC particle loop + //printf("End of MC particle loop \n"); + + +} + +//________________________________________________________________________ +void AliAnalysisTaskEfficiency::Terminate(Option_t *) +{ + +} + + +Bool_t AliAnalysisTaskEfficiency::SelectJouri(AliESDtrack* track) +{ + + Bool_t selected = kTRUE; + AliESDEvent* esdE = (AliESDEvent*) fInputEvent; + // > 50 Clusters + if (track->GetTPCclusters(0) < 50) selected = kFALSE; + + const AliESDVertex *vtx=esdE->GetPrimaryVertexSPD(); + if (!vtx->GetStatus()) selected = kFALSE; + + Double_t zv=vtx->GetZv(); + + + const AliExternalTrackParam *ct=track->GetTPCInnerParam(); + if (!ct) selected = kFALSE; + + Double_t xyz[3]; + ct->GetXYZ(xyz); + Float_t rv = TMath::Sqrt(xyz[0] * xyz[0] + xyz[1] * xyz[1]); + if (rv > 3.0) selected = kFALSE; + if (TMath::Abs(xyz[2] - zv)>2.5) selected = kFALSE; + return (selected); + +} + +Int_t AliAnalysisTaskEfficiency::GetIndexAC(AliESDtrack *track){ + if(!track)return -1; + + // Crude selection for A and C side + // just with eta + if (track->Eta() > 0) { + if (track->Charge() > 0) + return kPosA; + else + return kNegA; + } + else { // C side + if (track->Charge() > 0) + return kPosC; + else + return kNegC; + } + + return -1; + +} + + +Int_t AliAnalysisTaskEfficiency::ConvertePDG(Int_t pdg) +{ + + // function converts the pdg values for nuclei (d, t, 3He, + // 4He and their antiparticles) into indices for histo-array + // filled in UserExec + + if ( pdg == 1000010020 )return 0; + else if ( pdg == -1000010020 )return 1; + else if ( pdg == 1000010030 )return 2; + else if ( pdg == -1000010030 )return 3; + else if ( pdg == 1000020030 )return 4; + else if ( pdg == -1000020030 )return 5; + else if ( pdg == 1000020040 )return 6; + else if ( pdg == -1000020040 )return 7; + else return 9; + + +} diff --git a/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.h b/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.h new file mode 100644 index 00000000000..d5179e8c230 --- /dev/null +++ b/test/vmctest/scripts/efficiency/AliAnalysisTaskEfficiency.h @@ -0,0 +1,149 @@ +#ifndef AliAnalysisTaskEfficiency_cxx +#define AliAnalysisTaskEfficiency_cxx + +class TH1F; +class TH2F; +class TList; +class TProfile; + +class AliESDEvent; +class AliESDtrack; +class AliESDtrackCuts; + + +#include "AliAnalysisTaskSE.h" + +class AliAnalysisTaskEfficiency : public AliAnalysisTaskSE { + public: + AliAnalysisTaskEfficiency(const char *name = "AliAnalysisTaskEfficiency"); + virtual ~AliAnalysisTaskEfficiency() {} + + virtual void UserCreateOutputObjects(); + virtual void UserExec(Option_t *option); + virtual void Terminate(Option_t *); + virtual Bool_t SelectJouri(AliESDtrack* track); + virtual void SetTrackType(Int_t type) {fTrackType = type;} + virtual void SetCuts(AliESDtrackCuts* cuts){fCuts = cuts;} + virtual void SetFieldOn(Bool_t b = kTRUE) {fFieldOn = b;} + + private: + + // For crude division into A/C side and positive/negative + enum { kNegA = 0, kPosA, kNegC, kPosC, kTotalAC}; + enum { kTPC = 0, kSPD, kTotalVtx}; + Int_t GetIndexAC(AliESDtrack *track); + Int_t ConvertePDG(Int_t pdg); + + Bool_t fFieldOn; + + TList *fHists; // List of histos + TH1F *fHistRECpt; // pt of reconstructed tracks + TH1F *fHistMCpt; // pt of MC primaries + TH1F *fHistMCNRpt; // pt of not reconstructed primaries + TH1F *fHistFAKEpt; // pt of fake track + TH1F *fHistMULTpt; // pt of multiply reconstructed tracks + + TH1F *fHistRECeta; // eta of reconstructed tracks + TH1F *fHistMCeta; // eta of MC primaries + TH1F *fHistMCNReta; // eta of not reconstructed primaries + TH1F *fHistFAKEeta; // eta of fake track + TH1F *fHistMULTeta; // eta of multiply reconstructed tracks + + TH1F *fHistRECphi; // phi of reconstructed tracks + TH1F *fHistMCphi; // phi of MC primaries + TH1F *fHistMCNRphi; // phi of not reconstructed primaries + TH1F *fHistFAKEphi; // phi of fake track + TH1F *fHistMULTphi; // phi of multiply reconstructed tracks + + TH1F *fHistRECHPTeta; // eta of reconstructed high-pT tracks + TH1F *fHistMCHPTeta; // eta of MC primary high-pT tracks + TH1F *fHistMCNRHPTeta; // eta of not reconstructed primaries + TH1F *fHistFAKEHPTeta; // eta of fake high-pT tracks + + TH1F *fHistRECHPTphi; // phi of reconstructed high-pT tracks + TH1F *fHistMCHPTphi; // phi of MC primary high-pT tracks + TH1F *fHistMCNRHPTphi; // phi of not reconstructed primaries + TH1F *fHistFAKEHPTphi; // phi of fake high-pT tracks + + TH1F *fHistRecMult; // Reconstruction multiplicity + TH1F *fHistNCluster; // Number of clusters of suspicious tracks + + // Checks on ESD data only + TH1F *fh1VtxEff; // Vtx Reconstruction eff for TPC and SPD + TH2F *fh2PhiPadRow[kTotalAC]; // TPC pad row vs phi + TH2F *fh2PhiLayer[kTotalAC]; // ITS layer vs phi + TH2F *fh2MultSpdChips[2]; // Multiplicity vs fired chips + TH2F *fh2VtxTpcSpd; // Vtx correlation TPC <-> SPD + + // Correlation histos + TH2F* fh2VertexCorrelation[kTotalVtx]; // ESD z-vtx vs MC z-vtx + TH2F* fh2VertexCorrelationShift[kTotalVtx]; // (MC z-vtx - ESD z-vtx) vs MC z-vtx + TH1F* fh1VertexShift[kTotalVtx]; // (MC z-vtx - ESD z-vtx) + TH1F* fh1VertexShiftNorm[kTotalVtx]; // (MC z-vtx - ESD z-vtx) / (sigma_ESD-z-vtx) histogrammed + + TH2F* fh2EtaCorrelation; // ESD eta vs MC eta + TH2F* fh2EtaCorrelationShift; // (MC eta - ESD eta) vs MC eta + TH2F* fh2PhiCorrelation; // ESD phi vs MC phi + TH2F* fh2PhiCorrelationShift; // (MC phi - ESD phi) vs MC phi + TH2F* fh2PtCorrelation; // ESD pt vs MC pt + TH2F* fh2PtCorrelationShift; // (MC pt - ESD pt) vs MC pt + + //Tracked Deltaphi histos + TH1F* fHistDeltaphiprimaries; // Tracklet dPhi for primaries + TH1F* fHistDeltaphisecondaries; // Tracklet dPhi for secondaries + TH1F* fHistDeltaphireject; // Tracklet dPhi for rejected tracks + + //ITS-TPC matching + TH1F* fHistTPCITSdeltax; // TPC-ITS matching Delta_x + TH1F* fHistTPCITSdeltay; // TPC-ITS matching Delta_y + TH1F* fHistTPCITSdeltaz; // TPC-ITS matching Delta_z + TH1F* fHistTPCITSdeltar; // TPC-ITS matching Delta_r + //TRD-TPC matching + TH1F* fHistTPCTRDdeltax; // TPC-TRD matching Delta_x + TH1F* fHistTPCTRDdeltay; // TPC-TRD matching Delta_y + TH1F* fHistTPCTRDdeltaz; // TPC-TRD matching Delta_z + TH1F* fHistTPCTRDdeltar; // TPC-TRD matching Delta_r + TH1F* fHistTPCTRDdeltaphi; // TPC-TRD matching Delta_phi + // Pulls + TProfile* fPtPullsPos; // Pt pulls + TProfile* fPtPullsNeg; // Pt pulls + TProfile* fPtShiftPos; // Pt shift + TProfile* fPtShiftNeg; // Pt shift + // Others + Int_t fTrackType; + AliESDtrackCuts* fCuts; // List of cuts + // SPD Vertex + TH2F* fVertexRvsZ; + TH2F* fVertexRvsZC; + // Multi Fluctuations + TProfile* fEtaMultiFluc; // multiplicity fluctuations + TH1F* fEtaMulti; // multiplicity fluctuations + TH1F* fEtaMultiH; // multiplicity fluctuations + TProfile* fPhiMultiFluc; // multiplicity fluctuations + TH1F* fPhiMulti; // multiplicity fluctuations + TH1F* fPhiMultiH; // multiplicity fluctuations + + //for each nuclei type one histo (d, t, 3He, 4He) + TH1F *fHistRECptCharge[8]; // pt of reconstructed tracks + TH1F *fHistMCptCharge[8]; // pt of MC primaries + TH1F *fHistMCNRptCharge[8]; // pt of not reconstructed primaries + TH1F *fHistFAKEptCharge[8]; // pt of fake track + + TH1F *fHistRECetaCharge[8]; // eta of reconstructed tracks + TH1F *fHistMCetaCharge[8]; // eta of MC primaries + TH1F *fHistMCNRetaCharge[8]; // eta of not reconstructed primaries + TH1F *fHistFAKEetaCharge[8]; // eta of fake track + + TH1F *fHistRECphiCharge[8]; // phi of reconstructed tracks + TH1F *fHistMCphiCharge[8]; // phi of MC primaries + TH1F *fHistMCNRphiCharge[8]; // phi of not reconstructed primaries + TH1F *fHistFAKEphiCharge[8]; // phi of fake track + + + AliAnalysisTaskEfficiency(const AliAnalysisTaskEfficiency&); // not implemented + AliAnalysisTaskEfficiency& operator=(const AliAnalysisTaskEfficiency&); // not implemented + + ClassDef(AliAnalysisTaskEfficiency, 1); // example of analysis +}; + +#endif diff --git a/test/vmctest/scripts/efficiency/CreateAnalysisPlugin.C b/test/vmctest/scripts/efficiency/CreateAnalysisPlugin.C new file mode 100644 index 00000000000..2164e73b18f --- /dev/null +++ b/test/vmctest/scripts/efficiency/CreateAnalysisPlugin.C @@ -0,0 +1,145 @@ +AliAnalysisGrid* CreateAnalysisPlugin(TString analysisMode="full") +{ + AliAnalysisAlien *plugin = new AliAnalysisAlien(); + + // Overwrite all generated files, datasets and output results from a previous session + plugin->SetOverwriteMode(); + // Set the run mode (can be "full", "test", "offline", "submit" or "terminate") + + plugin->SetRunMode(TString (analysisMode.Data())); + + // Set versions of used packages + plugin->SetAPIVersion("V1.1x"); + plugin->SetROOTVersion("v5-28-00f"); + plugin->SetAliROOTVersion("v4-21-33-AN"); + // Declare input data to be processed. + + // Method 1: Create automatically XML collections using alien 'find' command. + // Define production directory LFN + // plugin->SetGridDataDir("/alice/sim/LHC10a6"); + //* plugin->SetGridDataDir("/alice/data/2010/LHC10b"); + // Set data search pattern + // plugin->SetDataPattern("*ESDs.root"); // simulated, tags not used + //* plugin->SetDataPattern("*ESDs/pass2/*ESDs.root"); // real data check reco pass and data base directory + //* plugin->SetRunPrefix("000"); // real data + // plugin->SetDataPattern("*tag.root"); // Use ESD tags (same applies for AOD's) + // ...then add run numbers to be considered + // plugin->AddRunNumber(125020); // simulated + + //for ESDs + plugin->SetGridDataDir("/alice/data/2011/LHC11a"); + plugin->SetDataPattern("*ESDs/pass1/*ESDs.root"); + plugin->SetRunPrefix("000"); + + //ESDs sim + // plugin->SetGridDataDir("/alice/sim/LHC10f9b"); + //plugin->SetDataPattern("*ESDs.root"); + + //for AODs + // plugin->SetGridDataDir("/alice/data/2010/LHC10c"); + // plugin->SetRunPrefix("000"); // real data + //plugin->SetDataPattern("*ESDs/pass2_recovery_900GeV/AOD017/*AOD.root"); + + //sim AODs + // plugin->SetGridDataDir("/alice/sim/LHC10d4a"); + //plugin->SetDataPattern("*AOD012/*AOD.root"); + + // TString runs ="120824:120823:120822:120821:120820:120758:120750:120741:120671:120617:120616:120505:120504:120503:120244:120079:120076:120073:120072:120069:120067:119862:119859:119856:119853:119849:119846:119845:119844:119842:119841:119163:119161:119159:119086:119085:119084:119079:119077:119067:119061:119047:119041:119037"; // dont forget last two runs + + //TString runs ="120829:120825"; + // TString runs="118506:118507:118512:118518:118556:118558:118560:118561"; + TString runs ="146801"; + + TObjArray* array = runs.Tokenize ( ":" ); + TObjString *str; + TString strr,strr2_1,strr2_2; + for ( Int_t i = 0;i < array->GetEntriesFast();i++ ) { + str = ( TObjString * ) array->At ( i ); + strr = str->GetString(); + if ( !strr.IsNull() ) { + plugin->AddRunNumber(strr.Atoi()); + } + } + + // Method 2: Declare existing data files (raw collections, xml collections, root file) + // If no path mentioned data is supposed to be in the work directory (see SetGridWorkingDir()) + // XML collections added via this method can be combined with the first method if + // the content is compatible (using or not tags) + // plugin->AddDataFile("tag.xml"); + // plugin->AddDataFile("/alice/data/2008/LHC08c/000057657/raw/Run57657.Merged.RAW.tag.root"); + + // Define alien work directory where all files will be copied. Relative to alien $HOME. + plugin->SetGridWorkingDir("146801"); + // Declare alien output directory. Relative to working directory. + plugin->SetGridOutputDir("output"); // In this case will be $HOME/work/output + // Declare the analysis source files names separated by blancs. To be compiled runtime + // using ACLiC on the worker nodes. + plugin->SetAnalysisSource("AliAnalysisTaskEfficiency.cxx"); + // plugin->SetAdditionalRootLibs("CORRFW PWG2resonances"); + // plugin->SetAdditionalRootLibs("PWG2resonances"); + // plugin->SetAdditionalRootLibs("PWG2resonances"); + // + plugin->SetAdditionalLibs("AliAnalysisTaskEfficiency.h AliAnalysisTaskEfficiency.cxx"); + // plugin->EnablePackage("PWG2resonances"); + // plugin->EnablePackage(""); + // plugin->EnablePackage(""); + // Declare all libraries (other than the default ones for the framework. These will be + // loaded by the generated analysis macro. Add all extra files (task .cxx/.h) here. + + // No need for output file names. Procedure is automatic. + // plugin->SetOutputFiles("Pt.ESD.1.root"); + // plugin->SetDefaultOutputs(); + // No need define the files to be archived. Note that this is handled automatically by the plugin. + // plugin->SetOutputArchive("log_archive.zip:stdout,stderr"); + // Set a name for the generated analysis macro (default MyAnalysis.C) Make this unique ! + plugin->SetAnalysisMacro("AnalysisTest.C"); + // Optionally set maximum number of input files/subjob (default 100, put 0 to ignore). The optimum for an analysis + // is correlated with the run time - count few hours TTL per job, not minutes ! + plugin->SetSplitMaxInputFileNumber(100); + // Optionally set number of failed jobs that will trigger killing waiting sub-jobs. + plugin->SetMaxInitFailed(5); + // Optionally resubmit threshold. + plugin->SetMasterResubmitThreshold(90); + // Optionally set time to live (default 30000 sec) + plugin->SetTTL(20000); + // Optionally set input format (default xml-single) + plugin->SetInputFormat("xml-single"); + // Optionally modify the name of the generated JDL (default analysis.jdl) + plugin->SetJDLName("TaskRsn.jdl"); + // Optionally modify job price (default 1) + plugin->SetPrice(1); + // Optionally modify split mode (default 'se') + plugin->SetSplitMode("se"); + + + //++++++++++++++++ PROOF ++++++++++++++++ + // Proof cluster + + plugin->SetProofCluster("alice-caf"); + //plugin->SetProofCluster("skaf.saske.sk"); + // plugin->SetProofCluster("skaf-test.saske.sk"); + // Dataset to be used + // plugin->SetProofDataSet("/alice/sim/LHC10a12_104316#esdTree"); + // plugin->SetProofDataSet("/alice/sim/LHC10a12_104157#esdTree"); + // plugin->SetProofDataSet("ds.txt"); + plugin->SetProofDataSet("g4g.txt"); + // May need to reset proof. Supported modes: 0-no reset, 1-soft, 2-hard + plugin->SetProofReset(0); + // May limit the number of workers per slave. If used with SetNproofWorkers, SetParallel(nproofworkers) will be called after connection + //plugin->SetNproofWorkersPerSlave(1); + // May request connection to alien upon connection to grid + // plugin->SetProofConnectGrid(kTRUE); + + // plugin->SetNproofWorkers(51); + // May use a specific version of root installed in proof + // plugin->SetRootVersionForProof("current"); + // May set the aliroot mode. Check http://aaf.cern.ch/node/83 + plugin->SetAliRootMode("default"); // Loads AF libs by default + // May request ClearPackages (individual ClearPackage not supported) + plugin->SetClearPackages(kFALSE); + // Plugin test mode works only providing a file containing test file locations + // plugin->SetFileForTestMode("AOD.txt"); + plugin->SetFileForTestMode("test.txt"); + //++++++++++++++ end PROOF ++++++++++++++++ + return plugin; +} diff --git a/test/vmctest/scripts/efficiency/SetupAnalysis.C b/test/vmctest/scripts/efficiency/SetupAnalysis.C new file mode 100644 index 00000000000..94648ff1163 --- /dev/null +++ b/test/vmctest/scripts/efficiency/SetupAnalysis.C @@ -0,0 +1,148 @@ +void SetupAnalysis(TString mode, + TString analysisMode="full", + Bool_t useMC = kFALSE, + Int_t nEvents=1.0*1e9, + Int_t nEventsSkip=0, + TString format="esd") +{ + + // ALICE stuff + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) mgr = new AliAnalysisManager("CAF train"); + + // Create and configure the alien handler plugin + gROOT->LoadMacro("CreateAnalysisPlugin.C"); + AliAnalysisGrid *alienHandler = CreateAnalysisPlugin(analysisMode); + if (!alienHandler) return; + mgr->SetGridHandler(alienHandler); + + // input handler for esd or AOD, real or MC data + InputHandlerSetup(format,useMC); + + // physics selection + if(!format.CompareTo("esd")){ + gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C"); + AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(kFALSE); + if(useMC) physSelTask->GetPhysicsSelection()->SetAnalyzeMC(); + AliPhysicsSelection* physSel = physSelTask->GetPhysicsSelection(); + physSel->AddBackgroundIdentification(new AliBackgroundSelection()); + } + + gROOT->ProcessLine(Form(".include %s/include",gSystem->ExpandPathName("$ALICE_ROOT"))); + + gROOT->LoadMacro("AliAnalysisTaskEfficiency.cxx+g"); + + // load and run AddTask macro + gROOT->LoadMacro("AddTaskEfficiency.C"); + + + AliAnalysisTaskSE* task1 = AddTaskEfficiency(-1); + if(!task1){ + Printf("AddTask could not be run."); + } + + // Run analysis + mgr->InitAnalysis(); + + if ((!mode.CompareTo("proof")) ||(!mode.CompareTo("local"))) { + mgr->StartAnalysis(mode.Data(),nEvents,nEventsSkip); + } + else { + mgr->StartAnalysis(mode.Data()); + + } + +} + + +TString GetFormatFromDataSet(TString dataset) { + + TString dsTreeName; + if (dataset.Contains("#")) { + Info("runSKAF.C",Form("Detecting format from dataset name '%s' ...",dataset.Data())); + dsTreeName=dataset(dataset.Last('#'),dataset.Length()); + } else { + Info("runSKAF.C",Form("Detecting format from dataset '%s' (may take while, depends on network connection) ...", + dataset.Data())); + TFileCollection *ds = gProof->GetDataSet(dataset.Data()); + if (!ds) { + Error(Form("Dataset %s doesn't exist on proof cluster!!!!",dataset.Data())); + return ""; + } + dsTreeName = ds->GetDefaultTreeName(); + } + + if (dsTreeName.Contains("esdTree")) { + Info("runSKAF.C","ESD input format detected ..."); + return "ESD"; + } else if (dsTreeName.Contains("aodTree")) { + Info("runSKAF.C","AOD input format detected ..."); + return "AOD"; + } else { + Error("runSKAF.C",Form("Tree %s is not supported !!!",dsTreeName.Data())); + Error("runSKAF.C",Form("Maybe set your DS to %s#esdTree or %s#aodTree", + dataset.Data(),dataset.Data())); + } + + return ""; +} + +Bool_t InputHandlerSetup(TString format = "esd", Bool_t useKine = kTRUE) +{ + format.ToLower(); + + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + + AliAnalysisDataContainer *cin = mgr->GetCommonInputContainer(); + + if (cin) return; + + if (!format.CompareTo("esd")) + { + AliESDInputHandler *esdInputHandler = + dynamic_cast(AliAnalysisManager::GetAnalysisManager()-> + GetInputEventHandler()); + + if (!esdInputHandler) + { + Info("CustomAnalysisTaskInputSetup", "Creating esdInputHandler ..."); + esdInputHandler = new AliESDInputHandler(); + mgr->SetInputEventHandler(esdInputHandler); + } + + if (useKine) + { + AliMCEventHandler* mcInputHandler = + dynamic_cast(AliAnalysisManager::GetAnalysisManager()-> + GetMCtruthEventHandler()); + + if (!mcInputHandler) + { + Info("CustomAnalysisTaskInputSetup", "Creating mcInputHandler ..."); + AliMCEventHandler* mcInputHandler = new AliMCEventHandler(); + mgr->SetMCtruthEventHandler(mcInputHandler); + } + } + + } + else if (!format.CompareTo("aod")) + { + AliAODInputHandler *aodInputHandler = + dynamic_cast(AliAnalysisManager::GetAnalysisManager()-> + GetInputEventHandler()); + + if (!aodInputHandler) + { + Info("CustomAnalysisTaskInputSetup", "Creating aodInputHandler ..."); + aodInputHandler = new AliAODInputHandler(); + mgr->SetInputEventHandler(aodInputHandler); + } + } + else + { + AliWarning("Wrong input format!!! Only ESD and AOD are supported. Skipping Task ..."); + return kFALSE; + } + + return kTRUE; +} diff --git a/test/vmctest/scripts/efficiency/doubleRatio.C b/test/vmctest/scripts/efficiency/doubleRatio.C new file mode 100644 index 00000000000..a0b15e4070a --- /dev/null +++ b/test/vmctest/scripts/efficiency/doubleRatio.C @@ -0,0 +1,310 @@ +// plot macro for double ratio of reconstruction efficiencies of different particles +// Author: Eva Sicking + +void doubleRatio(Int_t trackType=0, Int_t particle=3) +{ + + //particle names + TString partName[4]={"Deu", "Tri", "He3", "He4"}; + + // track type names + TString trackName[4]={"Global", "TPC", "ITS_SA", "ITS"}; + + //open file and lists + //output of AliAnalysisTaskEfficiency + TFile *f = new TFile("AnalysisResults.root"); + TList *list = (TList *)f->Get(Form("QAHists/QAHists_%s",trackName[trackType].Data())); + + //style + gROOT->SetStyle("Plain"); + gStyle->SetOptStat(0); + const Int_t NRGBs = 5; + const Int_t NCont = 500; + Double_t stops[NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 }; + Double_t red[NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 }; + Double_t green[NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 }; + Double_t blue[NRGBs] = { 0.51, 1.00, 0.12, 0.00, 0.00 }; + TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont); + gStyle->SetNumberContours(NCont); + + // Pt + TH1F* fHistRECpt = (TH1F*) list-> + FindObject(Form("fHistRECptChargeAnti%s",partName[particle].Data())); + TH1F* fHistMCpt = (TH1F*) list-> + FindObject(Form("fHistMCptChargeAnti%s",partName[particle].Data())); + TH1F* fHistFAKEpt = (TH1F*) list-> + FindObject(Form("fHistFAKEptChargeAnti%s",partName[particle].Data())); + TH1F* fHistMCNRpt = (TH1F*) list-> + FindObject(Form("fHistMCNRptChargeAnti%s",partName[particle].Data())); + + + c1 = new TCanvas(Form("pT_antinuclei_%s",trackName[trackType].Data()), + Form("pT_antinuclei_%s",trackName[trackType].Data()), + 100, 100, 1000,600 ); + c1->Divide(2,2); + c1->cd(1); + c1->GetPad(1)->SetLogy(); + + fHistMCpt->SetXTitle("p_{T} [GeV]"); + fHistMCpt->SetMinimum(1); + fHistMCpt->Draw(); + fHistMCpt->SetTitle(""); + fHistRECpt->SetLineColor(2); + fHistRECpt->Draw("same"); + fHistFAKEpt->SetLineColor(4); + fHistFAKEpt->Draw("same"); + fHistMCNRpt->SetLineColor(3); + fHistMCNRpt->Draw("same"); + + c1->cd(4); + fHistPtInEff = (TH1F*) fHistMCNRpt->Clone(); + fHistPtInEff->Divide(fHistMCpt); + fHistPtInEff->SetXTitle("p_{T} [GeV]"); + fHistPtInEff->SetTitle("Inefficiency from non-rec. MC tracks"); + fHistPtInEff->Draw(); + fHistPtInEff->SetMinimum(0); + fHistPtInEff->SetMaximum(1); + + c1->cd(3); + fHistPtEff = (TH1F*) fHistRECpt->Clone(); + fHistPtEff->Add(fHistFAKEpt, -1.); + fHistPtEff->Divide(fHistMCpt); + //fHistPtEff->Add(fHistPtInEff); + fHistPtEff->SetXTitle("p_{T} [GeV]"); + fHistPtEff->SetTitle("Efficiency"); + fHistPtEff->Draw(); + fHistPtEff->SetMinimum(0); + fHistPtEff->SetMaximum(1); + + c1->cd(2); + TLine *line = new TLine(0.07, 0.75, 0.2, 0.75); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.73, "MC"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.55, 0.2, 0.55); + line->SetLineColor(2); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.53, "Reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.35, 0.2, 0.35); + line->SetLineColor(4); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.33, "Fake Tracks"); + tex->SetLineWidth(2); + tex->Draw(); + + + line = new TLine(0.07, 0.15, 0.2, 0.15); + line->SetLineColor(3); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.13, "MC not reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + + + // Pt + TH1F* fHistRECpt2 = (TH1F*) list-> + FindObject(Form("fHistRECptCharge%s",partName[particle].Data())); + TH1F* fHistMCpt2 = (TH1F*) list-> + FindObject(Form("fHistMCptCharge%s",partName[particle].Data())); + TH1F* fHistFAKEpt2 = (TH1F*) list-> + FindObject(Form("fHistFAKEptCharge%s",partName[particle].Data())); + TH1F* fHistMCNRpt2 = (TH1F*) list-> + FindObject(Form("fHistMCNRptCharge%s",partName[particle].Data())); + + c2 = new TCanvas(Form("pT_nuclei_%s",trackName[trackType].Data()), + Form("pT_nuclei_%s",trackName[trackType].Data()), + 100, 100, 1000,600 ); + c2->Divide(2,2); + c2->cd(1); + c2->GetPad(1)->SetLogy(); + + fHistMCpt2->SetXTitle("p_{T} [GeV]"); + fHistMCpt2->SetMinimum(1); + fHistMCpt2->Draw(); + fHistMCpt2->SetTitle(""); + fHistRECpt2->SetLineColor(2); + fHistRECpt2->Draw("same"); + fHistFAKEpt2->SetLineColor(4); + fHistFAKEpt2->Draw("same"); + fHistMCNRpt2->SetLineColor(3); + fHistMCNRpt2->Draw("same"); + + c2->cd(4); + fHistPtInEff2 = (TH1F*) fHistMCNRpt2->Clone(); + fHistPtInEff2->Divide(fHistMCpt2); + fHistPtInEff2->SetXTitle("p_{T} [GeV]"); + fHistPtInEff2->SetTitle("Inefficiency from non-rec. MC tracks"); + fHistPtInEff2->Draw(); + fHistPtInEff2->SetMinimum(0); + fHistPtInEff2->SetMaximum(1); + + c2->cd(3); + fHistPtEff2 = (TH1F*) fHistRECpt2->Clone(); + fHistPtEff2->Add(fHistFAKEpt2, -1.); + fHistPtEff2->Divide(fHistMCpt2); + // fHistPtEff->Add(fHistPtInEff); + fHistPtEff2->SetXTitle("p_{T} [GeV]"); + fHistPtEff2->SetTitle("Efficiency"); + fHistPtEff2->Draw(); + fHistPtEff2->SetMinimum(0); + fHistPtEff2->SetMaximum(1); + + c2->cd(2); + TLine *line = new TLine(0.07, 0.75, 0.2, 0.75); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.73, "MC"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.55, 0.2, 0.55); + line->SetLineColor(2); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.53, "Reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.35, 0.2, 0.35); + line->SetLineColor(4); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.33, "Fake Tracks"); + tex->SetLineWidth(2); + tex->Draw(); + + + line = new TLine(0.07, 0.15, 0.2, 0.15); + line->SetLineColor(3); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.13, "MC not reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + + + + new TCanvas; + + fHistPtEffPosNeg = (TH1F*) fHistPtEff->Clone(); + fHistPtEffPosNeg->Divide(fHistPtEff2); + fHistPtEffPosNeg->SetLineWidth(2); + fHistPtEffPosNeg->Draw(); + fHistPtEffPosNeg->SetMinimum(0); + fHistPtEffPosNeg->SetMaximum(1.2); + fHistPtEff->Draw("same");//anti-particle + fHistPtEff->SetLineColor(kRed); + fHistPtEff->SetLineWidth(2); + fHistPtEff2->Draw("same"); + fHistPtEff2->SetLineColor(kBlack); + fHistPtEff2->SetLineWidth(2); + TF1 *func = new TF1("func", "1",-5,20); + func->SetLineWidth(1); + func->Draw("same"); + + + + + + + TLegend *legp= new TLegend(0.65,0.78,0.9,0.98); + legp->SetFillColor(kWhite); + legp->SetBorderSize(0); + Int_t number =4; + TF1 *fun[2]; + for(Int_t i=0;i<2;i++){ + fun[i]= new TF1(Form("fun%d",i),"gaus",-5.0,5.0); + fun[i]->SetLineColor(2**i); + fun[i]->SetLineStyle(1); + } + legp->AddEntry(fun[0],Form("%s",partName[particle].Data() ),"l"); + legp->AddEntry(fun[1],Form("Anti-%s",partName[particle].Data() ),"l"); + + + + + gROOT->SetStyle("Plain"); + gStyle->SetOptStat(0); + gStyle->SetPalette(1); + + + TCanvas *c = new TCanvas("c","PtAll" , 100, 100, 880, 680); + c->Draw(); + c ->cd(); + TPad *pad_spectrum = new TPad( "pad_spectrum","pad_spectrum", 0.0,0.39,1.0,0.95); + c ->cd(); + TPad *pad_ratio = new TPad( "pad_ratio","pad_ratio", 0.0,0.03,1.0,0.39); + + + c ->cd(); + pad_spectrum->SetTopMargin(0.); + pad_spectrum->SetBottomMargin(0.); + pad_spectrum->SetLeftMargin(0.12); + pad_spectrum->SetRightMargin(0.055); + pad_spectrum->Draw(); + pad_spectrum->cd(); + TH1F * histo = new TH1F("histo", "",150, -0.5,149.5 ); + histo->GetXaxis()->SetRangeUser(0,11); + // histo->GetXaxis()->SetRangeUser(25,35); + histo->GetYaxis()->SetRangeUser(0.01,1.3); + histo->SetXTitle("N_{charged}"); + histo->SetTitle(""); + histo->SetYTitle("efficiency "); + histo->Draw(); + histo->GetYaxis()->SetNdivisions(10); + histo->GetYaxis()->SetLabelSize(0.07); + histo->GetYaxis()->SetTitleSize(0.08); + histo->GetYaxis()->SetTitleOffset(0.6); + + fHistPtEff->Draw("same"); + fHistPtEff2->Draw("same"); + + legp->Draw(); + + + + c ->cd(); + pad_ratio->SetTopMargin(0.); + pad_ratio->SetBottomMargin(0.36); + pad_ratio->SetLeftMargin(0.12); + pad_ratio->SetRightMargin(0.055); + pad_ratio->Draw(); + pad_ratio->cd(); + TH1F * histo2 = new TH1F("histo2", "",150, -0.5,149.5 ); + histo2->GetXaxis()->SetRangeUser(0,90); + // histo->GetXaxis()->SetRangeUser(25,35); + histo2->SetXTitle("p_{T} (GeV/c)"); + histo2->SetTitle(""); + histo2->SetYTitle("Ratio Neg / Pos"); + histo2->Draw(); + histo2->GetXaxis()->SetRangeUser(0.0,11); + histo2->SetLabelSize(0.1,"X"); + histo2->SetLabelSize(0.1,"Y"); + histo2->SetTitleSize(0.09,"X"); + histo2->SetTitleSize(0.09,"Y"); + histo2->SetTitleOffset(0.3,"Y"); + histo2->GetYaxis()->SetRangeUser(0.59,1.16); + histo2->GetYaxis()->SetNdivisions(4); + histo2->GetYaxis()->SetLabelSize(0.12); + histo2->GetXaxis()->SetLabelSize(0.12); + histo2->GetXaxis()->SetTitleSize(0.15); + histo2->GetYaxis()->SetTitleSize(0.11); + + fHistPtEffPosNeg->Draw("same"); + func->Draw("same"); + func->SetLineStyle(2); + + c->Print(Form("%s.png",partName[particle].Data())); + +} diff --git a/test/vmctest/scripts/efficiency/macrosalida.C b/test/vmctest/scripts/efficiency/macrosalida.C new file mode 100644 index 00000000000..2e3bf1a8ce1 --- /dev/null +++ b/test/vmctest/scripts/efficiency/macrosalida.C @@ -0,0 +1,290 @@ +// monitor reconstruction efficiencies for different track types +// use output of AliAnalysisTaskEfficiency +// Authors: Veronica Canoa Roman, Eva Sicking + +void macrosalida(Int_t prod=0,Int_t trackType=0 ) +{ + + //postfix of analysis output files ("X") + TString prodName[6]={"A","B","C","D","E","F"}; + //names of physics lists using in analysed productions (LHC11d6x) + TString prodNameDetail[6]={"Geant4_QGSP_BERT_EMV_p02", + "Geant4_QGSP_BERT_CHIPS_p02", + "Geant4_QGSP_BERT_EMV_b01", + "Geant4_QGSP_BERT_CHIPS_b01", + "Geant4_QGSP_FTFP_BERT_b01", + "Geant3"}; + // name of track types + // ITS = left over tracks from global tracking, efficiency does not make sense + TString trackName[4]={"Global", "TPC", "ITS_SA", "ITS"}; + + //open file + // TFile *f = new TFile(Form("AnalysisResults%s.root",prodName[prod].Data())); + TFile *f = new TFile(Form("AnalysisResults.root",prodName[prod].Data())); + TList *list = (TList *)f->Get(Form("QAHists/QAHists_%s",trackName[trackType].Data())); + + //style + gROOT->SetStyle("Plain"); + gStyle->SetOptStat(0); + const Int_t NRGBs = 5; + const Int_t NCont = 500; + Double_t stops[NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 }; + Double_t red[NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 }; + Double_t green[NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 }; + Double_t blue[NRGBs] = { 0.51, 1.00, 0.12, 0.00, 0.00 }; + TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont); + gStyle->SetNumberContours(NCont); + + + + // Pt efficiency histos + TH1F* fHistRECpt = (TH1F*) list->FindObject("fHistRECpt"); + TH1F* fHistMCpt = (TH1F*) list->FindObject("fHistMCpt"); + TH1F* fHistFAKEpt = (TH1F*) list->FindObject("fHistFAKEpt"); + TH1F* fHistMCNRpt = (TH1F*) list->FindObject("fHistMCNRpt"); + + + c1 = new TCanvas(Form("pT_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + Form("pT_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + 100, 100, 1000,800 ); + c1->Divide(2,2); + c1->cd(1); + c1->GetPad(1)->SetLogy(); + + fHistMCpt->SetXTitle("p_{T} [GeV]"); + fHistMCpt->SetMinimum(1); + fHistMCpt->Draw(); + fHistMCpt->SetTitle(""); + fHistMCpt->SetLineWidth(2); + fHistRECpt->SetLineColor(2); + fHistRECpt->SetLineWidth(2); + fHistRECpt->Draw("same"); + fHistFAKEpt->SetLineColor(4); + fHistFAKEpt->SetLineWidth(2); + fHistFAKEpt->Draw("same"); + fHistMCNRpt->SetLineColor(3); + fHistMCNRpt->SetLineWidth(2); + fHistMCNRpt->Draw("same"); + + c1->cd(4); + fHistPtInEff = (TH1F*) fHistMCNRpt->Clone(); + fHistPtInEff->Divide(fHistMCpt); + fHistPtInEff->SetXTitle("p_{T} [GeV]"); + fHistPtInEff->SetTitle("Inefficiency from non-rec. MC tracks"); + fHistPtInEff->Draw(); + fHistPtInEff->SetMinimum(0); + fHistPtInEff->SetMaximum(1); + + c1->cd(3); + fHistPtEff = (TH1F*) fHistRECpt->Clone(); + fHistPtEff->Add(fHistFAKEpt, -1.); + fHistPtEff->Divide(fHistMCpt); + // fHistPtEff->Add(fHistPtInEff); + fHistPtEff->SetXTitle("p_{T} [GeV]"); + fHistPtEff->SetTitle("Efficiency"); + fHistPtEff->Draw(); + fHistPtEff->SetMinimum(0); + fHistPtEff->SetMaximum(1); + + c1->cd(2); + TLine *line = new TLine(0.07, 0.75, 0.2, 0.75); + line->SetLineWidth(3); + line->Draw(); + tex = new TLatex(0.25, 0.73, "MC"); + tex->SetLineWidth(2); + tex->SetTextSize(0.08); + tex->Draw(); + + line = new TLine(0.07, 0.55, 0.2, 0.55); + line->SetLineColor(2); + line->SetLineWidth(3); + line->Draw(); + tex = new TLatex(0.25, 0.53, "Reconstructed"); + tex->SetLineWidth(2); + tex->SetTextSize(0.08); + tex->Draw(); + + line = new TLine(0.07, 0.35, 0.2, 0.35); + line->SetLineColor(4); + line->SetLineWidth(3); + line->Draw(); + tex = new TLatex(0.25, 0.33, "Fake Tracks"); + tex->SetLineWidth(3); + tex->SetTextSize(0.08); + tex->Draw(); + + + line = new TLine(0.07, 0.15, 0.2, 0.15); + line->SetLineColor(3); + line->SetLineWidth(3); + line->Draw(); + tex = new TLatex(0.25, 0.13, "MC not reconstructed"); + tex->SetLineWidth(2); + tex->SetTextSize(0.08); + tex->Draw(); + + // c1->Print(Form("pT_%s_%s.png",prodNameDetail[prod].Data(),trackName[trackType].Data())); + + // Eta efficiency histos + c2 = new TCanvas(Form("eta_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + Form("eta_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + 100, 100, 1000,600 ); + c2->Divide(2,2); + c2->cd(1); + + TH1F* fHistRECeta = (TH1F*) list->FindObject("fHistRECeta"); + TH1F* fHistMCeta = (TH1F*) list->FindObject("fHistMCeta"); + TH1F* fHistFAKEeta = (TH1F*) list->FindObject("fHistFAKEeta"); + TH1F* fHistMCNReta = (TH1F*) list->FindObject("fHistMCNReta"); + fHistMCeta->SetMinimum(0.); + fHistMCeta->SetXTitle("#eta"); + fHistMCeta->Draw(); + fHistMCeta->SetTitle(""); + fHistRECeta->SetLineColor(2); + fHistRECeta->Draw("same"); + fHistFAKEeta->SetLineColor(4); + fHistFAKEeta->Draw("same"); + fHistMCNReta->SetLineColor(3); + fHistMCNReta->Draw("same"); + c2->cd(4); + fHistEtaInEff = (TH1F*) fHistMCNReta->Clone(); + fHistEtaInEff->Divide(fHistMCeta); + fHistEtaInEff->SetXTitle("#eta"); + fHistEtaInEff->SetTitle("Inefficiency from non-rec. MC tracks"); + fHistEtaInEff->Draw(); + fHistEtaInEff->SetMinimum(0); + fHistEtaInEff->SetMaximum(1); + + c2->cd(3); + fHistEtaEff = (TH1F*) fHistRECeta->Clone(); + fHistEtaEff->Add(fHistFAKEeta, -1.); + fHistEtaEff->Divide(fHistMCeta); + fHistEtaEff->SetXTitle("#eta"); + fHistEtaEff->SetTitle("Efficiency"); + fHistEtaEff->Draw(); + fHistEtaEff->SetMinimum(0); + fHistEtaEff->SetMaximum(1); + + c2->cd(2); + TLine *line = new TLine(0.07, 0.75, 0.2, 0.75); + line->Draw(); + line->SetLineWidth(2); + tex = new TLatex(0.25, 0.73, "MC"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.55, 0.2, 0.55); + line->SetLineColor(2); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.53, "Reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.35, 0.2, 0.35); + line->SetLineColor(4); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.33, "Fake Tracks"); + tex->SetLineWidth(2); + tex->Draw(); + + + line = new TLine(0.07, 0.15, 0.2, 0.15); + line->SetLineColor(3); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.13, "MC not reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + //c2->Print(Form("eta_%s_%s.png",prodNameDetail[prod].Data(),trackName[trackType].Data())); + + + // Phi efficiency histos + c3 = new TCanvas(Form("phi_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + Form("phi_%s_%s",prodNameDetail[prod].Data(), + trackName[trackType].Data()), + 100, 100, 1000,600 ); + c3->Divide(2,2); + c3->cd(1); + + TH1F* fHistRECphi = (TH1F*) list->FindObject("fHistRECphi"); + TH1F* fHistMCphi = (TH1F*) list->FindObject("fHistMCphi"); + TH1F* fHistFAKEphi = (TH1F*) list->FindObject("fHistFAKEphi"); + TH1F* fHistMCNRphi = (TH1F*) list->FindObject("fHistMCNRphi"); + + fHistMCphi->SetMinimum(0.); + fHistMCphi->SetXTitle("#phi"); + fHistMCphi->Draw(); + fHistMCphi->SetTitle(""); + fHistRECphi->SetLineColor(2); + fHistRECphi->Draw("same"); + fHistFAKEphi->SetLineColor(4); + fHistFAKEphi->Draw("same"); + fHistMCNRphi->SetLineColor(3); + fHistMCNRphi->Draw("same"); + + c3->cd(4); + + fHistPhiInEff = (TH1F*) fHistMCNRphi->Clone(); + fHistPhiInEff->Divide(fHistMCphi); + fHistPhiInEff->SetXTitle("#phi"); + fHistPhiInEff->SetTitle("Inefficiency from non-rec. MC tracks"); + fHistPhiInEff->Draw(); + fHistPhiInEff->SetMinimum(0); + fHistPhiInEff->SetMaximum(1); + + c3->cd(3); + fHistPhiEff = (TH1F*) fHistRECphi->Clone(); + fHistPhiEff->Add(fHistFAKEphi, -1.); + fHistPhiEff->Divide(fHistMCphi); + fHistPhiEff->SetXTitle("#phi"); + fHistPhiEff->SetTitle("Efficiency"); + fHistPhiEff->Draw(); + fHistPhiEff->SetMinimum(0); + fHistPhiEff->SetMaximum(1); + + + c3->cd(2); + TLine *line = new TLine(0.07, 0.75, 0.2, 0.75); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.73, "MC"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.55, 0.2, 0.55); + line->SetLineColor(2); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.53, "Reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + line = new TLine(0.07, 0.35, 0.2, 0.35); + line->SetLineColor(4); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.33, "Fake Tracks"); + tex->SetLineWidth(2); + tex->Draw(); + + + line = new TLine(0.07, 0.15, 0.2, 0.15); + line->SetLineColor(3); + line->SetLineWidth(2); + line->Draw(); + tex = new TLatex(0.25, 0.13, "MC not reconstructed"); + tex->SetLineWidth(2); + tex->Draw(); + + // c3->Print(Form("phi_%s_%s.png",prodNameDetail[prod].Data(),trackName[trackType].Data())); + + +} diff --git a/test/vmctest/scripts/efficiency/runALICE.C b/test/vmctest/scripts/efficiency/runALICE.C new file mode 100644 index 00000000000..e918087067a --- /dev/null +++ b/test/vmctest/scripts/efficiency/runALICE.C @@ -0,0 +1,65 @@ +void runALICE(TString mode="local", + TString analysisMode="full", + Bool_t useMC = kTRUE, + Int_t nEvents=1.0*1e6, + Int_t nEventsSkip=0, + TString format="esd") { + + gSystem->Load("libTree.so"); + gSystem->Load("libGeom.so"); + gSystem->Load("libVMC.so"); + gSystem->Load("libPhysics.so"); + gSystem->Load("libSTEERBase.so"); + gSystem->Load("libESD.so"); + gSystem->Load("libAOD.so"); + + gSystem->Load("libANALYSIS.so"); + gSystem->Load("libANALYSISalice.so"); + gSystem->Load("libCORRFW.so"); + + + // Create and configure the alien handler plugin + gROOT->LoadMacro("SetupAnalysis.C"); + SetupAnalysis(mode.Data(), + analysisMode.Data(), + useMC, + nEvents, + nEventsSkip, + format.Data()); + +} + +Int_t setupPar(const char* pararchivename) { + /////////////////// + // Setup PAR File// + /////////////////// + if (pararchivename) { + char processline[1024]; + sprintf(processline,".! tar xzf %s.par",pararchivename); + gROOT->ProcessLine(processline); + const char* ocwd = gSystem->WorkingDirectory(); + gSystem->ChangeDirectory(pararchivename); + + // check for BUILD.sh and execute + if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) { + // printf("*******************************\n"); + // printf("*** Building PAR archive ***\n"); + // printf("*******************************\n"); + + if (gSystem->Exec("PROOF-INF/BUILD.sh")) { + Error("runAnalysis","Cannot Build the PAR Archive! - Abort!"); + return -1; + } + } + // check for SETUP.C and execute + if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) { + // printf("*******************************\n"); + // printf("*** Setup PAR archive ***\n"); + // printf("*******************************\n"); + gROOT->Macro("PROOF-INF/SETUP.C"); + } + + gSystem->ChangeDirectory("../"); + } + return 1; +} -- 2.43.0