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