-#include "TChain.h"
-#include "TTree.h"
-#include "TCanvas.h"
-#include "TLorentzVector.h"
-#include "TGraphErrors.h"
-
-#include "AliAnalysisTask.h"
-#include "AliAnalysisManager.h"
-
-#include "AliESDEvent.h"
-#include "AliESDInputHandler.h"
-#include "AliAODEvent.h"
-#include "AliAODTrack.h"
-#include "AliAODInputHandler.h"
-#include "AliMCEventHandler.h"
-#include "AliMCEvent.h"
-#include "AliStack.h"
-
-#include "AliBalance.h"
-
-#include "AliAnalysisTaskBF.h"
-
-// Analysis task for the BF code
-// Authors: Panos Cristakoglou@cern.ch
-
-ClassImp(AliAnalysisTaskBF)
-
-//________________________________________________________________________
-AliAnalysisTaskBF::AliAnalysisTaskBF(const char *name)
- : AliAnalysisTask(name, ""), fESD(0), fAOD(0), fMC(0), fBalance(0) {
- // Constructor
-
- // Define input and output slots here
- // Input slot #0 works with a TChain
- DefineInput(0, TChain::Class());
- // Output slot #0 writes into a TH1 container
- DefineOutput(0, AliBalance::Class());
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskBF::ConnectInputData(Option_t *) {
- // Connect ESD or AOD here
- // Called once
- TString gAnalysisLevel = fBalance->GetAnalysisLevel();
-
- TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
- if (!tree) {
- Printf("ERROR: Could not read chain from input slot 0");
- }
- else {
- if(gAnalysisLevel == "ESD") {
- AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
-
- if (!esdH) {
- Printf("ERROR: Could not get ESDInputHandler");
- }
- else
- fESD = esdH->GetEvent();
- }//ESD
- else if(gAnalysisLevel == "AOD") {
- AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
-
- if (!aodH) {
- Printf("ERROR: Could not get AODInputHandler");
- } else
- fAOD = aodH->GetEvent();
- }//AOD
- else if(gAnalysisLevel == "MC") {
- AliMCEventHandler* mcH = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
- if (!mcH) {
- Printf("ERROR: Could not retrieve MC event handler");
- }
- else
- fMC = mcH->MCEvent();
- }//MC
- else
- Printf("Wrong analysis type: Only ESD, AOD and MC types are allowed!");
- }
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskBF::CreateOutputObjects() {
- // Create histograms
- // Called once
-
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskBF::Exec(Option_t *) {
- // Main loop
- // Called for each event
- TString gAnalysisLevel = fBalance->GetAnalysisLevel();
-
- TObjArray *array = new TObjArray();
-
- //ESD analysis
- if(gAnalysisLevel == "ESD") {
- if (!fESD) {
- Printf("ERROR: fESD not available");
- return;
- }
-
- Printf("There are %d tracks in this event", fESD->GetNumberOfTracks());
- for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
- AliESDtrack* track = fESD->GetTrack(iTracks);
- if (!track) {
- Printf("ERROR: Could not receive track %d", iTracks);
- continue;
- }
- array->Add(track);
- } //track loop
- }//ESD analysis
- //AOD analysis
- else if(gAnalysisLevel == "AOD") {
- if (!fAOD) {
- Printf("ERROR: fAOD not available");
- return;
- }
-
- Printf("There are %d tracks in this event", fAOD->GetNumberOfTracks());
- for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
- AliAODTrack* track = fAOD->GetTrack(iTracks);
- if (!track) {
- Printf("ERROR: Could not receive track %d", iTracks);
- continue;
- }
- array->Add(track);
- } //track loop
- }//AOD analysis
- //MC analysis
- else if(gAnalysisLevel == "MC") {
- if (!fMC) {
- Printf("ERROR: fMC not available");
- return;
- }
-
- Printf("There are %d tracks in this event", fMC->GetNumberOfPrimaries());
- for (Int_t iTracks = 0; iTracks < fMC->GetNumberOfPrimaries(); iTracks++) {
- AliMCParticle* track = dynamic_cast<AliMCParticle *>(fMC->GetTrack(iTracks));
- if (!track) {
- Printf("ERROR: Could not receive particle %d", iTracks);
- continue;
- }
- array->Add(track);
- } //track loop
- }//MC analysis
-
- fBalance->CalculateBalance(array);
-
- delete array;
-
- // Post output data.
- PostData(0, fBalance);
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskBF::Terminate(Option_t *) {
- // Draw result to the screen
- // Called once at the end of the query
-
- fBalance = dynamic_cast<AliBalance*> (GetOutputData(0));
- if (!fBalance) {
- Printf("ERROR: fBalance not available");
- return;
- }
-
- TGraphErrors *gr = fBalance->DrawBalance();
- gr->SetMarkerStyle(20);
- gr->Draw("AP");
-
- fBalance->PrintResults();
-}
+#include "TChain.h"\r
+#include "TList.h"\r
+#include "TCanvas.h"\r
+#include "TLorentzVector.h"\r
+#include "TGraphErrors.h"\r
+#include "TH1F.h"\r
+#include "TH2F.h"\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+#include "AliAnalysisManager.h"\r
+\r
+#include "AliESDVertex.h"\r
+#include "AliESDEvent.h"\r
+#include "AliESDInputHandler.h"\r
+#include "AliAODEvent.h"\r
+#include "AliAODTrack.h"\r
+#include "AliAODInputHandler.h"\r
+#include "AliMCEventHandler.h"\r
+#include "AliMCEvent.h"\r
+#include "AliStack.h"\r
+#include "AliESDtrackCuts.h"\r
+\r
+#include "AliAnalysisTaskBF.h"\r
+#include "AliBalance.h"\r
+\r
+\r
+// Analysis task for the BF code\r
+// Authors: Panos Cristakoglou@cern.ch\r
+\r
+ClassImp(AliAnalysisTaskBF)\r
+\r
+//________________________________________________________________________\r
+AliAnalysisTaskBF::AliAnalysisTaskBF(const char *name) \r
+: AliAnalysisTaskSE(name), \r
+ fBalance(0),\r
+ fRunShuffling(kFALSE),\r
+ fShuffledBalance(0),\r
+ fList(0),\r
+ fListBF(0),\r
+ fListBFS(0),\r
+ fHistEventStats(0),\r
+ fHistCentStats(0),\r
+ fHistTriggerStats(0),\r
+ fHistTrackStats(0),\r
+ fHistVx(0),\r
+ fHistVy(0),\r
+ fHistVz(0),\r
+ fHistClus(0),\r
+ fHistDCA(0),\r
+ fHistChi2(0),\r
+ fHistPt(0),\r
+ fHistEta(0),\r
+ fHistPhi(0),\r
+ fHistV0M(0),\r
+ fHistRefTracks(0),\r
+ fESDtrackCuts(0),\r
+ fCentralityEstimator("V0M"),\r
+ fUseCentrality(kFALSE),\r
+ fCentralityPercentileMin(0.), \r
+ fCentralityPercentileMax(5.),\r
+ fUseOfflineTrigger(kFALSE),\r
+ fVxMax(0.3),\r
+ fVyMax(0.3),\r
+ fVzMax(10.),\r
+ nAODtrackCutBit(128),\r
+ fPtMin(0.3),\r
+ fPtMax(1.5),\r
+ fEtaMin(-0.8),\r
+ fEtaMax(-0.8),\r
+ fDCAxyCut(-1),\r
+ fDCAzCut(-1),\r
+ fTPCchi2Cut(-1),\r
+ fNClustersTPCCut(-1){\r
+ // Constructor\r
+\r
+ // Define input and output slots here\r
+ // Input slot #0 works with a TChain\r
+ DefineInput(0, TChain::Class());\r
+ // Output slot #0 writes into a TH1 container\r
+ DefineOutput(1, TList::Class());\r
+ DefineOutput(2, TList::Class());\r
+ DefineOutput(3, TList::Class());\r
+}\r
+\r
+//________________________________________________________________________\r
+AliAnalysisTaskBF::~AliAnalysisTaskBF() {\r
+\r
+ // delete fBalance; \r
+ // delete fShuffledBalance; \r
+ // delete fList;\r
+ // delete fListBF; \r
+ // delete fListBFS;\r
+\r
+ // delete fHistEventStats; \r
+ // delete fHistTrackStats; \r
+ // delete fHistVx; \r
+ // delete fHistVy; \r
+ // delete fHistVz; \r
+\r
+ // delete fHistClus;\r
+ // delete fHistDCA;\r
+ // delete fHistChi2;\r
+ // delete fHistPt;\r
+ // delete fHistEta;\r
+ // delete fHistPhi;\r
+ // delete fHistV0M;\r
+}\r
+\r
+//________________________________________________________________________\r
+void AliAnalysisTaskBF::UserCreateOutputObjects() {\r
+ // Create histograms\r
+ // Called once\r
+ if(!fBalance) {\r
+ fBalance = new AliBalance();\r
+ fBalance->SetAnalysisLevel("ESD");\r
+ //fBalance->SetNumberOfBins(-1,16);\r
+ fBalance->SetInterval(-1,-0.8,0.8,16,0.,1.6);\r
+ }\r
+ if(fRunShuffling) {\r
+ if(!fShuffledBalance) {\r
+ fShuffledBalance = new AliBalance();\r
+ fShuffledBalance->SetAnalysisLevel("ESD");\r
+ //fShuffledBalance->SetNumberOfBins(-1,16);\r
+ fShuffledBalance->SetInterval(-1,-0.8,0.8,16,0.,1.6);\r
+ }\r
+ }\r
+\r
+ //QA list\r
+ fList = new TList();\r
+ fList->SetName("listQA");\r
+ fList->SetOwner();\r
+\r
+ //Balance Function list\r
+ fListBF = new TList();\r
+ fListBF->SetName("listBF");\r
+ fListBF->SetOwner();\r
+\r
+ if(fRunShuffling) {\r
+ fListBFS = new TList();\r
+ fListBFS->SetName("listBFShuffled");\r
+ fListBFS->SetOwner();\r
+ }\r
+\r
+ //Event stats.\r
+ TString gCutName[4] = {"Total","Offline trigger",\r
+ "Vertex","Analyzed"};\r
+ fHistEventStats = new TH1F("fHistEventStats",\r
+ "Event statistics;;N_{events}",\r
+ 4,0.5,4.5);\r
+ for(Int_t i = 1; i <= 4; i++)\r
+ fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());\r
+ fList->Add(fHistEventStats);\r
+\r
+ TString gCentName[9] = {"V0M","FMD","TRK","TKL","CL0","CL1","V0MvsFMD","TKLvsV0M","ZEMvsZDC"};\r
+ fHistCentStats = new TH2F("fHistCentStats",\r
+ "Centrality statistics;;Cent percentile",\r
+ 9,-0.5,8.5,220,-5,105);\r
+ for(Int_t i = 1; i <= 9; i++)\r
+ fHistCentStats->GetXaxis()->SetBinLabel(i,gCentName[i-1].Data());\r
+ fList->Add(fHistCentStats);\r
+\r
+ fHistTriggerStats = new TH1F("fHistTriggerStats","Trigger statistics;TriggerBit;N_{events}",130,0,130);\r
+ fList->Add(fHistTriggerStats);\r
+\r
+ fHistTrackStats = new TH1F("fHistTrackStats","Event statistics;TrackFilterBit;N_{events}",130,0,130);\r
+ fList->Add(fHistTrackStats);\r
+\r
+ // Vertex distributions\r
+ fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);\r
+ fList->Add(fHistVx);\r
+ fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);\r
+ fList->Add(fHistVy);\r
+ fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);\r
+ fList->Add(fHistVz);\r
+\r
+ // QA histograms\r
+ fHistClus = new TH2F("fHistClus","# Cluster (TPC vs. ITS)",10,0,10,200,0,200);\r
+ fList->Add(fHistClus);\r
+ fHistChi2 = new TH1F("fHistChi2","Chi2/NDF distribution",200,0,10);\r
+ fList->Add(fHistChi2);\r
+ fHistDCA = new TH2F("fHistDCA","DCA (xy vs. z)",400,-5,5,400,-5,5); \r
+ fList->Add(fHistDCA);\r
+ fHistPt = new TH1F("fHistPt","p_{T} distribution",200,0,10);\r
+ fList->Add(fHistPt);\r
+ fHistEta = new TH1F("fHistEta","#eta distribution",200,-2,2);\r
+ fList->Add(fHistEta);\r
+ fHistPhi = new TH1F("fHistPhi","#phi distribution",200,-20,380);\r
+ fList->Add(fHistPhi);\r
+ fHistV0M = new TH2F("fHistV0M","V0 Multiplicity C vs. A",500, 0, 20000, 500, 0, 20000);\r
+ fList->Add(fHistV0M);\r
+ TString gRefTrackName[6] = {"tracks","tracksPos","tracksNeg","tracksTPConly","clusITS0","clusITS1"};\r
+ fHistRefTracks = new TH2F("fHistRefTracks","Nr of Ref tracks/event vs. ref track estimator;;Nr of tracks",6, 0, 6, 400, 0, 20000);\r
+ for(Int_t i = 1; i <= 6; i++)\r
+ fHistRefTracks->GetXaxis()->SetBinLabel(i,gRefTrackName[i-1].Data());\r
+ fList->Add(fHistRefTracks);\r
+\r
+\r
+ // Balance function histograms\r
+\r
+ // Initialize histograms if not done yet\r
+ if(!fBalance->GetHistNp(0)){\r
+ AliWarning("Histograms not yet initialized! --> Will be done now");\r
+ AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
+ fBalance->InitHistograms();\r
+ }\r
+\r
+ if(fRunShuffling) {\r
+ if(!fShuffledBalance->GetHistNp(0)) {\r
+ AliWarning("Histograms (shuffling) not yet initialized! --> Will be done now");\r
+ AliWarning("--> Add 'gBalance->InitHistograms()' in your configBalanceFunction");\r
+ fShuffledBalance->InitHistograms();\r
+ }\r
+ }\r
+\r
+ for(Int_t a = 0; a < ANALYSIS_TYPES; a++){\r
+ fListBF->Add(fBalance->GetHistNp(a));\r
+ fListBF->Add(fBalance->GetHistNn(a));\r
+ fListBF->Add(fBalance->GetHistNpn(a));\r
+ fListBF->Add(fBalance->GetHistNnn(a));\r
+ fListBF->Add(fBalance->GetHistNpp(a));\r
+ fListBF->Add(fBalance->GetHistNnp(a));\r
+\r
+ if(fRunShuffling) {\r
+ fListBFS->Add(fShuffledBalance->GetHistNp(a));\r
+ fListBFS->Add(fShuffledBalance->GetHistNn(a));\r
+ fListBFS->Add(fShuffledBalance->GetHistNpn(a));\r
+ fListBFS->Add(fShuffledBalance->GetHistNnn(a));\r
+ fListBFS->Add(fShuffledBalance->GetHistNpp(a));\r
+ fListBFS->Add(fShuffledBalance->GetHistNnp(a));\r
+ } \r
+ }\r
+\r
+ if(fESDtrackCuts) fList->Add(fESDtrackCuts);\r
+\r
+ // Post output data.\r
+ PostData(1, fList);\r
+ PostData(2, fListBF);\r
+ if(fRunShuffling) PostData(3, fListBFS);\r
+}\r
+\r
+//________________________________________________________________________\r
+void AliAnalysisTaskBF::UserExec(Option_t *) {\r
+ // Main loop\r
+ // Called for each event\r
+ TString gAnalysisLevel = fBalance->GetAnalysisLevel();\r
+\r
+ TObjArray *array = new TObjArray();\r
+ AliESDtrack *track_TPC = NULL;\r
+\r
+\r
+ // vector holding the charges of all tracks\r
+ vector<Int_t> chargeVectorShuffle; // this will be shuffled\r
+ vector<Int_t> chargeVector; // original charge \r
+ \r
+ //ESD analysis\r
+ if(gAnalysisLevel == "ESD") {\r
+ AliESDEvent* gESD = dynamic_cast<AliESDEvent*>(InputEvent()); // from TaskSE\r
+ if (!gESD) {\r
+ Printf("ERROR: gESD not available");\r
+ return;\r
+ }\r
+\r
+ // store offline trigger bits\r
+ fHistTriggerStats->Fill(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected());\r
+\r
+ // event selection done in AliAnalysisTaskSE::Exec() --> this is not used\r
+ fHistEventStats->Fill(1); //all events\r
+ Bool_t isSelected = kTRUE;\r
+ if(fUseOfflineTrigger)\r
+ isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();\r
+ if(isSelected) {\r
+ fHistEventStats->Fill(2); //triggered events\r
+\r
+ if(fUseCentrality) {\r
+ //Centrality stuff\r
+ AliCentrality *centrality = gESD->GetCentrality();\r
+ //Int_t nCentrality = 0;\r
+ //nCentrality = centrality->GetCentralityClass5(fCentralityEstimator.Data());\r
+ //cout<<nCentrality<<" "<<centrality->IsEventInCentralityClass(fCentralityPercentileMin,fCentralityPercentileMax,fCentralityEstimator.Data())<<endl;\r
+ \r
+ // take only events inside centrality class\r
+ if(!centrality->IsEventInCentralityClass(fCentralityPercentileMin,\r
+ fCentralityPercentileMax,\r
+ fCentralityEstimator.Data()))\r
+ return;\r
+ \r
+ // centrality QA (V0M)\r
+ fHistV0M->Fill(gESD->GetVZEROData()->GetMTotV0A(), gESD->GetVZEROData()->GetMTotV0C());\r
+ }\r
+ \r
+ const AliESDVertex *vertex = gESD->GetPrimaryVertex();\r
+ if(vertex) {\r
+ if(vertex->GetNContributors() > 0) {\r
+ if(vertex->GetZRes() != 0) {\r
+ fHistEventStats->Fill(3); //events with a proper vertex\r
+ if(TMath::Abs(vertex->GetXv()) < fVxMax) {\r
+ if(TMath::Abs(vertex->GetYv()) < fVyMax) {\r
+ if(TMath::Abs(vertex->GetZv()) < fVzMax) {\r
+ fHistEventStats->Fill(4); //analayzed events\r
+ fHistVx->Fill(vertex->GetXv());\r
+ fHistVy->Fill(vertex->GetYv());\r
+ fHistVz->Fill(vertex->GetZv());\r
+ \r
+ //Printf("There are %d tracks in this event", gESD->GetNumberOfTracks());\r
+ for (Int_t iTracks = 0; iTracks < gESD->GetNumberOfTracks(); iTracks++) {\r
+ AliESDtrack* track = dynamic_cast<AliESDtrack *>(gESD->GetTrack(iTracks));\r
+ if (!track) {\r
+ Printf("ERROR: Could not receive track %d", iTracks);\r
+ continue;\r
+ } \r
+ \r
+ // take only TPC only tracks\r
+ track_TPC = new AliESDtrack();\r
+ if(!track->FillTPCOnlyTrack(*track_TPC)) continue;\r
+ \r
+ //ESD track cuts\r
+ if(fESDtrackCuts) \r
+ if(!fESDtrackCuts->AcceptTrack(track_TPC)) continue;\r
+ \r
+ // fill QA histograms\r
+ Float_t b[2];\r
+ Float_t bCov[3];\r
+ track_TPC->GetImpactParameters(b,bCov);\r
+ if (bCov[0]<=0 || bCov[2]<=0) {\r
+ AliDebug(1, "Estimated b resolution lower or equal zero!");\r
+ bCov[0]=0; bCov[2]=0;\r
+ }\r
+ \r
+ Int_t nClustersTPC = -1;\r
+ nClustersTPC = track_TPC->GetTPCNclsIter1(); // TPC standalone\r
+ //nClustersTPC = track->GetTPCclusters(0); // global track\r
+ Float_t chi2PerClusterTPC = -1;\r
+ if (nClustersTPC!=0) {\r
+ chi2PerClusterTPC = track_TPC->GetTPCchi2Iter1()/Float_t(nClustersTPC); // TPC standalone\r
+ //chi2PerClusterTPC = track->GetTPCchi2()/Float_t(nClustersTPC); // global track\r
+ }\r
+ \r
+ fHistClus->Fill(track_TPC->GetITSclusters(0),nClustersTPC);\r
+ fHistDCA->Fill(b[1],b[0]);\r
+ fHistChi2->Fill(chi2PerClusterTPC);\r
+ fHistPt->Fill(track_TPC->Pt());\r
+ fHistEta->Fill(track_TPC->Eta());\r
+ fHistPhi->Fill(track_TPC->Phi()*TMath::RadToDeg());\r
+ \r
+ // fill BF array\r
+ array->Add(track_TPC);\r
+ \r
+ // fill charge vector\r
+ chargeVector.push_back(track_TPC->Charge());\r
+ if(fRunShuffling){\r
+ chargeVectorShuffle.push_back(track_TPC->Charge());\r
+ }\r
+ \r
+ delete track_TPC;\r
+ \r
+ } //track loop\r
+ }//Vz cut\r
+ }//Vy cut\r
+ }//Vx cut\r
+ }//proper vertex resolution\r
+ }//proper number of contributors\r
+ }//vertex object valid\r
+ }//triggered event \r
+ }//ESD analysis\r
+ \r
+ //AOD analysis (vertex and track cuts also here!!!!)\r
+ else if(gAnalysisLevel == "AOD") {\r
+ AliAODEvent* gAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // from TaskSE\r
+ if(!gAOD) {\r
+ Printf("ERROR: gAOD not available");\r
+ return;\r
+ }\r
+\r
+ AliAODHeader *aodHeader = gAOD->GetHeader();\r
+\r
+ // store offline trigger bits\r
+ fHistTriggerStats->Fill(aodHeader->GetOfflineTrigger());\r
+\r
+ // event selection done in AliAnalysisTaskSE::Exec() --> this is not used\r
+ fHistEventStats->Fill(1); //all events\r
+ Bool_t isSelected = kTRUE;\r
+ if(fUseOfflineTrigger)\r
+ isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();\r
+ if(isSelected) {\r
+ fHistEventStats->Fill(2); //triggered events\r
+ \r
+ //Centrality stuff (centrality in AOD header)\r
+ if(fUseCentrality) {\r
+ Float_t fCentrality = aodHeader->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data());\r
+ // cout<<fCentralityEstimator.Data()<<" = "<<fCentrality<<" , others are V0M = "\r
+ // << aodHeader->GetCentralityP()->GetCentralityPercentile("V0M")\r
+ // <<" FMD = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("FMD")\r
+ // <<" TRK = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("TRK")\r
+ // <<" TKL = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("TKL")\r
+ // <<" CL0 ="<<aodHeader->GetCentralityP()->GetCentralityPercentile("CL0")\r
+ // <<" CL1 ="<<aodHeader->GetCentralityP()->GetCentralityPercentile("CL1")\r
+ // <<" V0MvsFMD = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("V0MvsFMD")\r
+ // <<" TKLvsV0M = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("TKLvsV0M")\r
+ // <<" ZEMvsZDC = "<<aodHeader->GetCentralityP()->GetCentralityPercentile("ZEMvsZDC")\r
+ // <<endl;\r
+ \r
+ // QA for centrality estimators\r
+ fHistCentStats->Fill(0.,aodHeader->GetCentralityP()->GetCentralityPercentile("V0M"));\r
+ fHistCentStats->Fill(1.,aodHeader->GetCentralityP()->GetCentralityPercentile("FMD"));\r
+ fHistCentStats->Fill(2.,aodHeader->GetCentralityP()->GetCentralityPercentile("TRK"));\r
+ fHistCentStats->Fill(3.,aodHeader->GetCentralityP()->GetCentralityPercentile("TKL"));\r
+ fHistCentStats->Fill(4.,aodHeader->GetCentralityP()->GetCentralityPercentile("CL0"));\r
+ fHistCentStats->Fill(5.,aodHeader->GetCentralityP()->GetCentralityPercentile("CL1"));\r
+ fHistCentStats->Fill(6.,aodHeader->GetCentralityP()->GetCentralityPercentile("V0MvsFMD"));\r
+ fHistCentStats->Fill(7.,aodHeader->GetCentralityP()->GetCentralityPercentile("TKLvsV0M"));\r
+ fHistCentStats->Fill(8.,aodHeader->GetCentralityP()->GetCentralityPercentile("ZEMvsZDC"));\r
+ \r
+ // take only events inside centrality class\r
+ if((fCentrality < fCentralityPercentileMin) || (fCentrality > fCentralityPercentileMax)) \r
+ return;\r
+ \r
+ // centrality QA (V0M)\r
+ fHistV0M->Fill(gAOD->GetVZEROData()->GetMTotV0A(), gAOD->GetVZEROData()->GetMTotV0C());\r
+ \r
+ // centrality QA (reference tracks)\r
+ fHistRefTracks->Fill(0.,aodHeader->GetRefMultiplicity());\r
+ fHistRefTracks->Fill(1.,aodHeader->GetRefMultiplicityPos());\r
+ fHistRefTracks->Fill(2.,aodHeader->GetRefMultiplicityNeg());\r
+ fHistRefTracks->Fill(3.,aodHeader->GetTPConlyRefMultiplicity());\r
+ fHistRefTracks->Fill(4.,aodHeader->GetNumberOfITSClusters(0));\r
+ fHistRefTracks->Fill(5.,aodHeader->GetNumberOfITSClusters(1));\r
+ fHistRefTracks->Fill(6.,aodHeader->GetNumberOfITSClusters(2));\r
+ fHistRefTracks->Fill(7.,aodHeader->GetNumberOfITSClusters(3));\r
+ fHistRefTracks->Fill(8.,aodHeader->GetNumberOfITSClusters(4));\r
+ }\r
+\r
+ const AliAODVertex *vertex = gAOD->GetPrimaryVertex();\r
+ \r
+ if(vertex) {\r
+ Double32_t fCov[6];\r
+ vertex->GetCovarianceMatrix(fCov);\r
+ \r
+ if(vertex->GetNContributors() > 0) {\r
+ if(fCov[5] != 0) {\r
+ fHistEventStats->Fill(3); //events with a proper vertex\r
+ if(TMath::Abs(vertex->GetX()) < fVxMax) {\r
+ if(TMath::Abs(vertex->GetY()) < fVyMax) {\r
+ if(TMath::Abs(vertex->GetZ()) < fVzMax) {\r
+ fHistEventStats->Fill(4); //analyzed events\r
+ fHistVx->Fill(vertex->GetX());\r
+ fHistVy->Fill(vertex->GetY());\r
+ fHistVz->Fill(vertex->GetZ());\r
+ \r
+ //Printf("There are %d tracks in this event", gAOD->GetNumberOfTracks());\r
+ for (Int_t iTracks = 0; iTracks < gAOD->GetNumberOfTracks(); iTracks++) {\r
+ AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(gAOD->GetTrack(iTracks));\r
+ if (!aodTrack) {\r
+ Printf("ERROR: Could not receive track %d", iTracks);\r
+ continue;\r
+ }\r
+ \r
+ // AOD track cuts\r
+ \r
+ // For ESD Filter Information: ANALYSIS/macros/AddTaskESDfilter.C\r
+ // take only TPC only tracks \r
+ fHistTrackStats->Fill(aodTrack->GetFilterMap());\r
+ if(!aodTrack->TestFilterBit(nAODtrackCutBit)) continue;\r
+ \r
+ Float_t pt = aodTrack->Pt();\r
+ Float_t eta = aodTrack->Eta();\r
+ \r
+ Float_t DCAxy = aodTrack->DCA(); // this is the DCA from global track (not exactly what is cut on)\r
+ Float_t DCAz = aodTrack->ZAtDCA(); // this is the DCA from global track (not exactly what is cut on)\r
+ \r
+ \r
+ // Kinematics cuts from ESD track cuts\r
+ if( pt < fPtMin || pt > fPtMax) continue;\r
+ if( eta < fEtaMin || eta > fEtaMax) continue;\r
+ \r
+ // Extra DCA cuts (for systematic studies [!= -1])\r
+ if( fDCAxyCut != -1 && fDCAxyCut != -1){\r
+ if(TMath::Sqrt((DCAxy*DCAxy)/(fDCAxyCut*fDCAxyCut)+(DCAz*DCAz)/(fDCAzCut*fDCAzCut)) > 1 ){\r
+ continue; // 2D cut\r
+ }\r
+ }\r
+ \r
+ // Extra TPC cuts (for systematic studies [!= -1])\r
+ if( fTPCchi2Cut != -1 && aodTrack->Chi2perNDF() > fTPCchi2Cut){\r
+ continue;\r
+ }\r
+ if( fNClustersTPCCut != -1 && aodTrack->GetTPCNcls() < fNClustersTPCCut){\r
+ continue;\r
+ }\r
+ \r
+ \r
+ // fill QA histograms\r
+ fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());\r
+ fHistDCA->Fill(DCAz,DCAxy);\r
+ fHistChi2->Fill(aodTrack->Chi2perNDF());\r
+ fHistPt->Fill(pt);\r
+ fHistEta->Fill(eta);\r
+ fHistPhi->Fill(aodTrack->Phi()*TMath::RadToDeg());\r
+ \r
+ // fill BF array\r
+ array->Add(aodTrack);\r
+ \r
+ // fill charge vector\r
+ chargeVector.push_back(aodTrack->Charge());\r
+ if(fRunShuffling) {\r
+ chargeVectorShuffle.push_back(aodTrack->Charge());\r
+ }\r
+ } //track loop\r
+ }//Vz cut\r
+ }//Vy cut\r
+ }//Vx cut\r
+ }//proper vertex resolution\r
+ }//proper number of contributors\r
+ }//vertex object valid\r
+ }//triggered event \r
+ }//AOD analysis\r
+\r
+ //MC analysis\r
+ else if(gAnalysisLevel == "MC") {\r
+ \r
+ AliMCEvent* mcEvent = MCEvent(); \r
+ if (!mcEvent) {\r
+ Printf("ERROR: mcEvent not available");\r
+ return;\r
+ }\r
+ \r
+ Printf("There are %d tracks in this event", mcEvent->GetNumberOfPrimaries());\r
+ for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfPrimaries(); iTracks++) {\r
+ AliMCParticle* track = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(iTracks));\r
+ if (!track) {\r
+ Printf("ERROR: Could not receive particle %d", iTracks);\r
+ continue;\r
+ }\r
+\r
+ if( track->Pt() < fPtMin || track->Pt() > fPtMax) continue;\r
+ if( track->Eta() < fEtaMin || track->Eta() > fEtaMax) continue;\r
+\r
+ array->Add(track);\r
+\r
+ // fill charge vector\r
+ chargeVector.push_back(track->Charge());\r
+ if(fRunShuffling) {\r
+ chargeVectorShuffle.push_back(track->Charge());\r
+ }\r
+\r
+ } //track loop\r
+ }//MC analysis\r
+ \r
+ // calculate balance function\r
+ fBalance->CalculateBalance(array,chargeVector);\r
+\r
+ if(fRunShuffling) {\r
+ // shuffle charges\r
+ random_shuffle( chargeVectorShuffle.begin(), chargeVectorShuffle.end() );\r
+ fShuffledBalance->CalculateBalance(array,chargeVectorShuffle);\r
+ }\r
+ \r
+ delete array;\r
+ \r
+\r
+} \r
+\r
+//________________________________________________________________________\r
+void AliAnalysisTaskBF::FinishTaskOutput(){\r
+ //Printf("END BF");\r
+\r
+ if (!fBalance) {\r
+ Printf("ERROR: fBalance not available");\r
+ return;\r
+ } \r
+ if(fRunShuffling) {\r
+ if (!fShuffledBalance) {\r
+ Printf("ERROR: fShuffledBalance not available");\r
+ return;\r
+ }\r
+ }\r
+\r
+}\r
+\r
+//________________________________________________________________________\r
+void AliAnalysisTaskBF::Terminate(Option_t *) {\r
+ // Draw result to the screen\r
+ // Called once at the end of the query\r
+\r
+ // not implemented ...\r
+\r
+}\r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff