]>
Commit | Line | Data |
---|---|---|
c7a4dac0 | 1 | /************************************************************************** |
20eef074 | 2 | * Copyright(c) 1998-1999, 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 | **************************************************************************/ | |
c7a4dac0 | 15 | |
5bfe44ce | 16 | /* |
85f5e9c2 | 17 | $Id$ |
5bfe44ce | 18 | */ |
f8d9a5b8 | 19 | |
5bfe44ce | 20 | #include <Riostream.h> |
f7a1cc68 | 21 | #include <TGeoGlobalMagField.h> |
b0f5e3fc | 22 | #include <TH1.h> |
9f033001 | 23 | #include <TString.h> |
e8189707 | 24 | #include "AliITS.h" |
e869281d | 25 | #include "AliITSdigitSPD.h" |
5bfe44ce | 26 | #include "AliITShit.h" |
9f033001 | 27 | #include "AliITSmodule.h" |
c7a4dac0 | 28 | #include "AliITSpList.h" |
5bfe44ce | 29 | #include "AliITSCalibrationSPD.h" |
30 | #include "AliITSsegmentationSPD.h" | |
b0f5e3fc | 31 | #include "AliITSsimulationSPD.h" |
f77f13c8 | 32 | #include "AliLog.h" |
5bfe44ce | 33 | #include "AliRun.h" |
a4005be7 | 34 | #include "AliMagF.h" |
bb7e41dd | 35 | #include "AliMathBase.h" |
9f033001 | 36 | |
bb6b39bf | 37 | //#define DEBUG |
e8189707 | 38 | |
b0f5e3fc | 39 | ClassImp(AliITSsimulationSPD) |
40 | //////////////////////////////////////////////////////////////////////// | |
5bfe44ce | 41 | // Version: 1 |
590d15ee | 42 | // Modified by D. Elia, G.E. Bruno, H. Tydesjo |
43 | // Fast diffusion code by Bjorn S. Nilsen | |
5bfe44ce | 44 | // March-April 2006 |
8e2f60b2 | 45 | // October 2007: GetCalibrationObjects() removed |
b0f5e3fc | 46 | // |
5bfe44ce | 47 | // Version: 0 |
48 | // Written by Boris Batyunya | |
49 | // December 20 1999 | |
21b825a4 | 50 | // |
5bfe44ce | 51 | // |
52 | // AliITSsimulationSPD is to do the simulation of SPDs. | |
53 | // | |
54 | //////////////////////////////////////////////////////////////////////// | |
55 | ||
c7a4dac0 | 56 | //______________________________________________________________________ |
5bfe44ce | 57 | AliITSsimulationSPD::AliITSsimulationSPD(): |
58 | AliITSsimulation(), | |
59 | fHis(0), | |
60 | fSPDname(), | |
a4005be7 | 61 | fCoupling(), |
62 | fLorentz(kFALSE), | |
20eef074 | 63 | fTanLorAng(0), |
64 | fStrobe(kTRUE), | |
65 | fStrobeLenght(4), | |
66 | fStrobePhase(-12.5e-9){ | |
67 | // Default constructor. | |
68 | // Inputs: | |
69 | // none. | |
70 | // Outputs: | |
71 | // none. | |
72 | // Return: | |
73 | // A default constructed AliITSsimulationSPD class. | |
74 | ||
75 | AliDebug(1,Form("Calling default constructor")); | |
5bfe44ce | 76 | // Init(); |
b0f5e3fc | 77 | } |
c7a4dac0 | 78 | //______________________________________________________________________ |
8ba39da9 | 79 | AliITSsimulationSPD::AliITSsimulationSPD(AliITSDetTypeSim *dettyp): |
80 | AliITSsimulation(dettyp), | |
5bfe44ce | 81 | fHis(0), |
82 | fSPDname(), | |
a4005be7 | 83 | fCoupling(), |
84 | fLorentz(kFALSE), | |
20eef074 | 85 | fTanLorAng(0), |
86 | fStrobe(kTRUE), | |
87 | fStrobeLenght(4), | |
88 | fStrobePhase(-12.5e-9){ | |
89 | // standard constructor | |
90 | // Inputs: | |
91 | // AliITSsegmentation *seg A pointer to the segmentation class | |
92 | // to be used for this simulation | |
93 | // AliITSCalibration *resp A pointer to the responce class to | |
94 | // be used for this simulation | |
95 | // Outputs: | |
96 | // none. | |
97 | // Return: | |
98 | // A default constructed AliITSsimulationSPD class. | |
99 | ||
100 | AliDebug(1,Form("Calling standard constructor ")); | |
101 | Init(); | |
aacedc3e | 102 | } |
103 | //______________________________________________________________________ | |
5bfe44ce | 104 | void AliITSsimulationSPD::Init(){ |
20eef074 | 105 | // Initilization |
106 | // Inputs: | |
107 | // none. | |
108 | // Outputs: | |
109 | // none. | |
110 | // Return: | |
111 | // none. | |
112 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
113 | ||
114 | SetModuleNumber(0); | |
115 | SetEventNumber(0); | |
116 | SetMap(new AliITSpList(GetNPixelsZ(),GetNPixelsX())); | |
117 | AliITSSimuParam* simpar = fDetType->GetSimuParam(); | |
118 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
119 | Double_t bias = simpar->GetSPDBiasVoltage(); | |
5bfe44ce | 120 | // cout << "Bias Voltage --> " << bias << endl; // dom |
20eef074 | 121 | simpar->SetDistanceOverVoltage(kmictocm*seg->Dy(),bias); |
5bfe44ce | 122 | // set kind of coupling ("old" or "new") |
20eef074 | 123 | char opt[20]; |
124 | simpar->GetSPDCouplingOption(opt); | |
125 | char *old = strstr(opt,"old"); | |
126 | if (old) { | |
127 | fCoupling=2; | |
128 | } else { | |
129 | fCoupling=1; | |
130 | } // end if | |
a0a6914c | 131 | SetLorentzDrift(simpar->GetSPDLorentzDrift()); |
132 | if (fLorentz) SetTanLorAngle(simpar->GetSPDLorentzHoleWeight()); | |
20eef074 | 133 | //SetStrobeGeneration(kFALSE); |
134 | if (fStrobe) GenerateStrobePhase(); | |
a4005be7 | 135 | } |
136 | //______________________________________________________________________ | |
137 | Bool_t AliITSsimulationSPD::SetTanLorAngle(Double_t WeightHole) { | |
20eef074 | 138 | // This function set the Tangent of the Lorentz angle. |
139 | // A weighted average is used for electrons and holes | |
140 | // Input: Double_t WeightHole: wheight for hole: it should be in the range [0,1] | |
141 | // output: Bool_t : kTRUE in case of success | |
142 | // | |
143 | if(!fDetType) { | |
144 | AliError("AliITSsimulationSPD::SetTanLorAngle: AliITSDetTypeSim* fDetType not set "); | |
145 | return kFALSE;} | |
146 | if(WeightHole<0) { | |
147 | WeightHole=0.; | |
148 | AliWarning("AliITSsimulationSPD::SetTanLorAngle: You have asked for negative Hole weight"); | |
149 | AliWarning("AliITSsimulationSPD::SetTanLorAngle: I'm going to use only electrons"); | |
150 | } | |
151 | if(WeightHole>1) { | |
152 | WeightHole=1.; | |
153 | AliWarning("AliITSsimulationSPD::SetTanLorAngle: You have asked for weight > 1"); | |
154 | AliWarning("AliITSsimulationSPD::SetTanLorAngle: I'm going to use only holes"); | |
155 | } | |
156 | Double_t WeightEle=1.-WeightHole; | |
157 | AliITSSimuParam* simpar = fDetType->GetSimuParam(); | |
a0a6914c | 158 | AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); |
159 | if (!fld) AliFatal("The field is not initialized"); | |
160 | Double_t bz = fld->SolenoidField(); | |
161 | fTanLorAng = TMath::Tan(WeightHole*simpar->LorentzAngleHole(bz) + | |
162 | WeightEle*simpar->LorentzAngleElectron(bz)); | |
20eef074 | 163 | return kTRUE; |
c7a4dac0 | 164 | } |
165 | //______________________________________________________________________ | |
5bfe44ce | 166 | AliITSsimulationSPD::~AliITSsimulationSPD(){ |
20eef074 | 167 | // destructor |
168 | // Inputs: | |
169 | // none. | |
170 | // Outputs: | |
171 | // none. | |
172 | // Return: | |
173 | // none. | |
174 | ||
175 | if (fHis) { | |
176 | fHis->Delete(); | |
177 | delete fHis; | |
178 | } // end if fHis | |
b0f5e3fc | 179 | } |
c7a4dac0 | 180 | //______________________________________________________________________ |
5bfe44ce | 181 | AliITSsimulationSPD::AliITSsimulationSPD(const |
182 | AliITSsimulationSPD | |
7537d03c | 183 | &s) : AliITSsimulation(s), |
184 | fHis(s.fHis), | |
185 | fSPDname(s.fSPDname), | |
a4005be7 | 186 | fCoupling(s.fCoupling), |
187 | fLorentz(s.fLorentz), | |
20eef074 | 188 | fTanLorAng(s.fTanLorAng), |
189 | fStrobe(s.fStrobe), | |
190 | fStrobeLenght(s.fStrobeLenght), | |
191 | fStrobePhase(s.fStrobePhase){ | |
192 | // Copy Constructor | |
193 | // Inputs: | |
194 | // AliITSsimulationSPD &s The original class for which | |
195 | // this class is a copy of | |
196 | // Outputs: | |
197 | // none. | |
198 | // Return: | |
c7a4dac0 | 199 | |
b0f5e3fc | 200 | } |
c7a4dac0 | 201 | //______________________________________________________________________ |
5bfe44ce | 202 | AliITSsimulationSPD& AliITSsimulationSPD::operator=(const |
20eef074 | 203 | AliITSsimulationSPD &s){ |
204 | // Assignment operator | |
205 | // Inputs: | |
206 | // AliITSsimulationSPD &s The original class for which | |
207 | // this class is a copy of | |
208 | // Outputs: | |
209 | // none. | |
210 | // Return: | |
211 | ||
212 | if(&s == this) return *this; | |
213 | this->fHis = s.fHis; | |
214 | fCoupling = s.fCoupling; | |
215 | fSPDname = s.fSPDname; | |
216 | fLorentz = s.fLorentz; | |
217 | fTanLorAng = s.fTanLorAng; | |
218 | fStrobe = s.fStrobe; | |
219 | fStrobeLenght = s.fStrobeLenght; | |
220 | fStrobePhase = s.fStrobePhase; | |
221 | return *this; | |
5bfe44ce | 222 | } |
85f5e9c2 | 223 | /* |
d2f55a22 | 224 | //______________________________________________________________________ |
5bfe44ce | 225 | AliITSsimulation& AliITSsimulationSPD::operator=(const |
20eef074 | 226 | AliITSsimulation &s){ |
227 | // Assignment operator | |
228 | // Inputs: | |
229 | // AliITSsimulationSPD &s The original class for which | |
230 | // this class is a copy of | |
231 | // Outputs: | |
232 | // none. | |
233 | // Return: | |
234 | ||
235 | if(&s == this) return *this; | |
236 | Error("AliITSsimulationSPD","Not allowed to make a = with " | |
237 | "AliITSsimulationSPD","Using default creater instead"); | |
238 | ||
239 | return *this; | |
5bfe44ce | 240 | } |
85f5e9c2 | 241 | */ |
5bfe44ce | 242 | //______________________________________________________________________ |
243 | void AliITSsimulationSPD::InitSimulationModule(Int_t module, Int_t event){ | |
20eef074 | 244 | // This function creates maps to build the list of tracks for each |
245 | // summable digit. Inputs defined by base class. | |
246 | // Inputs: | |
247 | // Int_t module // Module number to be simulated | |
248 | // Int_t event // Event number to be simulated | |
249 | // Outputs: | |
250 | // none | |
251 | // Returns: | |
252 | // none | |
253 | ||
254 | AliDebug(1,Form("(module=%d,event=%d)",module,event)); | |
255 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
256 | AliITSSimuParam* simpar = fDetType->GetSimuParam(); | |
257 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
258 | SetModuleNumber(module); | |
259 | SetEventNumber(event); | |
260 | simpar->SetDistanceOverVoltage(kmictocm*seg->Dy(),simpar->GetSPDBiasVoltage(module)); | |
261 | ClearMap(); | |
5bfe44ce | 262 | } |
263 | //_____________________________________________________________________ | |
264 | void AliITSsimulationSPD::SDigitiseModule(AliITSmodule *mod,Int_t, | |
20eef074 | 265 | Int_t event){ |
266 | // This function begins the work of creating S-Digits. Inputs defined | |
267 | // by base class. | |
268 | // Inputs: | |
269 | // AliITSmodule *mod // module | |
270 | // Int_t // not used | |
271 | // Int_t event // Event number | |
272 | // Outputs: | |
273 | // none | |
274 | // Return: | |
275 | // test // test returns kTRUE if the module contained hits | |
276 | // // test returns kFALSE if it did not contain hits | |
277 | ||
278 | AliDebug(1,Form("(mod=%p, ,event=%d)",mod,event)); | |
279 | if(!(mod->GetNhits())){ | |
280 | AliDebug(1,Form("In event %d module %d there are %d hits returning.", | |
5bfe44ce | 281 | event, mod->GetIndex(),mod->GetNhits())); |
20eef074 | 282 | return;// if module has no hits don't create Sdigits |
283 | } // end if | |
284 | SetModuleNumber(mod->GetIndex()); | |
285 | if (fStrobe) if(event != GetEventNumber()) GenerateStrobePhase(); | |
286 | SetEventNumber(event); | |
287 | InitSimulationModule( GetModuleNumber() , event ); | |
288 | // HitToSDigit(mod); | |
289 | HitToSDigitFast(mod); | |
ad7f2bfa | 290 | if (fDetType->GetSimuParam()->GetSPDAddNoisyFlag()) AddNoisyPixels(); |
291 | if (fDetType->GetSimuParam()->GetSPDRemoveDeadFlag()) RemoveDeadPixels(); | |
292 | ||
590d15ee | 293 | // cout << "After Remove in SDigitiseModule !!!!!" << endl; // dom |
294 | // cout << "Module " << mod->GetIndex() << " Event " << event << endl; // dom | |
20eef074 | 295 | WriteSDigits(); |
296 | ClearMap(); | |
5bfe44ce | 297 | } |
298 | //______________________________________________________________________ | |
299 | void AliITSsimulationSPD::WriteSDigits(){ | |
20eef074 | 300 | // This function adds each S-Digit to pList |
301 | // Inputs: | |
302 | // none. | |
303 | // Outputs: | |
304 | // none. | |
305 | // Return: | |
306 | // none | |
307 | Int_t ix, nix, iz, niz; | |
308 | static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS"); | |
309 | ||
310 | AliDebug(1,Form("Writing SDigits for module %d",GetModuleNumber())); | |
590d15ee | 311 | // cout << "WriteSDigits for module " << GetModuleNumber() << endl; // dom |
20eef074 | 312 | GetMap()->GetMaxMapIndex(niz, nix); |
313 | for(iz=0; iz<niz; iz++)for(ix=0; ix<nix; ix++){ | |
314 | if(GetMap()->GetSignalOnly(iz,ix)>0.0){ | |
590d15ee | 315 | // cout << " Signal gt 0 iz ix " << iz << ix << " Module " << GetModuleNumber() << endl; // dom |
20eef074 | 316 | aliITS->AddSumDigit(*(GetMap()->GetpListItem(iz,ix))); |
5bfe44ce | 317 | if(AliDebugLevel()>0) { |
318 | AliDebug(1,Form("%d, %d",iz,ix)); | |
319 | cout << *(GetMap()->GetpListItem(iz,ix)) << endl; | |
20eef074 | 320 | } // end if GetDebug |
321 | } // end if GetMap()->GetSignalOnly(iz,ix)>0.0 | |
322 | } // end for iz,ix | |
323 | return; | |
3a97c582 | 324 | } |
325 | //______________________________________________________________________ | |
326 | void AliITSsimulationSPD::FinishSDigitiseModule(){ | |
20eef074 | 327 | // This function calls SDigitsToDigits which creates Digits from SDigits |
328 | // Inputs: | |
329 | // none | |
330 | // Outputs: | |
331 | // none | |
332 | // Return | |
333 | // none | |
334 | ||
335 | AliDebug(1,"()"); | |
590d15ee | 336 | // cout << "FinishSDigitiseModule for module " << GetModuleNumber() << endl; // dom |
20eef074 | 337 | FrompListToDigits(); // Charge To Signal both adds noise and |
338 | ClearMap(); | |
339 | return; | |
5bfe44ce | 340 | } |
341 | //______________________________________________________________________ | |
342 | void AliITSsimulationSPD::DigitiseModule(AliITSmodule *mod,Int_t, | |
20eef074 | 343 | Int_t event){ |
344 | // This function creates Digits straight from the hits and then adds | |
345 | // electronic noise to the digits before adding them to pList | |
346 | // Each of the input variables is passed along to HitToSDigit | |
347 | // Inputs: | |
348 | // AliITSmodule *mod module | |
349 | // Int_t Dummy. | |
350 | // Int_t Dummy | |
351 | // Outputs: | |
352 | // none. | |
353 | // Return: | |
354 | // none. | |
355 | ||
356 | if (fStrobe) if(event != GetEventNumber()) GenerateStrobePhase(); | |
1a429a6b | 357 | AliDebug(1,Form("(mod=%p,,0)",mod)); |
20eef074 | 358 | // HitToSDigit(mod); |
359 | InitSimulationModule( mod->GetIndex(), event ); | |
360 | HitToSDigitFast(mod); | |
ad7f2bfa | 361 | |
362 | if (fDetType->GetSimuParam()->GetSPDAddNoisyFlag()) AddNoisyPixels(); | |
363 | if (fDetType->GetSimuParam()->GetSPDRemoveDeadFlag()) RemoveDeadPixels(); | |
590d15ee | 364 | // cout << "After Remove in DigitiseModule in module " << mod->GetIndex() << endl; // dom |
20eef074 | 365 | FrompListToDigits(); |
366 | ClearMap(); | |
c7a4dac0 | 367 | } |
368 | //______________________________________________________________________ | |
5bfe44ce | 369 | void AliITSsimulationSPD::HitToSDigit(AliITSmodule *mod){ |
20eef074 | 370 | // Does the charge distributions using Gaussian diffusion charge charing. |
371 | // Inputs: | |
372 | // AliITSmodule *mod Pointer to this module | |
373 | // Output: | |
374 | // none. | |
375 | // Return: | |
376 | // none. | |
377 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
378 | const Double_t kBunchLenght = 25e-9; // LHC clock | |
379 | TObjArray *hits = mod->GetHits(); | |
380 | Int_t nhits = hits->GetEntriesFast(); | |
381 | Int_t h,ix,iz,i; | |
382 | Int_t idtrack; | |
383 | Double_t x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0,de=0.0,ld=0.0; | |
384 | Double_t x,y,z,t,tp,st,dt=0.2,el,sig,sigx,sigz,fda; | |
385 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
386 | AliITSSimuParam *simpar = fDetType->GetSimuParam(); | |
387 | Double_t thick = 0.5*kmictocm*seg->Dy(); // Half Thickness | |
388 | simpar->GetSPDSigmaDiffusionAsymmetry(fda); // | |
389 | ||
390 | AliDebug(1,Form("(mod=%p) fCoupling=%d",mod,fCoupling)); | |
391 | if(nhits<=0) return; | |
392 | for(h=0;h<nhits;h++){ | |
393 | if(AliDebugLevel()>0) { | |
5bfe44ce | 394 | AliDebug(1,Form("Hits, %d", h)); |
395 | cout << *(mod->GetHit(h)) << endl; | |
20eef074 | 396 | } // end if GetDebug |
397 | // Check if the hit is inside readout window | |
398 | if (fStrobe) | |
399 | if ((mod->GetHit(h)->GetTOF() < fStrobePhase) || | |
400 | (mod->GetHit(h)->GetTOF() > (fStrobePhase+(Double_t)fStrobeLenght*kBunchLenght))) continue; | |
401 | if(!mod->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,idtrack)) continue; | |
402 | st = TMath::Sqrt(x1*x1+y1*y1+z1*z1); | |
403 | if(st>0.0){ | |
404 | st = (Double_t)((Int_t)(st/kmictocm)); // number of microns | |
405 | if(st<=1.0) st = 1.0; | |
406 | dt = 1.0/st; | |
407 | for(t=0.0;t<1.0;t+=dt){ // Integrate over t | |
408 | tp = t+0.5*dt; | |
409 | x = x0+x1*tp; | |
410 | y = y0+y1*tp; | |
411 | z = z0+z1*tp; | |
412 | if(!(seg->LocalToDet(x,z,ix,iz))) continue; // outside | |
413 | //el = res->GeVToCharge((Double_t)(dt*de)); | |
414 | el = dt * de / simpar->GetGeVToCharge(); | |
415 | if(GetDebug(1)){ | |
416 | if(el<=0.0) cout<<"el="<<el<<" dt="<<dt | |
417 | <<" de="<<de<<endl; | |
418 | } // end if GetDebug | |
419 | sig = simpar->SigmaDiffusion1D(TMath::Abs(thick + y)); | |
420 | sigx=sig; | |
421 | sigz=sig*fda; | |
422 | if (fLorentz) ld=(y+thick)*fTanLorAng; | |
423 | SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h); | |
424 | } // end for t | |
425 | } else { // st == 0.0 deposit it at this point | |
426 | x = x0; | |
427 | y = y0; | |
428 | z = z0; | |
429 | if(!(seg->LocalToDet(x,z,ix,iz))) continue; // outside | |
430 | //el = res->GeVToCharge((Double_t)de); | |
431 | el = de / simpar->GetGeVToCharge(); | |
432 | sig = simpar->SigmaDiffusion1D(TMath::Abs(thick + y)); | |
433 | sigx=sig; | |
434 | sigz=sig*fda; | |
435 | if (fLorentz) ld=(y+thick)*fTanLorAng; | |
436 | SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h); | |
437 | } // end if st>0.0 | |
438 | ||
439 | } // Loop over all hits h | |
440 | ||
441 | // Coupling | |
442 | switch (fCoupling) { | |
443 | default: | |
444 | break; | |
445 | case 1: //case 3: | |
446 | for(i=0;i<GetMap()->GetEntries();i++) | |
447 | if(GetMap()->GetpListItem(i)==0) continue; | |
448 | else{ | |
449 | GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); | |
450 | SetCoupling(iz,ix,idtrack,h); | |
451 | } // end for i | |
452 | break; | |
453 | case 2: // case 4: | |
454 | for(i=0;i<GetMap()->GetEntries();i++) | |
455 | if(GetMap()->GetpListItem(i)==0) continue; | |
456 | else{ | |
457 | GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); | |
458 | SetCouplingOld(iz,ix,idtrack,h); | |
459 | } // end for i | |
460 | break; | |
461 | } // end switch | |
462 | if(GetDebug(2))Info("HitToSDigit","Finished fCoupling=%d",fCoupling); | |
c7a4dac0 | 463 | } |
464 | //______________________________________________________________________ | |
5bfe44ce | 465 | void AliITSsimulationSPD::HitToSDigitFast(AliITSmodule *mod){ |
20eef074 | 466 | // Does the charge distributions using Gaussian diffusion charge charing. // Inputs: |
467 | // AliITSmodule *mod Pointer to this module | |
468 | // Output: | |
469 | // none. | |
470 | // Return: | |
471 | // none. | |
472 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
473 | const Int_t kn10=10; | |
474 | const Double_t kti[kn10]={7.443716945e-3,2.166976971e-1,3.397047841e-1, | |
475 | 4.325316833e-1,4.869532643e-1,5.130467358e-1, | |
476 | 5.674683167e-1,6.602952159e-1,7.833023029e-1, | |
477 | 9.255628306e-1}; | |
478 | const Double_t kwi[kn10]={1.477621124e-1,1.346333597e-1,1.095431813e-1, | |
479 | 7.472567455e-2,3.333567215e-2,3.333567215e-2, | |
480 | 7.472567455e-2,1.095431813e-1,1.346333597e-1, | |
481 | 1.477621124e-1}; | |
482 | const Double_t kBunchLenght = 25e-9; // LHC clock | |
483 | TObjArray *hits = mod->GetHits(); | |
484 | Int_t nhits = hits->GetEntriesFast(); | |
485 | Int_t h,ix,iz,i; | |
486 | Int_t idtrack; | |
487 | Double_t x0=0.0,x1=0.0,y0=0.0,y1=0.0,z0=0.0,z1=0.0,de=0.0,ld=0.0; | |
488 | Double_t x,y,z,t,st,el,sig,sigx,sigz,fda; | |
489 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
490 | AliITSSimuParam* simpar = fDetType->GetSimuParam(); | |
491 | Double_t thick = 0.5*kmictocm*seg->Dy(); // Half thickness | |
492 | simpar->GetSPDSigmaDiffusionAsymmetry(fda); | |
5bfe44ce | 493 | // cout << "Half Thickness " << thick << endl; // dom |
494 | // cout << "Diffusion asymm " << fda << endl; // dom | |
495 | ||
20eef074 | 496 | AliDebug(1,Form("(mod=%p) fCoupling=%d",mod,fCoupling)); |
497 | if(nhits<=0) return; | |
498 | for(h=0;h<nhits;h++){ | |
499 | if(AliDebugLevel()>0) { | |
500 | AliDebug(1,Form("Hits, %d", h)); | |
501 | cout << *(mod->GetHit(h)) << endl; | |
502 | } // end if GetDebug | |
503 | // Check if the hit is inside readout window | |
504 | if (fStrobe) | |
505 | if ((mod->GetHit(h)->GetTOF() < fStrobePhase) || | |
506 | (mod->GetHit(h)->GetTOF() > (fStrobePhase+(Double_t)fStrobeLenght*kBunchLenght))) continue; | |
507 | if(!mod->LineSegmentL(h,x0,x1,y0,y1,z0,z1,de,idtrack)) continue; | |
508 | st = TMath::Sqrt(x1*x1+y1*y1+z1*z1); | |
509 | if(st>0.0) for(i=0;i<kn10;i++){ // Integrate over t | |
510 | t = kti[i]; | |
511 | x = x0+x1*t; | |
512 | y = y0+y1*t; | |
513 | z = z0+z1*t; | |
514 | if(!(seg->LocalToDet(x,z,ix,iz))) continue; // outside | |
515 | el = kwi[i]*de/simpar->GetGeVToCharge(); | |
516 | if(GetDebug(1)){ | |
517 | if(el<=0.0) cout<<"el="<<el<<" kwi["<<i<<"]="<<kwi[i] | |
518 | <<" de="<<de<<endl; | |
519 | } // end if GetDebug | |
520 | sig = simpar->SigmaDiffusion1D(TMath::Abs(thick + y)); | |
521 | sigx=sig; | |
522 | sigz=sig*fda; | |
523 | if (fLorentz) ld=(y+thick)*fTanLorAng; | |
524 | SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h); | |
8e2f60b2 | 525 | // cout << "sigx sigz " << sigx << " " << sigz << endl; // dom |
20eef074 | 526 | } // end for i // End Integrate over t |
527 | else { // st == 0.0 deposit it at this point | |
528 | x = x0; | |
529 | y = y0; | |
530 | z = z0; | |
531 | if(!(seg->LocalToDet(x,z,ix,iz))) continue; // outside | |
532 | el = de / simpar->GetGeVToCharge(); | |
533 | sig = simpar->SigmaDiffusion1D(TMath::Abs(thick + y)); | |
534 | sigx=sig; | |
535 | sigz=sig*fda; | |
536 | if (fLorentz) ld=(y+thick)*fTanLorAng; | |
537 | SpreadChargeAsym(x,z,ix,iz,el,sigx,sigz,ld,idtrack,h); | |
538 | } // end if st>0.0 | |
539 | ||
540 | } // Loop over all hits h | |
541 | ||
542 | // Coupling | |
543 | switch (fCoupling) { | |
544 | default: | |
545 | break; | |
546 | case 1: // case 3: | |
547 | for(i=0;i<GetMap()->GetEntries();i++) | |
548 | if(GetMap()->GetpListItem(i)==0) continue; | |
549 | else{ | |
550 | GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); | |
551 | SetCoupling(iz,ix,idtrack,h); | |
552 | } // end for i | |
553 | break; | |
554 | case 2: // case 4: | |
555 | for(i=0;i<GetMap()->GetEntries();i++) | |
556 | if(GetMap()->GetpListItem(i)==0) continue; | |
557 | else{ | |
558 | GetMap()->GetMapIndex(GetMap()->GetpListItem(i)->GetIndex(),iz,ix); | |
559 | SetCouplingOld(iz,ix,idtrack,h); | |
560 | } // end for i | |
561 | break; | |
562 | } // end switch | |
563 | if(GetDebug(2))Info("HitToSDigit","Finished fCoupling=%d",fCoupling); | |
21b825a4 | 564 | } |
c7a4dac0 | 565 | //______________________________________________________________________ |
5bfe44ce | 566 | void AliITSsimulationSPD::SpreadCharge(Double_t x0,Double_t z0, |
20eef074 | 567 | Int_t ix0,Int_t iz0, |
a4005be7 | 568 | Double_t el,Double_t sig,Double_t ld, |
569 | Int_t t,Int_t hi){ | |
20eef074 | 570 | // Spreads the charge over neighboring cells. Assume charge is distributed |
571 | // as charge(x,z) = (el/2*pi*sig*sig)*exp(-arg) | |
572 | // arg=((x-x0)*(x-x0)/2*sig*sig)+((z-z0*z-z0)/2*sig*sig) | |
573 | // if fLorentz=kTRUE, then x0=x0+ld (Lorentz drift taken into account) | |
574 | // Defined this way, the integral over all x and z is el. | |
575 | // Inputs: | |
576 | // Double_t x0 x position of point where charge is liberated | |
577 | // Double_t z0 z position of point where charge is liberated | |
578 | // Int_t ix0 row of cell corresponding to point x0 | |
579 | // Int_t iz0 columb of cell corresponding to point z0 | |
580 | // Double_t el number of electrons liberated in this step | |
581 | // Double_t sig Sigma difusion for this step (y0 dependent) | |
582 | // Double_t ld lorentz drift in x for this step (y0 dependent) | |
583 | // Int_t t track number | |
584 | // Int_t ti hit track index number | |
585 | // Int_t hi hit "hit" index number | |
586 | // Outputs: | |
587 | // none. | |
588 | // Return: | |
589 | // none. | |
590 | const Int_t knx = 3,knz = 2; | |
591 | const Double_t kRoot2 = 1.414213562; // Sqrt(2). | |
592 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
593 | Int_t ix,iz,ixs,ixe,izs,ize; | |
594 | Float_t x,z; | |
595 | Double_t x1,x2,z1,z2,s,sp; | |
596 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
597 | ||
598 | ||
599 | if(GetDebug(4)) Info("SpreadCharge","(x0=%e,z0=%e,ix0=%d,iz0=%d,el=%e," | |
600 | "sig=%e,t=%d,i=%d)",x0,z0,ix0,iz0,el,sig,t,hi); | |
601 | if(sig<=0.0) { // if sig<=0 No diffusion to simulate. | |
602 | GetMap()->AddSignal(iz0,ix0,t,hi,GetModuleNumber(),el); | |
603 | if(GetDebug(2)){ | |
604 | cout << "sig<=0.0=" << sig << endl; | |
605 | } // end if GetDebug | |
606 | return; | |
607 | } // end if | |
608 | sp = 1.0/(sig*kRoot2); | |
609 | if(GetDebug(2)){ | |
610 | cout << "sig=" << sig << " sp=" << sp << endl; | |
611 | } // end if GetDebug | |
612 | ixs = TMath::Max(-knx+ix0,0); | |
613 | ixe = TMath::Min(knx+ix0,seg->Npx()-1); | |
614 | izs = TMath::Max(-knz+iz0,0); | |
615 | ize = TMath::Min(knz+iz0,seg->Npz()-1); | |
616 | for(ix=ixs;ix<=ixe;ix++) for(iz=izs;iz<=ize;iz++){ | |
617 | seg->DetToLocal(ix,iz,x,z); // pixel center | |
618 | x1 = x; | |
619 | z1 = z; | |
620 | x2 = x1 + 0.5*kmictocm*seg->Dpx(ix); // Upper | |
621 | x1 -= 0.5*kmictocm*seg->Dpx(ix); // Lower | |
622 | z2 = z1 + 0.5*kmictocm*seg->Dpz(iz); // Upper | |
623 | z1 -= 0.5*kmictocm*seg->Dpz(iz); // Lower | |
624 | x1 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) | |
625 | x2 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) | |
626 | z1 -= z0; // Distance from where track traveled | |
627 | z2 -= z0; // Distance from where track traveled | |
628 | s = 0.25; // Correction based on definision of Erfc | |
bb7e41dd | 629 | s *= AliMathBase::ErfcFast(sp*x1) - AliMathBase::ErfcFast(sp*x2); |
20eef074 | 630 | if(GetDebug(3)){ |
631 | cout <<"el="<<el<<" ix0="<<ix0<<" ix="<<ix<<" x0="<<x<< | |
632 | " iz0="<<iz0<<" iz="<<iz<<" z0="<<z<< | |
633 | " sp*x1="<<sp*x1<<" sp*x2="<<sp*x2<<" s="<<s; | |
634 | } // end if GetDebug | |
bb7e41dd | 635 | s *= AliMathBase::ErfcFast(sp*z1) - AliMathBase::ErfcFast(sp*z2); |
20eef074 | 636 | if(GetDebug(3)){ |
637 | cout<<" sp*z1="<<sp*z1<<" sp*z2="<<sp*z2<<" s="<<s<< endl; | |
638 | } // end if GetDebug | |
639 | GetMap()->AddSignal(iz,ix,t,hi,GetModuleNumber(),s*el); | |
640 | } // end for ix, iz | |
c7a4dac0 | 641 | } |
642 | //______________________________________________________________________ | |
5bfe44ce | 643 | void AliITSsimulationSPD::SpreadChargeAsym(Double_t x0,Double_t z0, |
20eef074 | 644 | Int_t ix0,Int_t iz0, |
645 | Double_t el,Double_t sigx,Double_t sigz, | |
646 | Double_t ld,Int_t t,Int_t hi){ | |
647 | // Spreads the charge over neighboring cells. Assume charge is distributed | |
648 | // as charge(x,z) = (el/2*pi*sigx*sigz)*exp(-arg) | |
649 | // arg=((x-x0)*(x-x0)/2*sigx*sigx)+((z-z0*z-z0)/2*sigz*sigz) | |
650 | // if fLorentz=kTRUE, then x0=x0+ld (Lorentz drift taken into account) | |
651 | // Defined this way, the integral over all x and z is el. | |
652 | // Inputs: | |
653 | // Double_t x0 x position of point where charge is liberated | |
654 | // Double_t z0 z position of point where charge is liberated | |
655 | // Int_t ix0 row of cell corresponding to point x0 | |
656 | // Int_t iz0 columb of cell corresponding to point z0 | |
657 | // Double_t el number of electrons liberated in this step | |
658 | // Double_t sigx Sigma difusion along x for this step (y0 dependent) | |
659 | // Double_t sigz Sigma difusion along z for this step (y0 dependent) | |
660 | // Double_t ld lorentz drift in x for this stip (y0 dependent) | |
661 | // Int_t t track number | |
662 | // Int_t ti hit track index number | |
663 | // Int_t hi hit "hit" index number | |
664 | // Outputs: | |
665 | // none. | |
666 | // Return: | |
667 | // none. | |
668 | const Int_t knx = 3,knz = 2; | |
669 | const Double_t kRoot2 = 1.414213562; // Sqrt(2). | |
670 | const Double_t kmictocm = 1.0e-4; // convert microns to cm. | |
671 | Int_t ix,iz,ixs,ixe,izs,ize; | |
672 | Float_t x,z; | |
673 | Double_t x1,x2,z1,z2,s,spx,spz; | |
674 | AliITSsegmentationSPD* seg = (AliITSsegmentationSPD*)GetSegmentationModel(0); | |
675 | ||
676 | ||
677 | if(GetDebug(4)) Info("SpreadChargeAsym","(x0=%e,z0=%e,ix0=%d,iz0=%d,el=%e," | |
a5a317a9 | 678 | "sigx=%e, sigz=%e, t=%d,i=%d)",x0,z0,ix0,iz0,el,sigx,sigz,t,hi); |
20eef074 | 679 | if(sigx<=0.0 || sigz<=0.0) { // if sig<=0 No diffusion to simulate. |
680 | GetMap()->AddSignal(iz0,ix0,t,hi,GetModuleNumber(),el); | |
681 | if(GetDebug(2)){ | |
682 | cout << "sigx<=0.0=" << sigx << endl; | |
683 | cout << "sigz<=0.0=" << sigz << endl; | |
684 | } // end if GetDebug | |
685 | return; | |
686 | } // end if | |
687 | spx = 1.0/(sigx*kRoot2); spz = 1.0/(sigz*kRoot2); | |
688 | if(GetDebug(2)){ | |
689 | cout << "sigx=" << sigx << " spx=" << spx << endl; | |
690 | cout << "sigz=" << sigz << " spz=" << spz << endl; | |
691 | } // end if GetDebug | |
692 | ixs = TMath::Max(-knx+ix0,0); | |
693 | ixe = TMath::Min(knx+ix0,seg->Npx()-1); | |
694 | izs = TMath::Max(-knz+iz0,0); | |
695 | ize = TMath::Min(knz+iz0,seg->Npz()-1); | |
696 | for(ix=ixs;ix<=ixe;ix++) for(iz=izs;iz<=ize;iz++){ | |
697 | seg->DetToLocal(ix,iz,x,z); // pixel center | |
698 | x1 = x; | |
699 | z1 = z; | |
700 | x2 = x1 + 0.5*kmictocm*seg->Dpx(ix); // Upper | |
701 | x1 -= 0.5*kmictocm*seg->Dpx(ix); // Lower | |
702 | z2 = z1 + 0.5*kmictocm*seg->Dpz(iz); // Upper | |
703 | z1 -= 0.5*kmictocm*seg->Dpz(iz); // Lower | |
704 | x1 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) | |
705 | x2 -= x0+ld; // Distance from where track traveled (taking into account the Lorentz drift) | |
706 | z1 -= z0; // Distance from where track traveled | |
707 | z2 -= z0; // Distance from where track traveled | |
708 | s = 0.25; // Correction based on definision of Erfc | |
bb7e41dd | 709 | s *= AliMathBase::ErfcFast(spx*x1) - AliMathBase::ErfcFast(spx*x2); |
20eef074 | 710 | if(GetDebug(3)){ |
711 | cout <<"el="<<el<<" ix0="<<ix0<<" ix="<<ix<<" x0="<<x<< | |
712 | " iz0="<<iz0<<" iz="<<iz<<" z0="<<z<< | |
713 | " spx*x1="<<spx*x1<<" spx*x2="<<spx*x2<<" s="<<s; | |
714 | } // end if GetDebug | |
bb7e41dd | 715 | s *= AliMathBase::ErfcFast(spz*z1) - AliMathBase::ErfcFast(spz*z2); |
20eef074 | 716 | if(GetDebug(3)){ |
717 | cout<<" spz*z1="<<spz*z1<<" spz*z2="<<spz*z2<<" s="<<s<< endl; | |
718 | } // end if GetDebug | |
719 | GetMap()->AddSignal(iz,ix,t,hi,GetModuleNumber(),s*el); | |
720 | } // end for ix, iz | |
c7a4dac0 | 721 | } |
722 | //______________________________________________________________________ | |
ad7f2bfa | 723 | void AliITSsimulationSPD::RemoveDeadPixels(){ |
724 | // Removes dead pixels on each module (ladder) | |
725 | // This should be called before going from sdigits to digits (FrompListToDigits) | |
726 | Int_t mod = GetModuleNumber(); | |
727 | AliITSCalibrationSPD* calObj = (AliITSCalibrationSPD*) fDetType->GetCalibrationModel(mod); | |
728 | ||
729 | Int_t nrDead = calObj->GetNrBad(); | |
730 | for (Int_t i=0; i<nrDead; i++) { | |
731 | GetMap()->DeleteHit(calObj->GetBadColAt(i), calObj->GetBadRowAt(i)); | |
732 | } | |
733 | } | |
734 | //______________________________________________________________________ | |
735 | void AliITSsimulationSPD::AddNoisyPixels() { | |
736 | // Adds noisy pixels on each module (ladder) | |
737 | // This should be called before going from sdigits to digits (FrompListToDigits) | |
738 | Int_t mod = GetModuleNumber(); | |
739 | AliITSCalibrationSPD* calObj = (AliITSCalibrationSPD*) fDetType->GetSPDNoisyModel(mod); | |
740 | ||
741 | Int_t nrNoisy = calObj->GetNrBad(); | |
742 | for (Int_t i=0; i<nrNoisy; i++) { | |
743 | // adding 10 times the threshold will for sure make this pixel fire... | |
744 | GetMap()->AddNoise(calObj->GetBadColAt(i), calObj->GetBadRowAt(i), mod, 10*GetThreshold()); | |
745 | } | |
5bfe44ce | 746 | } |
747 | //______________________________________________________________________ | |
590d15ee | 748 | void AliITSsimulationSPD::FrompListToDigits(){ |
20eef074 | 749 | // add noise and electronics, perform the zero suppression and add the |
750 | // digit to the list | |
751 | // Inputs: | |
752 | // none. | |
753 | // Outputs: | |
754 | // none. | |
755 | // Return: | |
756 | // none. | |
757 | static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS"); | |
758 | Int_t j,ix,iz; | |
759 | Double_t electronics; | |
760 | Double_t sig; | |
761 | const Int_t knmaxtrk=AliITSdigit::GetNTracks(); | |
762 | static AliITSdigitSPD dig; | |
763 | AliITSSimuParam *simpar = fDetType->GetSimuParam(); | |
764 | if(GetDebug(1)) Info("FrompListToDigits","()"); | |
765 | for(iz=0; iz<GetNPixelsZ(); iz++) for(ix=0; ix<GetNPixelsX(); ix++){ | |
5bfe44ce | 766 | // NEW (for the moment plugged by hand, in the future possibly read from Data Base) |
767 | // here parametrize the efficiency of the pixel along the row for the test columns (1,9,17,25) | |
768 | // if(iz==1 || iz == 9 || iz == 17 || iz == 25) { | |
769 | // Double_t eff,p1=0.,p2=0.; | |
770 | // Double_t x=ix; | |
771 | // switch (iz) { | |
772 | // case 1: p1=0.63460;p2=0.42438E-01;break; | |
773 | // case 9: p1=0.41090;p2=0.75914E-01;break; | |
774 | // case 17: p1=0.31883;p2=0.91502E-01;break; | |
775 | // case 25: p1=0.48828;p2=0.57975E-01;break; | |
776 | // } // end switch | |
777 | // eff=1.-p1*exp(-p2*x); | |
778 | // if (gRandom->Rndm() >= eff) continue; | |
779 | // } // end if | |
2ae37d58 | 780 | // END parametrize the efficiency |
781 | // | |
20eef074 | 782 | electronics = simpar->ApplySPDBaselineAndNoise(); |
783 | UpdateMapNoise(ix,iz,electronics); | |
784 | // | |
785 | // Apply Threshold and write Digits. | |
786 | sig = GetMap()->GetSignalOnly(iz,ix); | |
787 | FillHistograms(ix,iz,sig+electronics); | |
788 | if(GetDebug(3)){ | |
789 | cout<<sig<<"+"<<electronics<<">threshold("<<ix<<","<<iz | |
790 | <<")="<<GetThreshold() <<endl; | |
791 | } // end if GetDebug | |
ad7f2bfa | 792 | // if (sig+electronics <= GetThreshold()) continue; |
793 | if (GetMap()->GetSignal(iz,ix) <= GetThreshold()) continue; | |
20eef074 | 794 | dig.SetCoord1(iz); |
795 | dig.SetCoord2(ix); | |
796 | dig.SetSignal(1); | |
590d15ee | 797 | |
798 | // dig.SetSignalSPD((Int_t) GetMap()->GetSignal(iz,ix)); | |
20eef074 | 799 | Double_t aSignal = GetMap()->GetSignal(iz,ix); |
800 | if (TMath::Abs(aSignal)>2147483647.0) { | |
801 | //PH 2147483647 is the max. integer | |
802 | //PH This apparently is a problem which needs investigation | |
803 | AliWarning(Form("Too big or too small signal value %f",aSignal)); | |
804 | aSignal = TMath::Sign((Double_t)2147483647,aSignal); | |
805 | } | |
806 | dig.SetSignalSPD((Int_t)aSignal); | |
807 | ||
808 | for(j=0;j<knmaxtrk;j++){ | |
809 | if (j<GetMap()->GetNEntries()) { | |
810 | dig.SetTrack(j,GetMap()->GetTrack(iz,ix,j)); | |
811 | dig.SetHit(j,GetMap()->GetHit(iz,ix,j)); | |
812 | }else { // Default values | |
813 | dig.SetTrack(j,-3); | |
814 | dig.SetHit(j,-1); | |
815 | } // end if GetMap() | |
816 | } // end for j | |
817 | if(GetDebug(3)){ | |
818 | cout<<iz<<","<<ix<<","<<*(GetMap()->GetpListItem(iz,ix))<<endl; | |
819 | } // end if GetDebug | |
820 | aliITS->AddSimDigit(0,&dig); | |
ad7f2bfa | 821 | // simulate fo signal response for this pixel hit: |
822 | fDetType->ProcessSPDDigitForFastOr(fModule, dig.GetCoord1(), dig.GetCoord2()); | |
20eef074 | 823 | } // for ix/iz |
c7a4dac0 | 824 | } |
825 | //______________________________________________________________________ | |
5bfe44ce | 826 | void AliITSsimulationSPD::CreateHistograms(){ |
20eef074 | 827 | // create 1D histograms for tests |
828 | // Inputs: | |
829 | // none. | |
830 | // Outputs: | |
831 | // none. | |
832 | // Return: | |
833 | // none. | |
834 | ||
835 | if(GetDebug(1)) Info("CreateHistograms","create histograms"); | |
836 | ||
837 | fHis = new TObjArray(GetNPixelsZ()); | |
838 | fSPDname="spd_"; | |
839 | for(Int_t i=0;i<GetNPixelsZ();i++) { | |
840 | Char_t pixelz[4]; | |
ce189d7c | 841 | snprintf(pixelz,3,"%d",i); |
20eef074 | 842 | fSPDname.Append(pixelz); |
843 | fHis->AddAt(new TH1F(fSPDname.Data(),"SPD maps", | |
844 | GetNPixelsX(),0.,(Double_t)GetNPixelsX()),i); | |
845 | } // end for i | |
c7a4dac0 | 846 | } |
847 | //______________________________________________________________________ | |
5bfe44ce | 848 | void AliITSsimulationSPD::FillHistograms(Int_t ix,Int_t iz,Double_t v){ |
20eef074 | 849 | // Fill the histogram |
850 | // Inputs: | |
851 | // none. | |
852 | // Outputs: | |
853 | // none. | |
854 | // Return: | |
855 | // none. | |
856 | ||
857 | if(!GetHistArray()) return; // Only fill if setup. | |
858 | if(GetDebug(2)) Info("FillHistograms","fill histograms"); | |
859 | GetHistogram(iz)->Fill(ix,v); | |
c7a4dac0 | 860 | } |
861 | //______________________________________________________________________ | |
5bfe44ce | 862 | void AliITSsimulationSPD::ResetHistograms(){ |
20eef074 | 863 | // Reset histograms for this detector |
864 | // Inputs: | |
865 | // none. | |
866 | // Outputs: | |
867 | // none. | |
868 | // Return: | |
869 | // none. | |
870 | ||
871 | if(!GetHistArray()) return; // Only fill if setup. | |
872 | if(GetDebug(2)) Info("FillHistograms","fill histograms"); | |
873 | for ( int i=0;i<GetNPixelsZ();i++ ) { | |
874 | if (fHis->At(i)) ((TH1F*)fHis->At(i))->Reset(); | |
875 | } // end for i | |
21b825a4 | 876 | } |
5bfe44ce | 877 | |
c7a4dac0 | 878 | //______________________________________________________________________ |
590d15ee | 879 | void AliITSsimulationSPD::SetCoupling(Int_t col, Int_t row, Int_t ntrack, |
5bfe44ce | 880 | Int_t idhit) { |
20eef074 | 881 | // Take into account the coupling between adiacent pixels. |
882 | // The parameters probcol and probrow are the probability of the | |
883 | // signal in one pixel shared in the two adjacent pixels along | |
884 | // the column and row direction, respectively. | |
885 | // Note pList is goten via GetMap() and module is not need any more. | |
886 | // Otherwise it is identical to that coded by Tiziano Virgili (BSN). | |
887 | //Begin_Html | |
888 | /* | |
889 | <img src="picts/ITS/barimodel_3.gif"> | |
890 | </pre> | |
891 | <br clear=left> | |
892 | <font size=+2 color=red> | |
893 | <a href="mailto:tiziano.virgili@cern.ch"></a>. | |
894 | </font> | |
895 | <pre> | |
896 | */ | |
897 | //End_Html | |
898 | // Inputs: | |
899 | // Int_t col z cell index | |
900 | // Int_t row x cell index | |
901 | // Int_t ntrack track incex number | |
902 | // Int_t idhit hit index number | |
903 | // Outputs: | |
904 | // none. | |
905 | // Return: | |
906 | // none. | |
907 | Int_t j1,j2,flag=0; | |
908 | Double_t pulse1,pulse2; | |
909 | Double_t couplR=0.0,couplC=0.0; | |
910 | Double_t xr=0.; | |
911 | ||
912 | GetCouplings(couplC,couplR); | |
913 | if(GetDebug(3)) Info("SetCoupling","(col=%d,row=%d,ntrack=%d,idhit=%d) " | |
914 | "Calling SetCoupling couplC=%e couplR=%e", | |
915 | col,row,ntrack,idhit,couplC,couplR); | |
916 | j1 = col; | |
917 | j2 = row; | |
918 | pulse1 = GetMap()->GetSignalOnly(col,row); | |
919 | pulse2 = pulse1; | |
920 | for (Int_t isign=-1;isign<=1;isign+=2){// loop in col direction | |
921 | do{ | |
922 | j1 += isign; | |
923 | xr = gRandom->Rndm(); | |
924 | if ((j1<0) || (j1>GetNPixelsZ()-1) || (xr>couplC)){ | |
925 | j1 = col; | |
926 | flag = 1; | |
927 | }else{ | |
928 | UpdateMapSignal(row,j1,ntrack,idhit,pulse1); | |
929 | // flag = 0; | |
930 | flag = 1; // only first next!! | |
931 | } // end if | |
932 | } while(flag == 0); | |
933 | // loop in row direction | |
934 | do{ | |
935 | j2 += isign; | |
936 | xr = gRandom->Rndm(); | |
937 | if ((j2<0) || (j2>GetNPixelsX()-1) || (xr>couplR)){ | |
938 | j2 = row; | |
939 | flag = 1; | |
940 | }else{ | |
941 | UpdateMapSignal(j2,col,ntrack,idhit,pulse2); | |
942 | // flag = 0; | |
943 | flag = 1; // only first next!! | |
944 | } // end if | |
945 | } while(flag == 0); | |
946 | } // for isign | |
f8d9a5b8 | 947 | } |
948 | //______________________________________________________________________ | |
590d15ee | 949 | void AliITSsimulationSPD::SetCouplingOld(Int_t col, Int_t row, |
20eef074 | 950 | Int_t ntrack,Int_t idhit) { |
951 | // Take into account the coupling between adiacent pixels. | |
952 | // The parameters probcol and probrow are the fractions of the | |
953 | // signal in one pixel shared in the two adjacent pixels along | |
954 | // the column and row direction, respectively. | |
955 | //Begin_Html | |
956 | /* | |
957 | <img src="picts/ITS/barimodel_3.gif"> | |
958 | </pre> | |
959 | <br clear=left> | |
960 | <font size=+2 color=red> | |
961 | <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>. | |
962 | </font> | |
963 | <pre> | |
964 | */ | |
965 | //End_Html | |
966 | // Inputs: | |
967 | // Int_t col z cell index | |
968 | // Int_t row x cell index | |
969 | // Int_t ntrack track incex number | |
970 | // Int_t idhit hit index number | |
971 | // Int_t module module number | |
972 | // Outputs: | |
973 | // none. | |
974 | // Return: | |
975 | // none. | |
976 | Int_t j1,j2,flag=0; | |
977 | Double_t pulse1,pulse2; | |
978 | Double_t couplR=0.0,couplC=0.0; | |
979 | ||
980 | GetCouplings(couplC,couplR); | |
981 | ||
982 | // Debugging ... | |
5bfe44ce | 983 | // cout << "Threshold --> " << GetThreshold() << endl; // dom |
8e2f60b2 | 984 | // cout << "Couplings --> " << couplC << " " << couplR << endl; // dom |
5bfe44ce | 985 | |
20eef074 | 986 | if(GetDebug(3)) Info("SetCouplingOld","(col=%d,row=%d,ntrack=%d,idhit=%d) " |
987 | "Calling SetCoupling couplC=%e couplR=%e", | |
988 | col,row,ntrack,idhit,couplC,couplR); | |
989 | for (Int_t isign=-1;isign<=1;isign+=2){// loop in col direction | |
990 | pulse1 = GetMap()->GetSignalOnly(col,row); | |
991 | pulse2 = pulse1; | |
992 | j1 = col; | |
993 | j2 = row; | |
994 | do{ | |
995 | j1 += isign; | |
996 | pulse1 *= couplC; | |
997 | if ((j1<0)||(j1>GetNPixelsZ()-1)||(pulse1<GetThreshold())){ | |
998 | pulse1 = GetMap()->GetSignalOnly(col,row); | |
999 | j1 = col; | |
1000 | flag = 1; | |
1001 | }else{ | |
1002 | UpdateMapSignal(row,j1,ntrack,idhit,pulse1); | |
1003 | // flag = 0; | |
1004 | flag = 1; // only first next !! | |
1005 | } // end if | |
1006 | } while(flag == 0); | |
1007 | // loop in row direction | |
1008 | do{ | |
1009 | j2 += isign; | |
1010 | pulse2 *= couplR; | |
1011 | if((j2<0)||(j2>(GetNPixelsX()-1))||(pulse2<GetThreshold())){ | |
1012 | pulse2 = GetMap()->GetSignalOnly(col,row); | |
1013 | j2 = row; | |
1014 | flag = 1; | |
1015 | }else{ | |
1016 | UpdateMapSignal(j2,col,ntrack,idhit,pulse2); | |
1017 | // flag = 0; | |
1018 | flag = 1; // only first next!! | |
1019 | } // end if | |
1020 | } while(flag == 0); | |
1021 | } // for isign | |
21b825a4 | 1022 | } |
20eef074 | 1023 | //______________________________________________________________________ |
1024 | void AliITSsimulationSPD::GenerateStrobePhase() | |
1025 | { | |
1026 | // Generate randomly the strobe | |
1027 | // phase w.r.t to the LHC clock | |
1028 | // Done once per event | |
1029 | const Double_t kBunchLenght = 25e-9; // LHC clock | |
1030 | fStrobePhase = ((Double_t)gRandom->Integer(fStrobeLenght))*kBunchLenght- | |
1031 | (Double_t)fStrobeLenght*kBunchLenght+ | |
1032 | kBunchLenght/2; | |
1033 | } | |
1034 |