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