git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
fix macro control plots npart eccentricity
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTasks / AliAnalysisTaskQAflow.cxx
CommitLineData
2948ac5a 1#include "TMath.h"
2#include "TH1D.h"
3#include "TH2D.h"
4#include "TSeqCollection.h"
5#include "TObjArray.h"
6#include "TObjArray.h"
7#include "TChain.h"
8#include "TMCProcess.h"
9#include "TLorentzVector.h"
10#include "TDirectory.h"
11#include "TROOT.h"
076df7bf 12#include "TNtuple.h"
2948ac5a 13
14#include "AliLog.h"
15#include "AliVParticle.h"
16#include "AliMCParticle.h"
17#include "AliStack.h"
18#include "AliVEvent.h"
19#include "AliESDEvent.h"
20#include "AliMCEvent.h"
21#include "AliESDtrack.h"
22#include "AliFlowTrackCuts.h"
23#include "AliFlowEventCuts.h"
24#include "AliMultiplicity.h"
25#include "AliESDtrackCuts.h"
26#include "AliVertex.h"
076df7bf 27#include "AliFlowEventSimple.h"
28#include "AliFlowVector.h"
2948ac5a 29
30#include "AliAnalysisTask.h"
31#include "AliAnalysisManager.h"
32
33#include "AliAnalysisTaskQAflow.h"
34
35ClassImp(AliAnalysisTaskQAflow)
36
37//________________________________________________________________________
38AliAnalysisTaskQAflow::AliAnalysisTaskQAflow()
39 : AliAnalysisTaskSE(),
40 fOutput(NULL),
076df7bf 41 fFillNtuple(kFALSE),
42 fNtuple(NULL),
2948ac5a 43 fEventCuts(NULL),
44 fTrackCuts(NULL)
45{
46 // Default constructor
47}
48
49//________________________________________________________________________
50AliAnalysisTaskQAflow::AliAnalysisTaskQAflow(const char* name)
51 : AliAnalysisTaskSE(name),
52 fOutput(NULL),
076df7bf 53 fFillNtuple(kFALSE),
54 fNtuple(NULL),
2948ac5a 55 fEventCuts(NULL),
56 fTrackCuts(NULL)
57{
58 // Constructor
076df7bf 59 DefineInput(1, AliFlowEventSimple::Class());
2948ac5a 60 DefineOutput(1, TObjArray::Class());
61}
62
63//________________________________________________________________________
64void AliAnalysisTaskQAflow::UserCreateOutputObjects()
65{
66 // Called once at the beginning
67 fOutput=new TObjArray();
076df7bf 68 fNtuple = new TNtuple("flowQAtree","flowQAtree","meanpt:qx:qy:mult:refmult:trig:tpcvtxx:tpcvtxy:tpcvtxz:ntrackletsA:ntrackletsC");
2948ac5a 69 //TDirectory* before = gDirectory->mkdir("before cuts","before cuts");
70 //TDirectory* after = gDirectory->mkdir("after cuts","after cuts");
71 TObjArray* before = new TObjArray();
72 TObjArray* after = new TObjArray();
73 fOutput->Add(before);
74 fOutput->Add(after);
75
76 //define histograms
77
78 TH1* hist;
79 hist = new TH1I("tracklet_multiplicity","tracklet multiplicity",1000,0,20000);
80 before->Add(hist); after->Add(hist->Clone()); //0
81 hist = new TH1I("track_multiplicity", "standard TPC track multiplicity",1000,0,10000);
82 before->Add(hist); after->Add(hist->Clone()); //1
83 hist = new TH1I("refmult","refmult",1000, 0,10000);
84 before->Add(hist); after->Add(hist->Clone()); //2
85 hist = new TH1I("SPD_clusters","SPD clusters",1000,0,100000);
86 before->Add(hist); after->Add(hist->Clone()); //3
87 hist = new TH1D("primary_vertexZ","primary vertex z",100,-20,20);
88 before->Add(hist); after->Add(hist->Clone()); //4
89 hist = new TH1I("ITS_clusters_on_track", "ITS clusters on track", 8, 0, 8);
90 before->Add(hist); after->Add(hist->Clone()); //5
91 hist = new TH1I("TPC_clusters_on_track", "TPC clusters on track", 159, 1, 160);
92 before->Add(hist); after->Add(hist->Clone()); //6
93 hist = new TH1D("TPC_chi2_per_cluster","TPC #chi^{2}/cluster",100,0.0,5.0);
94 before->Add(hist); after->Add(hist->Clone()); //7
95 hist = new TH1D("DCA_xy","DCA xy", 1000, -5.0, 5.0 );
96 before->Add(hist); after->Add(hist->Clone()); //8
97 hist = new TH1D("DCA_z","DCA z", 1000, -5.0, 5.0 );
98 before->Add(hist); after->Add(hist->Clone()); //9
99 hist = new TH1D("phi_tracklets","#phi tracklets", 1000, 0.0, TMath::TwoPi() );
100 before->Add(hist); after->Add(hist->Clone()); //10
101 hist = new TH1D("phi_tracks","#phi tracks", 1000, 0.0, TMath::TwoPi() );
102 before->Add(hist); after->Add(hist->Clone()); //11
103 hist = new TH1D("eta_tracklets","#eta tracklets", 1000, -2.0, 2.0 );
104 before->Add(hist); after->Add(hist->Clone()); //12
105 hist = new TH1D("eta_tracks","#eta tracks", 1000, -2.0, 2.0 );
106 before->Add(hist); after->Add(hist->Clone()); //13
107 hist = new TH1D("TPC_vertex_z", "TPC vertex z", 100,-20.0,20.0);
108 before->Add(hist); after->Add(hist->Clone()); //14
74d0740e 109 hist = new TH1D("ptyield", "p_{t} spectrum", 10000,0.0,10.0);
110 before->Add(hist); after->Add(hist->Clone()); //15
2948ac5a 111
112 //post data here as it doesn't change anyway (the pointer to list anyway)
113
114 //restore dir add status
115 PostData(1, fOutput);
116}
117
118//________________________________________________________________________
119void AliAnalysisTaskQAflow::UserExec(Option_t *)
120{
121
122 //get teh input data
123 AliESDEvent* event = dynamic_cast<AliESDEvent*>(InputEvent());
124 if (!event)
125 {
126 AliFatal("no ESD event");
127 return;
128 }
129
130 //TObjArray* before = ((TDirectory*)fOutput->At(0))->GetList();
131 //TObjArray* after = ((TDirectory*)fOutput->At(1))->GetList();
132 TObjArray* before = (TObjArray*)fOutput->At(0);
133 TObjArray* after = (TObjArray*)fOutput->At(1);
134 TH1* htrackletmultB = dynamic_cast<TH1*>(before->At(0));
135 TH1* htrackletmultA = dynamic_cast<TH1*>(after->At(0));
136 TH1* htrackmultB = dynamic_cast<TH1*>(before->At(1));
137 TH1* htrackmultA = dynamic_cast<TH1*>(after->At(1));
138 TH1* hrefmultB = dynamic_cast<TH1*>(before->At(2));
139 TH1* hrefmultA = dynamic_cast<TH1*>(after->At(2));
140 TH1* hspdclustersB = dynamic_cast<TH1*>(before->At(3));
141 TH1* hspdclustersA = dynamic_cast<TH1*>(after->At(3));
142 TH1* hprimvtxzB = dynamic_cast<TH1*>(before->At(4));
143 TH1* hprimvtxzA = dynamic_cast<TH1*>(after->At(4));
144
145 TH1* hITSclsB = dynamic_cast<TH1*>(before->At(5));
146 TH1* hITSclsA = dynamic_cast<TH1*>(after->At(5));
147 TH1* hTPCclsB = dynamic_cast<TH1*>(before->At(6));
148 TH1* hTPCclsA = dynamic_cast<TH1*>(after->At(6));
149 TH1* hTPCchi2B = dynamic_cast<TH1*>(before->At(7));
150 TH1* hTPCchi2A = dynamic_cast<TH1*>(after->At(7));
151 TH1* hdcaxyB = dynamic_cast<TH1*>(before->At(8));
152 TH1* hdcaxyA = dynamic_cast<TH1*>(after->At(8));
153 TH1* hdcazB = dynamic_cast<TH1*>(before->At(9));
154 TH1* hdcazA = dynamic_cast<TH1*>(after->At(9));
155 TH1* hphitrackletsB = dynamic_cast<TH1*>(before->At(10));
156 TH1* hphitrackletsA = dynamic_cast<TH1*>(after->At(10));
157 TH1* hphitracksB = dynamic_cast<TH1*>(before->At(11));
158 TH1* hphitracksA = dynamic_cast<TH1*>(after->At(11));
159 TH1* hetatrackletsB = dynamic_cast<TH1*>(before->At(12));
160 TH1* hetatrackletsA = dynamic_cast<TH1*>(after->At(12));
161 TH1* hetatracksB = dynamic_cast<TH1*>(before->At(13));
162 TH1* hetatracksA = dynamic_cast<TH1*>(after->At(13));
163 TH1* hprimvtxzTPCB = dynamic_cast<TH1*>(before->At(14));
164 TH1* hprimvtxzTPCA = dynamic_cast<TH1*>(after->At(14));
74d0740e 165 TH1* hptyieldB = dynamic_cast<TH1*>(before->At(15));
166 TH1* hptyieldA = dynamic_cast<TH1*>(after->At(15));
2948ac5a 167
168 AliMultiplicity* tracklets = const_cast<AliMultiplicity*>(event->GetMultiplicity());
169 Int_t ntracklets=0;
170 Int_t nspdclusters=0;
076df7bf 171 Int_t nspd1clusters=0;
172 Int_t ntrackletsA=0;
173 Int_t ntrackletsC=0;
2948ac5a 174 if (tracklets)
175 {
176 ntracklets = tracklets->GetNumberOfTracklets();
177 nspdclusters = tracklets->GetNumberOfITSClusters(0,1);
076df7bf 178 nspd1clusters = tracklets->GetNumberOfITSClusters(1);
2948ac5a 179 for (Int_t i=0; i<tracklets->GetNumberOfTracklets(); i++)
180 {
181 Bool_t pass=fTrackCuts->IsSelected(tracklets,i);
182 Float_t phi=tracklets->GetPhi(i);
183 Float_t eta=tracklets->GetEta(i);
184 hphitrackletsB->Fill(phi); if (pass) hphitrackletsA->Fill(phi);
185 hetatrackletsB->Fill(eta); if (pass) hetatrackletsA->Fill(eta);
076df7bf 186 if (eta>0) ntrackletsC++;
187 else ntrackletsA++;
2948ac5a 188 }
189
190 }
0e3a89af 191 //Int_t trackmult = AliESDtrackCuts::GetReferenceMultiplicity(event,kTRUE);
192 AliFlowEventCuts newcuts;
193 newcuts.SetRefMultCuts(fTrackCuts);
194 Int_t trackmult = newcuts.RefMult(event);
2948ac5a 195 Int_t refmult = fEventCuts->RefMult(event);
196 Bool_t passevent = fEventCuts->IsSelected(event);
197 htrackletmultB->Fill(ntracklets); if (passevent) htrackletmultA->Fill(ntracklets);
198 htrackmultB->Fill(trackmult); if (passevent) htrackmultA->Fill( trackmult);
199 hrefmultB->Fill(refmult); if (passevent) hrefmultA->Fill( refmult);
200 hspdclustersB->Fill(nspdclusters); if (passevent) hspdclustersA->Fill( nspdclusters);
201 const AliVertex* vertex = event->GetPrimaryVertex();
076df7bf 202 Float_t vtxz=0.0;
203 Float_t vtxx=0.0;
204 Float_t vtxy=0.0;
2948ac5a 205 if (vertex)
206 {
207 vtxz = vertex->GetZ();
076df7bf 208 vtxx = vertex->GetX();
209 vtxy = vertex->GetY();
2948ac5a 210 hprimvtxzB->Fill(vtxz); if (passevent) hprimvtxzA->Fill(vtxz);
211 }
212 const AliVertex* vertextpc = event->GetPrimaryVertexTPC();
076df7bf 213 Float_t vtxTPCx=0.0;
214 Float_t vtxTPCy=0.0;
215 Float_t vtxTPCz=0.0;
2948ac5a 216 if (vertextpc)
217 {
076df7bf 218 vtxTPCx = vertextpc->GetX();
219 vtxTPCy = vertextpc->GetY();
220 vtxTPCz = vertextpc->GetZ();
221 hprimvtxzTPCB->Fill(vtxTPCz); if (passevent) hprimvtxzTPCA->Fill(vtxTPCz);
2948ac5a 222 }
223 fTrackCuts->SetEvent(event);
076df7bf 224 Float_t meanpt=0;
225 Int_t ntracks=fTrackCuts->GetNumberOfInputObjects();
226 for (Int_t i=0; i<ntracks; i++)
2948ac5a 227 {
228 TObject* obj = fTrackCuts->GetInputObject(i);
229 Bool_t pass = fTrackCuts->IsSelected(obj,i);
230 Float_t dcaxy=0.0;
231 Float_t dcaz=0.0;
232 Float_t tpcchi2=0.0;
233 Float_t tpcchi2percls=0.0;
234 Int_t ntpccls=0;
235 Int_t nitscls=0;
236 Float_t eta=0.0;
237 Float_t phi=0.0;
74d0740e 238 Float_t pt=0.0;
2948ac5a 239 AliESDtrack* track = dynamic_cast<AliESDtrack*>(fTrackCuts->GetTrack());
240 if (track)
241 {
242 track->GetImpactParameters(dcaxy,dcaz);
243 tpcchi2=track->GetTPCchi2();
244 ntpccls=track->GetTPCNcls();
245 eta=track->Eta();
246 phi=track->Phi();
74d0740e 247 pt=track->Pt();
248 meanpt+=pt;
2948ac5a 249 tpcchi2percls= (ntpccls==0)?0.0:tpcchi2/ntpccls;
250 nitscls = track->GetNcls(0);
251 hITSclsB->Fill(nitscls); if (pass) hITSclsA->Fill( nitscls);
252 hTPCclsB->Fill(ntpccls); if (pass) hTPCclsA->Fill( ntpccls);
253 hTPCchi2B->Fill(tpcchi2percls); if (pass) hTPCchi2A->Fill( tpcchi2percls);
254 hdcaxyB->Fill(dcaxy); if (pass) hdcaxyA->Fill( dcaxy);
255 hdcazB->Fill(dcaz); if (pass) hdcazA->Fill(dcaz);
256 hetatracksB->Fill(eta); if (pass) hetatracksA->Fill(eta);
257 hphitracksB->Fill(phi); if (pass) hphitracksA->Fill(phi);
74d0740e 258 hptyieldB->Fill(pt); if (pass) hptyieldA->Fill(pt);
2948ac5a 259 }
260 }
076df7bf 261 meanpt = meanpt/ntracks;
262
263 ///////////////////////////////////////////////////////////////////////
264 //the ntuple part/////////////
265 AliFlowEventSimple* flowevent = dynamic_cast<AliFlowEventSimple*>(GetInputData(1));
266 if (flowevent && fFillNtuple)
267 {
268 AliFlowVector qvec = flowevent->GetQ(2);
269 Double_t qx = qvec.X();
270 Double_t qy = qvec.Y();
271 TString triggers = event->GetFiredTriggerClasses();
272 Int_t trig=0;
273 if (triggers.Contains("CMBAC-B-NOPF-ALL")) trig+=1;
274 if (triggers.Contains("CMBS2C-B-NOPF-ALL")) trig+=10;
275 if (triggers.Contains("CMBS2A-B-NOPF-ALL")) trig+=100;
276 fNtuple->Fill(meanpt,qx,qy,trackmult,refmult,trig,vtxTPCx,vtxTPCy,vtxTPCz,ntrackletsA,ntrackletsC);
277 }//if flowevent
2948ac5a 278}
279
280//________________________________________________________________________
281void AliAnalysisTaskQAflow::Terminate(Option_t *)
282{
283
284}
285
286//________________________________________________________________________
287AliAnalysisTaskQAflow::~AliAnalysisTaskQAflow()
288{
289 //dtor
290 delete fTrackCuts;
291 delete fEventCuts;
292}