]>
Commit | Line | Data |
---|---|---|
df0d5480 | 1 | /************************************************************************** |
2 | * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | // Thil class computes background corrections for the FMD. The correction is computed | |
17 | // in eta,phi cells and the objects stored can be put into alien to use with analysis. | |
18 | // It is based on the AliFMDInput class that is used to loop over hits and primaries. | |
19 | // | |
20 | // Author: Hans Hjersing Dalsgaard, NBI, hans.dalsgaard@cern.ch | |
21 | // | |
22 | // | |
23 | ||
24 | #include "AliSimulation.h" | |
25 | #include "TStopwatch.h" | |
26 | #include "iostream" | |
27 | #include "TGrid.h" | |
28 | #include "AliRunLoader.h" | |
29 | #include "AliGeomManager.h" | |
30 | #include "AliFMDGeometry.h" | |
31 | #include "AliStack.h" | |
32 | #include "TParticle.h" | |
33 | #include "TDatabasePDG.h" | |
34 | #include "TParticlePDG.h" | |
35 | #include "AliRun.h" | |
36 | #include "AliFMDBackgroundCorrection.h" | |
37 | #include "TSystem.h" | |
38 | #include "AliCDBManager.h" | |
39 | #include "TTree.h" | |
40 | #include "TClonesArray.h" | |
41 | #include "TBranch.h" | |
42 | #include "AliFMDHit.h" | |
43 | #include "AliLoader.h" | |
44 | #include "AliFMD.h" | |
45 | #include "TH2F.h" | |
46 | #include "AliGenEventHeader.h" | |
47 | #include "AliHeader.h" | |
48 | #include "TFile.h" | |
49 | #include "TAxis.h" | |
50 | #include "AliCDBId.h" | |
51 | #include "AliCDBMetaData.h" | |
52 | #include "TROOT.h" | |
53 | #include "AliFMDParameters.h" | |
54 | #include "AliLog.h" | |
26c8b5ff | 55 | #include "TList.h" |
8dc823cc | 56 | #include "AliFMDAnaParameters.h" |
b82e76e0 | 57 | #include "AliFMDAnaCalibBackgroundCorrection.h" |
df0d5480 | 58 | |
59 | ClassImp(AliFMDBackgroundCorrection) | |
60 | //_____________________________________________________________________ | |
61 | AliFMDBackgroundCorrection::AliFMDBackgroundCorrection() : | |
62 | TNamed(), | |
ca3af5d9 | 63 | fCorrectionArray(), |
64 | fPrimaryList() | |
df0d5480 | 65 | {} |
66 | ||
67 | //_____________________________________________________________________ | |
68 | AliFMDBackgroundCorrection::AliFMDInputBG::AliFMDInputBG() : | |
69 | AliFMDInput(), | |
8dc823cc | 70 | fPrimaryArray(), |
71 | fHitArray(), | |
df0d5480 | 72 | fPrim(0), |
73 | fHits(0), | |
74 | fZvtxCut(0), | |
75 | fNvtxBins(0), | |
76 | fPrevTrack(-1), | |
77 | fPrevDetector(-1), | |
78 | fPrevRing('Q'), | |
79 | fNbinsEta(100) | |
80 | { | |
81 | AddLoad(kTracks); | |
82 | AddLoad(kHits); | |
83 | AddLoad(kKinematics); | |
84 | AddLoad(kHeader); | |
85 | } | |
86 | ||
87 | //_____________________________________________________________________ | |
88 | ||
ca3af5d9 | 89 | void AliFMDBackgroundCorrection::GenerateBackgroundCorrection(const Int_t nvtxbins, |
df0d5480 | 90 | Float_t zvtxcut, |
ca3af5d9 | 91 | const Int_t nBinsEta, |
df0d5480 | 92 | Bool_t storeInAlien, |
93 | Int_t runNo, | |
94 | Int_t endRunNo, | |
95 | const Char_t* filename, | |
96 | Bool_t simulate, | |
97 | Int_t nEvents) { | |
98 | ||
ca3af5d9 | 99 | //TGrid::Connect("alien:",0,0,"t"); |
df0d5480 | 100 | if(simulate) |
101 | Simulate(nEvents); | |
102 | else { | |
ca3af5d9 | 103 | //AliCDBManager::Instance()->SetDefaultStorage("alien://Folder=/alice/data/2008/LHC08d/OCDB/"); |
104 | AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT"); | |
df0d5480 | 105 | AliCDBManager::Instance()->SetRun(runNo); |
106 | ||
107 | #if defined(__CINT__) | |
108 | gSystem->Load("liblhapdf"); | |
109 | gSystem->Load("libEGPythia6"); | |
110 | gSystem->Load("libpythia6"); | |
111 | gSystem->Load("libAliPythia6"); | |
112 | gSystem->Load("libgeant321"); | |
113 | #endif | |
ca3af5d9 | 114 | // |
df0d5480 | 115 | } |
116 | ||
117 | //Setting up the geometry | |
118 | //----------------------------------------------- | |
119 | if (AliGeomManager::GetGeometry() == NULL) | |
120 | AliGeomManager::LoadGeometry(); | |
121 | ||
122 | AliFMDGeometry* geo = AliFMDGeometry::Instance(); | |
123 | geo->Init(); | |
124 | geo->InitTransformations(); | |
ca3af5d9 | 125 | |
df0d5480 | 126 | AliInfo("Processing hits and primaries "); |
127 | ||
128 | AliFMDInputBG input; | |
129 | input.SetVtxCutZ(zvtxcut); | |
130 | input.SetNvtxBins(nvtxbins); | |
131 | input.SetNbinsEta(nBinsEta); | |
132 | input.Run(); | |
133 | ||
134 | AliInfo(Form("Found %d primaries and %d hits.", input.GetNprim(),input.GetNhits())); | |
135 | ||
136 | TObjArray* hitArray = input.GetHits(); | |
137 | TObjArray* primaryArray = input.GetPrimaries(); | |
138 | fCorrectionArray.SetName("FMD_bg_correction"); | |
139 | fCorrectionArray.SetOwner(); | |
ca3af5d9 | 140 | |
141 | TList* primaryList = new TList(); | |
26c8b5ff | 142 | primaryList->SetName("primaries"); |
ca3af5d9 | 143 | |
26c8b5ff | 144 | TList* hitList = new TList(); |
145 | hitList->SetName("hits"); | |
146 | TList* corrList = new TList(); | |
147 | corrList->SetName("corrections"); | |
148 | ||
b82e76e0 | 149 | AliFMDAnaCalibBackgroundCorrection* background = new AliFMDAnaCalibBackgroundCorrection(); |
150 | ||
df0d5480 | 151 | for(Int_t det= 1; det <=3; det++) { |
152 | Int_t nRings = (det==1 ? 1 : 2); | |
153 | ||
154 | TObjArray* detArrayCorrection = new TObjArray(); | |
155 | detArrayCorrection->SetName(Form("FMD%d",det)); | |
156 | fCorrectionArray.AddAtAndExpand(detArrayCorrection,det); | |
157 | ||
158 | ||
159 | for(Int_t iring = 0; iring<nRings; iring++) { | |
160 | TObjArray* primRingArray = (TObjArray*)primaryArray->At(iring); | |
161 | Char_t ringChar = (iring == 0 ? 'I' : 'O'); | |
162 | TObjArray* vtxArrayCorrection = new TObjArray(); | |
163 | vtxArrayCorrection->SetName(Form("FMD%d%c",det,ringChar)); | |
164 | detArrayCorrection->AddAtAndExpand(vtxArrayCorrection,iring); | |
165 | ||
166 | for(Int_t vertexBin=0;vertexBin<nvtxbins;vertexBin++) { | |
167 | TObjArray* detArray = (TObjArray*)hitArray->At(det); | |
168 | TObjArray* vtxArray = (TObjArray*)detArray->At(iring); | |
169 | TH2F* hHits = (TH2F*)vtxArray->At(vertexBin); | |
26c8b5ff | 170 | hitList->Add(hHits); |
df0d5480 | 171 | TH2F* hPrimary = (TH2F*)primRingArray->At(vertexBin); |
26c8b5ff | 172 | primaryList->Add(hPrimary); |
df0d5480 | 173 | TH2F* hCorrection = (TH2F*)hHits->Clone(Form("FMD%d%c_vtxbin_%d_correction",det,ringChar,vertexBin)); |
ca3af5d9 | 174 | hCorrection->Divide(hPrimary); |
26c8b5ff | 175 | corrList->Add(hCorrection); |
8dc823cc | 176 | hCorrection->SetTitle(hCorrection->GetName()); |
df0d5480 | 177 | vtxArrayCorrection->AddAtAndExpand(hCorrection,vertexBin); |
b82e76e0 | 178 | background->SetBgCorrection(det,ringChar,vertexBin,hCorrection); |
df0d5480 | 179 | } |
180 | ||
181 | } | |
182 | } | |
183 | ||
184 | TAxis refAxis(nvtxbins,-1*zvtxcut,zvtxcut); | |
185 | ||
b82e76e0 | 186 | background->SetRefAxis(&refAxis); |
df0d5480 | 187 | |
8dc823cc | 188 | TFile* fout = new TFile(filename,"RECREATE"); |
df0d5480 | 189 | refAxis.Write("vertexbins"); |
26c8b5ff | 190 | |
191 | hitList->Write(); | |
192 | primaryList->Write(); | |
193 | corrList->Write(); | |
ca3af5d9 | 194 | |
8dc823cc | 195 | TObjArray* container = new TObjArray(); |
196 | container->SetOwner(); | |
197 | container->AddAtAndExpand(&refAxis,0); | |
198 | container->AddAtAndExpand(&fCorrectionArray,1); | |
199 | container->AddAtAndExpand(hitArray,2); | |
200 | container->AddAtAndExpand(primaryArray,3); | |
b82e76e0 | 201 | |
202 | ||
df0d5480 | 203 | if(storeInAlien) { |
204 | AliCDBManager* cdb = AliCDBManager::Instance(); | |
205 | cdb->SetDefaultStorage("local://$ALICE_ROOT"); | |
8dc823cc | 206 | AliCDBId id(AliFMDAnaParameters::GetBackgroundPath(),runNo,endRunNo); |
df0d5480 | 207 | |
208 | AliCDBMetaData* meta = new AliCDBMetaData; | |
209 | meta->SetResponsible(gSystem->GetUserInfo()->fRealName.Data()); | |
210 | meta->SetAliRootVersion(gROOT->GetVersion()); | |
211 | meta->SetBeamPeriod(1); | |
212 | meta->SetComment("Background Correction for FMD"); | |
b82e76e0 | 213 | meta->SetProperty("key1", background ); |
214 | cdb->Put(background, id, meta); | |
df0d5480 | 215 | |
216 | } | |
217 | ||
8dc823cc | 218 | fout->Close(); |
219 | ||
df0d5480 | 220 | |
ca3af5d9 | 221 | } |
222 | ||
df0d5480 | 223 | //_____________________________________________________________________ |
224 | void AliFMDBackgroundCorrection::Simulate(Int_t nEvents) { | |
225 | ||
226 | AliSimulation sim ; | |
227 | sim.SetRunNumber(0); | |
228 | TGrid::Connect("alien:",0,0,"t"); | |
229 | sim.SetDefaultStorage("alien://Folder=/alice/data/2008/LHC08d/OCDB/"); | |
230 | sim.SetConfigFile("Config.C"); | |
231 | sim.SetRunQA("FMD:"); | |
232 | TStopwatch timer; | |
233 | timer.Start(); | |
234 | sim.RunSimulation(nEvents); | |
235 | timer.Stop(); | |
236 | timer.Print(); | |
237 | ||
238 | } | |
239 | ||
240 | //_____________________________________________________________________ | |
ca3af5d9 | 241 | Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::ProcessHit( AliFMDHit* h, TParticle* p) { |
df0d5480 | 242 | |
df0d5480 | 243 | AliGenEventHeader* genHeader = fLoader->GetHeader()->GenEventHeader(); |
244 | TArrayF vertex; | |
245 | genHeader->PrimaryVertex(vertex); | |
246 | ||
247 | if(TMath::Abs(vertex.At(2)) > fZvtxCut) | |
248 | return kTRUE; | |
249 | ||
250 | Double_t delta = 2*fZvtxCut/fNvtxBins; | |
251 | Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta; | |
252 | Int_t vertexBin = (Int_t)vertexBinDouble; | |
253 | ||
ca3af5d9 | 254 | Int_t i = h->Track(); |
255 | ||
df0d5480 | 256 | if(h) |
257 | if(h->Q() != 0 && (i != fPrevTrack || h->Detector() != fPrevDetector || h->Ring() != fPrevRing)) { | |
258 | ||
259 | Int_t det = h->Detector(); | |
260 | Char_t ring = h->Ring(); | |
df0d5480 | 261 | Int_t iring = (ring == 'I' ? 0 : 1); |
262 | ||
263 | TObjArray* detArray = (TObjArray*)fHitArray.At(det); | |
264 | TObjArray* vtxArray = (TObjArray*)detArray->At(iring); | |
265 | TH2F* hHits = (TH2F*)vtxArray->At(vertexBin); | |
266 | ||
267 | Float_t phi = TMath::ATan2(h->Py(),h->Px()); | |
268 | if(phi<0) | |
269 | phi = phi+2*TMath::Pi(); | |
270 | Float_t theta = TMath::ATan2(TMath::Sqrt(TMath::Power(h->X(),2)+TMath::Power(h->Y(),2)),h->Z()+vertex.At(2)); | |
271 | Float_t eta = -1*TMath::Log(TMath::Tan(0.5*theta)); | |
272 | hHits->Fill(eta,phi); | |
273 | fHits++; | |
df0d5480 | 274 | fPrevDetector = det; |
275 | fPrevRing = ring; | |
276 | } | |
ca3af5d9 | 277 | |
df0d5480 | 278 | fPrevTrack = i; |
279 | return kTRUE; | |
280 | } | |
281 | ||
282 | //_____________________________________________________________________ | |
283 | Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Init() { | |
284 | ||
285 | fPrimaryArray.SetOwner(); | |
286 | fPrimaryArray.SetName("FMD_primary"); | |
287 | ||
288 | for(Int_t iring = 0; iring<2;iring++) { | |
289 | Char_t ringChar = (iring == 0 ? 'I' : 'O'); | |
290 | TObjArray* ringArray = new TObjArray(); | |
291 | ringArray->SetName(Form("FMD_%c",ringChar)); | |
292 | fPrimaryArray.AddAtAndExpand(ringArray,iring); | |
293 | Int_t nSec = (iring == 1 ? 40 : 20); | |
294 | for(Int_t v=0; v<fNvtxBins;v++) { | |
295 | ||
296 | TH2F* hPrimary = new TH2F(Form("hPrimary_FMD_%c_vtx%d",ringChar,v), | |
297 | Form("hPrimary_FMD_%c_vtx%d",ringChar,v), | |
298 | fNbinsEta, -6,6, nSec, 0,2*TMath::Pi()); | |
b82e76e0 | 299 | hPrimary->Sumw2(); |
df0d5480 | 300 | ringArray->AddAtAndExpand(hPrimary,v); |
301 | } | |
302 | } | |
303 | ||
304 | ||
305 | fHitArray.SetOwner(); | |
306 | fHitArray.SetName("FMD_hits"); | |
307 | ||
308 | for(Int_t det =1; det<=3;det++) | |
309 | { | |
310 | TObjArray* detArrayHits = new TObjArray(); | |
311 | detArrayHits->SetName(Form("FMD%d",det)); | |
312 | fHitArray.AddAtAndExpand(detArrayHits,det); | |
313 | Int_t nRings = (det==1 ? 1 : 2); | |
314 | for(Int_t ring = 0;ring<nRings;ring++) | |
315 | { | |
316 | Int_t nSec = (ring == 1 ? 40 : 20); | |
317 | Char_t ringChar = (ring == 0 ? 'I' : 'O'); | |
318 | TObjArray* vtxArrayHits = new TObjArray(); | |
319 | vtxArrayHits->SetName(Form("FMD%d%c",det,ringChar)); | |
320 | detArrayHits->AddAtAndExpand(vtxArrayHits,ring); | |
321 | for(Int_t v=0; v<fNvtxBins;v++) | |
322 | { | |
323 | ||
324 | TH2F* hHits = new TH2F(Form("hHits_FMD%d%c_vtx%d",det,ringChar,v), | |
325 | Form("hHits_FMD%d%c_vtx%d",det,ringChar,v), | |
326 | fNbinsEta, -6,6, nSec, 0, 2*TMath::Pi()); | |
b82e76e0 | 327 | hHits->Sumw2(); |
df0d5480 | 328 | vtxArrayHits->AddAtAndExpand(hHits,v); |
329 | ||
330 | } | |
331 | } | |
332 | ||
333 | } | |
334 | ||
335 | ||
336 | ||
337 | AliFMDInput::Init(); | |
338 | ||
339 | ||
340 | return kTRUE; | |
341 | } | |
ca3af5d9 | 342 | //_____________________________________________________________________ |
343 | ||
344 | Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Begin(Int_t event ) { | |
df0d5480 | 345 | |
ca3af5d9 | 346 | Bool_t retVal = AliFMDInput::Begin(event); |
347 | ||
348 | AliStack* partStack=fLoader->Stack(); | |
349 | ||
350 | Int_t nTracks=partStack->GetNtrack(); | |
351 | AliGenEventHeader* genHeader = fLoader->GetHeader()->GenEventHeader(); | |
352 | TArrayF vertex; | |
353 | genHeader->PrimaryVertex(vertex); | |
354 | ||
355 | if(TMath::Abs(vertex.At(2)) > fZvtxCut) | |
356 | return kTRUE; | |
357 | ||
358 | Double_t delta = 2*fZvtxCut/fNvtxBins; | |
359 | Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta; | |
360 | Int_t vertexBin = (Int_t)vertexBinDouble; | |
361 | ||
362 | for(Int_t j=0;j<nTracks;j++) | |
363 | { | |
364 | TParticle* p = partStack->Particle(j); | |
365 | TDatabasePDG* pdgDB = TDatabasePDG::Instance(); | |
366 | TParticlePDG* pdgPart = pdgDB->GetParticle(p->GetPdgCode()); | |
367 | Float_t charge = (pdgPart ? pdgPart->Charge() : 0); | |
368 | Float_t phi = TMath::ATan2(p->Py(),p->Px()); | |
369 | ||
370 | if(phi<0) | |
371 | phi = phi+2*TMath::Pi(); | |
372 | if(p->Theta() == 0) continue; | |
373 | Float_t eta = -1*TMath::Log(TMath::Tan(0.5*p->Theta())); | |
374 | ||
375 | Bool_t primary = (charge!=0)&&(TMath::Abs(p->Vx() - vertex.At(0))<0.01)&&(TMath::Abs(p->Vy() - vertex.At(1))<0.01)&&(TMath::Abs(p->Vz() - vertex.At(2))<0.01); | |
376 | if(primary) { | |
377 | fPrim++; | |
378 | TObjArray* innerArray = (TObjArray*)fPrimaryArray.At(0); | |
379 | TObjArray* outerArray = (TObjArray*)fPrimaryArray.At(1); | |
380 | ||
381 | TH2F* hPrimaryInner = (TH2F*)innerArray->At(vertexBin); | |
382 | TH2F* hPrimaryOuter = (TH2F*)outerArray->At(vertexBin); | |
383 | hPrimaryInner->Fill(eta,phi); | |
384 | hPrimaryOuter->Fill(eta,phi); | |
385 | } | |
386 | } | |
387 | ||
388 | return retVal; | |
389 | } | |
390 | //_____________________________________________________________________ | |
391 | // | |
392 | // EOF | |
393 | // | |
df0d5480 | 394 | |
395 |