e53090fc13f52465c8af406fc77cc7db624df722
[u/mrichter/AliRoot.git] / ITS / AliITSsimulationSPD.cxx
1 /**************************************************************************
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  **************************************************************************/
15
16 /* $Id$ */
17
18 #include <Riostream.h>
19 #include <TRandom.h>
20 #include <TH1.h>
21 #include <TMath.h>
22 #include <TString.h>
23 #include <TParticle.h>
24
25 #include "AliRun.h"
26 #include "AliITS.h"
27 #include "AliITShit.h"
28 #include "AliITSdigit.h"
29 #include "AliITSmodule.h"
30 #include "AliITSMapA2.h"
31 #include "AliITSpList.h"
32 #include "AliITSsimulationSPD.h"
33 #include "AliITSsegmentation.h"
34 #include "AliITSresponse.h"
35 #include "AliITSsegmentationSPD.h"
36 #include "AliITSresponseSPD.h"
37
38 //#define DEBUG
39
40 ClassImp(AliITSsimulationSPD)
41 ////////////////////////////////////////////////////////////////////////
42 // Version: 0
43 // Written by Rocco Caliandro
44 // from a model developed with T. Virgili and R.A. Fini
45 // June 15 2000
46 //
47 // AliITSsimulationSPD is the simulation of SPDs
48 //
49 //______________________________________________________________________
50 AliITSsimulationSPD::AliITSsimulationSPD(){
51     // Default constructor
52
53     fResponse     = 0;
54     fSegmentation = 0;
55     fHis          = 0;
56     fMapA2        = 0;
57
58 /*
59     fThresh       = 0.;
60     fSigma        = 0.;
61     fCouplCol     = 0.;
62     fCouplRow     = 0.; */
63 }
64 //______________________________________________________________________
65 AliITSsimulationSPD::AliITSsimulationSPD(AliITSsegmentation *seg,
66                                          AliITSresponse *resp) {
67     // Standard constructor
68
69     fResponse     = 0;
70     fSegmentation = 0;
71     fHis          = 0;
72     fMapA2        = 0;
73
74 /*
75     fThresh       = 0.;
76     fSigma        = 0.;
77     fCouplCol     = 0.;
78     fCouplRow     = 0.*/
79     Init((AliITSsegmentationSPD*)seg,(AliITSresponseSPD*)resp);
80 }
81 //______________________________________________________________________
82 void AliITSsimulationSPD::Init(AliITSsegmentationSPD *seg,
83                                AliITSresponseSPD *resp) {
84     // Initilizes the variables of AliITSsimulation SPD.
85
86     fHis = 0;
87     fResponse     = resp;
88     fSegmentation = seg;
89     fMapA2  = new AliITSMapA2(fSegmentation);
90     fpList  = new AliITSpList(GetNPixelsZ()+1,GetNPixelsX()+1);
91 /*
92     fResponse->Thresholds(fThresh,fSigma);
93     fResponse->GetNoiseParam(fCouplCol,fCouplRow);
94     fNPixelsZ = fSegmentation->Npz();
95     fNPixelsX = fSegmentation->Npx();
96 */
97 }
98 //______________________________________________________________________
99 AliITSsimulationSPD::~AliITSsimulationSPD() { 
100     // destructor
101
102     delete fMapA2;
103 //    delete fpList;
104
105     if (fHis) {
106         fHis->Delete(); 
107         delete fHis;     
108     } // end if
109 }
110 //______________________________________________________________________
111 AliITSsimulationSPD::AliITSsimulationSPD(const AliITSsimulationSPD &source) :
112     AliITSsimulation(source){
113     // Copy Constructor
114
115     if(&source == this) return;
116
117     this->fMapA2    = source.fMapA2;
118     this->fHis      = source.fHis;
119 /*
120     this->fThresh   = source.fThresh;
121     this->fSigma    = source.fSigma;
122     this->fCouplCol = source.fCouplCol;
123     this->fCouplRow = source.fCouplRow;
124     this->fNPixelsX = source.fNPixelsX;
125     this->fNPixelsZ = source.fNPixelsZ;
126 */
127     return;
128 }
129 //______________________________________________________________________
130 AliITSsimulationSPD& AliITSsimulationSPD::operator=(const AliITSsimulationSPD 
131                                                     &source) {
132     //    Assignment operator
133
134     if(&source == this) return *this;
135
136     this->fMapA2    = source.fMapA2;
137     this->fHis      = source.fHis;
138 /*
139     this->fThresh   = source.fThresh;
140     this->fSigma    = source.fSigma;
141     this->fCouplCol = source.fCouplCol;
142     this->fCouplRow = source.fCouplRow;
143     this->fNPixelsX = source.fNPixelsX;
144     this->fNPixelsZ = source.fNPixelsZ;
145 */
146     return *this;
147
148 //______________________________________________________________________
149 void AliITSsimulationSPD::InitSimulationModule(Int_t module,Int_t event){
150     // Creates maps to build the list of tracks for each sumable digit
151     // Inputs:
152     //   Int_t module    // Module number to be simulated
153     //   Int_t event     // Event number to be simulated
154     // Outputs:
155     //   none.
156     // Return
157     //    none.
158  
159     fModule = module;
160     fEvent  = event;
161     fMapA2->ClearMap();
162     fpList->ClearMap();
163 }
164 //______________________________________________________________________
165 void AliITSsimulationSPD::FinishSDigitiseModule(){
166     // Does the Sdigits to Digits work
167     // Inputs:
168     //   none.
169     // Outputs:
170     //   none.
171     // Return:
172     //   none.
173
174     SDigitsToDigits(fModule,fpList);
175 }
176 //______________________________________________________________________
177 void AliITSsimulationSPD::SDigitiseModule(AliITSmodule *mod, Int_t dummy0,
178                                              Int_t dummy1) {
179     // Sum digitize module
180     if (!(mod->GetNhits())) return; // if module has no hits then no Sdigits.
181     Int_t    number     = 10000;
182     Int_t    *frowpixel = new Int_t[number];
183     Int_t    *fcolpixel = new Int_t[number];
184     Double_t *fenepixel = new Double_t[number];
185
186     dummy0 = dummy1; // remove unsued variable warning.
187     fModule = mod->GetIndex();
188
189     // Array of pointers to store the track index of the digits
190     // leave +1, otherwise pList crashes when col=256, row=192
191
192     HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,fpList);
193
194     WriteSDigits(fpList);
195
196     // clean memory
197     delete[] frowpixel;
198     delete[] fcolpixel;
199     delete[] fenepixel;
200     fMapA2->ClearMap();
201     fpList->ClearMap();
202 }
203 //______________________________________________________________________
204 void AliITSsimulationSPD::DigitiseModule(AliITSmodule *mod, Int_t dummy0,
205                                              Int_t dummy1) {
206     // digitize module. Also need to digitize modules with only noise.
207
208     Int_t    number     = 10000;
209     Int_t    *frowpixel = new Int_t[number];
210     Int_t    *fcolpixel = new Int_t[number];
211     Double_t *fenepixel = new Double_t[number];
212
213     dummy0 = dummy1; // remove unsued variable warning.
214     // Array of pointers to store the track index of the digits
215     // leave +1, otherwise pList crashes when col=256, row=192
216     fModule = mod->GetIndex();
217     // noise setting
218     SetFluctuations(fpList,fModule);
219
220     HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,fpList);
221
222     // apply mask to SPD module
223     SetMask();
224
225     CreateDigit(fModule,fpList);
226
227     // clean memory
228     delete[] frowpixel;
229     delete[] fcolpixel;
230     delete[] fenepixel;
231     fMapA2->ClearMap();
232     fpList->ClearMap();
233 }
234 //______________________________________________________________________
235 void AliITSsimulationSPD::SDigitsToDigits(Int_t module,AliITSpList *pList) {
236     // sum digits to Digits.
237
238 #ifdef DEBUG
239     cout << "Entering AliITSsimulatinSPD::SDigitsToDigits for module=";
240     cout << module << endl;
241 #endif
242     fModule = module;
243
244     // noise setting
245     SetFluctuations(pList,module);
246
247     fMapA2->ClearMap(); // since noise is in pList aready. Zero Map so that
248     // noise is not doubled when calling FillMapFrompList.
249
250     FillMapFrompList(pList);
251
252     // apply mask to SPD module
253     SetMask();
254
255     CreateDigit(module,pList);
256
257     fMapA2->ClearMap();
258     pList->ClearMap();
259 }
260 //______________________________________________________________________
261 void AliITSsimulationSPD::UpdateMapSignal(Int_t row,Int_t col,Int_t trk,
262                                           Int_t hit,Int_t mod,Double_t ene,
263                                           AliITSpList *pList) {
264     // updates the Map of signal, adding the energy  (ene) released by
265     // the current track
266
267     fMapA2->AddSignal(row,col,ene);
268     pList->AddSignal(row,col,trk,hit,mod,ene);
269 }
270 //______________________________________________________________________
271 void AliITSsimulationSPD::UpdateMapNoise(Int_t row,Int_t col,Int_t mod,
272                                          Double_t ene,AliITSpList *pList) {
273     // updates the Map of noise, adding the energy  (ene) give my noise
274
275     fMapA2->AddSignal(row,col,ene);
276     pList->AddNoise(row,col,mod,ene);
277 }
278 //______________________________________________________________________
279 void AliITSsimulationSPD::HitsToAnalogDigits(AliITSmodule *mod,
280                                              Int_t *frowpixel,Int_t *fcolpixel,
281                                              Double_t *fenepixel,
282                                              AliITSpList *pList) {
283     // Loops over all hits to produce Analog/floting point digits. This
284     // is also the first task in producing standard digits.
285     
286     // loop over hits in the module
287     Int_t hitpos,nhits = mod->GetNhits();
288     for (hitpos=0;hitpos<nhits;hitpos++) {
289         HitToDigit(mod,hitpos,frowpixel,fcolpixel,fenepixel,pList);
290     }// end loop over digits
291 }
292 //______________________________________________________________________
293 void AliITSsimulationSPD::HitToDigit(AliITSmodule *mod,Int_t hitpos,
294                                      Int_t *frowpixel,Int_t *fcolpixel,
295                                      Double_t *fenepixel,AliITSpList *pList) {
296     //  Steering function to determine the digits associated to a given
297     // hit (hitpos)
298     // The digits are created by charge sharing (ChargeSharing) and by
299     // capacitive coupling (SetCoupling). At all the created digits is
300     // associated the track number of the hit (ntrack)
301     Double_t x1l=0.0,y1l=0.0,z1l=0.0,x2l=0.0,y2l=0.0,z2l=0.0;
302     Int_t r1,r2,c1,c2,row,col,npixel = 0;
303     Int_t ntrack;
304     Double_t ene=0.0,etot=0.0;
305     const Float_t kconv = 10000.;     // cm -> microns
306     const Float_t kconv1= 0.277e9;    // GeV -> electrons equivalent
307
308     if(!(mod->LineSegmentL(hitpos,x1l,x2l,y1l,y2l,z1l,z2l,etot,ntrack)))return;
309
310     x2l += x1l; y2l += y1l; z2l += z1l; // Convert to ending coordinate.
311     // positions shifted and converted in microns
312     x1l   = x1l*kconv + fSegmentation->Dx()/2.;
313     z1l   = z1l*kconv + fSegmentation->Dz()/2.;
314     // positions  shifted and converted in microns
315     x2l   = x2l*kconv + fSegmentation->Dx()/2.;
316     z2l   = z2l*kconv + fSegmentation->Dz()/2.;
317     etot *= kconv1; // convert from GeV to electrons equivalent.
318     Int_t module = mod->GetIndex();
319
320     // to account for the effective sensitive area
321     // introduced in geometry 
322     if (z1l<0 || z1l>fSegmentation->Dz()) return;
323     if (z2l<0 || z2l>fSegmentation->Dz()) return;
324     if (x1l<0 || x1l>fSegmentation->Dx()) return;
325     if (x2l<0 || x2l>fSegmentation->Dx()) return;
326
327     //Get the col and row number starting from 1
328     // the x direction is not inverted for the second layer!!!
329     fSegmentation->GetPadIxz(x1l, z1l, c1, r1); 
330     fSegmentation->GetPadIxz(x2l, z2l, c2, r2);
331
332     // to account for unexpected equal entrance and 
333     // exit coordinates
334     if (x1l==x2l) x2l=x2l+x2l*0.1;
335     if (z1l==z2l) z2l=z2l+z2l*0.1;
336
337     if ((r1==r2) && (c1==c2)){
338         // no charge sharing
339         npixel = 1;              
340         frowpixel[npixel-1] = r1;
341         fcolpixel[npixel-1] = c1;
342         fenepixel[npixel-1] = etot;
343     } else {
344         // charge sharing
345         ChargeSharing(x1l,z1l,x2l,z2l,c1,r1,c2,r2,etot,
346                       npixel,frowpixel,fcolpixel,fenepixel);
347     } // end if r1==r2 && c1==c2.
348
349     for (Int_t npix=0;npix<npixel;npix++){
350         row = frowpixel[npix];
351         col = fcolpixel[npix];
352         ene = fenepixel[npix];
353         UpdateMapSignal(row,col,ntrack,hitpos,module,ene,pList); 
354         // Starting capacitive coupling effect
355         SetCoupling(row,col,ntrack,hitpos,module,pList); 
356     } // end for npix
357 }
358 //______________________________________________________________________
359 void AliITSsimulationSPD::ChargeSharing(Float_t x1l,Float_t z1l,Float_t x2l,
360                                         Float_t z2l,Int_t c1,Int_t r1,Int_t c2,
361                                         Int_t r2,Float_t etot,
362                                         Int_t &npixel,Int_t *frowpixel,
363                                         Int_t *fcolpixel,Double_t *fenepixel){
364     //  Take into account the geometrical charge sharing when the track
365     //  crosses more than one pixel.
366     //
367     //Begin_Html
368     /*
369       <img src="picts/ITS/barimodel_2.gif">
370       </pre>
371       <br clear=left>
372       <font size=+2 color=red>
373       <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
374       </font>
375       <pre>
376     */
377     //End_Html
378     //Float_t dm;
379     Float_t xa,za,xb,zb,dx,dz,dtot,refr,refm,refc;
380     Float_t refn=0.;
381     Float_t arefm, arefr, arefn, arefc, azb, az2l, axb, ax2l;
382     Int_t   dirx,dirz,rb,cb;
383     Int_t flag,flagrow,flagcol;
384     Double_t epar;
385
386     npixel = 0;
387     xa     = x1l;
388     za     = z1l;
389 //    dx     = x1l-x2l;
390 //    dz     = z1l-z2l;
391     dx     = x2l-x1l;
392     dz     = z2l-z1l;
393     dtot   = TMath::Sqrt((dx*dx)+(dz*dz));   
394     if (dtot==0.0) dtot = 0.01;
395     dirx   = (Int_t) TMath::Sign((Float_t)1,dx);
396     dirz   = (Int_t) TMath::Sign((Float_t)1,dz);
397
398     // calculate the x coordinate of  the pixel in the next column    
399     // and the z coordinate of  the pixel in the next row
400     Float_t xpos, zpos;
401
402     fSegmentation->GetPadCxz(c1, r1-1, xpos, zpos); 
403
404     Float_t xsize = fSegmentation->Dpx(0);
405     Float_t zsize = fSegmentation->Dpz(r1-1);
406     
407     if (dirx == 1) refr = xpos+xsize/2.;
408     else refr = xpos-xsize/2.;
409
410     if (dirz == 1) refn = zpos+zsize/2.;
411     else refn = zpos-zsize/2.;
412
413     flag = 0;
414     flagrow = 0;
415     flagcol = 0;
416     do{
417         // calculate the x coordinate of the intersection with the pixel
418         // in the next cell in row  direction
419       if(dz!=0)
420         refm = dx*((refn - z1l)/dz) + x1l;
421       else
422         refm = refr+dirx*xsize;
423    
424         // calculate the z coordinate of the intersection with the pixel
425         // in the next cell in column direction
426       if (dx!=0)
427         refc = dz*((refr - x1l)/dx) + z1l;
428       else
429         refc = refn+dirz*zsize;
430
431         arefm = refm * dirx;
432         arefr = refr * dirx;
433         arefn = refn * dirz;
434         arefc = refc * dirz;
435
436         if ((arefm < arefr) && (arefn < arefc)){
437             // the track goes in the pixel in the next cell in row direction
438             xb = refm;
439             zb = refn;
440             cb = c1;
441             rb = r1 + dirz;
442             azb = zb * dirz;
443             az2l = z2l * dirz;
444             if (rb == r2) flagrow=1;
445             if (azb > az2l) {
446                 zb = z2l;
447                 xb = x2l;
448             } // end if
449             // shift to the pixel in the next cell in row direction
450             Float_t zsizeNext = fSegmentation->Dpz(rb-1);
451             //to account for cell at the borders of the detector
452             if(zsizeNext==0) zsizeNext = zsize;
453             refn += zsizeNext*dirz;
454         }else {
455             // the track goes in the pixel in the next cell in column direction
456             xb = refr;
457             zb = refc;
458             cb = c1 + dirx;
459             rb = r1;
460             axb = xb * dirx;
461             ax2l = x2l * dirx;
462             if (cb == c2) flagcol=1;
463             if (axb > ax2l) {
464                 zb = z2l;
465                 xb = x2l;
466             } // end ifaxb > ax2l
467
468             // shift to the pixel in the next cell in column direction
469             Float_t xsizeNext = fSegmentation->Dpx(cb-1);
470             //to account for cell at the borders of the detector
471             if(xsizeNext==0) xsizeNext = xsize;
472             refr += xsizeNext*dirx;
473         } // end if (arefm < arefr) && (arefn < arefc)
474
475         //calculate the energy lost in the crossed pixel      
476         epar = TMath::Sqrt((xb-xa)*(xb-xa)+(zb-za)*(zb-za)); 
477         epar = etot*(epar/dtot);
478
479         //store row, column and energy lost in the crossed pixel
480         frowpixel[npixel] = r1;
481         fcolpixel[npixel] = c1;
482         fenepixel[npixel] = epar;
483         npixel++;
484  
485         // the exit point of the track is reached
486         if (epar == 0) flag = 1;
487         if ((r1 == r2) && (c1 == c2)) flag = 1;
488         if (flag!=1) {
489             r1 = rb;
490             c1 = cb;
491             xa = xb;
492             za = zb;
493         } // end if flag!=1
494     } while (flag==0);
495 }
496 //______________________________________________________________________
497 void AliITSsimulationSPD::SetCoupling(Int_t row, Int_t col, Int_t ntrack,
498                                       Int_t idhit,Int_t module,
499                                       AliITSpList *pList) {
500     //  Take into account the coupling between adiacent pixels.
501     //  The parameters probcol and probrow are the probability of the
502     //  signal in one pixel shared in the two adjacent pixels along
503     //  the column and row direction, respectively.
504     //
505     //Begin_Html
506     /*
507       <img src="picts/ITS/barimodel_3.gif">
508       </pre>
509       <br clear=left>
510       <font size=+2 color=red>
511       <a href="mailto:tiziano.virgili@cern.ch"></a>.
512       </font>
513       <pre>
514     */
515     //End_Html
516     Int_t j1,j2,flag=0;
517     Double_t pulse1,pulse2;
518     Float_t couplR=0.0,couplC=0.0;
519     Double_t xr=0.;
520
521     GetCouplings(couplR,couplC);
522     j1 = row;
523     j2 = col;
524     pulse1 = fMapA2->GetSignal(row,col);
525     pulse2 = pulse1;
526     for (Int_t isign=-1;isign<=1;isign+=2){// loop in row direction
527         do{
528             j1 += isign;
529             //   pulse1 *= couplR; 
530           xr = gRandom->Rndm();
531           //   if ((j1<0) || (j1>GetNPixelsZ()-1) || (pulse1<GetThreshold())){
532             if ((j1<0) || (j1>GetNPixelsZ()-1) || (xr>couplR)){
533                 j1 = row;
534                 flag = 1;
535             }else{
536                 UpdateMapSignal(j1,col,ntrack,idhit,module,pulse1,pList);
537               //  flag = 0;
538               flag = 1; // only first next!!
539             } // end if
540         } while(flag == 0);
541         // loop in column direction
542       do{
543           j2 += isign;
544           // pulse2 *= couplC; 
545           xr = gRandom->Rndm();
546           //  if ((j2<0) || (j2>(GetNPixelsX()-1)) || (pulse2<GetThreshold())){
547             if ((j2<0) || (j2>GetNPixelsX()-1) || (xr>couplC)){
548               j2 = col;
549               flag = 1;
550           }else{
551               UpdateMapSignal(row,j2,ntrack,idhit,module,pulse2,pList);
552               //  flag = 0;
553               flag = 1; // only first next!!
554           } // end if
555       } while(flag == 0);
556     } // for isign
557 }
558 //______________________________________________________________________
559 void AliITSsimulationSPD::SetCouplingOld(Int_t row, Int_t col, Int_t ntrack,
560                                          Int_t idhit,Int_t module,
561                                          AliITSpList *pList) {
562     //  Take into account the coupling between adiacent pixels.
563     //  The parameters probcol and probrow are the fractions of the
564     //  signal in one pixel shared in the two adjacent pixels along
565     //  the column and row direction, respectively.
566     //
567     //Begin_Html
568     /*
569       <img src="picts/ITS/barimodel_3.gif">
570       </pre>
571       <br clear=left>
572       <font size=+2 color=red>
573       <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
574       </font>
575       <pre>
576     */
577     //End_Html
578     Int_t j1,j2,flag=0;
579     Double_t pulse1,pulse2;
580     Float_t couplR=0.0,couplC=0.0;
581
582     GetCouplings(couplR,couplC);
583     j1 = row;
584     j2 = col;
585     pulse1 = fMapA2->GetSignal(row,col);
586     pulse2 = pulse1;
587     for (Int_t isign=-1;isign<=1;isign+=2){// loop in row direction
588         do{
589             j1 += isign;
590             pulse1 *= couplR;
591             if ((j1<0) || (j1>GetNPixelsZ()-1) || (pulse1<GetThreshold())){
592                 pulse1 = fMapA2->GetSignal(row,col);
593                 j1 = row;
594                 flag = 1;
595             }else{
596                 UpdateMapSignal(j1,col,ntrack,idhit,module,pulse1,pList);
597                 flag = 0;
598             } // end if
599         } while(flag == 0);
600         // loop in column direction
601       do{
602           j2 += isign;
603           pulse2 *= couplC;
604           if ((j2<0) || (j2>(GetNPixelsX()-1)) || (pulse2<GetThreshold())){
605               pulse2 = fMapA2->GetSignal(row,col);
606               j2 = col;
607               flag = 1;
608           }else{
609               UpdateMapSignal(row,j2,ntrack,idhit,module,pulse2,pList);
610               flag = 0;
611           } // end if
612       } while(flag == 0);
613     } // for isign
614 }
615 //______________________________________________________________________
616 void AliITSsimulationSPD::CreateDigit(Int_t module,AliITSpList *pList) {
617     // The pixels are fired if the energy deposited inside them is above
618     // the threshold parameter ethr. Fired pixed are interpreted as digits
619     // and stored in the file digitfilename. One also needs to write out
620     // cases when there is only noise (nhits==0).
621
622     static AliITS *aliITS  = (AliITS*)gAlice->GetModule("ITS");
623
624     Int_t size = AliITSdigitSPD::GetNTracks();
625     Int_t * digits = new Int_t[size];
626     Int_t * tracks = new Int_t[size];
627     Int_t * hits = new Int_t[size];
628     Float_t * charges = new Float_t[size]; 
629     Int_t j1;
630
631     module=0; // remove unused variable warning.
632     for(j1=0;j1<size;j1++){tracks[j1]=-3;hits[j1]=-1;charges[j1]=0.0;}
633     for (Int_t r=1;r<=GetNPixelsZ();r++) {
634         for (Int_t c=1;c<=GetNPixelsX();c++) {
635             // check if the deposited energy in a pixel is above the
636             // threshold 
637             Float_t signal = (Float_t) fMapA2->GetSignal(r,c);
638             if ( signal > GetThreshold()) {
639                 digits[0] = r-1;  // digits starts from 0
640                 digits[1] = c-1;  // digits starts from 0
641                 //digits[2] = 1;  
642                 digits[2] =  (Int_t) signal;  // the signal is stored in
643                                               //  electrons
644                 for(j1=0;j1<size;j1++){
645                     if(j1<pList->GetNEnteries()){
646                         tracks[j1] = pList->GetTrack(r,c,j1);
647                         hits[j1]   = pList->GetHit(r,c,j1);
648                         //}else{
649                         //tracks[j1] = -3;
650                         //hits[j1]   = -1;
651                     } // end if
652                     //charges[j1] = 0;
653                 } // end for j1
654                 Float_t phys = 0;
655                 aliITS->AddSimDigit(0,phys,digits,tracks,hits,charges);
656 #ifdef DEBUG
657                 cout << " CreateSPDDigit mod=" << fModule << " r,c=" << r;
658                 cout <<","<<c<< " sig=" << fpList->GetSignalOnly(r,c);
659                 cout << " noise=" << fpList->GetNoise(r,c);
660                 cout << " Msig="<< signal << " Thres=" << GetThreshold()<<endl;
661 #endif
662             } // end if of threshold condition
663         } // for c
664     }// end do on pixels
665     delete [] digits;
666     delete [] tracks;
667     delete [] hits;
668     delete [] charges;
669
670 }
671 //______________________________________________________________________
672 void AliITSsimulationSPD::SetFluctuations(AliITSpList *pList,Int_t module) {
673     //  Set the electronic noise and threshold non-uniformities to all the
674     //  pixels in a detector.
675     //  The parameter fSigma is the squared sum of the sigma due to noise
676     //  and the sigma of the threshold distribution among pixels.
677     //
678     //Begin_Html
679     /*
680       <img src="picts/ITS/barimodel_1.gif">
681       </pre>
682       <br clear=left>
683       <font size=+2 color=red>
684       <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
685       </font>
686       <pre>
687     */
688     //End_Html
689     Float_t  thr=0.0,sigm=0.0;
690     Double_t signal,sigma;
691     Int_t iz,ix;
692
693     GetThresholds(thr,sigm);
694     sigma = (Double_t) sigm;
695     for(iz=1;iz<=GetNPixelsZ();iz++){
696         for(ix=1;ix<=GetNPixelsX();ix++){
697             signal = sigma*gRandom->Gaus(); 
698             fMapA2->SetHit(iz,ix,signal);
699             // insert in the label-signal-hit list the pixels fired
700             // only by noise
701             pList->AddNoise(iz,ix,module,signal);
702         } // end of loop on pixels
703     } // end of loop on pixels
704 }
705 //______________________________________________________________________
706 void AliITSsimulationSPD::SetMask() {
707     //  Apply a mask to the SPD module. 1% of the pixel channels are
708     //  masked. When the database will be ready, the masked pixels
709     //  should be read from it.
710     Double_t signal;
711     Int_t iz,ix,im;
712     Float_t totMask;
713     Float_t perc = ((AliITSresponseSPD*)fResponse)->GetFractionDead();
714     // in this way we get the same set of random numbers for all runs.
715     // This is a cluge for now.
716     static TRandom *rnd = new TRandom();
717
718     totMask= perc*GetNPixelsZ()*GetNPixelsX();
719     for(im=1;im<totMask;im++){
720         do{
721             ix=(Int_t)(rnd->Rndm()*(GetNPixelsX()-1.) + 1.);
722         } while(ix<=0 || ix>GetNPixelsX());
723         do{
724             iz=(Int_t)(rnd->Rndm()*(GetNPixelsZ()-1.) + 1.);
725         } while(iz<=0 || iz>GetNPixelsZ());
726         signal = -1.;
727         fMapA2->SetHit(iz,ix,signal);
728     } // end loop on masked pixels
729 }
730 //______________________________________________________________________
731 void AliITSsimulationSPD::CreateHistograms() {
732     // Create Histograms
733     Int_t i;
734
735     fHis=new TObjArray(GetNPixelsZ());
736     for(i=0;i<GetNPixelsZ();i++) {
737         TString spdname("spd_");
738         Char_t candnum[4];
739         sprintf(candnum,"%d",i+1);
740         spdname.Append(candnum);
741         (*fHis)[i] = new TH1F(spdname.Data(),"SPD maps",
742                               GetNPixelsX(),0.,(Float_t) GetNPixelsX());
743     } // end for i
744 }
745 //______________________________________________________________________
746 void AliITSsimulationSPD::ResetHistograms() {
747     // Reset histograms for this detector
748     Int_t i;
749
750     for(i=0;i<GetNPixelsZ();i++ ) {
751         if ((*fHis)[i])    ((TH1F*)(*fHis)[i])->Reset();
752     } // end for i
753 }
754 //______________________________________________________________________
755 void AliITSsimulationSPD::WriteSDigits(AliITSpList *pList){
756     // Fills the Summable digits Tree
757     Int_t i,ni,j,nj;
758     static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
759
760     pList->GetMaxMapIndex(ni,nj);
761     for(i=0;i<ni;i++)for(j=0;j<nj;j++){
762         if(pList->GetSignalOnly(i,j)>0.0){
763             aliITS->AddSumDigit(*(pList->GetpListItem(i,j)));
764 #ifdef DEBUG
765             cout << "pListSPD: " << *(pList->GetpListItem(i,j)) << endl;
766             cout << " CreateSPDSDigit mod=" << fModule << " r,c=";
767             cout << i  <<","<< j << " sig=" << fpList->GetSignalOnly(i,j);
768             cout << " noise=" << fpList->GetNoise(i,j) <<endl;
769 #endif
770         } // end if
771     } // end for i,j
772     return;
773 }
774 //______________________________________________________________________
775 void AliITSsimulationSPD::FillMapFrompList(AliITSpList *pList){
776     // Fills fMap2A from the pList of Summable digits
777     Int_t ix,iz;
778
779     for(iz=0;iz<GetNPixelsZ();iz++)for(ix=0;ix<GetNPixelsX();ix++) 
780         fMapA2->AddSignal(iz,ix,pList->GetSignal(iz,ix));
781     return;
782 }