7c3e5162 |
1 | |
2 | #include <TROOT.h> |
3 | #include <TSystem.h> |
4 | #include <TInterpreter.h> |
5 | #include <TChain.h> |
6 | #include <TFile.h> |
7 | #include <TList.h> |
8 | #include <iostream> |
bb8a464f |
9 | #include "TH1F.h" |
7c3e5162 |
10 | #include "TH2F.h" |
11 | #include "AliFMDAnalysisTaskDndeta.h" |
12 | #include "AliAnalysisManager.h" |
13 | #include "AliESDFMD.h" |
14 | #include "AliESDEvent.h" |
15 | #include "AliAODEvent.h" |
16 | #include "AliAODHandler.h" |
17 | #include "AliMCEventHandler.h" |
18 | #include "AliStack.h" |
19 | #include "AliLog.h" |
20 | #include "AliESDVertex.h" |
21 | #include "TMath.h" |
22 | #include "AliFMDAnaParameters.h" |
78f6f750 |
23 | //#include "AliFMDGeometry.h" |
bb8a464f |
24 | #include "AliGenEventHeader.h" |
cc066cb9 |
25 | #include "AliHeader.h" |
78f6f750 |
26 | //#include "TDatabasePDG.h" |
27 | //#include "TParticlePDG.h" |
cc066cb9 |
28 | #include "AliFMDStripIndex.h" |
7c3e5162 |
29 | ClassImp(AliFMDAnalysisTaskDndeta) |
30 | |
31 | |
32 | AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta() |
33 | : fDebug(0), |
34 | fOutputList(0), |
35 | fInputList(0), |
36 | fArray(0), |
37 | fInputArray(0), |
38 | fVertexString(0x0), |
39 | fNevents(), |
bb8a464f |
40 | fNMCevents(), |
41 | fStandalone(kTRUE), |
cc066cb9 |
42 | fMCevent(0), |
b83cd7a7 |
43 | fLastTrackByStrip(0), |
b3546e91 |
44 | fPrimary(kTRUE) |
7c3e5162 |
45 | { |
46 | // Default constructor |
47 | DefineInput (0, TList::Class()); |
48 | DefineOutput(0, TList::Class()); |
49 | } |
50 | //_____________________________________________________________________ |
51 | AliFMDAnalysisTaskDndeta::AliFMDAnalysisTaskDndeta(const char* name, Bool_t SE): |
52 | AliAnalysisTask(name, "Density"), |
53 | fDebug(0), |
54 | fOutputList(0), |
55 | fInputList(0), |
56 | fArray(), |
57 | fInputArray(0), |
58 | fVertexString(0x0), |
59 | fNevents(), |
bb8a464f |
60 | fNMCevents(), |
61 | fStandalone(kTRUE), |
cc066cb9 |
62 | fMCevent(0), |
b83cd7a7 |
63 | fLastTrackByStrip(0), |
b3546e91 |
64 | fPrimary(kTRUE) |
7c3e5162 |
65 | { |
66 | fStandalone = SE; |
67 | if(fStandalone) { |
68 | DefineInput (0, TList::Class()); |
69 | DefineInput(1, TObjString::Class()); |
70 | DefineOutput(0, TList::Class()); |
71 | |
72 | } |
73 | } |
74 | //_____________________________________________________________________ |
75 | void AliFMDAnalysisTaskDndeta::CreateOutputObjects() |
76 | { |
77 | AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance(); |
78 | |
79 | fArray.SetName("FMD"); |
80 | fArray.SetOwner(); |
81 | |
82 | if(!fOutputList) |
83 | fOutputList = new TList(); |
84 | fOutputList->SetName("BackgroundCorrected"); |
85 | |
86 | |
87 | TH2F* hMult = 0; |
cc066cb9 |
88 | TH1F* hHits = 0; |
bb8a464f |
89 | TH1F* hPrimVertexBin = 0; |
7c3e5162 |
90 | |
bb8a464f |
91 | |
5ea1e0a9 |
92 | TH2F* hBgTmp = pars->GetBackgroundCorrection(1, 'I', 0); |
bb8a464f |
93 | TH1F* hPrimary = new TH1F("hMultvsEta","hMultvsEta", |
5ea1e0a9 |
94 | hBgTmp->GetNbinsX(), |
95 | hBgTmp->GetXaxis()->GetXmin(), |
96 | hBgTmp->GetXaxis()->GetXmax()); |
bb8a464f |
97 | hPrimary->Sumw2(); |
98 | fOutputList->Add(hPrimary); |
7c3e5162 |
99 | Int_t nVtxbins = pars->GetNvtxBins(); |
100 | |
bb8a464f |
101 | |
7c3e5162 |
102 | for(Int_t det =1; det<=3;det++) |
103 | { |
104 | TObjArray* detArray = new TObjArray(); |
105 | detArray->SetName(Form("FMD%d",det)); |
106 | fArray.AddAtAndExpand(detArray,det); |
107 | Int_t nRings = (det==1 ? 1 : 2); |
108 | for(Int_t ring = 0;ring<nRings;ring++) |
109 | { |
110 | Char_t ringChar = (ring == 0 ? 'I' : 'O'); |
111 | Int_t nSec = (ring == 0 ? 20 : 40); |
112 | |
113 | TObjArray* vtxArray = new TObjArray(); |
114 | vtxArray->SetName(Form("FMD%d%c",det,ringChar)); |
115 | detArray->AddAtAndExpand(vtxArray,ring); |
116 | for(Int_t i = 0; i< nVtxbins; i++) { |
117 | TH2F* hBg = pars->GetBackgroundCorrection(det, ringChar, i); |
118 | hMult = new TH2F(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i), |
119 | hBg->GetNbinsX(), |
120 | hBg->GetXaxis()->GetXmin(), |
121 | hBg->GetXaxis()->GetXmax(), |
122 | nSec, 0, 2*TMath::Pi()); |
cc066cb9 |
123 | |
124 | hHits = new TH1F(Form("hHits_FMD%d%c_vtxbin%d",det,ringChar,i),Form("hHits_FMD%d%c_vtxbin%d",det,ringChar,i), |
125 | hBg->GetNbinsX(), |
126 | hBg->GetXaxis()->GetXmin(), |
127 | hBg->GetXaxis()->GetXmax()); |
128 | |
129 | |
130 | |
7c3e5162 |
131 | hMult->Sumw2(); |
cc066cb9 |
132 | hHits->Sumw2(); |
7c3e5162 |
133 | fOutputList->Add(hMult); |
cc066cb9 |
134 | fOutputList->Add(hHits); |
7c3e5162 |
135 | vtxArray->AddAtAndExpand(hMult,i); |
136 | |
137 | } |
138 | } |
139 | } |
140 | |
bb8a464f |
141 | for(Int_t i = 0; i< nVtxbins; i++) { |
142 | |
143 | hPrimVertexBin = new TH1F(Form("primmult_vtxbin%d",i), |
144 | Form("primmult_vtxbin%d",i), |
5ea1e0a9 |
145 | hBgTmp->GetNbinsX(), |
146 | hBgTmp->GetXaxis()->GetXmin(), |
147 | hBgTmp->GetXaxis()->GetXmax()); |
bb8a464f |
148 | hPrimVertexBin->Sumw2(); |
149 | fOutputList->Add(hPrimVertexBin); |
150 | |
151 | } |
152 | |
7c3e5162 |
153 | fNevents.SetBins(nVtxbins,0,nVtxbins); |
154 | fNevents.SetName("nEvents"); |
bb8a464f |
155 | fNMCevents.SetBins(nVtxbins,0,nVtxbins); |
156 | fNMCevents.SetName("nMCEvents"); |
7c3e5162 |
157 | fOutputList->Add(&fNevents); |
bb8a464f |
158 | fOutputList->Add(&fNMCevents); |
7c3e5162 |
159 | |
160 | } |
161 | //_____________________________________________________________________ |
162 | void AliFMDAnalysisTaskDndeta::ConnectInputData(Option_t */*option*/) |
163 | { |
164 | if(fStandalone) { |
165 | fInputList = (TList*)GetInputData(0); |
166 | fVertexString = (TObjString*)GetInputData(1); |
167 | } |
168 | } |
169 | //_____________________________________________________________________ |
170 | void AliFMDAnalysisTaskDndeta::Exec(Option_t */*option*/) |
171 | { |
172 | Int_t vtxbin = fVertexString->GetString().Atoi(); |
173 | fNevents.Fill(vtxbin); |
7c3e5162 |
174 | for(UShort_t det=1;det<=3;det++) { |
175 | //TObjArray* detInputArray = (TObjArray*)fInputArray->At(det); |
176 | TObjArray* detArray = (TObjArray*)fArray.At(det); |
177 | Int_t nRings = (det==1 ? 1 : 2); |
178 | for (UShort_t ir = 0; ir < nRings; ir++) { |
179 | Char_t ringChar = (ir == 0 ? 'I' : 'O'); |
180 | //TObjArray* vtxInputArray = (TObjArray*)detInputArray->At(ir); |
181 | TObjArray* vtxArray = (TObjArray*)detArray->At(ir); |
182 | TH2F* hMultTotal = (TH2F*)vtxArray->At(vtxbin); |
183 | |
184 | |
185 | TH2F* hMultInput = (TH2F*)fInputList->FindObject(Form("mult_FMD%d%c_vtxbin%d",det,ringChar,vtxbin)); |
1282ce49 |
186 | |
7c3e5162 |
187 | hMultTotal->Add(hMultInput); |
188 | |
189 | |
190 | } |
191 | } |
bb8a464f |
192 | |
b3546e91 |
193 | if(fMCevent && fPrimary) |
bb8a464f |
194 | ProcessPrimary(); |
195 | |
7c3e5162 |
196 | if(fStandalone) { |
197 | PostData(0, fOutputList); |
198 | } |
199 | |
200 | } |
201 | //_____________________________________________________________________ |
202 | void AliFMDAnalysisTaskDndeta::Terminate(Option_t */*option*/) { |
203 | |
1282ce49 |
204 | |
7c3e5162 |
205 | AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance(); |
206 | |
207 | Int_t nVtxbins = pars->GetNvtxBins(); |
208 | |
209 | for(UShort_t det=1;det<=3;det++) { |
210 | TObjArray* detArray = (TObjArray*)fArray.At(det); |
211 | Int_t nRings = (det==1 ? 1 : 2); |
212 | for (UShort_t ir = 0; ir < nRings; ir++) { |
213 | TObjArray* vtxArray = (TObjArray*)detArray->At(ir); |
1282ce49 |
214 | Char_t ringChar = (ir == 0 ? 'I' : 'O'); |
7c3e5162 |
215 | for(Int_t i =0; i<nVtxbins; i++) { |
216 | TH2F* hMultTotal = (TH2F*)vtxArray->At(i); |
1282ce49 |
217 | TH1D* hMultProj = hMultTotal->ProjectionX(Form("dNdeta_FMD%d%c_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY()); |
218 | fOutputList->Add(hMultProj); |
7c3e5162 |
219 | } |
220 | } |
221 | } |
1282ce49 |
222 | |
bb8a464f |
223 | } |
224 | //_____________________________________________________________________ |
225 | void AliFMDAnalysisTaskDndeta::ProcessPrimary() { |
cc066cb9 |
226 | |
227 | fLastTrackByStrip.Reset(-1); |
228 | |
bb8a464f |
229 | AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance(); |
230 | |
231 | AliMCParticle* particle = 0; |
232 | AliStack* stack = fMCevent->Stack(); |
233 | |
234 | TH1F* hPrimary = (TH1F*)fOutputList->FindObject("hMultvsEta"); |
cc066cb9 |
235 | AliHeader* header = fMCevent->Header(); |
236 | AliGenEventHeader* genHeader = header->GenEventHeader(); |
237 | |
238 | TArrayF vertex; |
239 | genHeader->PrimaryVertex(vertex); |
240 | if(TMath::Abs(vertex.At(2)) > pars->GetVtxCutZ()) |
241 | return; |
242 | Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins(); |
243 | Double_t vertexBinDouble = (vertex.At(2) + pars->GetVtxCutZ()) / delta; |
244 | Int_t vertexBin = (Int_t)vertexBinDouble; |
bb8a464f |
245 | |
246 | Bool_t firstTrack = kTRUE; |
247 | Int_t nTracks = fMCevent->GetNumberOfTracks(); |
248 | for(Int_t i = 0 ;i<nTracks;i++) { |
249 | particle = fMCevent->GetTrack(i); |
250 | if(!particle) |
251 | continue; |
6289b3e8 |
252 | |
bb8a464f |
253 | if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) { |
254 | hPrimary->Fill(particle->Eta()); |
cc066cb9 |
255 | |
6289b3e8 |
256 | |
bb8a464f |
257 | TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin)); |
258 | hPrimVtxBin->Fill(particle->Eta()); |
259 | if(firstTrack) { |
260 | fNMCevents.Fill(vertexBin); |
261 | firstTrack = kFALSE; |
262 | } |
263 | } |
b3546e91 |
264 | |
cc066cb9 |
265 | for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) { |
bb8a464f |
266 | |
cc066cb9 |
267 | AliTrackReference* ref = particle->GetTrackReference(j); |
268 | UShort_t det,sec,strip; |
269 | Char_t ring; |
270 | if(ref->DetectorId() != AliTrackReference::kFMD) |
271 | continue; |
272 | AliFMDStripIndex::Unpack(ref->UserId(),det,ring,sec,strip); |
b83cd7a7 |
273 | Float_t thisStripTrack = fLastTrackByStrip(det,ring,sec,strip); |
cc066cb9 |
274 | if(particle->Charge() != 0 && i != thisStripTrack ) { |
78f6f750 |
275 | //Double_t x,y,z; |
cc066cb9 |
276 | |
f58a4769 |
277 | Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));//-1*TMath::Log(TMath::Tan(0.5*theta)); |
cc066cb9 |
278 | TH1F* hHits = (TH1F*)fOutputList->FindObject(Form("hHits_FMD%d%c_vtxbin%d",det,ring,vertexBin)); |
279 | hHits->Fill(eta); |
280 | Float_t nstrips = (ring =='O' ? 256 : 512); |
281 | |
f58a4769 |
282 | //if(det == 1 && ring == 'I') |
283 | // std::cout<<"hit in "<<det<<" "<<ring<<" "<<sec<<" "<<strip<<" "<<std::endl; |
b83cd7a7 |
284 | fLastTrackByStrip(det,ring,sec,strip) = (Float_t)i; |
cc066cb9 |
285 | |
286 | if(strip >0) |
b83cd7a7 |
287 | fLastTrackByStrip(det,ring,sec,strip-1) = (Float_t)i; |
cc066cb9 |
288 | if(strip < (nstrips - 1)) |
b83cd7a7 |
289 | fLastTrackByStrip(det,ring,sec,strip+1) = (Float_t)i; |
cc066cb9 |
290 | |
291 | |
292 | } |
b3546e91 |
293 | } |
cc066cb9 |
294 | |
295 | |
b3546e91 |
296 | } |
bb8a464f |
297 | |
7c3e5162 |
298 | } |
299 | //_____________________________________________________________________ |
300 | // |
301 | // |
302 | // EOF |