9 #include "TLorentzVector.h"
11 #include "AliAnalysisTask.h"
12 #include "AliAnalysisManager.h"
14 #include "AliESDEvent.h"
15 #include "AliESDInputHandler.h"
16 #include "AliAODEvent.h"
17 #include "AliAODInputHandler.h"
18 #include "AliMCEventHandler.h"
19 #include "AliMCEvent.h"
21 #include "AliCFContainer.h"
22 #include "AliVertexerTracks.h"
23 #include "AliESDVertex.h"
25 #include "AliProtonAnalysis.h"
26 #include "AliAnalysisTaskProtons.h"
28 // Analysis task creating a the 2d y-p_t spectrum of p and antip
29 // Author: Panos Cristakoglou
31 ClassImp(AliAnalysisTaskProtons)
33 //________________________________________________________________________
34 AliAnalysisTaskProtons::AliAnalysisTaskProtons()
35 : AliAnalysisTask(), fESD(0), fAOD(0), fMC(0), fAnalysisType("ESD"),
36 fList(0), fProtonAnalysis(0),
37 fElectronFunction(0), fMuonFunction(0),
38 fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
39 fFunctionUsed(kFALSE) {
44 //________________________________________________________________________
45 AliAnalysisTaskProtons::AliAnalysisTaskProtons(const char *name)
46 : AliAnalysisTask(name, ""), fESD(0), fAOD(0), fMC(0), fAnalysisType("ESD"),
47 fList(0), fProtonAnalysis(0),
48 fElectronFunction(0), fMuonFunction(0),
49 fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
50 fFunctionUsed(kFALSE),
51 fTriggerMode(kMB2), fProtonAnalysisMode(kTPC),
52 fVxMax(0), fVyMax(0), fVzMax(0) {
55 // Define input and output slots here
56 // Input slot #0 works with a TChain
57 DefineInput(0, TChain::Class());
58 // Output slot #0 writes into a TList container
59 DefineOutput(0, TList::Class());
62 //________________________________________________________________________
63 void AliAnalysisTaskProtons::ConnectInputData(Option_t *) {
64 // Connect ESD or AOD here
67 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
69 Printf("ERROR: Could not read chain from input slot 0");
71 if(fAnalysisType == "ESD") {
72 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
75 Printf("ERROR: Could not get ESDInputHandler");
77 fESD = esdH->GetEvent();
79 else if(fAnalysisType == "AOD") {
80 AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
83 Printf("ERROR: Could not get AODInputHandler");
85 fAOD = aodH->GetEvent();
87 else if(fAnalysisType == "MC") {
88 AliMCEventHandler* mcH = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
90 Printf("ERROR: Could not retrieve MC event handler");
96 Printf("Wrong analysis type: Only ESD, AOD and MC types are allowed!");
100 //________________________________________________________________________
101 void AliAnalysisTaskProtons::CreateOutputObjects() {
104 //Prior probabilities
105 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
107 //proton analysis object
108 fProtonAnalysis = new AliProtonAnalysis();
110 if(fAnalysisType == "ESD") {
111 //Use of TPConly tracks
112 if(fProtonAnalysisMode == kTPC) {
113 fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
114 fProtonAnalysis->UseTPCOnly();
115 //fProtonAnalysis->SetTPCpid();
116 fProtonAnalysis->SetMinTPCClusters(100);
117 fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
118 fProtonAnalysis->SetMaxCov11(0.5);
119 fProtonAnalysis->SetMaxCov22(0.5);
120 fProtonAnalysis->SetMaxCov33(0.5);
121 fProtonAnalysis->SetMaxCov44(0.5);
122 fProtonAnalysis->SetMaxCov55(0.5);
123 fProtonAnalysis->SetMaxSigmaToVertexTPC(2.0);
124 //fProtonAnalysis->SetMaxDCAXYTPC(1.5);
125 //fProtonAnalysis->SetMaxDCAZTPC(1.5);
127 //Use of HybridTPC tracks
128 else if(fProtonAnalysisMode == kHybrid) {
129 fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
130 fProtonAnalysis->UseHybridTPC();
131 fProtonAnalysis->SetTPCpid();
132 fProtonAnalysis->SetMinTPCClusters(110);
133 fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
134 fProtonAnalysis->SetMaxCov11(0.5);
135 fProtonAnalysis->SetMaxCov22(0.5);
136 fProtonAnalysis->SetMaxCov33(0.5);
137 fProtonAnalysis->SetMaxCov44(0.5);
138 fProtonAnalysis->SetMaxCov55(0.5);
139 fProtonAnalysis->SetMaxSigmaToVertex(2.0);
140 /*fProtonAnalysis->SetMaxDCAXY(1.5);
141 fProtonAnalysis->SetMaxDCAZ(1.5);*/
142 fProtonAnalysis->SetPointOnITSLayer6();
143 fProtonAnalysis->SetPointOnITSLayer5();
144 //fProtonAnalysis->SetPointOnITSLayer4();
145 //fProtonAnalysis->SetPointOnITSLayer3();
146 fProtonAnalysis->SetPointOnITSLayer2();
147 fProtonAnalysis->SetPointOnITSLayer1();
148 fProtonAnalysis->SetMinITSClusters(5);
151 else if(fProtonAnalysisMode == kGlobal) {
152 fProtonAnalysis->InitAnalysisHistograms(20, -1.0, 1.0, 48, 0.3, 1.5);
153 fProtonAnalysis->SetMinTPCClusters(110);
154 fProtonAnalysis->SetMaxChi2PerTPCCluster(2.2);
155 fProtonAnalysis->SetMaxCov11(0.5);
156 fProtonAnalysis->SetMaxCov22(0.5);
157 fProtonAnalysis->SetMaxCov33(0.5);
158 fProtonAnalysis->SetMaxCov44(0.5);
159 fProtonAnalysis->SetMaxCov55(0.5);
160 fProtonAnalysis->SetMaxSigmaToVertex(2.0);
161 //fProtonAnalysis->SetMaxDCAXY(2.0);
162 //fProtonAnalysis->SetMaxDCAZ(2.0);
163 fProtonAnalysis->SetTPCRefit();
164 fProtonAnalysis->SetPointOnITSLayer1();
165 fProtonAnalysis->SetPointOnITSLayer2();
166 //fProtonAnalysis->SetPointOnITSLayer3();
167 //fProtonAnalysis->SetPointOnITSLayer4();
168 fProtonAnalysis->SetPointOnITSLayer5();
169 fProtonAnalysis->SetPointOnITSLayer6();
170 fProtonAnalysis->SetMinITSClusters(5);
171 fProtonAnalysis->SetITSRefit();
172 fProtonAnalysis->SetESDpid();
176 fProtonAnalysis->SetPriorProbabilityFunctions(fElectronFunction,
182 fProtonAnalysis->SetPriorProbabilities(partFrac);
184 else if(fAnalysisType == "MC")
185 fProtonAnalysis->InitAnalysisHistograms(10, -0.5, 0.5, 16, 0.5, 0.9);
188 fList->Add(fProtonAnalysis->GetProtonYPtHistogram());
189 fList->Add(fProtonAnalysis->GetAntiProtonYPtHistogram());
190 fList->Add(fProtonAnalysis->GetEventHistogram());
191 fList->Add(fProtonAnalysis->GetProtonContainer());
192 fList->Add(fProtonAnalysis->GetAntiProtonContainer());
195 //________________________________________________________________________
196 void AliAnalysisTaskProtons::Exec(Option_t *) {
198 // Called for each event
200 if(fAnalysisType == "ESD") {
202 Printf("ERROR: fESD not available");
206 if(IsEventTriggered(fESD,fTriggerMode)) {
207 const AliESDVertex *vertex = GetVertex(fESD,fProtonAnalysisMode,
208 fVxMax,fVyMax,fVzMax);
210 Printf("Proton ESD analysis task: There are %d tracks in this event", fESD->GetNumberOfTracks());
211 fProtonAnalysis->Analyze(fESD,vertex);
212 }//reconstructed vertex
216 else if(fAnalysisType == "AOD") {
218 Printf("ERROR: fAOD not available");
222 Printf("Proton AOD analysis task: There are %d tracks in this event", fAOD->GetNumberOfTracks());
223 fProtonAnalysis->Analyze(fAOD);
226 else if(fAnalysisType == "MC") {
228 Printf("ERROR: Could not retrieve MC event");
232 AliStack* stack = fMC->Stack();
234 Printf("ERROR: Could not retrieve the stack");
237 Printf("Proton MC analysis task: There are %d primaries in this event", stack->GetNprimary());
238 fProtonAnalysis->Analyze(stack,kFALSE);//kTRUE in case of inclusive measurement
245 //________________________________________________________________________
246 void AliAnalysisTaskProtons::Terminate(Option_t *) {
247 // Draw result to the screen
248 // Called once at the end of the query
250 fList = dynamic_cast<TList*> (GetOutputData(0));
252 Printf("ERROR: fList not available");
256 TH2F *fHistYPtProtons = (TH2F *)fList->At(0);
257 TH2F *fHistYPtAntiProtons = (TH2F *)fList->At(1);
259 TCanvas *c1 = new TCanvas("c1","p-\bar{p}",200,0,800,400);
260 c1->SetFillColor(10); c1->SetHighLightColor(10); c1->Divide(2,1);
262 c1->cd(1)->SetLeftMargin(0.15); c1->cd(1)->SetBottomMargin(0.15);
263 if (fHistYPtProtons) fHistYPtProtons->DrawCopy("colz");
264 c1->cd(2)->SetLeftMargin(0.15); c1->cd(2)->SetBottomMargin(0.15);
265 if (fHistYPtAntiProtons) fHistYPtAntiProtons->DrawCopy("colz");
268 //________________________________________________________________________
269 Bool_t AliAnalysisTaskProtons::IsEventTriggered(const AliESDEvent *esd,
270 TriggerMode trigger) {
271 // check if the event was triggered
272 ULong64_t triggerMask = esd->GetTriggerMask();
274 // definitions from p-p.cfg
275 ULong64_t spdFO = (1 << 14);
276 ULong64_t v0left = (1 << 11);
277 ULong64_t v0right = (1 << 12);
281 if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
286 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
291 if (triggerMask & spdFO)
300 //________________________________________________________________________
301 const AliESDVertex* AliAnalysisTaskProtons::GetVertex(AliESDEvent* esd,
306 // Get the vertex from the ESD and returns it if the vertex is valid
307 // Second argument decides which vertex is used (this selects
308 // also the quality criteria that are applied)
309 const AliESDVertex* vertex = 0;
311 vertex = esd->GetPrimaryVertexSPD();
312 else if(mode == kTPC){
313 Double_t kBz = esd->GetMagneticField();
314 AliVertexerTracks vertexer(kBz);
315 vertexer.SetTPCMode();
316 AliESDVertex *vTPC = vertexer.FindPrimaryVertex(esd);
317 esd->SetPrimaryVertexTPC(vTPC);
318 for (Int_t i=0; i<esd->GetNumberOfTracks(); i++) {
319 AliESDtrack *t = esd->GetTrack(i);
320 t->RelateToVertexTPC(vTPC, kBz, kVeryBig);
323 vertex = esd->GetPrimaryVertexTPC();
325 else if(mode == kGlobal)
326 vertex = esd->GetPrimaryVertex();
328 Printf("GetVertex: ERROR: Invalid second argument %d", mode);
330 if(!vertex) return 0;
332 // check Ncontributors
333 if(vertex->GetNContributors() <= 0) return 0;
336 Double_t zRes = vertex->GetZRes();
337 if(zRes == 0) return 0;
340 if(TMath::Abs(vertex->GetXv()) > gVxMax) return 0;
341 if(TMath::Abs(vertex->GetYv()) > gVyMax) return 0;
342 if(TMath::Abs(vertex->GetZv()) > gVzMax) return 0;