]>
Commit | Line | Data |
---|---|---|
762e8424 | 1 | // $Id$ |
2 | // | |
3 | // Jet modelling task. | |
4 | // | |
5 | // Author: Salvatore Aiola, Constantin Loizides | |
6 | ||
7 | #include <TClonesArray.h> | |
8 | #include <TLorentzVector.h> | |
9 | #include <TRandom3.h> | |
10 | ||
11 | #include "AliAnalysisManager.h" | |
12 | #include "AliVEvent.h" | |
13 | #include "AliVCluster.h" | |
14 | #include "AliEMCALDigit.h" | |
15 | #include "AliEMCALRecPoint.h" | |
16 | #include "AliPicoTrack.h" | |
17 | #include "AliEMCALGeometry.h" | |
18 | #include "AliLog.h" | |
19 | ||
20 | #include "AliJetModelBaseTask.h" | |
21 | ||
22 | ClassImp(AliJetModelBaseTask) | |
23 | ||
24 | //________________________________________________________________________ | |
25 | AliJetModelBaseTask::AliJetModelBaseTask() : | |
26 | AliAnalysisTaskSE("AliJetModelBaseTask"), | |
27 | fGeomName(), | |
28 | fTracksName(), | |
29 | fOutTracksName(), | |
30 | fCaloName(), | |
31 | fOutCaloName(), | |
32 | fSuffix(), | |
33 | fEtaMin(-1), | |
34 | fEtaMax(1), | |
35 | fPhiMin(0), | |
36 | fPhiMax(TMath::Pi() * 2), | |
37 | fPtMin(0), | |
38 | fPtMax(0), | |
39 | fCopyArray(kTRUE), | |
40 | fNClusters(0), | |
41 | fNTracks(0), | |
42 | fGeom(0), | |
43 | fClusters(0), | |
44 | fOutClusters(0), | |
45 | fTracks(0), | |
46 | fOutTracks(0) | |
47 | { | |
48 | // Default constructor. | |
49 | } | |
50 | ||
51 | //________________________________________________________________________ | |
52 | AliJetModelBaseTask::AliJetModelBaseTask(const char *name) : | |
53 | AliAnalysisTaskSE(name), | |
54 | fGeomName(""), | |
55 | fTracksName("PicoTracks"), | |
56 | fOutTracksName("PicoTracksEmbedded"), | |
57 | fCaloName("CaloClustersCorr"), | |
58 | fOutCaloName("CaloClustersCorrEmbedded"), | |
59 | fSuffix("Processed"), | |
60 | fEtaMin(-1), | |
61 | fEtaMax(1), | |
62 | fPhiMin(0), | |
63 | fPhiMax(TMath::Pi() * 2), | |
64 | fPtMin(50), | |
65 | fPtMax(60), | |
66 | fCopyArray(kTRUE), | |
67 | fNClusters(0), | |
68 | fNTracks(1), | |
69 | fGeom(0), | |
70 | fClusters(0), | |
71 | fOutClusters(0), | |
72 | fTracks(0), | |
73 | fOutTracks(0) | |
74 | { | |
75 | // Standard constructor. | |
76 | ||
77 | } | |
78 | ||
79 | //________________________________________________________________________ | |
80 | AliJetModelBaseTask::~AliJetModelBaseTask() | |
81 | { | |
82 | // Destructor | |
83 | } | |
84 | ||
85 | //________________________________________________________________________ | |
86 | void AliJetModelBaseTask::Init() | |
87 | { | |
88 | // Init task. | |
89 | ||
90 | if (fPtMax < fPtMin) { | |
91 | AliWarning (Form("PtMax (%f) < PtMin (%f), setting PtMax = PtMin = %f", fPtMax, fPtMin, fPtMin)); | |
92 | fPtMax = fPtMin; | |
93 | } | |
94 | ||
95 | if (fEtaMax < fEtaMin) { | |
96 | AliWarning (Form("EtaMax (%f) < EtaMin (%f), setting EtaMax = EtaMin = %f", fEtaMax, fEtaMin, fEtaMin)); | |
97 | fEtaMax = fEtaMin; | |
98 | } | |
99 | ||
100 | if (fPhiMax < fPhiMin) { | |
101 | AliWarning (Form("PhiMax (%f) < PhiMin (%f), setting PhiMax = PhiMin = %f", fPhiMax, fPhiMin, fPhiMin)); | |
102 | fPhiMax = fPhiMin; | |
103 | } | |
104 | ||
105 | if (fNTracks > 0) { | |
106 | fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName)); | |
107 | ||
108 | if (strcmp(fTracks->GetClass()->GetName(), "AliPicoTrack")) { | |
109 | AliError("Can only embed PicoTracks!"); | |
110 | return; | |
111 | } | |
112 | ||
113 | if (!fTracks) { | |
114 | AliError(Form("Couldn't retrieve tracks with name %s!", fTracksName.Data())); | |
115 | return; | |
116 | } | |
117 | ||
118 | if (!fOutTracks) { | |
119 | fOutTracksName = fTracksName; | |
120 | if (fCopyArray) { | |
121 | fOutTracksName += fSuffix; | |
122 | fOutTracks = new TClonesArray(*fTracks); | |
123 | fOutTracks->SetName(fOutTracksName); | |
124 | } | |
125 | else { | |
126 | fOutTracks = fTracks; | |
127 | } | |
128 | } | |
129 | ||
130 | if (fCopyArray) { | |
131 | if (!(InputEvent()->FindListObject(fOutTracksName))) | |
132 | InputEvent()->AddObject(fOutTracks); | |
133 | fOutTracks->Clear(); | |
134 | } | |
135 | } | |
136 | ||
137 | if (fNClusters > 0) { | |
138 | fClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName)); | |
139 | ||
140 | if (!fClusters) { | |
141 | AliError(Form("Couldn't retrieve clusters with name %s!", fCaloName.Data())); | |
142 | return; | |
143 | } | |
144 | ||
145 | if (!fOutClusters) { | |
146 | fOutCaloName = fCaloName; | |
147 | if (fCopyArray) { | |
148 | fOutCaloName += fSuffix; | |
149 | fOutClusters = new TClonesArray(*fClusters); | |
150 | fOutClusters->SetName(fOutCaloName); | |
151 | } | |
152 | else { | |
153 | fOutClusters = fClusters; | |
154 | } | |
155 | } | |
156 | ||
157 | if (fCopyArray) { | |
158 | if (!(InputEvent()->FindListObject(fOutCaloName))) | |
159 | InputEvent()->AddObject(fOutClusters); | |
160 | fOutClusters->Clear(); | |
161 | } | |
762e8424 | 162 | |
bb838b50 | 163 | if (!fGeom) { |
164 | if (fGeomName.Length() > 0) { | |
165 | fGeom = AliEMCALGeometry::GetInstance(fGeomName); | |
166 | if (!fGeom) | |
167 | AliError(Form("Could not get geometry with name %s!", fGeomName.Data())); | |
168 | } else { | |
169 | fGeom = AliEMCALGeometry::GetInstance(); | |
170 | if (!fGeom) | |
171 | AliError("Could not get default geometry!"); | |
172 | } | |
762e8424 | 173 | } |
bb838b50 | 174 | |
3d0b9d27 | 175 | const Double_t EmcalMinEta = fGeom->GetArm1EtaMin(); |
176 | const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax(); | |
177 | const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad(); | |
178 | const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad(); | |
179 | ||
180 | if (fEtaMax > EmcalMaxEta) fEtaMax = EmcalMaxEta; | |
181 | if (fEtaMax < EmcalMinEta) fEtaMax = EmcalMinEta; | |
182 | if (fEtaMin > EmcalMaxEta) fEtaMin = EmcalMaxEta; | |
183 | if (fEtaMin < EmcalMinEta) fEtaMin = EmcalMinEta; | |
762e8424 | 184 | |
3d0b9d27 | 185 | if (fPhiMax > EmcalMaxPhi) fPhiMax = EmcalMaxPhi; |
186 | if (fPhiMax < EmcalMinPhi) fPhiMax = EmcalMinPhi; | |
187 | if (fPhiMin > EmcalMaxPhi) fPhiMin = EmcalMaxPhi; | |
188 | if (fPhiMin < EmcalMinPhi) fPhiMin = EmcalMinPhi; | |
762e8424 | 189 | } |
190 | } | |
191 | ||
192 | //________________________________________________________________________ | |
193 | void AliJetModelBaseTask::GetRandomCell(Double_t &eta, Double_t &phi, Int_t &absId) | |
194 | { | |
195 | Int_t repeats = 0; | |
05debff6 | 196 | Double_t rndEta = eta; |
197 | Double_t rndPhi = phi; | |
762e8424 | 198 | do { |
05debff6 | 199 | if (eta < -100) |
200 | rndEta = GetRandomEta(); | |
201 | if (phi < 0) | |
202 | rndPhi = GetRandomPhi(); | |
203 | fGeom->GetAbsCellIdFromEtaPhi(rndEta, rndPhi, absId); | |
762e8424 | 204 | repeats++; |
205 | } while (absId == -1 && repeats < 100); | |
206 | ||
207 | if (!(absId > -1)) { | |
208 | AliWarning(Form("Could not extract random cluster! Random eta-phi extracted more than 100 times!\n" | |
209 | "eta [%f, %f], phi [%f, %f]\n", fEtaMin, fEtaMax, fPhiMin, fPhiMax)); | |
210 | } | |
05debff6 | 211 | else { |
212 | eta = rndEta; | |
213 | phi = rndPhi; | |
214 | } | |
762e8424 | 215 | } |
216 | ||
217 | //________________________________________________________________________ | |
218 | Double_t AliJetModelBaseTask::GetRandomEta() | |
219 | { | |
220 | return gRandom->Rndm() * (fEtaMax - fEtaMin) + fEtaMin; | |
221 | } | |
222 | ||
223 | //________________________________________________________________________ | |
224 | Double_t AliJetModelBaseTask::GetRandomPhi() | |
225 | { | |
226 | return gRandom->Rndm() * (fPhiMax - fPhiMin) + fPhiMin; | |
227 | } | |
228 | ||
229 | //________________________________________________________________________ | |
230 | Double_t AliJetModelBaseTask::GetRandomPt() | |
231 | { | |
232 | return gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin; | |
233 | } | |
234 | ||
235 | //________________________________________________________________________ | |
236 | AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Double_t eta, Double_t phi) | |
237 | { | |
238 | Int_t absId = 0; | |
05debff6 | 239 | if (eta < -100 || phi < 0) { |
762e8424 | 240 | GetRandomCell(eta, phi, absId); |
241 | } | |
242 | else { | |
243 | fGeom->EtaPhiFromIndex(absId, eta, phi); | |
244 | } | |
245 | ||
246 | if (absId == -1) { | |
247 | AliWarning(Form("Unable to embed cluster in eta = %f, phi = %f!" | |
248 | " Maybe the eta-phi range is not inside the EMCal acceptance (eta = [%f, %f], phi = [%f, %f])", | |
249 | eta, phi, fEtaMin, fEtaMax, fPhiMin, fPhiMax)); | |
250 | return 0; | |
251 | } | |
252 | ||
253 | if (e < 0) { | |
254 | Double_t pt = GetRandomPt(); | |
255 | TLorentzVector nPart; | |
256 | nPart.SetPtEtaPhiM(pt, eta, phi, 0); | |
257 | e = nPart.E(); | |
258 | } | |
259 | ||
260 | return AddCluster(e, absId); | |
261 | } | |
262 | ||
263 | //________________________________________________________________________ | |
264 | AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Int_t absId) | |
265 | { | |
266 | const Int_t nClusters = fOutClusters->GetEntriesFast(); | |
267 | ||
268 | TClonesArray digits("AliEMCALDigit", 1); | |
269 | ||
270 | AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(digits.New(0)); | |
271 | digit->SetId(absId); | |
272 | digit->SetIndexInList(0); | |
273 | digit->SetType(AliEMCALDigit::kHG); | |
274 | digit->SetAmplitude(e); | |
275 | ||
276 | AliEMCALRecPoint recPoint(""); | |
277 | recPoint.AddDigit(*digit, e, kFALSE); | |
278 | recPoint.EvalGlobalPosition(0, &digits); | |
279 | ||
280 | TVector3 gpos; | |
281 | recPoint.GetGlobalPosition(gpos); | |
282 | Float_t g[3]; | |
283 | gpos.GetXYZ(g); | |
284 | ||
285 | AliVCluster *cluster = static_cast<AliVCluster*>(fOutClusters->New(nClusters)); | |
286 | cluster->SetType(AliVCluster::kEMCALClusterv1); | |
287 | cluster->SetE(recPoint.GetEnergy()); | |
288 | cluster->SetPosition(g); | |
289 | cluster->SetNCells(1); | |
290 | UShort_t shortAbsId = absId; | |
291 | cluster->SetCellsAbsId(&shortAbsId); | |
292 | Double32_t fract = 1; | |
293 | cluster->SetCellsAmplitudeFraction(&fract); | |
294 | cluster->SetID(nClusters); | |
295 | cluster->SetEmcCpvDistance(-1); | |
296 | cluster->SetChi2(100); // MC flag! | |
297 | ||
298 | return cluster; | |
299 | } | |
300 | ||
301 | //________________________________________________________________________ | |
302 | AliPicoTrack* AliJetModelBaseTask::AddTrack(Double_t pt, Double_t eta, Double_t phi) | |
303 | { | |
304 | const Int_t nTracks = fOutTracks->GetEntriesFast(); | |
305 | ||
306 | if (pt < 0) | |
307 | pt = GetRandomPt(); | |
05debff6 | 308 | if (eta < -100) |
762e8424 | 309 | eta = GetRandomEta(); |
310 | if (phi < 0) | |
311 | phi = GetRandomPhi(); | |
3d0b9d27 | 312 | |
762e8424 | 313 | AliPicoTrack *track = new ((*fOutTracks)[nTracks]) AliPicoTrack(pt, |
314 | eta, | |
315 | phi, | |
316 | 1, | |
317 | 100, // MC flag! | |
318 | 0, | |
319 | 0, | |
320 | kFALSE); | |
321 | return track; | |
322 | } | |
323 | ||
324 | //________________________________________________________________________ | |
325 | void AliJetModelBaseTask::Run() | |
326 | { | |
327 | ||
328 | } | |
329 | ||
330 | //________________________________________________________________________ | |
331 | void AliJetModelBaseTask::UserExec(Option_t *) | |
332 | { | |
333 | // Execute per event. | |
334 | ||
335 | Init(); | |
336 | ||
337 | Run(); | |
338 | ||
339 | } |