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 **************************************************************************/
18 #include <Riostream.h>
30 #include "AliITSMapA2.h"
31 #include "AliITSRawData.h"
32 #include "AliITSdigitSPD.h"
33 #include "AliITSetfSDD.h"
34 #include "AliITSmodule.h"
35 #include "AliITSpList.h"
36 #include "AliITSresponseSDD.h"
37 #include "AliITSCalibrationSDD.h"
38 #include "AliITSsegmentationSDD.h"
39 #include "AliITSsimulationSDD.h"
43 ClassImp(AliITSsimulationSDD)
44 ////////////////////////////////////////////////////////////////////////
46 // Written by Piergiorgio Cerello //
47 // November 23 1999 //
49 // AliITSsimulationSDD is the simulation of SDDs. //
50 ////////////////////////////////////////////////////////////////////////
52 //______________________________________________________________________
53 Int_t power(Int_t b, Int_t e) {
54 // compute b to the e power, where both b and e are Int_ts.
57 for(i=0; i<e; i++) power *= b;
60 //______________________________________________________________________
61 void FastFourierTransform(AliITSetfSDD *alisddetf,Double_t *real,
62 Double_t *imag,Int_t direction) {
63 // Do a Fast Fourier Transform
65 Int_t samples = alisddetf->GetSamples();
66 Int_t l = (Int_t) ((log((Float_t) samples)/log(2.))+0.5);
69 Int_t m2 = samples/m1;
72 for(j=0; j<samples; j += m1) {
74 for(k=j; k<= j+m-1; k++) {
75 Double_t wsr = alisddetf->GetWeightReal(p);
76 Double_t wsi = alisddetf->GetWeightImag(p);
77 if(direction == -1) wsi = -wsi;
78 Double_t xr = *(real+k+m);
79 Double_t xi = *(imag+k+m);
80 *(real+k+m) = wsr*(*(real+k)-xr) - wsi*(*(imag+k)-xi);
81 *(imag+k+m) = wsr*(*(imag+k)-xi) + wsi*(*(real+k)-xr);
92 for(j=0; j<samples; j++) {
96 for(i1=1; i1<=l; i1++) {
99 p = p + p + j2 - j1 - j1;
102 Double_t xr = *(real+j);
103 Double_t xi = *(imag+j);
104 *(real+j) = *(real+p);
105 *(imag+j) = *(imag+p);
110 if(direction == -1) {
111 for(i=0; i<samples; i++) {
112 *(real+i) /= samples;
113 *(imag+i) /= samples;
115 } // end if direction == -1
118 //______________________________________________________________________
119 AliITSsimulationSDD::AliITSsimulationSDD():
142 fCrosstalkFlag(kFALSE),
147 // Default constructor
149 SetPerpendTracksFlag();
154 //______________________________________________________________________
155 AliITSsimulationSDD::AliITSsimulationSDD(const AliITSsimulationSDD &source) :
156 AliITSsimulation(source),
158 fHitMap2(source.fHitMap2),
159 fHitSigMap2(source.fHitSigMap2),
160 fHitNoiMap2(source.fHitNoiMap2),
161 fStream(source.fStream),
162 fElectronics(source.fElectronics),
165 fOutZR(source.fOutZR),
166 fOutZI(source.fOutZI),
167 fAnodeFire(source.fAnodeFire),
173 fTreeB(source.fTreeB),
174 fParam(source.fParam),
175 fFileName(source.fFileName),
177 fCheckNoise(source.fCheckNoise),
178 fCrosstalkFlag(source.fCrosstalkFlag),
179 fDoFFT(source.fDoFFT),
180 fNofMaps(source.fNofMaps),
181 fMaxNofSamples(source.fMaxNofSamples),
182 fScaleSize(source.fScaleSize){
183 // Copy constructor to satify Coding roules only.
186 //______________________________________________________________________
187 AliITSsimulationSDD& AliITSsimulationSDD::operator=(const AliITSsimulationSDD &src){
188 // Assignment operator to satify Coding roules only.
190 if(this==&src) return *this;
191 Error("AliITSsimulationSDD","Not allowed to make a = with "
192 "AliITSsimulationSDD Using default creater instead");
195 //______________________________________________________________________
196 AliITSsimulation& AliITSsimulationSDD::operator=(const AliITSsimulation &src){
197 // Assignment operator to satify Coding roules only.
199 if(this==&src) return *this;
200 Error("AliITSsimulationSSD","Not allowed to make a = with "
201 "AliITSsimulationSDD Using default creater instead");
205 //______________________________________________________________________
206 AliITSsimulationSDD::AliITSsimulationSDD(AliITSDetTypeSim* dettyp):
207 AliITSsimulation(dettyp),
229 fCrosstalkFlag(kFALSE),
234 // Default Constructor
237 //______________________________________________________________________
238 void AliITSsimulationSDD::Init(){
239 // Standard Constructor
242 SetPerpendTracksFlag();
247 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
249 AliITSresponseSDD* res = (AliITSresponseSDD*)fDetType->GetResponse(1);
250 fpList = new AliITSpList( seg->Npz(),
251 fScaleSize*seg->Npx() );
252 fHitSigMap2 = new AliITSMapA2(seg,fScaleSize,1);
253 fHitNoiMap2 = new AliITSMapA2(seg,fScaleSize,1);
254 fHitMap2 = fHitSigMap2;
256 fNofMaps = seg->Npz();
257 fMaxNofSamples = seg->Npx();
258 fAnodeFire = new Bool_t [fNofMaps];
260 Float_t sddWidth = seg->Dz();
262 Float_t anodePitch = seg->Dpz(dummy);
263 Double_t timeStep = (Double_t)seg->Dpx(dummy);
265 if(anodePitch*(fNofMaps/2) > sddWidth) {
266 Warning("AliITSsimulationSDD",
267 "Too many anodes %d or too big pitch %f \n",
268 fNofMaps/2,anodePitch);
272 fElectronics = new AliITSetfSDD(timeStep/fScaleSize,
275 char opt1[20], opt2[20];
276 res->ParamOptions(opt1,opt2);
279 const char *kopt=res->ZeroSuppOption();
285 Bool_t write = res->OutputOption();
286 if(write && strstr(kopt,"2D")) MakeTreeB();
288 fITS = (AliITS*)gAlice->GetModule("ITS");
289 Int_t size = fNofMaps*fMaxNofSamples;
290 fStream = new AliITSInStream(size);
292 fInZR = new Double_t [fScaleSize*fMaxNofSamples];
293 fInZI = new Double_t [fScaleSize*fMaxNofSamples];
294 fOutZR = new Double_t [fScaleSize*fMaxNofSamples];
295 fOutZI = new Double_t [fScaleSize*fMaxNofSamples];
297 //______________________________________________________________________
298 AliITSsimulationSDD::~AliITSsimulationSDD() {
313 if(fTreeB) delete fTreeB;
314 if(fInZR) delete [] fInZR;
315 if(fInZI) delete [] fInZI;
316 if(fOutZR) delete [] fOutZR;
317 if(fOutZI) delete [] fOutZI;
318 if(fAnodeFire) delete [] fAnodeFire;
320 //______________________________________________________________________
321 void AliITSsimulationSDD::InitSimulationModule( Int_t module, Int_t event ) {
322 // create maps to build the lists of tracks for each summable digit
326 memset(fAnodeFire,0,sizeof(Bool_t)*fNofMaps);
328 //______________________________________________________________________
329 void AliITSsimulationSDD::ClearMaps() {
332 fHitSigMap2->ClearMap();
333 fHitNoiMap2->ClearMap();
335 //______________________________________________________________________
336 void AliITSsimulationSDD::SDigitiseModule(AliITSmodule *mod,Int_t md,Int_t ev){
337 // digitize module using the "slow" detector simulator creating
340 TObjArray *fHits = mod->GetHits();
341 Int_t nhits = fHits->GetEntriesFast();
344 InitSimulationModule( md, ev );
345 HitsToAnalogDigits( mod ); // fills fHitMap2 which is = fHitSigmap2
346 ChargeToSignal( fModule,kFALSE,kTRUE ); // - Process signal adding gain without adding noise
347 fHitMap2 = fHitNoiMap2; // - Swap to noise map
348 ChargeToSignal( fModule,kTRUE,kFALSE ); // - Process only noise
349 fHitMap2 = fHitSigMap2; // - Return to signal map
353 //______________________________________________________________________
354 Bool_t AliITSsimulationSDD::AddSDigitsToModule(TClonesArray *pItemArray,
356 // Add Summable digits to module maps.
357 AliITSresponseSDD* res = (AliITSresponseSDD*)fDetType->GetResponse(1);
358 Int_t nItems = pItemArray->GetEntries();
359 Double_t maxadc = res->MaxAdc();
362 // cout << "Adding "<< nItems <<" SDigits to module " << fModule << endl;
363 for( Int_t i=0; i<nItems; i++ ) {
364 AliITSpListItem * pItem = (AliITSpListItem *)(pItemArray->At( i ));
365 if( pItem->GetModule() != fModule ) {
366 Error( "AliITSsimulationSDD","Error reading, SDigits module "
367 "%d != current module %d: exit",
368 pItem->GetModule(), fModule );
372 if(pItem->GetSignal()>0.0 ) sig = kTRUE;
374 fpList->AddItemTo( mask, pItem ); // Add SignalAfterElect + noise
375 AliITSpListItem * pItem2 = fpList->GetpListItem( pItem->GetIndex() );
376 Double_t sigAE = pItem2->GetSignalAfterElect();
377 if( sigAE >= maxadc ) sigAE = maxadc-1; // avoid overflow signal
380 fpList->GetMapIndex( pItem->GetIndex(), ia, it );
381 fHitMap2->SetHit( ia, it, sigAE );
382 fAnodeFire[ia] = kTRUE;
386 //______________________________________________________________________
387 void AliITSsimulationSDD::FinishSDigitiseModule() {
388 // digitize module using the "slow" detector simulator from
389 // the sum of summable digits.
393 //______________________________________________________________________
394 void AliITSsimulationSDD::DigitiseModule(AliITSmodule *mod,Int_t md,Int_t ev){
395 // create maps to build the lists of tracks for each digit
397 TObjArray *fHits = mod->GetHits();
398 Int_t nhits = fHits->GetEntriesFast();
400 InitSimulationModule( md, ev );
402 if( !nhits && fCheckNoise ) {
403 ChargeToSignal( fModule,kTRUE,kFALSE ); // process noise
410 HitsToAnalogDigits( mod );
411 ChargeToSignal( fModule,kTRUE,kTRUE ); // process signal + noise
413 for( Int_t i=0; i<fNofMaps; i++ ) {
414 for( Int_t j=0; j<fMaxNofSamples; j++ ) {
415 Int_t jdx = j*fScaleSize;
416 Int_t index = fpList->GetHitIndex( i, j );
417 AliITSpListItem pItemTmp2( fModule, index, 0. );
418 // put the fScaleSize analog digits in only one
419 for( Int_t ik=0; ik<fScaleSize; ik++ ) {
420 AliITSpListItem *pItemTmp = fpList->GetpListItem( i, jdx+ik );
421 if( pItemTmp == 0 ) continue;
422 pItemTmp2.Add( pItemTmp );
424 fpList->DeleteHit( i, j );
425 fpList->AddItemTo( 0, &pItemTmp2 );
431 //______________________________________________________________________
432 void AliITSsimulationSDD::FinishDigits() {
433 // introduce the electronics effects and do zero-suppression if required
435 if( fCrosstalkFlag ) ApplyCrosstalk(fModule);
437 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
438 const char *kopt = res->GetZeroSuppOption();
439 ZeroSuppression( kopt );
441 //______________________________________________________________________
442 void AliITSsimulationSDD::HitsToAnalogDigits( AliITSmodule *mod ) {
443 // create maps to build the lists of tracks for each digit
445 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
446 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
447 TObjArray *hits = mod->GetHits();
448 Int_t nhits = hits->GetEntriesFast();
450 // Int_t arg[6] = {0,0,0,0,0,0};
452 Int_t nofAnodes = fNofMaps/2;
453 Double_t sddLength = seg->Dx();
454 Double_t sddWidth = seg->Dz();
455 Double_t anodePitch = seg->Dpz(dummy);
456 Double_t timeStep = seg->Dpx(dummy);
457 Double_t driftSpeed ; // drift velocity (anode dependent)
458 //Float_t maxadc = res->GetMaxAdc();
459 //Float_t topValue = res->GetDynamicRange();
460 Double_t norm = res->GetMaxAdc()/res->GetDynamicRange(); // maxadc/topValue;
461 Double_t cHloss = res->GetChargeLoss();
462 Float_t dfCoeff, s1; res->DiffCoeff(dfCoeff,s1); // Signal 2d Shape
463 Double_t eVpairs = res->GetGeVToCharge()*1.0E9; // 3.6 eV by def.
464 Double_t nsigma = res->GetNSigmaIntegration(); //
465 Int_t nlookups = res->GetGausNLookUp(); //
466 Float_t jitter = res->GetJitterError(); //
468 // Piergiorgio's part (apart for few variables which I made float
469 // when i thought that can be done
470 // Fill detector maps with GEANT hits
471 // loop over hits in the module
473 const Float_t kconv = 1.0e+6; // GeV->KeV
475 Int_t hitDetector; // detector number (lay,lad,hitDetector)
476 Int_t iWing; // which detector wing/side.
477 Int_t detector; // 2*(detector-1)+iWing
478 Int_t ii,kk,ka,kt; // loop indexs
479 Int_t ia,it,index; // sub-pixel integration indexies
480 Int_t iAnode; // anode number.
481 Int_t timeSample; // time buckett.
482 Int_t anodeWindow; // anode direction charge integration width
483 Int_t timeWindow; // time direction charge integration width
484 Int_t jamin,jamax; // anode charge integration window
485 Int_t jtmin,jtmax; // time charge integration window
486 Int_t ndiv; // Anode window division factor.
487 Int_t nsplit; // the number of splits in anode and time windows==1.
488 Int_t nOfSplits; // number of times track length is split into
489 Float_t nOfSplitsF; // Floating point version of nOfSplits.
490 Float_t kkF; // Floating point version of loop index kk.
491 Double_t pathInSDD; // Track length in SDD.
492 Double_t drPath; // average position of track in detector. in microns
493 Double_t drTime; // Drift time
494 Double_t nmul; // drift time window multiplication factor.
495 Double_t avDrft; // x position of path length segment in cm.
496 Double_t avAnode; // Anode for path length segment in Anode number (float)
497 Double_t xAnode; // Floating point anode number.
498 Double_t driftPath; // avDrft in microns.
499 Double_t width; // width of signal at anodes.
500 Double_t depEnergy; // Energy deposited in this GEANT step.
501 Double_t xL[3],dxL[3]; // local hit coordinates and diff.
502 Double_t sigA; // sigma of signal at anode.
503 Double_t sigT; // sigma in time/drift direction for track segment
504 Double_t aStep,aConst; // sub-pixel size and offset anode
505 Double_t tStep,tConst; // sub-pixel size and offset time
506 Double_t amplitude; // signal amplitude for track segment in nanoAmpere
507 Double_t chargeloss; // charge loss for track segment.
508 Double_t anodeAmplitude; // signal amplitude in anode direction
509 Double_t aExpo; // exponent of Gaussian anode direction
510 Double_t timeAmplitude; // signal amplitude in time direction
511 Double_t tExpo; // exponent of Gaussian time direction
512 // Double_t tof; // Time of flight in ns of this step.
514 for(ii=0; ii<nhits; ii++) {
515 if(!mod->LineSegmentL(ii,xL[0],dxL[0],xL[1],dxL[1],xL[2],dxL[2],
516 depEnergy,itrack)) continue;
517 xL[0] += 0.0001*gRandom->Gaus( 0, jitter ); //
518 xAnode=10000.*(xL[2]+0.5*dxL[2])/anodePitch + nofAnodes/2;;
519 driftSpeed = res->GetDriftSpeedAtAnode(xAnode);
520 if(timeStep*fMaxNofSamples < sddLength/driftSpeed) {
521 Warning("AliITSsimulationSDD",
522 "Time Interval > Allowed Time Interval\n");
525 hitDetector = mod->GetDet();
526 //tof = 1.E+09*(mod->GetHit(ii)->GetTOF()); // tof in ns.
527 //if(tof>sddLength/driftSpeed) continue; // hit happed too late.
529 // scale path to simulate a perpendicular track
530 // continue if the particle did not lose energy
531 // passing through detector
534 Form("fTrack = %d hit=%d module=%d This particle has passed without losing energy!",
535 itrack,ii,mod->GetIndex()));
537 } // end if !depEnergy
539 pathInSDD = TMath::Sqrt(dxL[0]*dxL[0]+dxL[1]*dxL[1]+dxL[2]*dxL[2]);
541 if (fFlag && pathInSDD) { depEnergy *= (0.03/pathInSDD); }
542 drPath = 10000.*(dxL[0]+2.*xL[0])*0.5;
543 if(drPath < 0) drPath = -drPath;
544 drPath = sddLength-drPath;
546 AliDebug(1, // this should be fixed at geometry level
547 Form("negative drift path drPath=%e sddLength=%e dxL[0]=%e xL[0]=%e",
548 drPath,sddLength,dxL[0],xL[0]));
550 } // end if drPath < 0
552 // Compute number of segments to brake step path into
553 drTime = drPath/driftSpeed; // Drift Time
554 sigA = TMath::Sqrt(2.*dfCoeff*drTime+s1*s1);// Sigma along the anodes
555 // calcuate the number of time the path length should be split into.
556 nOfSplits = (Int_t) (1. + 10000.*pathInSDD/sigA);
557 if(fFlag) nOfSplits = 1;
559 // loop over path segments, init. some variables.
560 depEnergy /= nOfSplits;
561 nOfSplitsF = (Float_t) nOfSplits;
562 for(kk=0;kk<nOfSplits;kk++) { // loop over path segments
563 kkF = (Float_t) kk + 0.5;
564 avDrft = xL[0]+dxL[0]*kkF/nOfSplitsF;
565 avAnode = xL[2]+dxL[2]*kkF/nOfSplitsF;
566 driftPath = 10000.*avDrft;
568 iWing = 2; // Assume wing is 2
569 if(driftPath < 0) { // if wing is not 2 it is 1.
571 driftPath = -driftPath;
572 } // end if driftPath < 0
573 driftPath = sddLength-driftPath;
574 detector = 2*(hitDetector-1) + iWing;
576 AliDebug(1, // this should be fixed at geometry level
577 Form("negative drift path driftPath=%e sddLength=%e avDrft=%e dxL[0]=%e xL[0]=%e",
578 driftPath,sddLength,avDrft,dxL[0],xL[0]));
580 } // end if driftPath < 0
583 drTime = driftPath/driftSpeed; // drift time for segment.
584 timeSample = (Int_t) (fScaleSize*drTime/timeStep + 1);
585 // compute time Sample including tof information. The tof only
586 // effects the time of the signal is recoreded and not the
588 // timeSample = (Int_t) (fScaleSize*(drTime+tof)/timeStep + 1);
589 if(timeSample > fScaleSize*fMaxNofSamples) {
590 Warning("HitsToAnalogDigits","Wrong Time Sample: %e",
593 } // end if timeSample > fScaleSize*fMaxNoofSamples
596 xAnode = 10000.*(avAnode)/anodePitch + nofAnodes/2; // +1?
597 if(xAnode*anodePitch > sddWidth || xAnode*anodePitch < 0.)
598 Warning("HitsToAnalogDigits",
599 "Exceedubg sddWidth=%e Z = %e",
600 sddWidth,xAnode*anodePitch);
601 iAnode = (Int_t) (1.+xAnode); // xAnode?
602 if(iAnode < 1 || iAnode > nofAnodes) {
603 Warning("HitToAnalogDigits","Wrong iAnode: 1<%d>%d (xanode=%e)",
604 iAnode,nofAnodes, xAnode);
606 } // end if iAnode < 1 || iAnode > nofAnodes
608 // store straight away the particle position in the array
609 // of particles and take idhit=ii only when part is entering (this
610 // requires FillModules() in the macro for analysis) :
612 // Sigma along the anodes for track segment.
613 sigA = TMath::Sqrt(2.*dfCoeff*drTime+s1*s1);
614 sigT = sigA/driftSpeed;
615 // Peak amplitude in nanoAmpere
616 amplitude = fScaleSize*160.*depEnergy/
617 (timeStep*eVpairs*2.*acos(-1.)*sigT*sigA);
618 amplitude *= timeStep/25.; // WARNING!!!!! Amplitude scaling to
619 // account for clock variations
620 // (reference value: 40 MHz)
621 chargeloss = 1.-cHloss*driftPath/1000;
622 amplitude *= chargeloss;
623 width = 2.*nsigma/(nlookups-1);
631 } // end if drTime > 1200.
633 nsplit = 4; // hard-wired //nsplit=4;nsplit = (nsplit+1)/2*2;
634 // Sub-pixel size see computation of aExpo and tExpo.
635 aStep = anodePitch/(nsplit*fScaleSize*sigA);
636 aConst = xAnode*anodePitch/sigA;
637 tStep = timeStep/(nsplit*fScaleSize*sigT);
638 tConst = drTime/sigT;
639 // Define SDD window corresponding to the hit
640 anodeWindow = (Int_t)(fScaleSize*nsigma*sigA/anodePitch+1);
641 timeWindow = (Int_t) (fScaleSize*nsigma*sigT/timeStep+1.);
642 jamin = (iAnode - anodeWindow/ndiv - 1)*fScaleSize*nsplit +1;
643 jamax = (iAnode + anodeWindow/ndiv)*fScaleSize*nsplit;
644 if(jamin <= 0) jamin = 1;
645 if(jamax > fScaleSize*nofAnodes*nsplit)
646 jamax = fScaleSize*nofAnodes*nsplit;
647 // jtmin and jtmax are Hard-wired
648 jtmin = (Int_t)(timeSample-timeWindow*nmul-1)*nsplit+1;
649 jtmax = (Int_t)(timeSample+timeWindow*nmul)*nsplit;
650 if(jtmin <= 0) jtmin = 1;
651 if(jtmax > fScaleSize*fMaxNofSamples*nsplit)
652 jtmax = fScaleSize*fMaxNofSamples*nsplit;
653 // Spread the charge in the anode-time window
654 for(ka=jamin; ka <=jamax; ka++) {
655 ia = (ka-1)/(fScaleSize*nsplit) + 1;
657 Warning("HitsToAnalogDigits","ia < 1: ");
660 if(ia > nofAnodes) ia = nofAnodes;
661 aExpo = (aStep*(ka-0.5)-aConst);
662 if(TMath::Abs(aExpo) > nsigma) anodeAmplitude = 0.;
664 dummy = (Int_t) ((aExpo+nsigma)/width);
665 anodeAmplitude = amplitude*res->GetGausLookUp(dummy);
666 } // end if TMath::Abs(aEspo) > nsigma
667 // index starts from 0
668 index = ((detector+1)%2)*nofAnodes+ia-1;
669 if(anodeAmplitude) for(kt=jtmin; kt<=jtmax; kt++) {
670 it = (kt-1)/nsplit+1; // it starts from 1
672 Warning("HitsToAnalogDigits","it < 1:");
675 if(it>fScaleSize*fMaxNofSamples)
676 it = fScaleSize*fMaxNofSamples;
677 tExpo = (tStep*(kt-0.5)-tConst);
678 if(TMath::Abs(tExpo) > nsigma) timeAmplitude = 0.;
680 dummy = (Int_t) ((tExpo+nsigma)/width);
681 timeAmplitude = anodeAmplitude*
682 res->GetGausLookUp(dummy);
683 } // end if TMath::Abs(tExpo) > nsigma
684 // build the list of Sdigits for this module
687 // arg[2] = itrack; // track number
688 // arg[3] = ii-1; // hit number.
689 timeAmplitude *= norm;
691 // ListOfFiredCells(arg,timeAmplitude,alst,padr);
692 Double_t charge = timeAmplitude;
693 charge += fHitMap2->GetSignal(index,it-1);
694 fHitMap2->SetHit(index, it-1, charge);
695 fpList->AddSignal(index,it-1,itrack,ii-1,
696 mod->GetIndex(),timeAmplitude);
697 fAnodeFire[index] = kTRUE;
698 } // end if anodeAmplitude and loop over time in window
699 } // loop over anodes in window
700 } // end loop over "sub-hits"
701 } // end loop over hits
704 //____________________________________________
705 void AliITSsimulationSDD::AddDigit( Int_t i, Int_t j, Int_t signal ) {
707 Int_t size = AliITSdigit::GetNTracks();
710 Int_t * tracks = new Int_t[size];
711 Int_t * hits = new Int_t[size];
713 Float_t * charges = new Float_t[size];
719 AliITSpListItem *pItem = fpList->GetpListItem( i, j );
722 for( Int_t l=0; l<size; l++ ) {
728 Int_t idtrack = pItem->GetTrack( 0 );
729 if( idtrack >= 0 ) phys = pItem->GetSignal();
732 for( Int_t l=0; l<size; l++ ) if(l<pItem->GetMaxKept()) {
733 tracks[l] = pItem->GetTrack( l );
734 hits[l] = pItem->GetHit( l );
735 charges[l] = pItem->GetSignal( l );
743 fITS->AddSimDigit( 1, phys, digits, tracks, hits, charges );
748 //______________________________________________________________________
749 void AliITSsimulationSDD::ChargeToSignal(Int_t mod,Bool_t bAddNoise, Bool_t bAddGain) {
750 // add baseline, noise, gain, electronics and ADC saturation effects
751 // apply dead channels
753 char opt1[20], opt2[20];
754 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(mod);
755 res->GetParamOptions(opt1,opt2);
761 Float_t maxadc = res->GetMaxAdc();
763 for (i=0;i<fNofMaps;i++) {
764 if( !fAnodeFire[i] ) continue;
765 baseline = res->GetBaseline(i);
766 noise = res->GetNoise(i);
767 gain = res->GetChannelGain(i);
768 if(res->IsDead()) gain=0;
769 for(k=0; k<fScaleSize*fMaxNofSamples; k++) {
770 fInZR[k] = fHitMap2->GetSignal(i,k);
771 if(bAddGain) fInZR[k]*=gain;
773 contrib = (baseline + noise*gRandom->Gaus());
779 for(k=0; k<fMaxNofSamples; k++) {
780 Double_t newcont = 0.;
781 Double_t maxcont = 0.;
782 for(kk=0;kk<fScaleSize;kk++) {
783 newcont = fInZR[fScaleSize*k+kk];
784 if(newcont > maxcont) maxcont = newcont;
787 if (newcont >= maxadc) newcont = maxadc -1;
788 if(newcont >= baseline){
789 Warning("","newcont=%d>=baseline=%d",newcont,baseline);
792 fHitMap2->SetHit(i,k,newcont);
795 FastFourierTransform(fElectronics,&fInZR[0],&fInZI[0],1);
796 for(k=0; k<fScaleSize*fMaxNofSamples; k++) {
797 Double_t rw = fElectronics->GetTraFunReal(k);
798 Double_t iw = fElectronics->GetTraFunImag(k);
799 fOutZR[k] = fInZR[k]*rw - fInZI[k]*iw;
800 fOutZI[k] = fInZR[k]*iw + fInZI[k]*rw;
802 FastFourierTransform(fElectronics,&fOutZR[0],&fOutZI[0],-1);
803 for(k=0; k<fMaxNofSamples; k++) {
804 Double_t newcont1 = 0.;
805 Double_t maxcont1 = 0.;
806 for(kk=0;kk<fScaleSize;kk++) {
807 newcont1 = fOutZR[fScaleSize*k+kk];
808 if(newcont1 > maxcont1) maxcont1 = newcont1;
811 if (newcont1 >= maxadc) newcont1 = maxadc -1;
812 fHitMap2->SetHit(i,k,newcont1);
815 } // end for i loop over anodes
819 //______________________________________________________________________
820 void AliITSsimulationSDD::ApplyCrosstalk(Int_t mod) {
821 // function add the crosstalk effect to signal
822 // temporal function, should be checked...!!!
823 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
825 Int_t fNofMaps = seg->Npz();
826 Int_t fMaxNofSamples = seg->Npx();
828 // create and inizialice crosstalk map
829 Float_t* ctk = new Float_t[fNofMaps*fMaxNofSamples+1];
831 Error( "ApplyCrosstalk", "no memory for temporal map: exit \n" );
834 memset( ctk, 0, sizeof(Float_t)*(fNofMaps*fMaxNofSamples+1) );
835 AliITSCalibrationSDD* calibr = (AliITSCalibrationSDD*)GetCalibrationModel(mod);
836 for( Int_t z=0; z<fNofMaps; z++ ) {
837 Double_t baseline = calibr->GetBaseline(z);
843 for( Int_t l=0; l<fMaxNofSamples; l++ ) {
844 Float_t fadc = (Float_t)fHitMap2->GetSignal( z, l );
845 if( fadc > baseline ) {
846 if( on == kFALSE && l<fMaxNofSamples-4 ) {
847 Float_t fadc1 = (Float_t)fHitMap2->GetSignal( z, l+1 );
848 if( fadc1 < fadc ) continue;
855 else { // end fadc > baseline
859 // make smooth derivative
860 Float_t* dev = new Float_t[fMaxNofSamples+1];
861 memset( dev, 0, sizeof(Float_t)*(fMaxNofSamples+1) );
863 Error( "ApplyCrosstalk",
864 "no memory for temporal array: exit \n" );
867 for( Int_t i=tstart; i<tstop; i++ ) {
868 if( i > 2 && i < fMaxNofSamples-2 )
869 dev[i] = -0.2*fHitMap2->GetSignal( z,i-2 )
870 -0.1*fHitMap2->GetSignal( z,i-1 )
871 +0.1*fHitMap2->GetSignal( z,i+1 )
872 +0.2*fHitMap2->GetSignal( z,i+2 );
875 // add crosstalk contribution to neibourg anodes
876 for( Int_t i=tstart; i<tstop; i++ ) {
878 Int_t i1 = (Int_t)((i-tstart)*.61+tstart+0.5); //
879 Float_t ctktmp = -dev[i1] * 0.25;
881 ctk[anode*fMaxNofSamples+i] += ctktmp;
884 if( anode < fNofMaps ) {
885 ctk[anode*fMaxNofSamples+i] += ctktmp;
890 } // if( nTsteps > 2 )
892 } // if( on == kTRUE )
897 for( Int_t a=0; a<fNofMaps; a++ )
898 for( Int_t t=0; t<fMaxNofSamples; t++ ) {
899 Float_t signal = fHitMap2->GetSignal(a,t)+ctk[a*fMaxNofSamples+t];
900 fHitMap2->SetHit( a, t, signal );
906 //______________________________________________________________________
907 void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl,
909 // Returns the compression alogirthm parameters
914 //______________________________________________________________________
915 void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl) const{
916 // returns the compression alogirthm parameters
922 //______________________________________________________________________
923 void AliITSsimulationSDD::SetCompressParam(){
924 // Sets the compression alogirthm parameters
925 AliITSCalibrationSDD* calibr = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
926 for(Int_t ian = 0; ian<fNofMaps;ian++){
927 fD[ian] = (Int_t)(calibr->GetBaseline(ian));
928 fT1[ian] = (Int_t)(2.*calibr->GetNoiseAfterElectronics(ian)+0.5);
929 fT2[ian] = 0; // used by 2D clustering - not defined yet
930 fTol[ian] = 0; // used by 2D clustering - not defined yet
933 //______________________________________________________________________
934 Int_t AliITSsimulationSDD::Convert10to8(Int_t signal) const {
935 // To the 10 to 8 bit lossive compression.
936 // code from Davide C. and Albert W.
938 if (signal < 128) return signal;
939 if (signal < 256) return (128+((signal-128)>>1));
940 if (signal < 512) return (192+((signal-256)>>3));
941 if (signal < 1024) return (224+((signal-512)>>4));
945 //______________________________________________________________________
946 AliITSMap* AliITSsimulationSDD::HitMap(Int_t i){
947 //Return the correct map.
949 return ((i==0)? fHitMap1 : fHitMap2);
952 //______________________________________________________________________
953 void AliITSsimulationSDD::ZeroSuppression(const char *option) {
954 // perform the zero suppresion
955 if (strstr(option,"2D")) {
956 //Init2D(); // activate if param change module by module
958 } else if (strstr(option,"1D")) {
959 //Init1D(); // activate if param change module by module
961 } else StoreAllDigits();
963 //______________________________________________________________________
964 void AliITSsimulationSDD::Init2D(){
965 // read in and prepare arrays: fD, fT1, fT2
966 // savemu[nanodes], savesigma[nanodes]
967 // read baseline and noise from file - either a .root file and in this
968 // case data should be organised in a tree with one entry for each
969 // module => reading should be done accordingly
970 // or a classic file and do smth. like this ( code from Davide C. and
972 // Read 2D zero-suppression parameters for SDD
974 if (!strstr(fParam.Data(),"file")) return;
978 Float_t *savemu = new Float_t [fNofMaps];
979 Float_t *savesigma = new Float_t [fNofMaps];
980 char input[100],basel[100],par[100];
983 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
985 res->Thresholds(tmp1,tmp2);
986 Int_t minval = static_cast<Int_t>(tmp1);
988 res->Filenames(input,basel,par);
991 filtmp = gSystem->ExpandPathName(fFileName.Data());
992 FILE *param = fopen(filtmp,"r");
996 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
998 Error("Init2D","Anode number not in increasing order!",filtmp);
1000 } // end if pos != na+1
1002 savesigma[na] = sigma;
1003 if ((2.*sigma) < mu) {
1004 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
1007 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
1008 if (tempTh < 0) tempTh=0;
1010 tempTh = (Int_t)floor(mu+3.0*sigma+0.5) - minval;
1011 if (tempTh < 0) tempTh=0;
1016 Error("Init2D","THE FILE %s DOES NOT EXIST !",filtmp);
1023 delete [] savesigma;
1025 //______________________________________________________________________
1026 void AliITSsimulationSDD::Compress2D(){
1027 // simple ITS cluster finder -- online zero-suppression conditions
1031 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1033 res->Thresholds(tmp1,tmp2);
1034 Int_t minval = static_cast<Int_t>(tmp1);
1035 Bool_t write = res->OutputOption();
1036 Bool_t do10to8 = res->Do10to8();
1037 Int_t nz, nl, nh, low, i, j;
1039 for (i=0; i<fNofMaps; i++) {
1040 CompressionParam(i,db,tl,th);
1045 for (j=0; j<fMaxNofSamples; j++) {
1046 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1047 signal -= db; // if baseline eq. is done here
1048 if (signal <= 0) {nz++; continue;}
1049 if ((signal - tl) < minval) low++;
1050 if ((signal - th) >= minval) {
1053 FindCluster(i,j,signal,minval,cond);
1055 ((TMath::Abs(fHitMap2->GetSignal(i,j-1))-th)>=minval)){
1056 if(do10to8) signal = Convert10to8(signal);
1057 AddDigit(i,j,signal);
1058 } // end if cond&&j&&()
1059 } else if ((signal - tl) >= minval) nl++;
1060 } // end for j loop time samples
1061 if (write) TreeB()->Fill(nz,nl,nh,low,i+1);
1062 } //end for i loop anodes
1066 sprintf(hname,"TNtuple%d_%d",fModule,fEvent);
1067 TreeB()->Write(hname);
1072 //______________________________________________________________________
1073 void AliITSsimulationSDD::FindCluster(Int_t i,Int_t j,Int_t signal,
1074 Int_t minval,Bool_t &cond){
1075 // Find clusters according to the online 2D zero-suppression algorithm
1076 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1077 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
1079 Bool_t do10to8 = res->Do10to8();
1080 Bool_t high = kFALSE;
1082 fHitMap2->FlagHit(i,j);
1084 // check the online zero-suppression conditions
1086 const Int_t kMaxNeighbours = 4;
1089 Int_t xList[kMaxNeighbours], yList[kMaxNeighbours];
1090 seg->Neighbours(i,j,&nn,xList,yList);
1092 for (in=0; in<nn; in++) {
1095 if (fHitMap2->TestHit(ix,iy)==kUnused) {
1096 CompressionParam(ix,dbx,tlx,thx);
1097 Int_t qn = (Int_t)(fHitMap2->GetSignal(ix,iy));
1098 qn -= dbx; // if baseline eq. is done here
1099 if ((qn-tlx) < minval) {
1100 fHitMap2->FlagHit(ix,iy);
1103 if ((qn - thx) >= minval) high=kTRUE;
1105 if(do10to8) signal = Convert10to8(signal);
1106 AddDigit(i,j,signal);
1108 if(do10to8) qns = Convert10to8(qn);
1110 if (!high) AddDigit(ix,iy,qns);
1112 if(!high) fHitMap2->FlagHit(ix,iy);
1113 } // end if qn-tlx < minval
1115 } // end for in loop over neighbours
1117 //______________________________________________________________________
1118 void AliITSsimulationSDD::Init1D(){
1119 // this is just a copy-paste of input taken from 2D algo
1120 // Torino people should give input
1121 // Read 1D zero-suppression parameters for SDD
1123 if (!strstr(fParam.Data(),"file")) return;
1125 Int_t na,pos,tempTh;
1127 Float_t *savemu = new Float_t [fNofMaps];
1128 Float_t *savesigma = new Float_t [fNofMaps];
1129 char input[100],basel[100],par[100];
1132 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1134 res->Thresholds(tmp1,tmp2);
1135 Int_t minval = static_cast<Int_t>(tmp1);
1137 res->Filenames(input,basel,par);
1140 // set first the disable and tol param
1143 filtmp = gSystem->ExpandPathName(fFileName.Data());
1144 FILE *param = fopen(filtmp,"r");
1148 fscanf(param,"%d %d %d %d ", &fT2[0], &fT2[1], &fTol[0], &fTol[1]);
1149 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
1151 Error("Init1D","Anode number not in increasing order!",filtmp);
1153 } // end if pos != na+1
1155 savesigma[na]=sigma;
1156 if ((2.*sigma) < mu) {
1157 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
1160 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
1161 if (tempTh < 0) tempTh=0;
1166 Error("Init1D","THE FILE %s DOES NOT EXIST !",filtmp);
1173 delete [] savesigma;
1175 //______________________________________________________________________
1176 void AliITSsimulationSDD::Compress1D(){
1177 // 1D zero-suppression algorithm (from Gianluca A.)
1178 Int_t dis,tol,thres,decr,diff;
1179 UChar_t *str=fStream->Stream();
1181 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1183 Bool_t do10to8=res->Do10to8();
1187 for (k=0; k<2; k++) {
1190 for (i=0; i<fNofMaps/2; i++) {
1191 Bool_t firstSignal=kTRUE;
1192 Int_t idx=i+k*fNofMaps/2;
1193 if( !fAnodeFire[idx] ) continue;
1194 CompressionParam(idx,decr,thres);
1196 for (j=0; j<fMaxNofSamples; j++) {
1197 Int_t signal=(Int_t)(fHitMap2->GetSignal(idx,j));
1198 signal -= decr; // if baseline eq.
1199 if(do10to8) signal = Convert10to8(signal);
1200 if (signal <= thres) {
1204 // write diff in the buffer for HuffT
1205 str[counter]=(UChar_t)diff;
1208 } // end if signal <= thres
1210 if (diff > 127) diff=127;
1211 if (diff < -128) diff=-128;
1213 // tol has changed to 8 possible cases ? - one can write
1214 // this if(TMath::Abs(diff)<tol) ... else ...
1215 if(TMath::Abs(diff)<tol) diff=0;
1216 // or keep it as it was before
1217 AddDigit(idx,j,last+diff);
1219 AddDigit(idx,j,signal);
1220 } // end if singal < dis
1222 // write diff in the buffer used to compute Huffman tables
1223 if (firstSignal) str[counter]=(UChar_t)signal;
1224 else str[counter]=(UChar_t)diff;
1228 } // end for j loop time samples
1229 } // end for i loop anodes one half of detector
1233 fStream->CheckCount(counter);
1235 // open file and write out the stream of diff's
1236 static Bool_t open=kTRUE;
1237 static TFile *outFile;
1238 Bool_t write = res->OutputOption();
1239 TDirectory *savedir = gDirectory;
1243 SetFileName("stream.root");
1244 cout<<"filename "<<fFileName<<endl;
1245 outFile=new TFile(fFileName,"recreate");
1246 cout<<"I have opened "<<fFileName<<" file "<<endl;
1253 fStream->ClearStream();
1255 // back to galice.root file
1256 if(savedir) savedir->cd();
1258 //______________________________________________________________________
1259 void AliITSsimulationSDD::StoreAllDigits(){
1260 // if non-zero-suppressed data
1261 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1263 Bool_t do10to8 = res->Do10to8();
1264 Int_t i, j, digits[3];
1266 for (i=0; i<fNofMaps; i++) {
1267 for (j=0; j<fMaxNofSamples; j++) {
1268 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1269 if(do10to8) signal = Convert10to8(signal);
1273 fITS->AddRealDigit(1,digits);
1277 //______________________________________________________________________
1278 void AliITSsimulationSDD::CreateHistograms(Int_t scale){
1279 // Creates histograms of maps for debugging
1282 fHis=new TObjArray(fNofMaps);
1283 for (i=0;i<fNofMaps;i++) {
1284 TString sddName("sdd_");
1286 sprintf(candNum,"%d",i+1);
1287 sddName.Append(candNum);
1288 fHis->AddAt(new TH1F(sddName.Data(),"SDD maps",scale*fMaxNofSamples,
1289 0.,(Float_t) scale*fMaxNofSamples), i);
1292 //______________________________________________________________________
1293 void AliITSsimulationSDD::FillHistograms(){
1294 // fill 1D histograms from map
1298 for( Int_t i=0; i<fNofMaps; i++) {
1299 TH1F *hist =(TH1F *)fHis->UncheckedAt(i);
1300 Int_t nsamples = hist->GetNbinsX();
1301 for( Int_t j=0; j<nsamples; j++) {
1302 Double_t signal=fHitMap2->GetSignal(i,j);
1303 hist->Fill((Float_t)j,signal);
1307 //______________________________________________________________________
1308 void AliITSsimulationSDD::ResetHistograms(){
1309 // Reset histograms for this detector
1312 for (i=0;i<fNofMaps;i++ ) {
1313 if (fHis->At(i)) ((TH1F*)fHis->At(i))->Reset();
1316 //______________________________________________________________________
1317 TH1F *AliITSsimulationSDD::GetAnode(Int_t wing, Int_t anode) {
1318 // Fills a histogram from a give anode.
1320 if (!fHis) return 0;
1322 if(wing <=0 || wing > 2) {
1323 Warning("GetAnode","Wrong wing number: %d",wing);
1325 } // end if wing <=0 || wing >2
1326 if(anode <=0 || anode > fNofMaps/2) {
1327 Warning("GetAnode","Wrong anode number: %d",anode);
1329 } // end if ampde <=0 || andoe > fNofMaps/2
1331 Int_t index = (wing-1)*fNofMaps/2 + anode-1;
1332 return (TH1F*)(fHis->At(index));
1334 //______________________________________________________________________
1335 void AliITSsimulationSDD::WriteToFile(TFile *hfile) {
1336 // Writes the histograms to a file
1342 for(i=0; i<fNofMaps; i++) fHis->At(i)->Write(); //fAdcs[i]->Write();
1345 //______________________________________________________________________
1346 Float_t AliITSsimulationSDD::GetNoise() {
1347 // Returns the noise value
1348 //Bool_t do10to8=GetResp()->Do10to8();
1349 //noise will always be in the liniar part of the signal
1351 Int_t threshold = fT1[0];
1352 char opt1[20], opt2[20];
1353 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1355 res->GetParamOptions(opt1,opt2);
1357 Double_t noise,baseline;
1358 //GetBaseline(fModule);
1360 TCanvas *c2 = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("c2");
1361 if(c2) delete c2->GetPrimitive("noisehist");
1362 if(c2) delete c2->GetPrimitive("anode");
1363 else c2=new TCanvas("c2");
1365 c2->SetFillColor(0);
1367 TH1F *noisehist = new TH1F("noisehist","noise",100,0.,(float)2*threshold);
1368 TH1F *anode = new TH1F("anode","Anode Projection",fMaxNofSamples,0.,
1369 (float)fMaxNofSamples);
1371 for (i=0;i<fNofMaps;i++) {
1372 CompressionParam(i,decr,threshold);
1373 baseline = res->GetBaseline(i);
1374 noise = res->GetNoise(i);
1376 for (k=0;k<fMaxNofSamples;k++) {
1377 Float_t signal=(Float_t)fHitMap2->GetSignal(i,k);
1378 //if (signal <= (float)threshold) noisehist->Fill(signal-baseline);
1379 if (signal <= (float)(threshold+decr)) noisehist->Fill(signal);
1380 anode->Fill((float)k,signal);
1385 TF1 *gnoise = new TF1("gnoise","gaus",0.,threshold);
1386 noisehist->Fit("gnoise","RQ");
1389 Float_t mnoise = gnoise->GetParameter(1);
1390 cout << "mnoise : " << mnoise << endl;
1391 Float_t rnoise = gnoise->GetParameter(2);
1392 cout << "rnoise : " << rnoise << endl;
1396 //______________________________________________________________________
1397 void AliITSsimulationSDD::WriteSDigits(){
1398 // Fills the Summable digits Tree
1399 static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
1401 for( Int_t i=0; i<fNofMaps; i++ ) {
1402 if( !fAnodeFire[i] ) continue;
1403 for( Int_t j=0; j<fMaxNofSamples; j++ ) {
1404 Double_t sig = fHitMap2->GetSignal( i, j );
1406 Int_t jdx = j*fScaleSize;
1407 Int_t index = fpList->GetHitIndex( i, j );
1408 AliITSpListItem pItemTmp2( fModule, index, 0. );
1409 // put the fScaleSize analog digits in only one
1410 for( Int_t ik=0; ik<fScaleSize; ik++ ) {
1411 AliITSpListItem *pItemTmp = fpList->GetpListItem(i,jdx+ik);
1412 if( pItemTmp == 0 ) continue;
1413 pItemTmp2.Add( pItemTmp );
1415 pItemTmp2.AddSignalAfterElect( fModule, index, sig );
1416 pItemTmp2.AddNoise(fModule,index,fHitNoiMap2->GetSignal(i,j));
1417 aliITS->AddSumDigit( pItemTmp2 );
1418 } // end if (sig > 0.2)
1423 //______________________________________________________________________
1424 void AliITSsimulationSDD::PrintStatus() const {
1425 // Print SDD simulation Parameters
1427 cout << "**************************************************" << endl;
1428 cout << " Silicon Drift Detector Simulation Parameters " << endl;
1429 cout << "**************************************************" << endl;
1430 cout << "Flag for Perpendicular tracks: " << (Int_t) fFlag << endl;
1431 cout << "Flag for noise checking: " << (Int_t) fCheckNoise << endl;
1432 cout << "Flag to switch off electronics: " << (Int_t) fDoFFT << endl;
1433 cout << "Number pf Anodes used: " << fNofMaps << endl;
1434 cout << "Number of Time Samples: " << fMaxNofSamples << endl;
1435 cout << "Scale size factor: " << fScaleSize << endl;
1436 cout << "**************************************************" << endl;