| 409f8c84 |
1 | #include <iostream.h> |
| 2 | #include <TRandom.h> |
| 3 | #include <TH1.h> |
| 4 | #include <TMath.h> |
| 5 | #include <TString.h> |
| 6 | #include <TParticle.h> |
| 7 | |
| 8 | |
| 9 | #include "AliRun.h" |
| 10 | #include "AliITS.h" |
| 11 | #include "AliITShit.h" |
| 12 | #include "AliITSdigit.h" |
| 13 | #include "AliITSmodule.h" |
| 14 | #include "AliITSMapA2.h" |
| 2cc6b29a |
15 | #include "AliITSpList.h" |
| 409f8c84 |
16 | #include "AliITSsimulationSPDdubna.h" |
| 17 | #include "AliITSsegmentation.h" |
| 18 | #include "AliITSresponse.h" |
| 19 | |
| 20 | |
| 21 | |
| 22 | |
| 23 | ClassImp(AliITSsimulationSPDdubna) |
| 24 | //////////////////////////////////////////////////////////////////////// |
| 25 | // Version: 0 |
| 26 | // Written by Boris Batyunya |
| 27 | // December 20 1999 |
| 28 | // |
| 29 | // AliITSsimulationSPDdubna is the simulation of SPDs |
| 2cc6b29a |
30 | //______________________________________________________________________ |
| 31 | |
| 32 | |
| 33 | AliITSsimulationSPDdubna::AliITSsimulationSPDdubna(){ |
| 34 | // constructor |
| 35 | |
| 36 | fResponse = 0; |
| 37 | fSegmentation = 0; |
| 38 | fMapA2 = 0; |
| 39 | fpList = 0; |
| 40 | fModule = 0; |
| 41 | fEvent = 0; |
| 42 | fHis = 0; |
| 43 | fNoise = 0.; |
| 44 | fBaseline = 0.; |
| 45 | fNPixelsZ = 0; |
| 46 | fNPixelsX = 0; |
| 409f8c84 |
47 | } |
| 2cc6b29a |
48 | //______________________________________________________________________ |
| 49 | AliITSsimulationSPDdubna::AliITSsimulationSPDdubna(AliITSsegmentation *seg, |
| 50 | AliITSresponse *resp){ |
| 51 | // standard constructor |
| 409f8c84 |
52 | |
| 2cc6b29a |
53 | fHis = 0; |
| 54 | fResponse = resp; |
| 55 | fSegmentation = seg; |
| 56 | fModule = 0; |
| 57 | fEvent = 0; |
| 409f8c84 |
58 | |
| 2cc6b29a |
59 | fNPixelsZ=fSegmentation->Npz(); |
| 60 | fNPixelsX=fSegmentation->Npx(); |
| 409f8c84 |
61 | |
| 2cc6b29a |
62 | fResponse->GetNoiseParam(fNoise,fBaseline); |
| 409f8c84 |
63 | |
| 2cc6b29a |
64 | fMapA2 = new AliITSMapA2(fSegmentation); |
| 409f8c84 |
65 | |
| 2cc6b29a |
66 | fpList = new AliITSpList(fNPixelsZ+1,fNPixelsX+1); |
| 409f8c84 |
67 | |
| 2cc6b29a |
68 | } |
| 69 | //______________________________________________________________________ |
| 70 | AliITSsimulationSPDdubna::~AliITSsimulationSPDdubna(){ |
| 71 | // destructor |
| 409f8c84 |
72 | |
| 2cc6b29a |
73 | delete fMapA2; |
| 409f8c84 |
74 | |
| 2cc6b29a |
75 | if (fHis) { |
| 76 | fHis->Delete(); |
| 77 | delete fHis; |
| 78 | } // end if fHis |
| 409f8c84 |
79 | } |
| 2cc6b29a |
80 | //______________________________________________________________________ |
| 81 | AliITSsimulationSPDdubna::AliITSsimulationSPDdubna(const |
| 82 | AliITSsimulationSPDdubna |
| 83 | &source){ |
| 84 | // Copy Constructor |
| 85 | if(&source == this) return; |
| 86 | this->fMapA2 = source.fMapA2; |
| 87 | this->fNoise = source.fNoise; |
| 88 | this->fBaseline = source.fBaseline; |
| 89 | this->fNPixelsX = source.fNPixelsX; |
| 90 | this->fNPixelsZ = source.fNPixelsZ; |
| 91 | this->fHis = source.fHis; |
| 92 | return; |
| 93 | } |
| 94 | //______________________________________________________________________ |
| 95 | AliITSsimulationSPDdubna& AliITSsimulationSPDdubna::operator=(const |
| 96 | AliITSsimulationSPDdubna &source){ |
| 97 | // Assignment operator |
| 98 | if(&source == this) return *this; |
| 99 | this->fMapA2 = source.fMapA2; |
| 100 | this->fNoise = source.fNoise; |
| 101 | this->fBaseline = source.fBaseline; |
| 102 | this->fNPixelsX = source.fNPixelsX; |
| 103 | this->fNPixelsZ = source.fNPixelsZ; |
| 104 | this->fHis = source.fHis; |
| 105 | return *this; |
| 106 | } |
| 107 | //______________________________________________________________________ |
| 108 | void AliITSsimulationSPDdubna::InitSimulationModule(Int_t module, Int_t event){ |
| 109 | // This function creates maps to build the list of tracks for each |
| 110 | // summable digit. |
| 111 | // |
| 112 | // Inputs: |
| 113 | // Int_t module // Module number to be simulated |
| 114 | // Int_t event // Event number to be simulated |
| 115 | // |
| 116 | // Outputs: |
| 117 | // none |
| 118 | // |
| 119 | // Returns: |
| 120 | // none |
| 409f8c84 |
121 | |
| 2cc6b29a |
122 | fModule = module; |
| 123 | fEvent = event; |
| 124 | fMapA2->ClearMap(); |
| 125 | fpList->ClearMap(); |
| 126 | } |
| 127 | //_____________________________________________________________________ |
| 128 | void AliITSsimulationSPDdubna::SDigitiseModule(AliITSmodule *mod, Int_t mask, |
| 129 | Int_t event){ |
| 130 | // This function begins the work of creating S-Digits |
| 131 | // |
| 132 | // Inputs: |
| 133 | // AliITSmodule *mod // module |
| 134 | // Int_t mask // mask to be applied to the module |
| 135 | // |
| 136 | // Outputs: |
| 137 | // none |
| 138 | // |
| 139 | // Return: |
| 140 | // test // test returns kTRUE if the module contained hits |
| 141 | // // test returns kFALSE if it did not contain hits |
| 409f8c84 |
142 | |
| 2cc6b29a |
143 | Int_t module = 0; |
| 409f8c84 |
144 | |
| 2cc6b29a |
145 | if(!(mod->GetNhits())) return;// if module has no hits don't create Sdigits |
| 146 | fModule = mod->GetIndex(); |
| 147 | HitToSDigit(mod, module, mask, fpList); |
| 148 | WriteSDigits(fpList); |
| 149 | fMapA2->ClearMap(); |
| 150 | fpList->ClearMap(); |
| 409f8c84 |
151 | } |
| 2cc6b29a |
152 | //______________________________________________________________________ |
| 153 | void AliITSsimulationSPDdubna::WriteSDigits(AliITSpList *pList){ |
| 154 | // This function adds each S-Digit to pList |
| 155 | // |
| 156 | // Inputs: |
| 157 | // AliITSpList *pList |
| 158 | // |
| 159 | // Outputs: |
| 160 | // none |
| 161 | // |
| 162 | // Return: |
| 163 | // none |
| 164 | Int_t i, ni, j, nj; |
| 165 | static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS"); |
| 166 | |
| 167 | pList->GetMaxMapIndex(ni, nj); |
| 168 | for(i=0; i<ni; i++)for(j=0; j<nj; j++){ |
| 169 | if(pList->GetSignalOnly(i, j)>0.0){ |
| 170 | aliITS->AddSumDigit(*(pList->GetpListItem(i, j))); |
| 171 | } // end if pList |
| 172 | } // end for i,j |
| 173 | return; |
| 174 | } |
| 175 | //______________________________________________________________________ |
| 176 | void AliITSsimulationSPDdubna::FinishSDigitiseModule(){ |
| 177 | // This function calls SDigitsToDigits which creates Digits from SDigits |
| 178 | // |
| 179 | // Inputs: |
| 180 | // none |
| 181 | // |
| 182 | // Outputs: |
| 183 | // none |
| 184 | // Return |
| 185 | // none |
| 409f8c84 |
186 | |
| 2cc6b29a |
187 | SDigitsToDigits(fModule, fpList); |
| 188 | return; |
| 189 | } |
| 190 | //______________________________________________________________________ |
| 191 | void AliITSsimulationSPDdubna::SDigitsToDigits(Int_t module, |
| 192 | AliITSpList *pList){ |
| 193 | // This function adds electronic noise to the S-Digits and then adds them |
| 194 | // to a new pList |
| 195 | // |
| 196 | // Inputs: |
| 197 | // Int_t module // module number |
| 198 | // AliITSpList *pList // pList |
| 199 | // |
| 200 | // Outputs: |
| 201 | // pList is passed along to the functions ChargeToSignal and GetList |
| 202 | // |
| 203 | // Return: |
| 204 | // none |
| 409f8c84 |
205 | |
| 2cc6b29a |
206 | fModule = module; |
| 207 | ChargeToSignal(pList); // Charge To Signal both adds noise and |
| 208 | fMapA2->ClearMap(); |
| 209 | pList->ClearMap(); |
| 409f8c84 |
210 | } |
| 2cc6b29a |
211 | //______________________________________________________________________ |
| 212 | void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, |
| 213 | Int_t dummy){ |
| 214 | // This function creates Digits straight from the hits and then adds |
| 215 | // electronic noise to the digits before adding them to pList |
| 216 | // |
| 217 | // Inputs: |
| 218 | // AliITSmodule *mod // module |
| 219 | // Int_t module // module number Dummy. |
| 220 | // Int_t dummy |
| 221 | // |
| 222 | // Outputs: |
| 223 | // Each of the input variables is passed along to HitToSDigit |
| 224 | // |
| 225 | // Return: |
| 226 | // none |
| 409f8c84 |
227 | |
| 2cc6b29a |
228 | fModule = mod->GetIndex(); //This calls the module for HitToSDigit |
| 229 | HitToSDigit(mod,fModule, dummy, fpList); |
| 230 | ChargeToSignal(fpList); |
| 231 | fMapA2->ClearMap(); |
| 232 | fpList->ClearMap(); |
| 233 | } |
| 234 | //______________________________________________________________________ |
| 235 | void AliITSsimulationSPDdubna::UpdateMapSignal(Int_t i, Int_t j, Int_t trk, |
| 236 | Int_t ht, Int_t module, |
| 237 | Double_t signal, |
| 238 | AliITSpList *pList){ |
| 239 | // This function adds a signal to the pList from the pList class |
| 240 | // |
| 241 | // Inputs: |
| 242 | // Int_t i // row number |
| 243 | // Int_t j // column number |
| 244 | // Int_t trk // track number |
| 245 | // Int_t ht // hit number |
| 246 | // Double_t signal // signal strength |
| 247 | // AliITSpList *pList // pList |
| 248 | // |
| 249 | // Outputs: |
| 250 | // All of the inputs are passed to AliITSpList::AddSignal |
| 251 | // Int_t ix // row number |
| 252 | // Int_t iz // column number |
| 253 | // Double_t sig // signal strength |
| 254 | // // These three variables are defined to preserve the |
| 255 | // // assignments used in the function AliITSMapA2::AddSignal |
| 256 | // |
| 257 | // Return: |
| 258 | // none |
| 259 | Int_t iz = j; |
| 260 | Int_t ix = i; |
| 261 | Double_t sig = signal; |
| 262 | |
| 263 | fMapA2->AddSignal(iz, ix, sig); |
| 264 | pList->AddSignal(i, j, trk, ht, fModule, signal); |
| 265 | } |
| 266 | //______________________________________________________________________ |
| 267 | void AliITSsimulationSPDdubna::UpdateMapNoise(Int_t i, Int_t j, Int_t ix, |
| 268 | Int_t iz, Int_t fModule, |
| 269 | Double_t sig, Float_t noise, |
| 270 | AliITSpList *pList){ |
| 271 | // This function adds noise to data in the MapA2 as well as the pList |
| 272 | // |
| 273 | // Inputs: |
| 274 | // Int_t i == ix // row number |
| 275 | // Int_t j == iz // column number |
| 276 | // Int_t mod // module number |
| 277 | // Double_t sig // signal strength |
| 278 | // Double_t noise // electronic noise generated by ChargeToSignal |
| 279 | // AliITSpList *pList // pList |
| 280 | // |
| 281 | // Outputs: |
| 282 | // All of the inputs are passed to AliITSMapA2::AddSignal or |
| 283 | // AliITSpList::AddNoise |
| 284 | // |
| 285 | // Return: |
| 286 | // none |
| 409f8c84 |
287 | |
| 2cc6b29a |
288 | fMapA2->AddSignal(iz, ix, sig); |
| 289 | pList->AddNoise(i, j, fModule, noise); |
| 290 | } |
| 291 | //______________________________________________________________________ |
| 292 | void AliITSsimulationSPDdubna::HitToDigit(AliITSmodule *mod, Int_t module, |
| 293 | Int_t dummy){ |
| 294 | DigitiseModule(mod, module, dummy); |
| 295 | } |
| 296 | //______________________________________________________________________ |
| 297 | void AliITSsimulationSPDdubna::HitToSDigit(AliITSmodule *mod, Int_t module, |
| 298 | Int_t dummy, AliITSpList *pList){ |
| 299 | // digitize module |
| 409f8c84 |
300 | const Float_t kEnToEl = 2.778e+8; // GeV->charge in electrons |
| 301 | // for 3.6 eV/pair |
| 302 | const Float_t kconv = 10000.; // cm -> microns |
| 303 | |
| 304 | Float_t spdLength = fSegmentation->Dz(); |
| 305 | Float_t spdWidth = fSegmentation->Dx(); |
| 338e4f06 |
306 | Float_t spdThickness = fSegmentation->Dy(); |
| 409f8c84 |
307 | Float_t difCoef, dum; |
| 308 | fResponse->DiffCoeff(difCoef,dum); |
| 338e4f06 |
309 | if(spdThickness > 290) difCoef = 0.00613; |
| 409f8c84 |
310 | |
| 311 | Float_t zPix0 = 1e+6; |
| 312 | Float_t xPix0 = 1e+6; |
| 313 | Float_t yPrev = 1e+6; |
| 314 | |
| 315 | Float_t zPitch = fSegmentation->Dpz(0); |
| 316 | Float_t xPitch = fSegmentation->Dpx(0); |
| 317 | |
| 318 | TObjArray *fHits = mod->GetHits(); |
| 2cc6b29a |
319 | module = mod->GetIndex(); |
| 409f8c84 |
320 | Int_t nhits = fHits->GetEntriesFast(); |
| 321 | if (!nhits) return; |
| 322 | |
| 2cc6b29a |
323 | cout<<"len,wid,thickness,nx,nz,pitchx,pitchz,difcoef ="<<spdLength<<"," |
| 324 | <<spdWidth<<","<<spdThickness<<","<<fNPixelsX<<","<<fNPixelsZ<<"," |
| 325 | <<xPitch<<","<<zPitch<<","<<difCoef<<endl; |
| 326 | // Array of pointers to the label-signal list |
| 409f8c84 |
327 | Int_t indexRange[4] = {0,0,0,0}; |
| 328 | |
| 329 | // Fill detector maps with GEANT hits |
| 330 | // loop over hits in the module |
| 331 | static Bool_t first; |
| 332 | Int_t lasttrack=-2; |
| 333 | Int_t hit, iZi, jz, jx; |
| 338e4f06 |
334 | Int_t idhit=-1; //! |
| 335 | cout<<"SPDdubna: module,nhits ="<<module<<","<<nhits<<endl; |
| 409f8c84 |
336 | for (hit=0;hit<nhits;hit++) { |
| 337 | AliITShit *iHit = (AliITShit*) fHits->At(hit); |
| 338e4f06 |
338 | //Int_t layer = iHit->GetLayer(); |
| 339 | Float_t yPix0 = -spdThickness/2; |
| 340 | |
| 341 | // work with the idtrack=entry number in the TreeH |
| 342 | //Int_t idhit,idtrack; //! |
| 343 | //mod->GetHitTrackAndHitIndex(hit,idtrack,idhit); //! |
| 344 | //Int_t idtrack=mod->GetHitTrackIndex(hit); |
| 345 | // or store straight away the particle position in the array |
| 346 | // of particles : |
| 409f8c84 |
347 | if(iHit->StatusEntering()) idhit=hit; |
| 348 | Int_t itrack = iHit->GetTrack(); |
| 349 | Int_t dray = 0; |
| 350 | |
| 351 | if (lasttrack != itrack || hit==(nhits-1)) first = kTRUE; |
| 352 | |
| 338e4f06 |
353 | //Int_t parent = iHit->GetParticle()->GetFirstMother(); |
| 409f8c84 |
354 | Int_t partcode = iHit->GetParticle()->GetPdgCode(); |
| 355 | |
| 2cc6b29a |
356 | // partcode (pdgCode): 11 - e-, 13 - mu-, 22 - gamma, 111 - pi0, |
| 357 | // 211 - pi+, 310 - K0s, 321 - K+, 2112 - n, 2212 - p, 3122 - lambda |
| 409f8c84 |
358 | |
| 359 | Float_t pmod = iHit->GetParticle()->P(); // total momentum at the |
| 360 | // vertex |
| 361 | pmod *= 1000; |
| 362 | |
| 409f8c84 |
363 | if(partcode == 11 && pmod < 6) dray = 1; // delta ray is e- |
| 364 | // at p < 6 MeV/c |
| 365 | |
| 409f8c84 |
366 | // Get hit z and x(r*phi) cordinates for each module (detector) |
| 367 | // in local system. |
| 368 | |
| 369 | Float_t zPix = kconv*iHit->GetZL(); |
| 370 | Float_t xPix = kconv*iHit->GetXL(); |
| 371 | Float_t yPix = kconv*iHit->GetYL(); |
| 372 | |
| 373 | // Get track status |
| 374 | Int_t status = iHit->GetTrackStatus(); |
| 409f8c84 |
375 | |
| 376 | // Check boundaries |
| 377 | if(zPix > spdLength/2) { |
| 2cc6b29a |
378 | //cout<<"!!! SPD: z outside ="<<zPix<<endl; |
| 379 | zPix = spdLength/2 - 10; |
| 409f8c84 |
380 | } |
| 381 | if(zPix < 0 && zPix < -spdLength/2) { |
| 2cc6b29a |
382 | //cout<<"!!! SPD: z outside ="<<zPix<<endl; |
| 383 | zPix = -spdLength/2 + 10; |
| 409f8c84 |
384 | } |
| 385 | if(xPix > spdWidth/2) { |
| 2cc6b29a |
386 | //cout<<"!!! SPD: x outside ="<<xPix<<endl; |
| 387 | xPix = spdWidth/2 - 10; |
| 409f8c84 |
388 | } |
| 389 | if(xPix < 0 && xPix < -spdWidth/2) { |
| 2cc6b29a |
390 | //cout<<"!!! SPD: x outside ="<<xPix<<endl; |
| 391 | xPix = -spdWidth/2 + 10; |
| 409f8c84 |
392 | } |
| 393 | Int_t trdown = 0; |
| 394 | |
| 395 | // enter Si or after event in Si |
| 396 | if (status == 66 ) { |
| 2cc6b29a |
397 | zPix0 = zPix; |
| 398 | xPix0 = xPix; |
| 399 | yPrev = yPix; |
| 400 | } // end if status == 66 |
| 409f8c84 |
401 | |
| 402 | Float_t depEnergy = iHit->GetIonization(); |
| 403 | // skip if the input point to Si |
| 404 | |
| 405 | if(depEnergy <= 0.) continue; |
| 406 | |
| 407 | // if track returns to the opposite direction: |
| 408 | if (yPix < yPrev) { |
| 409 | trdown = 1; |
| 2cc6b29a |
410 | } // end if yPix < yPrev |
| 409f8c84 |
411 | |
| 412 | // take into account the holes diffusion inside the Silicon |
| 413 | // the straight line between the entrance and exit points in Si is |
| 414 | // divided into the several steps; the diffusion is considered |
| 415 | // for each end point of step and charge |
| 416 | // is distributed between the pixels through the diffusion. |
| 409f8c84 |
417 | |
| 418 | // ---------- the diffusion in Z (beam) direction ------- |
| 409f8c84 |
419 | Float_t charge = depEnergy*kEnToEl; // charge in e- |
| 420 | Float_t drPath = 0.; |
| 421 | Float_t tang = 0.; |
| 422 | Float_t sigmaDif = 0.; |
| 423 | Float_t zdif = zPix - zPix0; |
| 424 | Float_t xdif = xPix - xPix0; |
| 425 | Float_t ydif = TMath::Abs(yPix - yPrev); |
| 426 | Float_t ydif0 = TMath::Abs(yPrev - yPix0); |
| 427 | |
| 428 | if(ydif < 1) continue; // ydif is not zero |
| 429 | |
| 430 | Float_t projDif = sqrt(xdif*xdif + zdif*zdif); |
| 431 | |
| 432 | Int_t ndZ = (Int_t)TMath::Abs(zdif/zPitch) + 1; |
| 433 | Int_t ndX = (Int_t)TMath::Abs(xdif/xPitch) + 1; |
| 434 | |
| 435 | // number of the steps along the track: |
| 436 | Int_t nsteps = ndZ; |
| 437 | if(ndX > ndZ) nsteps = ndX; |
| 338e4f06 |
438 | if(nsteps < 20) nsteps = 20; // minimum number of the steps |
| 409f8c84 |
439 | |
| 440 | if (projDif < 5 ) { |
| 2cc6b29a |
441 | drPath = (yPix-yPix0)*1.e-4; |
| 442 | drPath = TMath::Abs(drPath); // drift path in cm |
| 443 | sigmaDif = difCoef*sqrt(drPath); // sigma diffusion in cm |
| 444 | sigmaDif = sigmaDif*kconv; // sigma diffusion in microns |
| 445 | nsteps = 1; |
| 446 | } // end if projDif < 5 |
| 409f8c84 |
447 | |
| 448 | if(projDif > 5) tang = ydif/projDif; |
| 449 | Float_t dCharge = charge/nsteps; // charge in e- for one step |
| 450 | Float_t dZ = zdif/nsteps; |
| 451 | Float_t dX = xdif/nsteps; |
| 452 | |
| 2cc6b29a |
453 | for (iZi = 1; iZi <= nsteps;iZi++) { |
| 454 | Float_t dZn = iZi*dZ; |
| 409f8c84 |
455 | Float_t dXn = iZi*dX; |
| 456 | Float_t zPixn = zPix0 + dZn; |
| 457 | Float_t xPixn = xPix0 + dXn; |
| 458 | |
| 459 | if(projDif >= 5) { |
| 2cc6b29a |
460 | Float_t dProjn = sqrt(dZn*dZn+dXn*dXn); |
| 461 | drPath = dProjn*tang*1.e-4; // drift path for iZi+1 step in cm |
| 462 | if(trdown == 0) { |
| 463 | drPath = TMath::Abs(drPath) + ydif0*1.e-4; |
| 464 | }// end if trdow ==0 |
| 465 | if(trdown == 1) { |
| 466 | drPath = ydif0*1.e-4 - TMath::Abs(drPath); |
| 467 | drPath = TMath::Abs(drPath); |
| 468 | } // end if trdown == 1 |
| 469 | sigmaDif = difCoef*sqrt(drPath); |
| 470 | sigmaDif = sigmaDif*kconv; // sigma diffusion in microns |
| 471 | } // end if projdif >= 5 |
| 409f8c84 |
472 | |
| 473 | zPixn = (zPixn + spdLength/2.); |
| 474 | xPixn = (xPixn + spdWidth/2.); |
| 475 | Int_t nZpix, nXpix; |
| 476 | fSegmentation->GetPadIxz(xPixn,zPixn,nXpix,nZpix); |
| 477 | zPitch = fSegmentation->Dpz(nZpix); |
| 478 | fSegmentation->GetPadTxz(xPixn,zPixn); |
| 479 | // set the window for the integration |
| 480 | Int_t jzmin = 1; |
| 481 | Int_t jzmax = 3; |
| 482 | if(nZpix == 1) jzmin =2; |
| 483 | if(nZpix == fNPixelsZ) jzmax = 2; |
| 484 | |
| 485 | Int_t jxmin = 1; |
| 486 | Int_t jxmax = 3; |
| 487 | if(nXpix == 1) jxmin =2; |
| 488 | if(nXpix == fNPixelsX) jxmax = 2; |
| 489 | |
| 490 | Float_t zpix = nZpix; |
| 491 | Float_t dZright = zPitch*(zpix - zPixn); |
| 492 | Float_t dZleft = zPitch - dZright; |
| 493 | |
| 494 | Float_t xpix = nXpix; |
| 495 | Float_t dXright = xPitch*(xpix - xPixn); |
| 496 | Float_t dXleft = xPitch - dXright; |
| 497 | |
| 498 | Float_t dZprev = 0.; |
| 499 | Float_t dZnext = 0.; |
| 500 | Float_t dXprev = 0.; |
| 501 | Float_t dXnext = 0.; |
| 502 | |
| 503 | for(jz=jzmin; jz <=jzmax; jz++) { |
| 504 | if(jz == 1) { |
| 2cc6b29a |
505 | dZprev = -zPitch - dZleft; |
| 506 | dZnext = -dZleft; |
| 507 | } else if(jz == 2) { |
| 508 | dZprev = -dZleft; |
| 509 | dZnext = dZright; |
| 510 | } else if(jz == 3) { |
| 511 | dZprev = dZright; |
| 512 | dZnext = dZright + zPitch; |
| 513 | } // end if jz |
| 409f8c84 |
514 | // kz changes from 1 to the fNofPixels(270) |
| 515 | Int_t kz = nZpix + jz -2; |
| 516 | |
| 517 | Float_t zArg1 = dZprev/sigmaDif; |
| 518 | Float_t zArg2 = dZnext/sigmaDif; |
| 519 | Float_t zProb1 = TMath::Erfc(zArg1); |
| 520 | Float_t zProb2 = TMath::Erfc(zArg2); |
| 521 | Float_t dZCharge =0.5*(zProb1-zProb2)*dCharge; |
| 522 | |
| 523 | |
| 524 | // ----------- holes diffusion in X(r*phi) direction -------- |
| 525 | |
| 526 | if(dZCharge > 1.) { |
| 2cc6b29a |
527 | for(jx=jxmin; jx <=jxmax; jx++) { |
| 528 | if(jx == 1) { |
| 529 | dXprev = -xPitch - dXleft; |
| 530 | dXnext = -dXleft; |
| 531 | } else if(jx == 2) { |
| 532 | dXprev = -dXleft; |
| 533 | dXnext = dXright; |
| 534 | } else if(jx == 3) { |
| 535 | dXprev = dXright; |
| 536 | dXnext = dXright + xPitch; |
| 537 | } // end if jx |
| 538 | Int_t kx = nXpix + jx -2; |
| 539 | Float_t xArg1 = dXprev/sigmaDif; |
| 540 | Float_t xArg2 = dXnext/sigmaDif; |
| 541 | Float_t xProb1 = TMath::Erfc(xArg1); |
| 542 | Float_t xProb2 = TMath::Erfc(xArg2); |
| 543 | Float_t dXCharge =0.5*(xProb1-xProb2)*dZCharge; |
| 544 | |
| 545 | if(dXCharge > 1.) { |
| 546 | Int_t index = kz-1; |
| 547 | if (first) { |
| 548 | indexRange[0]=indexRange[1]=index; |
| 549 | indexRange[2]=indexRange[3]=kx-1; |
| 550 | first=kFALSE; |
| 551 | } // end if first |
| 552 | indexRange[0]=TMath::Min(indexRange[0],kz-1); |
| 553 | indexRange[1]=TMath::Max(indexRange[1],kz-1); |
| 554 | indexRange[2]=TMath::Min(indexRange[2],kx-1); |
| 555 | indexRange[3]=TMath::Max(indexRange[3],kx-1); |
| 556 | |
| 557 | // build the list of digits for this module |
| 558 | Double_t signal=fMapA2->GetSignal(index,kx-1); |
| 559 | signal+=dXCharge; |
| 560 | fMapA2->SetHit(index,kx-1,(double)signal); |
| 561 | |
| 562 | // The calling sequence for UpdateMapSignal was |
| 563 | // moved into the (dx > 1 e-) loop because it |
| 564 | // needs to call signal which is defined inside |
| 565 | // this loop |
| 566 | Int_t i = kx-1; |
| 567 | Int_t j = kz-1; |
| 568 | Int_t trk = mod->GetHitTrackIndex(hit); |
| 569 | Int_t ht = hit; |
| 570 | fModule = module;//Defined because functions |
| 571 | // called by UpdateMapSignal |
| 572 | // expect module to be an |
| 573 | // integer |
| 574 | UpdateMapSignal(j,i,trk,ht,fModule,signal,pList); |
| 575 | } // dXCharge > 1 e- |
| 576 | } // jx loop |
| 409f8c84 |
577 | } // dZCharge > 1 e- |
| 578 | } // jz loop |
| 579 | } // iZi loop |
| 409f8c84 |
580 | if (status == 65) { // the step is inside of Si |
| 2cc6b29a |
581 | zPix0 = zPix; |
| 582 | xPix0 = xPix; |
| 583 | } // end if status == 65 |
| 584 | yPrev = yPix; |
| 409f8c84 |
585 | } // hit loop inside the module |
| 409f8c84 |
586 | } |
| 2cc6b29a |
587 | //______________________________________________________________________ |
| 588 | void AliITSsimulationSPDdubna::ChargeToSignal(AliITSpList *pList){ |
| 589 | // add noise and electronics, perform the zero suppression and add the |
| 590 | // digit to the list |
| 591 | |
| 592 | AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS"); |
| 593 | |
| 594 | Float_t threshold = (float)fResponse->MinVal(); |
| 595 | |
| 596 | Int_t digits[3], tracks[3], hits[3], gi, j1; |
| 597 | Float_t charges[3]; |
| 598 | Float_t electronics; |
| 599 | Float_t signal; |
| 600 | Float_t phys; |
| 601 | Double_t sig; |
| 602 | Int_t module = 0; |
| 603 | for(Int_t iz=0; iz<fNPixelsZ; iz++){ |
| 604 | for(Int_t ix=0; ix<fNPixelsX; ix++){ |
| 605 | electronics = fBaseline + fNoise*gRandom->Gaus(); |
| 606 | signal = (float)pList->GetSignalOnly(ix, iz); |
| 607 | sig = Double_t (signal); // sig will be passed along to |
| 608 | // UpdateMapNoise this is necessary so |
| 609 | // that a signal without electronic |
| 610 | // noise is passed along |
| 611 | signal += electronics; |
| 612 | gi =iz*fNPixelsX+ix; // global index |
| 613 | if (signal > threshold) { |
| 614 | digits[0]=iz; |
| 615 | digits[1]=ix; |
| 616 | digits[2]=1; |
| 617 | for(j1=0;j1<3;j1++){ |
| 618 | if (pList->GetTrack(ix, iz, gi)) { |
| 619 | //b.b. tracks[j1]=-3; |
| 620 | tracks[j1] = (Int_t)(pList->GetTrack(ix, iz, j1)+j1); |
| 621 | hits[j1] = (Int_t)(pList->GetHit(ix, iz, j1)+j1+6); |
| 622 | }else { |
| 623 | tracks[j1]=-2; //noise |
| 624 | hits[j1] = -1; |
| 625 | } // end if pList |
| 626 | charges[j1] = 0; |
| 627 | } // end for j1 |
| 628 | |
| 629 | if(tracks[0] == tracks[1] && tracks[0] == tracks[2]) { |
| 630 | tracks[1] = -3; |
| 631 | hits[1] = -1; |
| 632 | tracks[2] = -3; |
| 633 | hits[2] = -1; |
| 634 | } else if(tracks[0] == tracks[1] && tracks[0] != tracks[2]) { |
| 635 | tracks[1] = -3; |
| 636 | hits[1] = -1; |
| 637 | } else if(tracks[0] == tracks[2] && tracks[0] != tracks[1]) { |
| 638 | tracks[2] = -3; |
| 639 | hits[2] = -1; |
| 640 | } else if(tracks[1] == tracks[2] && tracks[0] != tracks[1]) { |
| 641 | tracks[2] = -3; |
| 642 | hits[2] = -1; |
| 643 | } // end if |
| 644 | |
| 645 | phys = 0; |
| 646 | |
| 647 | Int_t i = ix; // These variables are declared so to be |
| 648 | Int_t j = iz; // passed along to UpdateMapNoise and |
| 649 | Float_t noise = electronics; // in that function |
| 650 | UpdateMapNoise(j, i, ix, iz, fModule, sig, noise, pList); |
| 651 | aliITS->AddSimDigit(0, phys, digits, tracks, hits, charges); |
| 652 | } // |
| 653 | } // |
| 654 | } // |
| 409f8c84 |
655 | } |
| 2cc6b29a |
656 | //______________________________________________________________________ |
| 657 | void AliITSsimulationSPDdubna::CreateHistograms(){ |
| 658 | // create 1D histograms for tests |
| 659 | |
| 660 | printf("SPD - create histograms\n"); |
| 661 | |
| 662 | fHis=new TObjArray(fNPixelsZ); |
| 663 | TString spdName("spd_"); |
| 664 | for (Int_t i=0;i<fNPixelsZ;i++) { |
| 665 | Char_t pixelz[4]; |
| 666 | sprintf(pixelz,"%d",i+1); |
| 667 | spdName.Append(pixelz); |
| 668 | //PH (*fHis)[i] = new TH1F(spdName.Data(),"SPD maps", |
| 669 | //PH fNPixelsX,0.,(Float_t) fNPixelsX); |
| 670 | fHis->AddAt(new TH1F(spdName.Data(),"SPD maps", |
| 671 | fNPixelsX,0.,(Float_t) fNPixelsX), i); |
| 672 | } // end for i |
| 409f8c84 |
673 | } |
| 2cc6b29a |
674 | //______________________________________________________________________ |
| 675 | void AliITSsimulationSPDdubna::ResetHistograms(){ |
| 409f8c84 |
676 | // |
| 677 | // Reset histograms for this detector |
| 678 | // |
| 679 | |
| 680 | for ( int i=0;i<fNPixelsZ;i++ ) { |
| 2cc6b29a |
681 | //PH if ((*fHis)[i]) ((TH1F*)(*fHis)[i])->Reset(); |
| 2682e810 |
682 | if (fHis->At(i)) ((TH1F*)fHis->At(i))->Reset(); |
| 2cc6b29a |
683 | } // end for i |
| 409f8c84 |
684 | } |
git://git.uio.no / u / mrichter / AliRoot.git / blame