Partical update of SetDefault.
[u/mrichter/AliRoot.git] / ITS / AliITSsimulationSSD.cxx
CommitLineData
b0f5e3fc 1#include <stdio.h>
8a00af9a 2#include <stdlib.h>
fd61217e 3#include <iostream.h>
b0f5e3fc 4#include <TObjArray.h>
1ca7869b 5#include <TRandom.h>
fd61217e 6#include <TMath.h>
b0f5e3fc 7
fd61217e 8#include "AliITSmodule.h"
9#include "AliITSMapA2.h"
9e8d6423 10#include "AliITSsegmentationSSD.h"
11#include "AliITSresponseSSD.h"
b0f5e3fc 12#include "AliITSsimulationSSD.h"
fd61217e 13//#include "AliITSdictSSD.h"
b0f5e3fc 14#include "AliITSdcsSSD.h"
15#include "AliITS.h"
16#include "AliRun.h"
17
18
19ClassImp(AliITSsimulationSSD);
fd61217e 20////////////////////////////////////////////////////////////////////////
21// Version: 0
22// Written by Enrico Fragiacomo
23// July 2000
24//
25// AliITSsimulationSSD is the simulation of SSDs.
26
b0f5e3fc 27//------------------------------------------------------------
28AliITSsimulationSSD::AliITSsimulationSSD(AliITSsegmentation *seg,
29 AliITSresponse *resp){
30 // Constructor
31
32 fSegmentation = seg;
33 fResponse = resp;
fd61217e 34 Float_t noise[2] = {0.,0.};
35 fResponse->GetNoiseParam(noise[0],noise[1]); // retrieves noise parameters
36 cout<<"nois1,2 ="<<noise[0]<<","<<noise[1]<<endl;
37 fDCS = new AliITSdcsSSD(seg,resp);
b0f5e3fc 38
39 fNstrips = fSegmentation->Npx();
40 fPitch = fSegmentation->Dpx(0);
fd61217e 41 cout<<" Dx,Dz ="<<fSegmentation->Dx()<<","<<fSegmentation->Dz()<<endl;
42 cout<<"fNstrips="<<fNstrips<<" fPitch="<<fPitch<<endl;
43
44 //fP = new TArrayF(fNstrips+1);
45 //fN = new TArrayF(fNstrips+1);
e8189707 46
fd61217e 47 fMapA2 = new AliITSMapA2(fSegmentation);
48
49 //fTracksP = new AliITSdictSSD[fNstrips+1];
50 //fTracksN = new AliITSdictSSD[fNstrips+1];
b0f5e3fc 51
fd61217e 52 fSteps = 100; // still hard-wired - set in SetDetParam and get it via
b0f5e3fc 53 // fDCS together with the others eventually
b0f5e3fc 54}
fd61217e 55
b0f5e3fc 56//___________________________________________________________________________
57AliITSsimulationSSD& AliITSsimulationSSD::operator=(AliITSsimulationSSD
58 &source){
59// Operator =
fd61217e 60
b0f5e3fc 61 if(this==&source) return *this;
62
63 this->fDCS = new AliITSdcsSSD(*(source.fDCS));
fd61217e 64 //this->fN = new TArrayF(*(source.fN));
65 //this->fP = new TArrayF(*(source.fP));
66 this->fMapA2 = source.fMapA2;
67 //this->fTracksP = new AliITSdictSSD(*(source.fTracksP));
68 //this->fTracksN = new AliITSdictSSD(*(source.fTracksN));
b0f5e3fc 69 this->fNstrips = source.fNstrips;
70 this->fPitch = source.fPitch;
71 this->fSteps = source.fSteps;
72 return *this;
73}
74//_____________________________________________________________
75AliITSsimulationSSD::AliITSsimulationSSD(AliITSsimulationSSD &source){
76 // copy constructor
fd61217e 77
b0f5e3fc 78 *this = source;
79}
80//____________________________________________________________________________
81AliITSsimulationSSD::~AliITSsimulationSSD() {
fd61217e 82 // destructor
83 //if(fP) delete fP;
84 //if(fN) delete fN;
85 delete fMapA2;
86 //if(fTracksP) delete [] fTracksP;
87 //if(fTracksN) delete [] fTracksN;
b0f5e3fc 88 delete fDCS;
89}
fd61217e 90
b0f5e3fc 91//_______________________________________________________________
b0f5e3fc 92void AliITSsimulationSSD::DigitiseModule(AliITSmodule *mod,Int_t module,
93 Int_t dummy) {
fd61217e 94 // Digitizes hits for one SSD module
95
96 TObjArray *hits = mod->GetHits();
97 Int_t nhits = hits->GetEntriesFast();
98 if (!nhits) return;
99 //cout<<"!! module, nhits ="<<module<<","<<nhits<<endl; //b.b.
100
101 Double_t x0=0.0, y0=0.0, z0=0.0;
102 Double_t x1=0.0, y1=0.0, z1=0.0;
103 Double_t de=0.0;
104 Int_t maxNdigits = 2*fNstrips;
105 Float_t **pList = new Float_t* [maxNdigits];
106 memset(pList,0,sizeof(Float_t*)*maxNdigits);
107 Int_t indexRange[4] = {0,0,0,0};
108 static Bool_t first=kTRUE;
109 Int_t lasttrack = -2;
110 Int_t idtrack = -2;
b0f5e3fc 111
fd61217e 112 for(Int_t i=0; i<nhits; i++) {
113 // LineSegmentL returns 0 if the hit is entering
114 // If hits is exiting returns positions of entering and exiting hits
115 // Returns also energy loss
116 if (mod->LineSegmentL(i, x0, x1, y0, y1, z0, z1, de, idtrack)) {
117 //cout<<"!! befor HitToDigit: hit ="<<i<<endl; //b.b.
118 HitToDigit(module, x0, y0, z0, x1, y1, z1, de, indexRange, first);
119
120 if (lasttrack != idtrack || i==(nhits-1)) {
121 GetList(idtrack,pList,indexRange);
122 first=kTRUE;
123 }
124 lasttrack=idtrack;
b0f5e3fc 125 }
fd61217e 126 } // end loop over hits
127
128 ApplyNoise();
129 ApplyCoupling();
e8189707 130
fd61217e 131 ChargeToSignal(pList);
b0f5e3fc 132
fd61217e 133 fMapA2->ClearMap();
b0f5e3fc 134}
135
136//---------------------------------------------------------------
fd61217e 137void AliITSsimulationSSD::HitToDigit(Int_t module, Double_t x0, Double_t y0,
138 Double_t z0, Double_t x1, Double_t y1,
139 Double_t z1, Double_t de,
140 Int_t *indexRange, Bool_t first) {
141 // Turns hits in SSD module into one or more digits.
142
143 AliITSsegmentationSSD *seg = new AliITSsegmentationSSD();
144 // AliITSresponseSSD *res = new AliITSresponseSSD();
145
146 Float_t tang[2] = {0.0,0.0};
147 seg->Angles(tang[0], tang[1]); // stereo<< -> tan(stereo)~=stereo
148 //fSegmentation->Angles(tang[0], tang[1]); // stereo<< -> tan(stereo)~=stereo
149 Double_t x, y, z;
150 Double_t dex=0.0, dey=0.0, dez=0.0;
151 Double_t pairs;
152 Double_t ionE = 3.62E-9; // ionization energy of Si (GeV)
153 //ionE = (Double_t) res->GetIonE();
b0f5e3fc 154
fd61217e 155 Double_t sigma[2] = {0.,0.}; // standard deviation of the diffusion gaussian
156
157 Double_t D[2] = {11.,30.}; // diffusion constant {h,e} (cm**2/sec)
158 //D[0] = (Double_t) res->GetDiffusionConstantP();
159 //D[1] = (Double_t) res->GetDiffusionConstantN();
160
161 Double_t tdrift[2] = {0.,0.}; // time of drift
162 Double_t vdrift[2] = {0.86E6,2.28E6}; // drift velocity (cm/sec)
163 //vdrift[0] = (Double_t) res->GetDriftVelocityP();
164 //vdrift[1] = (Double_t) res->GetDriftVelocityN();
165
166 Double_t w;
167 Double_t inf[2], sup[2], par0[2];
168 // Steps in the module are determined "manually" (i.e. No Geant)
169 // NumOfSteps divide path between entering and exiting hits in steps
170 Int_t numOfSteps = NumOfSteps(x1, y1, z1, dex, dey, dez);
171
172 // Enery loss is equally distributed among steps
173 de = de/numOfSteps;
174 pairs = de/ionE; // e-h pairs generated
175
176 //cout<<"Dy ="<<seg->Dy()<<endl;
177 //cout<<"numOfSteps ="<<numOfSteps<<endl;
178 //cout<<"dex,dey,dez ="<<dex<<","<<dey<<","<<dez<<endl;
179 //cout<<"y0,y1 ="<<y0<<","<<y1<<endl;
180 for(Int_t j=0; j<numOfSteps; j++) { // stepping
181 // cout<<"step number ="<<j<<endl;
182 x = x0 + (j+0.5)*dex;
183 y = y0 + (j+0.5)*dey;
184 if ( y > (seg->Dy()/2 +10)*1.0E-4 ) {
185 //if ( y > (seg->Dy()*1.0E-4/2) ) {
186 //if ( y > (fSegmentation->Dy()*1.0E-4/2) ) {
187 // check if particle is within the detector
188 cout<<"AliITSsimulationSSD::HitToDigit: Warning: hit out of detector y0,y,dey,j ="<<y0<<","<<y<<","<<dey<<","<<j<<endl;
189 return;
190 };
191 z = z0 + (j+0.5)*dez;
192
193 // calculate drift time
194 tdrift[0] = (y+(seg->Dy()*1.0E-4)/2) / vdrift[0]; // y is the minimum path
195 tdrift[1] = ((seg->Dy()*1.0E-4)/2-y) / vdrift[1]; // y is the minimum path
196 //tdrift[0] = (y+(fSegmentation->Dy()*1.0E-4)/2) / vdrift[0]; // y is the minimum path
197 //tdrift[1] = ((fSegmentation->Dy()*1.0E-4)/2-y) / vdrift[1]; // y is the minimum path
198
199 for(Int_t k=0; k<2; k++) { // both sides remember: 0=Pside 1=Nside
200
201 tang[k]=TMath::Tan(tang[k]);
202
203 // w is the coord. perpendicular to the strips
204 if(k==0) {
205 w = (x+(seg->Dx()*1.0E-4)/2) - (z+(seg->Dz()*1.0E-4)/2)*tang[k];
206 //w = (x+(fSegmentation->Dx()*1.0E-4)/2) - (z+(fSegmentation->Dz()*1.0E-4)/2)*tang[k];
207 //cout<<"k,x,z,w ="<<k<<","<<x<<","<<z<<","<<w<<endl;
208 }
209 else {
210 w = (x+(seg->Dx()*1.0E-4)/2) + (z-(seg->Dz()*1.0E-4)/2)*tang[k];
211 //w = (x+(fSegmentation->Dx()*1.0E-4)/2) + (z-(fSegmentation->Dz()*1.0E-4)/2)*tang[k];
212 //cout<<"k,x,z,w ="<<k<<","<<x<<","<<z<<","<<w<<endl;
213 }
214 w = w / (fPitch*1.0E-4); // w is converted in units of pitch
215
216 if((w<(-0.5)) || (w>(fNstrips-0.5))) {
217 // this check rejects hits in regions not covered by strips
218 // 0.5 takes into account boundaries
219 if(k==0) cout<<"AliITSsimulationSSD::HitToDigit: Warning: no strip in this region of P side"<<endl;
220 else cout<<"AliITSsimulationSSD::HitToDigit: Warning: no strip in this region of N side"<<endl;
221 return;
222 }
957d57f9 223
fd61217e 224 // sigma is the standard deviation of the diffusion gaussian
b0f5e3fc 225
fd61217e 226 if(tdrift[k]<0) return;
b0f5e3fc 227
fd61217e 228 sigma[k] = TMath::Sqrt(2*D[k]*tdrift[k]);
229 sigma[k] = sigma[k] /(fPitch*1.0E-4); //units of Pitch
230 if(sigma[k]==0.0) {
231 cout<<"AliITSsimulationSSD::DigitiseModule: Error: sigma=0"<<endl;
232 exit(0);
233 }
b0f5e3fc 234
fd61217e 235 par0[k] = pairs;
236 // we integrate the diffusion gaussian from -3sigma to 3sigma
237 inf[k] = w - 3*sigma[k]; // 3 sigma from the gaussian average
238 sup[k] = w + 3*sigma[k]; // 3 sigma from the gaussian average
239 // IntegrateGaussian does the actual integration of diffusion gaussian
240 IntegrateGaussian(k, par0[k], w, sigma[k], inf[k], sup[k],
241 indexRange, first);
242 } // end for loop over side (0=Pside, 1=Nside)
243 } // end stepping
244 delete seg;
b0f5e3fc 245}
246
b0f5e3fc 247//____________________________________________________________________
b0f5e3fc 248
249void AliITSsimulationSSD::ApplyNoise() {
250 // Apply Noise.
fd61217e 251
252 Float_t signal;
253 Float_t noise[2] = {0.,0.};
254 fResponse->GetNoiseParam(noise[0],noise[1]); // retrieves noise parameters
255
256 for(Int_t k=0;k<2;k++){ // both sides (0=Pside, 1=Nside)
257 for(Int_t ix=0;ix<fNstrips;ix++){ // loop over strips
258 signal = (Float_t) fMapA2->GetSignal(k,ix);// retrieves signal from map
259
260 signal += gRandom->Gaus(0,noise[k]); // add noise to signal
261 if(signal<0.) signal=0.0; // in case noise is negative...
262
263 fMapA2->SetHit(k,ix,(Double_t)signal); // give back signal to map
264 } // loop over strip
265 } // loop over k (P or N side)
b0f5e3fc 266}
267
268//_________________________________________________________________________
269
270void AliITSsimulationSSD::ApplyCoupling() {
fd61217e 271 // Apply the effect of electronic coupling between channels
272 Float_t signal, signalLeft=0, signalRight=0;
273
274 for(Int_t ix=0;ix<fNstrips;ix++){
275 if(ix>0.) signalLeft = (Float_t) fMapA2->GetSignal(0,ix-1)*fDCS->GetCouplingPL();
276 else signalLeft = 0.0;
277 if(ix<(fNstrips-1)) signalRight = (Float_t) fMapA2->GetSignal(0,ix+1)*fDCS->GetCouplingPR();
278 else signalRight = 0.0;
279 signal = (Float_t) fMapA2->GetSignal(0,ix);
280 signal += signalLeft + signalRight;
281 fMapA2->SetHit(0,ix,(Double_t)signal);
282
283 if(ix>0.) signalLeft = (Float_t) fMapA2->GetSignal(1,ix-1)*fDCS->GetCouplingNL();
284 else signalLeft = 0.0;
285 if(ix<(fNstrips-1)) signalRight = (Float_t) fMapA2->GetSignal(1,ix+1)*fDCS->GetCouplingNR();
286 else signalRight = 0.0;
287 signal = (Float_t) fMapA2->GetSignal(1,ix);
288 signal += signalLeft + signalRight;
289 fMapA2->SetHit(1,ix,(Double_t)signal);
290 } // loop over strips
b0f5e3fc 291}
292
fd61217e 293//____________________________________________________________________________
294Float_t AliITSsimulationSSD::F(Float_t av, Float_t x, Float_t s) {
295 // Computes the integral of a gaussian using Error Function
296 Float_t sqrt2 = TMath::Sqrt(2.0);
297 Float_t sigm2 = sqrt2*s;
298 Float_t integral;
299
300 integral = 0.5 * TMath::Erf( (x - av) / sigm2);
301 return integral;
302}
b0f5e3fc 303
fd61217e 304//_________________________________________________________________________
305void AliITSsimulationSSD::IntegrateGaussian(Int_t k,Double_t par, Double_t w,
306 Double_t sigma,
307 Double_t inf, Double_t sup,
308 Int_t *indexRange, Bool_t first) {
309 // integrate the diffusion gaussian
310 // remind: inf and sup are w-3sigma and w+3sigma
311 // we could define them here instead of passing them
312 // this way we are free to introduce asimmetry
313
314 Double_t a=0.0, b=0.0;
315 Double_t signal = 0.0, dXCharge1 = 0.0, dXCharge2 = 0.0;
316 // dXCharge1 and 2 are the charge to two neighbouring strips
317 // Watch that we only involve at least two strips
318 // Numbers greater than 2 of strips in a cluster depend on
319 // geometry of the track and delta rays, not charge diffusion!
320
321 Double_t strip = TMath::Floor(w); // clostest strip on the left
322
323 if ( TMath::Abs((strip - w)) < 0.5) {
324 // gaussian mean is closer to strip on the left
325 a = inf; // integration starting point
326 if((strip+0.5)<=sup) {
327 // this means that the tail of the gaussian goes beyond
328 // the middle point between strips ---> part of the signal
329 // is given to the strip on the right
330 b = strip + 0.5; // integration stopping point
331 dXCharge1 = F( w, b, sigma) - F(w, a, sigma);
332 dXCharge2 = F( w, sup, sigma) - F(w ,b, sigma);
333 }
334 else {
335 // this means that all the charge is given to the strip on the left
336 b = sup;
337 dXCharge1 = 0.9973; // gaussian integral at 3 sigmas
338 dXCharge2 = 0.0;
339 }
b0f5e3fc 340
fd61217e 341 dXCharge1 = par * dXCharge1; // normalize by mean of number of carriers
342 dXCharge2 = par * dXCharge2;
343
344 // for the time being, signal is the charge
345 // in ChargeToSignal signal is converted in ADC channel
346 signal = fMapA2->GetSignal(k,strip);
347 signal += dXCharge1;
348
349 fMapA2->SetHit(k,strip,(Double_t)signal);
350 if(((Int_t) strip) < (fNstrips-1)) {
351 // strip doesn't have to be the last (remind: last=fNstrips-1)
352 // otherwise part of the charge is lost
353 signal = fMapA2->GetSignal(k,(strip+1));
354 signal += dXCharge2;
355 fMapA2->SetHit(k,(strip+1),(Double_t)signal);
356 }
357
358 if(dXCharge1 > 1.) {
359 if (first) {
360 indexRange[k*2+0]=indexRange[k*2+1]=(Int_t) strip;
361 first=kFALSE;
362 }
b0f5e3fc 363
fd61217e 364 indexRange[k*2+0]=TMath::Min(indexRange[k*2+0],(Int_t) strip);
365 indexRange[k*2+1]=TMath::Max(indexRange[k*2+1],(Int_t) strip);
366 } // dXCharge > 1 e-
367
368 }
369 else {
370 // gaussian mean is closer to strip on the right
371 strip++; // move to strip on the rigth
372 b = sup; // now you know where to stop integrating
373 if((strip-0.5)>=inf) {
374 // tail of diffusion gaussian on the left goes left of
375 // middle point between strips
376 a = strip - 0.5; // integration starting point
377 dXCharge1 = F(w, b, sigma) - F(w, a, sigma);
378 dXCharge2 = F(w, a, sigma) - F(w, inf, sigma);
379 }
380 else {
381 a = inf;
382 dXCharge1 = 0.9973; // gaussian integral at 3 sigmas
383 dXCharge2 = 0.0;
b0f5e3fc 384 }
fd61217e 385
386 dXCharge1 = par * dXCharge1; // normalize by means of carriers
387 dXCharge2 = par * dXCharge2;
388
389 // for the time being, signal is the charge
390 // in ChargeToSignal signal is converted in ADC channel
391 signal = fMapA2->GetSignal(k,strip);
392 signal += dXCharge1;
393 fMapA2->SetHit(k,strip,(Double_t)signal);
394 if(((Int_t) strip) > 0) {
395 // strip doesn't have to be the first
396 // otherwise part of the charge is lost
397 signal = fMapA2->GetSignal(k,(strip-1));
398 signal += dXCharge2;
399 fMapA2->SetHit(k,(strip-1),(Double_t)signal);
400 }
401
402 if(dXCharge1 > 1.) {
403 if (first) {
404 indexRange[k*2+0]=indexRange[k*2+1]=(Int_t) strip;
405 first=kFALSE;
406 }
407
408 indexRange[k*2+0]=TMath::Min(indexRange[k*2+0],(Int_t) strip);
409 indexRange[k*2+1]=TMath::Max(indexRange[k*2+1],(Int_t) strip);
410 } // dXCharge > 1 e-
411 }
b0f5e3fc 412}
413
b0f5e3fc 414
fd61217e 415 //_________________________________________________________________________
416
417Int_t AliITSsimulationSSD::NumOfSteps(Double_t x, Double_t y, Double_t z,
418 Double_t & dex,Double_t & dey,Double_t & dez) {
419 // number of steps
420 // it also returns steps for each coord
421 //AliITSsegmentationSSD *seg = new AliITSsegmentationSSD();
422
423 Double_t step = 25E-4;
424 //step = (Double_t) seg->GetStepSize(); // step size (cm)
425 Int_t numOfSteps = (Int_t) (TMath::Sqrt(x*x+y*y+z*z)/step);
426
427 if (numOfSteps < 1) numOfSteps = 1; // one step, at least
428
429 // we could condition the stepping depending on the incident angle
430 // of the track
431 dex = x/numOfSteps;
432 dey = y/numOfSteps;
433 dez = z/numOfSteps;
b0f5e3fc 434
fd61217e 435 return numOfSteps;
b0f5e3fc 436
fd61217e 437}
438
439//---------------------------------------------
440void AliITSsimulationSSD::GetList(Int_t label,Float_t **pList,Int_t *indexRange) {
441 // loop over nonzero digits
442 Int_t ix,globalIndex;
443 Float_t signal=0.;
444 Float_t highest,middle,lowest;
445 // printf("SPD-GetList: indexRange[0] indexRange[1] indexRange[2] indexRange[3] %d %d %d %d\n",indexRange[0], indexRange[1], indexRange[2], indexRange[3]);
446
447 for(Int_t k=0; k<2; k++) {
448 for(ix=indexRange[k*2+0];ix<indexRange[k*2+1]+1;ix++){
449 if(indexRange[k*2+0]<indexRange[k*2+1])
450 signal=fMapA2->GetSignal(k,ix);
451
452 globalIndex = k*fNstrips+ix; // globalIndex starts from 0!
453 if(!pList[globalIndex]){
454
455 //
456 // Create new list (6 elements - 3 signals and 3 tracks + total sig)
457 //
458
459 pList[globalIndex] = new Float_t [6];
460 // set list to -1
461
462 *pList[globalIndex] = -2.;
463 *(pList[globalIndex]+1) = -2.;
464 *(pList[globalIndex]+2) = -2.;
465 *(pList[globalIndex]+3) = 0.;
466 *(pList[globalIndex]+4) = 0.;
467 *(pList[globalIndex]+5) = 0.;
468
469 *pList[globalIndex] = (float)label;
470 *(pList[globalIndex]+3) = signal;
e8189707 471
e8189707 472 }
fd61217e 473 else{
474
475 // check the signal magnitude
476
477 highest = *(pList[globalIndex]+3);
478 middle = *(pList[globalIndex]+4);
479 lowest = *(pList[globalIndex]+5);
480
481 signal -= (highest+middle+lowest);
482
483 //
484 // compare the new signal with already existing list
485 //
486
487 if(signal<lowest) continue; // neglect this track
488
489 if (signal>highest){
490 *(pList[globalIndex]+5) = middle;
491 *(pList[globalIndex]+4) = highest;
492 *(pList[globalIndex]+3) = signal;
493
494 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
495 *(pList[globalIndex]+1) = *pList[globalIndex];
496 *pList[globalIndex] = label;
497 }
498 else if (signal>middle){
499 *(pList[globalIndex]+5) = middle;
500 *(pList[globalIndex]+4) = signal;
501
502 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
503 *(pList[globalIndex]+1) = label;
504 }
505 else{
506 *(pList[globalIndex]+5) = signal;
507 *(pList[globalIndex]+2) = label;
508 }
509 }
510 } // end of loop pixels in x
511 } // end of loop over pixels in z
512}
513
514//---------------------------------------------
515void AliITSsimulationSSD::ChargeToSignal(Float_t **pList) {
516 // charge to signal
517
518 AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
519
520 Float_t threshold = 0.;
521
522 Int_t digits[3], tracks[3],hits[3],gi,j1;
523 Float_t charges[3];
524 Float_t signal,phys;
525 Float_t noise[2] = {0.,0.};
526 fResponse->GetNoiseParam(noise[0],noise[1]);
527
528 for(Int_t k=0;k<2;k++){ // both sides (0=Pside, 1=Nside)
529
530 // Threshold for zero-suppression
531 // It can be defined in AliITSresponseSSD
532 // threshold = (Float_t)fResponse->MinVal(k);
533 // I prefer to think adjusting the threshold "manually", looking
534 // at the scope, and considering noise standard deviation
535 threshold = 4.0*noise[k]; // 4 times noise is a choice
536 //cout<<"SSD: k,thresh ="<<k<<","<<threshold<<endl;
537 for(Int_t ix=0;ix<fNstrips;ix++){ // loop over strips
538
539 signal = (Float_t) fMapA2->GetSignal(k,ix);
540
541 gi =k*fNstrips+ix; // global index
542 if (signal > threshold) {
543 digits[0]=k;
544 digits[1]=ix;
545
546 // convert to ADC signal
547 // conversion factor are rather arbitrary (need tuning)
548 // minimum ionizing particle --> ~30000 pairs --> ADC channel 50
549 signal = signal*50.0/30000.0;
550 //cout<<"SSD: 1 signal ="<<signal<<endl;
551 if(signal>1000.) signal = 1000.0; // if exceeding, accumulate last one
552 //cout<<"SSD: 2 signal ="<<signal<<endl;
553 digits[2]=(Int_t) signal;
554
555 //gi =k*fNstrips+ix; // global index
556 for(j1=0;j1<3;j1++){
557 if (pList[gi]) {
558 tracks[j1] = (Int_t)(*(pList[gi]+j1));
559 }
560 else {
561 tracks[j1]=-2; //noise
562 }
563 charges[j1] = 0;
564 }
565
566 phys=0;
567
568 hits[0]=0;
569 hits[1]=0;
570 hits[2]=0;
571 aliITS->AddSimDigit(2,phys,digits,tracks,hits,charges); // finally add digit
572
573 //if(pList[gi]) delete [] pList[gi];
574 }
575 if(pList[gi]) delete [] pList[gi];
b0f5e3fc 576 }
fd61217e 577 }
578 delete [] pList;
b0f5e3fc 579}
580
581
b0f5e3fc 582
9e8d6423 583
b0f5e3fc 584
b0f5e3fc 585
b0f5e3fc 586
b0f5e3fc 587
b0f5e3fc 588
589
fd61217e 590
591
592
593