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