1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
24 #include <TParticle.h>
28 #include "AliITShit.h"
29 #include "AliITSdigit.h"
30 #include "AliITSmodule.h"
31 #include "AliITSMapA2.h"
32 #include "AliITSpList.h"
33 #include "AliITSsimulationSPD.h"
34 #include "AliITSsegmentation.h"
35 #include "AliITSresponse.h"
36 #include "AliITSsegmentationSPD.h"
37 #include "AliITSresponseSPD.h"
40 ClassImp(AliITSsimulationSPD)
41 ////////////////////////////////////////////////////////////////////////
43 // Written by Rocco Caliandro
44 // from a model developed with T. Virgili and R.A. Fini
47 // AliITSsimulationSPD is the simulation of SPDs
49 //______________________________________________________________________
50 AliITSsimulationSPD::AliITSsimulationSPD(){
51 // Default constructor
63 //______________________________________________________________________
64 AliITSsimulationSPD::AliITSsimulationSPD(AliITSsegmentation *seg,
65 AliITSresponse *resp) {
66 // Standard constructor
77 Init((AliITSsegmentationSPD*)seg,(AliITSresponseSPD*)resp);
79 //______________________________________________________________________
80 void AliITSsimulationSPD::Init(AliITSsegmentationSPD *seg,
81 AliITSresponseSPD *resp) {
82 // Initilizes the variables of AliITSsimulation SPD.
87 fMapA2 = new AliITSMapA2(fSegmentation);
89 fResponse->Thresholds(fThresh,fSigma);
90 fResponse->GetNoiseParam(fCouplCol,fCouplRow);
91 fNPixelsZ = fSegmentation->Npz();
92 fNPixelsX = fSegmentation->Npx();
95 //______________________________________________________________________
96 AliITSsimulationSPD::~AliITSsimulationSPD() {
106 //______________________________________________________________________
107 AliITSsimulationSPD::AliITSsimulationSPD(const AliITSsimulationSPD &source){
110 if(&source == this) return;
112 this->fMapA2 = source.fMapA2;
113 this->fHis = source.fHis;
115 this->fThresh = source.fThresh;
116 this->fSigma = source.fSigma;
117 this->fCouplCol = source.fCouplCol;
118 this->fCouplRow = source.fCouplRow;
119 this->fNPixelsX = source.fNPixelsX;
120 this->fNPixelsZ = source.fNPixelsZ;
124 //______________________________________________________________________
125 AliITSsimulationSPD& AliITSsimulationSPD::operator=(const AliITSsimulationSPD
127 // Assignment operator
129 if(&source == this) return *this;
131 this->fMapA2 = source.fMapA2;
132 this->fHis = source.fHis;
134 this->fThresh = source.fThresh;
135 this->fSigma = source.fSigma;
136 this->fCouplCol = source.fCouplCol;
137 this->fCouplRow = source.fCouplRow;
138 this->fNPixelsX = source.fNPixelsX;
139 this->fNPixelsZ = source.fNPixelsZ;
143 //______________________________________________________________________
144 void AliITSsimulationSPD::SDigitiseModule(AliITSmodule *mod, Int_t dummy0,
146 // Sum digitize module
147 if (!(mod->GetNhits())) return; // if module has no hits then no Sdigits.
148 Int_t number = 10000;
149 Int_t *frowpixel = new Int_t[number];
150 Int_t *fcolpixel = new Int_t[number];
151 Double_t *fenepixel = new Double_t[number];
153 // Array of pointers to store the track index of the digits
154 // leave +1, otherwise pList crashes when col=256, row=192
155 AliITSpList *pList = new AliITSpList(GetNPixelsZ()+1,GetNPixelsX()+1);
157 HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,pList);
168 //______________________________________________________________________
169 void AliITSsimulationSPD::DigitiseModule(AliITSmodule *mod, Int_t dummy0,
171 // digitize module. Also need to digitize modules with only noise.
173 Int_t number = 10000;
174 Int_t *frowpixel = new Int_t[number];
175 Int_t *fcolpixel = new Int_t[number];
176 Double_t *fenepixel = new Double_t[number];
177 Int_t module = mod->GetIndex();
179 // Array of pointers to store the track index of the digits
180 // leave +1, otherwise pList crashes when col=256, row=192
181 AliITSpList *pList = new AliITSpList(GetNPixelsZ()+1,GetNPixelsX()+1);
184 SetFluctuations(pList,module);
186 HitsToAnalogDigits(mod,frowpixel,fcolpixel,fenepixel,pList);
188 // apply mask to SPD module
191 CreateDigit(module,pList);
200 //______________________________________________________________________
201 void AliITSsimulationSPD::SDigitsToDigits(Int_t module,AliITSpList *pList) {
202 // sum digits to Digits.
204 FillMapFrompList(pList);
207 SetFluctuations(pList,module);
209 // apply mask to SPD module
212 CreateDigit(module,pList);
214 //______________________________________________________________________
215 void AliITSsimulationSPD::UpdateMapSignal(Int_t row,Int_t col,Int_t trk,
216 Int_t hit,Int_t mod,Double_t ene,
217 AliITSpList *pList) {
218 // updates the Map of signal, adding the energy (ene) released by
221 fMapA2->AddSignal(row,col,ene);
222 pList->AddSignal(row,col,trk,hit,mod,ene);
224 //______________________________________________________________________
225 void AliITSsimulationSPD::UpdateMapNoise(Int_t row,Int_t col,Int_t mod,
226 Double_t ene,AliITSpList *pList) {
227 // updates the Map of noise, adding the energy (ene) give my noise
229 fMapA2->AddSignal(row,col,ene);
230 pList->AddNoise(row,col,mod,ene);
232 //______________________________________________________________________
233 void AliITSsimulationSPD::HitsToAnalogDigits(AliITSmodule *mod,
234 Int_t *frowpixel,Int_t *fcolpixel,
236 AliITSpList *pList) {
237 // Loops over all hits to produce Analog/floting point digits. This
238 // is also the first task in producing standard digits.
240 // loop over hits in the module
241 Int_t hitpos,nhits = mod->GetNhits();
242 for (hitpos=0;hitpos<nhits;hitpos++) {
243 HitToDigit(mod,hitpos,frowpixel,fcolpixel,fenepixel,pList);
244 }// end loop over digits
246 //______________________________________________________________________
247 void AliITSsimulationSPD::HitToDigit(AliITSmodule *mod,Int_t hitpos,
248 Int_t *frowpixel,Int_t *fcolpixel,
249 Double_t *fenepixel,AliITSpList *pList) {
250 // Steering function to determine the digits associated to a given
252 // The digits are created by charge sharing (ChargeSharing) and by
253 // capacitive coupling (SetCoupling). At all the created digits is
254 // associated the track number of the hit (ntrack)
255 Double_t x1l=0.0,y1l=0.0,z1l=0.0,x2l=0.0,y2l=0.0,z2l=0.0;
256 Int_t r1,r2,c1,c2,row,col,npixel = 0;
258 Double_t ene=0.0,etot=0.0;
259 const Float_t kconv = 10000.; // cm -> microns
260 const Float_t kconv1= 0.277e9; // GeV -> electrons equivalent
262 if(!(mod->LineSegmentL(hitpos,x1l,x2l,y1l,y2l,z1l,z2l,etot,ntrack)))return;
264 x2l += x1l; y2l += y1l; z2l += z1l; // Convert to ending coordinate.
265 // positions shifted and converted in microns
266 x1l = x1l*kconv + fSegmentation->Dx()/2.;
267 z1l = z1l*kconv + fSegmentation->Dz()/2.;
268 // positions shifted and converted in microns
269 x2l = x2l*kconv + fSegmentation->Dx()/2.;
270 z2l = z2l*kconv + fSegmentation->Dz()/2.;
271 etot *= kconv1; // convert from GeV to electrons equivalent.
272 Int_t module = mod->GetIndex();
274 // to account for the effective sensitive area
275 // introduced in geometry
276 if (z1l<0 || z1l>fSegmentation->Dz()) return;
277 if (z2l<0 || z2l>fSegmentation->Dz()) return;
278 if (x1l<0 || x1l>fSegmentation->Dx()) return;
279 if (x2l<0 || x2l>fSegmentation->Dx()) return;
281 //Get the col and row number starting from 1
282 // the x direction is not inverted for the second layer!!!
283 fSegmentation->GetPadIxz(x1l, z1l, c1, r1);
284 fSegmentation->GetPadIxz(x2l, z2l, c2, r2);
286 // to account for unexpected equal entrance and
288 if (x1l==x2l) x2l=x2l+x2l*0.000001;
289 if (z1l==z2l) z2l=z2l+z2l*0.000001;
291 if ((r1==r2) && (c1==c2)){
294 frowpixel[npixel-1] = r1;
295 fcolpixel[npixel-1] = c1;
296 fenepixel[npixel-1] = etot;
299 ChargeSharing(x1l,z1l,x2l,z2l,c1,r1,c2,r2,etot,
300 npixel,frowpixel,fcolpixel,fenepixel);
301 } // end if r1==r2 && c1==c2.
303 for (Int_t npix=0;npix<npixel;npix++){
304 row = frowpixel[npix];
305 col = fcolpixel[npix];
306 ene = fenepixel[npix];
307 UpdateMapSignal(row,col,ntrack,hitpos,module,ene,pList);
308 // Starting capacitive coupling effect
309 SetCoupling(row,col,ntrack,hitpos,module,pList);
312 //______________________________________________________________________
313 void AliITSsimulationSPD::ChargeSharing(Float_t x1l,Float_t z1l,Float_t x2l,
314 Float_t z2l,Int_t c1,Int_t r1,Int_t c2,
315 Int_t r2,Float_t etot,
316 Int_t &npixel,Int_t *frowpixel,
317 Int_t *fcolpixel,Double_t *fenepixel){
318 // Take into account the geometrical charge sharing when the track
319 // crosses more than one pixel.
323 <img src="picts/ITS/barimodel_2.gif">
326 <font size=+2 color=red>
327 <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
332 Float_t xa,za,xb,zb,dx,dz,dtot,dm,refr,refm,refc;
334 Float_t arefm, arefr, arefn, arefc, azb, az2l, axb, ax2l;
335 Int_t dirx,dirz,rb,cb;
336 Int_t flag,flagrow,flagcol;
342 dx = TMath::Abs(x1l-x2l);
343 dz = TMath::Abs(z1l-z2l);
344 dtot = TMath::Sqrt((dx*dx)+(dz*dz));
345 dm = (x2l - x1l) / (z2l - z1l);
346 dirx = (Int_t) ((x2l - x1l) / dx);
347 dirz = (Int_t) ((z2l - z1l) / dz);
349 // calculate the x coordinate of the pixel in the next column
350 // and the z coordinate of the pixel in the next row
353 fSegmentation->GetPadCxz(c1, r1-1, xpos, zpos);
355 Float_t xsize = fSegmentation->Dpx(0);
356 Float_t zsize = fSegmentation->Dpz(r1-1);
358 if (dirx == 1) refr = xpos+xsize/2.;
359 else refr = xpos-xsize/2.;
361 if (dirz == 1) refn = zpos+zsize/2.;
362 else refn = zpos-zsize/2.;
368 // calculate the x coordinate of the intersection with the pixel
369 // in the next cell in row direction
370 refm = (refn - z1l)*dm + x1l;
372 // calculate the z coordinate of the intersection with the pixel
373 // in the next cell in column direction
374 refc = (refr - x1l)/dm + z1l;
381 if ((arefm < arefr) && (arefn < arefc)){
382 // the track goes in the pixel in the next cell in row direction
389 if (rb == r2) flagrow=1;
394 // shift to the pixel in the next cell in row direction
395 Float_t zsizeNext = fSegmentation->Dpz(rb-1);
396 //to account for cell at the borders of the detector
397 if(zsizeNext==0) zsizeNext = zsize;
398 refn += zsizeNext*dirz;
400 // the track goes in the pixel in the next cell in column direction
407 if (cb == c2) flagcol=1;
411 } // end ifaxb > ax2l
413 // shift to the pixel in the next cell in column direction
414 Float_t xsizeNext = fSegmentation->Dpx(cb-1);
415 //to account for cell at the borders of the detector
416 if(xsizeNext==0) xsizeNext = xsize;
417 refr += xsizeNext*dirx;
418 } // end if (arefm < arefr) && (arefn < arefc)
420 //calculate the energy lost in the crossed pixel
421 epar = TMath::Sqrt((xb-xa)*(xb-xa)+(zb-za)*(zb-za));
422 epar = etot*(epar/dtot);
424 //store row, column and energy lost in the crossed pixel
425 frowpixel[npixel] = r1;
426 fcolpixel[npixel] = c1;
427 fenepixel[npixel] = epar;
430 // the exit point of the track is reached
431 if (epar == 0) flag = 1;
432 if ((r1 == r2) && (c1 == c2)) flag = 1;
441 //______________________________________________________________________
442 void AliITSsimulationSPD::SetCoupling(Int_t row, Int_t col, Int_t ntrack,
443 Int_t idhit,Int_t module,
444 AliITSpList *pList) {
445 // Take into account the coupling between adiacent pixels.
446 // The parameters probcol and probrow are the fractions of the
447 // signal in one pixel shared in the two adjacent pixels along
448 // the column and row direction, respectively.
452 <img src="picts/ITS/barimodel_3.gif">
455 <font size=+2 color=red>
456 <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
462 Double_t pulse1,pulse2;
463 Float_t couplR=0.0,couplC=0.0;
465 GetCouplings(couplR,couplC);
468 pulse1 = fMapA2->GetSignal(row,col);
470 for (Int_t isign=-1;isign<=1;isign+=2){// loop in row direction
474 if ((j1<0) || (j1>GetNPixelsZ()-1) || (pulse1<GetThreshold())){
475 pulse1 = fMapA2->GetSignal(row,col);
479 UpdateMapSignal(j1,col,ntrack,idhit,module,pulse1,pList);
483 // loop in column direction
487 if ((j2<0) || (j2>(GetNPixelsX()-1)) || (pulse2<GetThreshold())){
488 pulse2 = fMapA2->GetSignal(row,col);
492 UpdateMapSignal(row,j2,ntrack,idhit,module,pulse2,pList);
498 //______________________________________________________________________
499 void AliITSsimulationSPD::CreateDigit(Int_t module,AliITSpList *pList) {
500 // The pixels are fired if the energy deposited inside them is above
501 // the threshold parameter ethr. Fired pixed are interpreted as digits
502 // and stored in the file digitfilename. One also needs to write out
503 // cases when there is only noise (nhits==0).
505 static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
513 for (Int_t r=1;r<=GetNPixelsZ();r++) {
514 for (Int_t c=1;c<=GetNPixelsX();c++) {
515 // check if the deposited energy in a pixel is above the
517 Float_t signal = (Float_t) fMapA2->GetSignal(r,c);
518 if ( signal > GetThreshold()) {
519 digits[0] = r-1; // digits starts from 0
520 digits[1] = c-1; // digits starts from 0
522 digits[2] = (Int_t) signal; // the signal is stored in
525 tracks[j1] = pList->GetTrack(r,c,j1);
526 hits[j1] = pList->GetHit(r,c,j1);
530 aliITS->AddSimDigit(0,phys,digits,tracks,hits,charges);
531 } // end if of threshold condition
535 //______________________________________________________________________
536 void AliITSsimulationSPD::SetFluctuations(AliITSpList *pList,Int_t module) {
537 // Set the electronic noise and threshold non-uniformities to all the
538 // pixels in a detector.
539 // The parameter fSigma is the squared sum of the sigma due to noise
540 // and the sigma of the threshold distribution among pixels.
544 <img src="picts/ITS/barimodel_1.gif">
547 <font size=+2 color=red>
548 <a href="mailto:Rocco.Caliandro@ba.infn.it"></a>.
553 Float_t thr=0.0,sigm=0.0;
554 Double_t signal,sigma;
557 GetThresholds(thr,sigm);
558 sigma = (Double_t) sigm;
559 for(iz=1;iz<=GetNPixelsZ();iz++){
560 for(ix=1;ix<=GetNPixelsX();ix++){
561 signal = sigma*gRandom->Gaus();
562 fMapA2->SetHit(iz,ix,signal);
563 // insert in the label-signal-hit list the pixels fired
565 pList->AddNoise(iz,ix,module,signal);
566 } // end of loop on pixels
567 } // end of loop on pixels
569 //______________________________________________________________________
570 void AliITSsimulationSPD::SetMask() {
571 // Apply a mask to the SPD module. 1% of the pixel channels are
572 // masked. When the database will be ready, the masked pixels
573 // should be read from it.
577 Float_t perc = ((AliITSresponseSPD*)fResponse)->GetFractionDead();
578 // in this way we get the same set of random numbers for all runs.
579 // This is a cluge for now.
580 static TRandom *rnd = new TRandom();
582 totMask= perc*GetNPixelsZ()*GetNPixelsX();
583 for(im=1;im<totMask;im++){
585 ix=(Int_t)(rnd->Rndm()*(GetNPixelsX()-1.) + 1.);
586 } while(ix<=0 || ix>GetNPixelsX());
588 iz=(Int_t)(rnd->Rndm()*(GetNPixelsZ()-1.) + 1.);
589 } while(iz<=0 || iz>GetNPixelsZ());
591 fMapA2->SetHit(iz,ix,signal);
592 } // end loop on masked pixels
594 //______________________________________________________________________
595 void AliITSsimulationSPD::CreateHistograms() {
599 fHis=new TObjArray(GetNPixelsZ());
600 for(i=0;i<GetNPixelsZ();i++) {
601 TString spdname("spd_");
603 sprintf(candnum,"%d",i+1);
604 spdname.Append(candnum);
605 (*fHis)[i] = new TH1F(spdname.Data(),"SPD maps",
606 GetNPixelsX(),0.,(Float_t) GetNPixelsX());
609 //______________________________________________________________________
610 void AliITSsimulationSPD::ResetHistograms() {
611 // Reset histograms for this detector
614 for(i=0;i<GetNPixelsZ();i++ ) {
615 if ((*fHis)[i]) ((TH1F*)(*fHis)[i])->Reset();
618 //______________________________________________________________________
619 void AliITSsimulationSPD::WriteSDigits(AliITSpList *pList){
620 // Fills the Summable digits Tree
622 static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
624 pList->GetMaxMapIndex(ni,nj);
625 for(i=0;i<ni;i++)for(j=0;j<nj;j++){
626 if(pList->GetSignalOnly(i,j)>0.0){
627 aliITS->AddSumDigit(*(pList->GetpListItem(i,j)));
628 // cout << "pListSPD: " << *(pList->GetpListItem(i,j)) << endl;
633 //______________________________________________________________________
634 void AliITSsimulationSPD::FillMapFrompList(AliITSpList *pList){
635 // Fills fMap2A from the pList of Summable digits
638 for(k=0;k<GetNPixelsZ();k++)for(ix=0;ix<GetNPixelsX();ix++)
639 fMapA2->AddSignal(k,ix,pList->GetSignal(k,ix));