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