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();
261 Float_t anodePitch = seg->Dpz(0);
262 Double_t timeStep = (Double_t)seg->Dpx(0);
264 if(anodePitch*(fNofMaps/2) > sddWidth) {
265 Warning("AliITSsimulationSDD",
266 "Too many anodes %d or too big pitch %f \n",
267 fNofMaps/2,anodePitch);
271 fElectronics = new AliITSetfSDD(timeStep/fScaleSize,
274 char opt1[20], opt2[20];
275 res->ParamOptions(opt1,opt2);
278 const char *kopt=res->ZeroSuppOption();
284 Bool_t write = res->OutputOption();
285 if(write && strstr(kopt,"2D")) MakeTreeB();
287 fITS = (AliITS*)gAlice->GetModule("ITS");
288 Int_t size = fNofMaps*fMaxNofSamples;
289 fStream = new AliITSInStream(size);
291 fInZR = new Double_t [fScaleSize*fMaxNofSamples];
292 fInZI = new Double_t [fScaleSize*fMaxNofSamples];
293 fOutZR = new Double_t [fScaleSize*fMaxNofSamples];
294 fOutZI = new Double_t [fScaleSize*fMaxNofSamples];
296 //______________________________________________________________________
297 AliITSsimulationSDD::~AliITSsimulationSDD() {
312 if(fTreeB) delete fTreeB;
313 if(fInZR) delete [] fInZR;
314 if(fInZI) delete [] fInZI;
315 if(fOutZR) delete [] fOutZR;
316 if(fOutZI) delete [] fOutZI;
317 if(fAnodeFire) delete [] fAnodeFire;
319 //______________________________________________________________________
320 void AliITSsimulationSDD::InitSimulationModule( Int_t module, Int_t event ) {
321 // create maps to build the lists of tracks for each summable digit
325 memset(fAnodeFire,0,sizeof(Bool_t)*fNofMaps);
327 //______________________________________________________________________
328 void AliITSsimulationSDD::ClearMaps() {
331 fHitSigMap2->ClearMap();
332 fHitNoiMap2->ClearMap();
334 //______________________________________________________________________
335 void AliITSsimulationSDD::SDigitiseModule(AliITSmodule *mod,Int_t md,Int_t ev){
336 // digitize module using the "slow" detector simulator creating
339 TObjArray *fHits = mod->GetHits();
340 Int_t nhits = fHits->GetEntriesFast();
343 InitSimulationModule( md, ev );
344 HitsToAnalogDigits( mod ); // fills fHitMap2 which is = fHitSigmap2
345 ChargeToSignal( fModule,kFALSE,kTRUE ); // - Process signal adding gain without adding noise
346 fHitMap2 = fHitNoiMap2; // - Swap to noise map
347 ChargeToSignal( fModule,kTRUE,kFALSE ); // - Process only noise
348 fHitMap2 = fHitSigMap2; // - Return to signal map
352 //______________________________________________________________________
353 Bool_t AliITSsimulationSDD::AddSDigitsToModule(TClonesArray *pItemArray,
355 // Add Summable digits to module maps.
356 AliITSresponseSDD* res = (AliITSresponseSDD*)fDetType->GetResponse(1);
357 Int_t nItems = pItemArray->GetEntries();
358 Double_t maxadc = res->MaxAdc();
361 // cout << "Adding "<< nItems <<" SDigits to module " << fModule << endl;
362 for( Int_t i=0; i<nItems; i++ ) {
363 AliITSpListItem * pItem = (AliITSpListItem *)(pItemArray->At( i ));
364 if( pItem->GetModule() != fModule ) {
365 Error( "AliITSsimulationSDD","Error reading, SDigits module "
366 "%d != current module %d: exit",
367 pItem->GetModule(), fModule );
371 if(pItem->GetSignal()>0.0 ) sig = kTRUE;
373 fpList->AddItemTo( mask, pItem ); // Add SignalAfterElect + noise
374 AliITSpListItem * pItem2 = fpList->GetpListItem( pItem->GetIndex() );
375 Double_t sigAE = pItem2->GetSignalAfterElect();
376 if( sigAE >= maxadc ) sigAE = maxadc-1; // avoid overflow signal
379 fpList->GetMapIndex( pItem->GetIndex(), ia, it );
380 fHitMap2->SetHit( ia, it, sigAE );
381 fAnodeFire[ia] = kTRUE;
385 //______________________________________________________________________
386 void AliITSsimulationSDD::FinishSDigitiseModule() {
387 // digitize module using the "slow" detector simulator from
388 // the sum of summable digits.
392 //______________________________________________________________________
393 void AliITSsimulationSDD::DigitiseModule(AliITSmodule *mod,Int_t md,Int_t ev){
394 // create maps to build the lists of tracks for each digit
396 TObjArray *fHits = mod->GetHits();
397 Int_t nhits = fHits->GetEntriesFast();
399 InitSimulationModule( md, ev );
401 if( !nhits && fCheckNoise ) {
402 ChargeToSignal( fModule,kTRUE,kFALSE ); // process noise
409 HitsToAnalogDigits( mod );
410 ChargeToSignal( fModule,kTRUE,kTRUE ); // process signal + noise
412 for( Int_t i=0; i<fNofMaps; i++ ) {
413 for( Int_t j=0; j<fMaxNofSamples; j++ ) {
414 Int_t jdx = j*fScaleSize;
415 Int_t index = fpList->GetHitIndex( i, j );
416 AliITSpListItem pItemTmp2( fModule, index, 0. );
417 // put the fScaleSize analog digits in only one
418 for( Int_t ik=0; ik<fScaleSize; ik++ ) {
419 AliITSpListItem *pItemTmp = fpList->GetpListItem( i, jdx+ik );
420 if( pItemTmp == 0 ) continue;
421 pItemTmp2.Add( pItemTmp );
423 fpList->DeleteHit( i, j );
424 fpList->AddItemTo( 0, &pItemTmp2 );
430 //______________________________________________________________________
431 void AliITSsimulationSDD::FinishDigits() {
432 // introduce the electronics effects and do zero-suppression if required
434 if( fCrosstalkFlag ) ApplyCrosstalk(fModule);
436 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
437 const char *kopt = res->GetZeroSuppOption();
438 ZeroSuppression( kopt );
440 //______________________________________________________________________
441 void AliITSsimulationSDD::HitsToAnalogDigits( AliITSmodule *mod ) {
442 // create maps to build the lists of tracks for each digit
443 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
444 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
445 TObjArray *hits = mod->GetHits();
446 Int_t nhits = hits->GetEntriesFast();
448 // Int_t arg[6] = {0,0,0,0,0,0};
449 Int_t nofAnodes = fNofMaps/2;
450 Double_t sddLength = seg->Dx();
451 Double_t sddWidth = seg->Dz();
452 Double_t anodePitch = seg->Dpz(0);
453 Double_t timeStep = seg->Dpx(0);
454 Double_t driftSpeed ; // drift velocity (anode dependent)
455 //Float_t maxadc = res->GetMaxAdc();
456 //Float_t topValue = res->GetDynamicRange();
457 Double_t norm = res->GetMaxAdc()/res->GetDynamicRange(); // maxadc/topValue;
458 Double_t cHloss = res->GetChargeLoss();
459 Float_t dfCoeff, s1; res->DiffCoeff(dfCoeff,s1); // Signal 2d Shape
460 Double_t eVpairs = res->GetGeVToCharge()*1.0E9; // 3.6 eV by def.
461 Double_t nsigma = res->GetNSigmaIntegration(); //
462 Int_t nlookups = res->GetGausNLookUp(); //
463 Float_t jitter = res->GetJitterError(); //
465 // Piergiorgio's part (apart for few variables which I made float
466 // when i thought that can be done
467 // Fill detector maps with GEANT hits
468 // loop over hits in the module
470 const Float_t kconv = 1.0e+6; // GeV->KeV
472 Int_t iWing; // which detector wing/side.
473 Int_t ii,kk,ka,kt; // loop indexs
474 Int_t ia,it,index; // sub-pixel integration indexies
475 Int_t iAnode; // anode number.
476 Int_t timeSample; // time buckett.
477 Int_t anodeWindow; // anode direction charge integration width
478 Int_t timeWindow; // time direction charge integration width
479 Int_t jamin,jamax; // anode charge integration window
480 Int_t jtmin,jtmax; // time charge integration window
481 Int_t ndiv; // Anode window division factor.
482 Int_t nsplit; // the number of splits in anode and time windows==1.
483 Int_t nOfSplits; // number of times track length is split into
484 Float_t nOfSplitsF; // Floating point version of nOfSplits.
485 Float_t kkF; // Floating point version of loop index kk.
486 Double_t pathInSDD; // Track length in SDD.
487 Double_t drPath; // average position of track in detector. in microns
488 Double_t drTime; // Drift time
489 Double_t nmul; // drift time window multiplication factor.
490 Double_t avDrft; // x position of path length segment in cm.
491 Double_t avAnode; // Anode for path length segment in Anode number (float)
492 Double_t zAnode; // Floating point anode number.
493 Double_t driftPath; // avDrft in microns.
494 Double_t width; // width of signal at anodes.
495 Double_t depEnergy; // Energy deposited in this GEANT step.
496 Double_t xL[3],dxL[3]; // local hit coordinates and diff.
497 Double_t sigA; // sigma of signal at anode.
498 Double_t sigT; // sigma in time/drift direction for track segment
499 Double_t aStep,aConst; // sub-pixel size and offset anode
500 Double_t tStep,tConst; // sub-pixel size and offset time
501 Double_t amplitude; // signal amplitude for track segment in nanoAmpere
502 Double_t chargeloss; // charge loss for track segment.
503 Double_t anodeAmplitude; // signal amplitude in anode direction
504 Double_t aExpo; // exponent of Gaussian anode direction
505 Double_t timeAmplitude; // signal amplitude in time direction
506 Double_t tExpo; // exponent of Gaussian time direction
507 // Double_t tof; // Time of flight in ns of this step.
509 for(ii=0; ii<nhits; ii++) {
510 if(!mod->LineSegmentL(ii,xL[0],dxL[0],xL[1],dxL[1],xL[2],dxL[2],
511 depEnergy,itrack)) continue;
513 if(xloc>0) iWing=0; // left side, carlos channel 0
514 else iWing=1; // right side
516 Float_t zloc=xL[2]+0.5*dxL[2];
517 zAnode=seg->GetAnodeFromLocal(xloc,zloc); // anode number in the range 0.-511.
518 driftSpeed = res->GetDriftSpeedAtAnode(zAnode);
519 if(timeStep*fMaxNofSamples < sddLength/driftSpeed) {
520 AliWarning("Time Interval > Allowed Time Interval\n");
524 // scale path to simulate a perpendicular track
525 // continue if the particle did not lose energy
526 // passing through detector
529 Form("fTrack = %d hit=%d module=%d This particle has passed without losing energy!",
530 itrack,ii,mod->GetIndex()));
532 } // end if !depEnergy
534 xL[0] += 0.0001*gRandom->Gaus( 0, jitter ); //
535 pathInSDD = TMath::Sqrt(dxL[0]*dxL[0]+dxL[1]*dxL[1]+dxL[2]*dxL[2]);
537 if (fFlag && pathInSDD) { depEnergy *= (0.03/pathInSDD); }
538 drPath = TMath::Abs(10000.*(dxL[0]+2.*xL[0])*0.5);
539 drPath = sddLength-drPath;
541 AliDebug(1, // this should be fixed at geometry level
542 Form("negative drift path drPath=%e sddLength=%e dxL[0]=%e xL[0]=%e",
543 drPath,sddLength,dxL[0],xL[0]));
545 } // end if drPath < 0
547 // Compute number of segments to brake step path into
548 drTime = drPath/driftSpeed; // Drift Time
549 sigA = TMath::Sqrt(2.*dfCoeff*drTime+s1*s1);// Sigma along the anodes
550 // calcuate the number of time the path length should be split into.
551 nOfSplits = (Int_t) (1. + 10000.*pathInSDD/sigA);
552 if(fFlag) nOfSplits = 1;
554 // loop over path segments, init. some variables.
555 depEnergy /= nOfSplits;
556 nOfSplitsF = (Float_t) nOfSplits;
557 Float_t theAverage=0.,theSteps=0.;
558 for(kk=0;kk<nOfSplits;kk++) { // loop over path segments
559 kkF = (Float_t) kk + 0.5;
560 avDrft = xL[0]+dxL[0]*kkF/nOfSplitsF;
561 avAnode = xL[2]+dxL[2]*kkF/nOfSplitsF;
564 zAnode = seg->GetAnodeFromLocal(avDrft,avAnode);
565 driftSpeed = res->GetDriftSpeedAtAnode(zAnode);
566 driftPath = TMath::Abs(10000.*avDrft);
567 driftPath = sddLength-driftPath;
569 AliDebug(1, // this should be fixed at geometry level
570 Form("negative drift path driftPath=%e sddLength=%e avDrft=%e dxL[0]=%e xL[0]=%e",
571 driftPath,sddLength,avDrft,dxL[0],xL[0]));
573 } // end if driftPath < 0
574 drTime = driftPath/driftSpeed; // drift time for segment.
575 timeSample = (Int_t) (fScaleSize*drTime/timeStep + 1); // time bin in range 1-256 !!!
576 if(timeSample > fScaleSize*fMaxNofSamples) {
577 AliWarning(Form("Wrong Time Sample: %e",timeSample));
579 } // end if timeSample > fScaleSize*fMaxNoofSamples
581 if(zAnode>nofAnodes) zAnode-=nofAnodes; // to have the anode number between 0. and 256.
582 if(zAnode*anodePitch > sddWidth || zAnode*anodePitch < 0.)
583 AliWarning(Form("Exceeding sddWidth=%e Z = %e",sddWidth,zAnode*anodePitch));
584 iAnode = (Int_t) (1.+zAnode); // iAnode in range 1-256 !!!!
585 if(iAnode < 1 || iAnode > nofAnodes) {
586 AliWarning(Form("Wrong iAnode: 1<%d>%d (xanode=%e)",iAnode,nofAnodes, zAnode));
588 } // end if iAnode < 1 || iAnode > nofAnodes
590 // store straight away the particle position in the array
591 // of particles and take idhit=ii only when part is entering (this
592 // requires FillModules() in the macro for analysis) :
594 // Sigma along the anodes for track segment.
595 sigA = TMath::Sqrt(2.*dfCoeff*drTime+s1*s1);
596 sigT = sigA/driftSpeed;
597 // Peak amplitude in nanoAmpere
598 amplitude = fScaleSize*160.*depEnergy/
599 (timeStep*eVpairs*2.*acos(-1.)*sigT*sigA);
600 amplitude *= timeStep/25.; // WARNING!!!!! Amplitude scaling to
601 // account for clock variations
602 // (reference value: 40 MHz)
603 chargeloss = 1.-cHloss*driftPath/1000.;
604 amplitude *= chargeloss;
605 width = 2.*nsigma/(nlookups-1);
613 } // end if drTime > 1200.
615 nsplit = 4; // hard-wired //nsplit=4;nsplit = (nsplit+1)/2*2;
616 // Sub-pixel size see computation of aExpo and tExpo.
617 aStep = anodePitch/(nsplit*fScaleSize*sigA);
618 aConst = zAnode*anodePitch/sigA;
619 tStep = timeStep/(nsplit*fScaleSize*sigT);
620 tConst = drTime/sigT;
621 // Define SDD window corresponding to the hit
622 anodeWindow = (Int_t)(fScaleSize*nsigma*sigA/anodePitch+1);
623 timeWindow = (Int_t) (fScaleSize*nsigma*sigT/timeStep+1.);
624 jamin = (iAnode - anodeWindow/ndiv - 2)*fScaleSize*nsplit +1;
625 jamax = (iAnode + anodeWindow/ndiv + 1)*fScaleSize*nsplit;
626 if(jamin <= 0) jamin = 1;
627 if(jamax > fScaleSize*nofAnodes*nsplit)
628 jamax = fScaleSize*nofAnodes*nsplit;
629 // jtmin and jtmax are Hard-wired
630 jtmin = (Int_t)(timeSample-timeWindow*nmul-1)*nsplit+1;
631 jtmax = (Int_t)(timeSample+timeWindow*nmul)*nsplit;
632 if(jtmin <= 0) jtmin = 1;
633 if(jtmax > fScaleSize*fMaxNofSamples*nsplit)
634 jtmax = fScaleSize*fMaxNofSamples*nsplit;
635 // Spread the charge in the anode-time window
636 for(ka=jamin; ka <=jamax; ka++) {
637 ia = (ka-1)/(fScaleSize*nsplit) + 1;
639 Warning("HitsToAnalogDigits","ia < 1: ");
642 if(ia > nofAnodes) ia = nofAnodes;
643 aExpo = (aStep*(ka-0.5)-aConst);
644 if(TMath::Abs(aExpo) > nsigma) anodeAmplitude = 0.;
646 Int_t theBin = (Int_t) ((aExpo+nsigma)/width+0.5);
647 anodeAmplitude = amplitude*res->GetGausLookUp(theBin);
648 } // end if TMath::Abs(aEspo) > nsigma
649 // index starts from 0
650 index = iWing*nofAnodes+ia-1;
652 for(kt=jtmin; kt<=jtmax; kt++) {
653 it = (kt-1)/nsplit+1; // it starts from 1
655 Warning("HitsToAnalogDigits","it < 1:");
658 if(it>fScaleSize*fMaxNofSamples)
659 it = fScaleSize*fMaxNofSamples;
660 tExpo = (tStep*(kt-0.5)-tConst);
661 if(TMath::Abs(tExpo) > nsigma) timeAmplitude = 0.;
663 Int_t theBin = (Int_t) ((tExpo+nsigma)/width+0.5);
664 timeAmplitude = anodeAmplitude*res->GetGausLookUp(theBin);
665 } // end if TMath::Abs(tExpo) > nsigma
666 // build the list of Sdigits for this module
669 // arg[2] = itrack; // track number
670 // arg[3] = ii-1; // hit number.
671 timeAmplitude *= norm;
673 // ListOfFiredCells(arg,timeAmplitude,alst,padr);
674 Double_t charge = timeAmplitude;
675 charge += fHitMap2->GetSignal(index,it-1);
676 fHitMap2->SetHit(index, it-1, charge);
677 fpList->AddSignal(index,it-1,itrack,ii-1,
678 mod->GetIndex(),timeAmplitude);
679 fAnodeFire[index] = kTRUE;
680 } // end loop over time in window
681 } // end if anodeAmplitude
682 } // loop over anodes in window
683 } // end loop over "sub-hits"
684 } // end loop over hits
687 //____________________________________________
688 void AliITSsimulationSDD::AddDigit( Int_t i, Int_t j, Int_t signal ) {
690 Int_t size = AliITSdigit::GetNTracks();
693 Int_t * tracks = new Int_t[size];
694 Int_t * hits = new Int_t[size];
696 Float_t * charges = new Float_t[size];
702 AliITSpListItem *pItem = fpList->GetpListItem( i, j );
705 for( Int_t l=0; l<size; l++ ) {
711 Int_t idtrack = pItem->GetTrack( 0 );
712 if( idtrack >= 0 ) phys = pItem->GetSignal();
715 for( Int_t l=0; l<size; l++ ) if(l<pItem->GetMaxKept()) {
716 tracks[l] = pItem->GetTrack( l );
717 hits[l] = pItem->GetHit( l );
718 charges[l] = pItem->GetSignal( l );
726 fITS->AddSimDigit( 1, phys, digits, tracks, hits, charges );
731 //______________________________________________________________________
732 void AliITSsimulationSDD::ChargeToSignal(Int_t mod,Bool_t bAddNoise, Bool_t bAddGain) {
733 // add baseline, noise, gain, electronics and ADC saturation effects
734 // apply dead channels
736 char opt1[20], opt2[20];
737 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(mod);
738 res->GetParamOptions(opt1,opt2);
744 Float_t maxadc = res->GetMaxAdc();
746 for (i=0;i<fNofMaps;i++) {
747 if( !fAnodeFire[i] ) continue;
748 baseline = res->GetBaseline(i);
749 noise = res->GetNoise(i);
750 gain = res->GetChannelGain(i);
751 if(res->IsBad()) gain=0.;
752 if( res->IsChipBad(res->GetChip(i)) ){
753 printf("Chip bad mod %d chip %d anode %d\n",mod,res->GetChip(i),i);
756 for(k=0; k<fScaleSize*fMaxNofSamples; k++) {
757 fInZR[k] = fHitMap2->GetSignal(i,k);
758 if(bAddGain) fInZR[k]*=gain;
760 contrib = (baseline + noise*gRandom->Gaus());
766 for(k=0; k<fMaxNofSamples; k++) {
767 Double_t newcont = 0.;
768 Double_t maxcont = 0.;
769 for(kk=0;kk<fScaleSize;kk++) {
770 newcont = fInZR[fScaleSize*k+kk];
771 if(newcont > maxcont) maxcont = newcont;
774 if (newcont >= maxadc) newcont = maxadc -1;
775 if(newcont >= baseline){
776 Warning("","newcont=%d>=baseline=%d",newcont,baseline);
779 fHitMap2->SetHit(i,k,newcont);
782 FastFourierTransform(fElectronics,&fInZR[0],&fInZI[0],1);
783 for(k=0; k<fScaleSize*fMaxNofSamples; k++) {
784 Double_t rw = fElectronics->GetTraFunReal(k);
785 Double_t iw = fElectronics->GetTraFunImag(k);
786 fOutZR[k] = fInZR[k]*rw - fInZI[k]*iw;
787 fOutZI[k] = fInZR[k]*iw + fInZI[k]*rw;
789 FastFourierTransform(fElectronics,&fOutZR[0],&fOutZI[0],-1);
790 for(k=0; k<fMaxNofSamples; k++) {
791 Double_t newcont1 = 0.;
792 Double_t maxcont1 = 0.;
793 for(kk=0;kk<fScaleSize;kk++) {
794 newcont1 = fOutZR[fScaleSize*k+kk];
795 if(newcont1 > maxcont1) maxcont1 = newcont1;
798 if (newcont1 >= maxadc) newcont1 = maxadc -1;
799 fHitMap2->SetHit(i,k,newcont1);
802 } // end for i loop over anodes
806 //______________________________________________________________________
807 void AliITSsimulationSDD::ApplyCrosstalk(Int_t mod) {
808 // function add the crosstalk effect to signal
809 // temporal function, should be checked...!!!
810 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
812 Int_t fNofMaps = seg->Npz();
813 Int_t fMaxNofSamples = seg->Npx();
815 // create and inizialice crosstalk map
816 Float_t* ctk = new Float_t[fNofMaps*fMaxNofSamples+1];
818 Error( "ApplyCrosstalk", "no memory for temporal map: exit \n" );
821 memset( ctk, 0, sizeof(Float_t)*(fNofMaps*fMaxNofSamples+1) );
822 AliITSCalibrationSDD* calibr = (AliITSCalibrationSDD*)GetCalibrationModel(mod);
823 for( Int_t z=0; z<fNofMaps; z++ ) {
824 Double_t baseline = calibr->GetBaseline(z);
830 for( Int_t l=0; l<fMaxNofSamples; l++ ) {
831 Float_t fadc = (Float_t)fHitMap2->GetSignal( z, l );
832 if( fadc > baseline ) {
833 if( on == kFALSE && l<fMaxNofSamples-4 ) {
834 Float_t fadc1 = (Float_t)fHitMap2->GetSignal( z, l+1 );
835 if( fadc1 < fadc ) continue;
842 else { // end fadc > baseline
846 // make smooth derivative
847 Float_t* dev = new Float_t[fMaxNofSamples+1];
848 memset( dev, 0, sizeof(Float_t)*(fMaxNofSamples+1) );
850 Error( "ApplyCrosstalk",
851 "no memory for temporal array: exit \n" );
854 for( Int_t i=tstart; i<tstop; i++ ) {
855 if( i > 2 && i < fMaxNofSamples-2 )
856 dev[i] = -0.2*fHitMap2->GetSignal( z,i-2 )
857 -0.1*fHitMap2->GetSignal( z,i-1 )
858 +0.1*fHitMap2->GetSignal( z,i+1 )
859 +0.2*fHitMap2->GetSignal( z,i+2 );
862 // add crosstalk contribution to neibourg anodes
863 for( Int_t i=tstart; i<tstop; i++ ) {
865 Int_t i1 = (Int_t)((i-tstart)*.61+tstart+0.5); //
866 Float_t ctktmp = -dev[i1] * 0.25;
868 ctk[anode*fMaxNofSamples+i] += ctktmp;
871 if( anode < fNofMaps ) {
872 ctk[anode*fMaxNofSamples+i] += ctktmp;
877 } // if( nTsteps > 2 )
879 } // if( on == kTRUE )
884 for( Int_t a=0; a<fNofMaps; a++ )
885 for( Int_t t=0; t<fMaxNofSamples; t++ ) {
886 Float_t signal = fHitMap2->GetSignal(a,t)+ctk[a*fMaxNofSamples+t];
887 fHitMap2->SetHit( a, t, signal );
893 //______________________________________________________________________
894 void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl,
896 // Returns the compression alogirthm parameters
901 //______________________________________________________________________
902 void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl) const{
903 // returns the compression alogirthm parameters
909 //______________________________________________________________________
910 void AliITSsimulationSDD::SetCompressParam(){
911 // Sets the compression alogirthm parameters
912 AliITSCalibrationSDD* calibr = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
913 for(Int_t ian = 0; ian<fNofMaps;ian++){
914 fD[ian] = (Int_t)(calibr->GetBaseline(ian));
915 fT1[ian] = (Int_t)(2.*calibr->GetNoiseAfterElectronics(ian)+0.5);
916 fT2[ian] = 0; // used by 2D clustering - not defined yet
917 fTol[ian] = 0; // used by 2D clustering - not defined yet
920 //______________________________________________________________________
921 Int_t AliITSsimulationSDD::Convert10to8(Int_t signal) const {
922 // To the 10 to 8 bit lossive compression.
923 // code from Davide C. and Albert W.
925 if (signal < 128) return signal;
926 if (signal < 256) return (128+((signal-128)>>1));
927 if (signal < 512) return (192+((signal-256)>>3));
928 if (signal < 1024) return (224+((signal-512)>>4));
932 //______________________________________________________________________
933 AliITSMap* AliITSsimulationSDD::HitMap(Int_t i){
934 //Return the correct map.
936 return ((i==0)? fHitMap1 : fHitMap2);
939 //______________________________________________________________________
940 void AliITSsimulationSDD::ZeroSuppression(const char *option) {
941 // perform the zero suppresion
942 if (strstr(option,"2D")) {
943 //Init2D(); // activate if param change module by module
945 } else if (strstr(option,"1D")) {
946 //Init1D(); // activate if param change module by module
948 } else StoreAllDigits();
950 //______________________________________________________________________
951 void AliITSsimulationSDD::Init2D(){
952 // read in and prepare arrays: fD, fT1, fT2
953 // savemu[nanodes], savesigma[nanodes]
954 // read baseline and noise from file - either a .root file and in this
955 // case data should be organised in a tree with one entry for each
956 // module => reading should be done accordingly
957 // or a classic file and do smth. like this ( code from Davide C. and
959 // Read 2D zero-suppression parameters for SDD
961 if (!strstr(fParam.Data(),"file")) return;
965 Float_t *savemu = new Float_t [fNofMaps];
966 Float_t *savesigma = new Float_t [fNofMaps];
967 char input[100],basel[100],par[100];
970 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
972 res->Thresholds(tmp1,tmp2);
973 Int_t minval = static_cast<Int_t>(tmp1);
975 res->Filenames(input,basel,par);
978 filtmp = gSystem->ExpandPathName(fFileName.Data());
979 FILE *param = fopen(filtmp,"r");
983 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
985 Error("Init2D","Anode number not in increasing order!",filtmp);
987 } // end if pos != na+1
989 savesigma[na] = sigma;
990 if ((2.*sigma) < mu) {
991 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
994 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
995 if (tempTh < 0) tempTh=0;
997 tempTh = (Int_t)floor(mu+3.0*sigma+0.5) - minval;
998 if (tempTh < 0) tempTh=0;
1003 Error("Init2D","THE FILE %s DOES NOT EXIST !",filtmp);
1010 delete [] savesigma;
1012 //______________________________________________________________________
1013 void AliITSsimulationSDD::Compress2D(){
1014 // simple ITS cluster finder -- online zero-suppression conditions
1018 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1020 res->Thresholds(tmp1,tmp2);
1021 Int_t minval = static_cast<Int_t>(tmp1);
1022 Bool_t write = res->OutputOption();
1023 Bool_t do10to8 = res->Do10to8();
1024 Int_t nz, nl, nh, low, i, j;
1026 for (i=0; i<fNofMaps; i++) {
1027 CompressionParam(i,db,tl,th);
1032 for (j=0; j<fMaxNofSamples; j++) {
1033 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1034 signal -= db; // if baseline eq. is done here
1035 if (signal <= 0) {nz++; continue;}
1036 if ((signal - tl) < minval) low++;
1037 if ((signal - th) >= minval) {
1040 FindCluster(i,j,signal,minval,cond);
1042 ((TMath::Abs(fHitMap2->GetSignal(i,j-1))-th)>=minval)){
1043 if(do10to8) signal = Convert10to8(signal);
1044 AddDigit(i,j,signal);
1045 } // end if cond&&j&&()
1046 } else if ((signal - tl) >= minval) nl++;
1047 } // end for j loop time samples
1048 if (write) TreeB()->Fill(nz,nl,nh,low,i+1);
1049 } //end for i loop anodes
1053 sprintf(hname,"TNtuple%d_%d",fModule,fEvent);
1054 TreeB()->Write(hname);
1059 //______________________________________________________________________
1060 void AliITSsimulationSDD::FindCluster(Int_t i,Int_t j,Int_t signal,
1061 Int_t minval,Bool_t &cond){
1062 // Find clusters according to the online 2D zero-suppression algorithm
1063 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1064 AliITSsegmentationSDD* seg = (AliITSsegmentationSDD*)GetSegmentationModel(1);
1066 Bool_t do10to8 = res->Do10to8();
1067 Bool_t high = kFALSE;
1069 fHitMap2->FlagHit(i,j);
1071 // check the online zero-suppression conditions
1073 const Int_t kMaxNeighbours = 4;
1076 Int_t xList[kMaxNeighbours], yList[kMaxNeighbours];
1077 seg->Neighbours(i,j,&nn,xList,yList);
1079 for (in=0; in<nn; in++) {
1082 if (fHitMap2->TestHit(ix,iy)==kUnused) {
1083 CompressionParam(ix,dbx,tlx,thx);
1084 Int_t qn = (Int_t)(fHitMap2->GetSignal(ix,iy));
1085 qn -= dbx; // if baseline eq. is done here
1086 if ((qn-tlx) < minval) {
1087 fHitMap2->FlagHit(ix,iy);
1090 if ((qn - thx) >= minval) high=kTRUE;
1092 if(do10to8) signal = Convert10to8(signal);
1093 AddDigit(i,j,signal);
1095 if(do10to8) qns = Convert10to8(qn);
1097 if (!high) AddDigit(ix,iy,qns);
1099 if(!high) fHitMap2->FlagHit(ix,iy);
1100 } // end if qn-tlx < minval
1102 } // end for in loop over neighbours
1104 //______________________________________________________________________
1105 void AliITSsimulationSDD::Init1D(){
1106 // this is just a copy-paste of input taken from 2D algo
1107 // Torino people should give input
1108 // Read 1D zero-suppression parameters for SDD
1110 if (!strstr(fParam.Data(),"file")) return;
1112 Int_t na,pos,tempTh;
1114 Float_t *savemu = new Float_t [fNofMaps];
1115 Float_t *savesigma = new Float_t [fNofMaps];
1116 char input[100],basel[100],par[100];
1119 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1121 res->Thresholds(tmp1,tmp2);
1122 Int_t minval = static_cast<Int_t>(tmp1);
1124 res->Filenames(input,basel,par);
1127 // set first the disable and tol param
1130 filtmp = gSystem->ExpandPathName(fFileName.Data());
1131 FILE *param = fopen(filtmp,"r");
1135 fscanf(param,"%d %d %d %d ", &fT2[0], &fT2[1], &fTol[0], &fTol[1]);
1136 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
1138 Error("Init1D","Anode number not in increasing order!",filtmp);
1140 } // end if pos != na+1
1142 savesigma[na]=sigma;
1143 if ((2.*sigma) < mu) {
1144 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
1147 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
1148 if (tempTh < 0) tempTh=0;
1153 Error("Init1D","THE FILE %s DOES NOT EXIST !",filtmp);
1160 delete [] savesigma;
1162 //______________________________________________________________________
1163 void AliITSsimulationSDD::Compress1D(){
1164 // 1D zero-suppression algorithm (from Gianluca A.)
1165 Int_t dis,tol,thres,decr,diff;
1166 UChar_t *str=fStream->Stream();
1168 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1170 Bool_t do10to8=res->Do10to8();
1174 for (k=0; k<2; k++) {
1177 for (i=0; i<fNofMaps/2; i++) {
1178 Bool_t firstSignal=kTRUE;
1179 Int_t idx=i+k*fNofMaps/2;
1180 if( !fAnodeFire[idx] ) continue;
1181 CompressionParam(idx,decr,thres);
1183 for (j=0; j<fMaxNofSamples; j++) {
1184 Int_t signal=(Int_t)(fHitMap2->GetSignal(idx,j));
1185 signal -= decr; // if baseline eq.
1186 if(do10to8) signal = Convert10to8(signal);
1187 if (signal <= thres) {
1191 // write diff in the buffer for HuffT
1192 str[counter]=(UChar_t)diff;
1195 } // end if signal <= thres
1197 if (diff > 127) diff=127;
1198 if (diff < -128) diff=-128;
1200 // tol has changed to 8 possible cases ? - one can write
1201 // this if(TMath::Abs(diff)<tol) ... else ...
1202 if(TMath::Abs(diff)<tol) diff=0;
1203 // or keep it as it was before
1204 AddDigit(idx,j,last+diff);
1206 AddDigit(idx,j,signal);
1207 } // end if singal < dis
1209 // write diff in the buffer used to compute Huffman tables
1210 if (firstSignal) str[counter]=(UChar_t)signal;
1211 else str[counter]=(UChar_t)diff;
1215 } // end for j loop time samples
1216 } // end for i loop anodes one half of detector
1220 fStream->CheckCount(counter);
1222 // open file and write out the stream of diff's
1223 static Bool_t open=kTRUE;
1224 static TFile *outFile;
1225 Bool_t write = res->OutputOption();
1226 TDirectory *savedir = gDirectory;
1230 SetFileName("stream.root");
1231 cout<<"filename "<<fFileName<<endl;
1232 outFile=new TFile(fFileName,"recreate");
1233 cout<<"I have opened "<<fFileName<<" file "<<endl;
1240 fStream->ClearStream();
1242 // back to galice.root file
1243 if(savedir) savedir->cd();
1245 //______________________________________________________________________
1246 void AliITSsimulationSDD::StoreAllDigits(){
1247 // if non-zero-suppressed data
1248 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1250 Bool_t do10to8 = res->Do10to8();
1251 Int_t i, j, digits[3];
1253 for (i=0; i<fNofMaps; i++) {
1254 for (j=0; j<fMaxNofSamples; j++) {
1255 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1256 if(do10to8) signal = Convert10to8(signal);
1260 fITS->AddRealDigit(1,digits);
1264 //______________________________________________________________________
1265 void AliITSsimulationSDD::CreateHistograms(Int_t scale){
1266 // Creates histograms of maps for debugging
1269 fHis=new TObjArray(fNofMaps);
1270 for (i=0;i<fNofMaps;i++) {
1271 TString sddName("sdd_");
1273 sprintf(candNum,"%d",i+1);
1274 sddName.Append(candNum);
1275 fHis->AddAt(new TH1F(sddName.Data(),"SDD maps",scale*fMaxNofSamples,
1276 0.,(Float_t) scale*fMaxNofSamples), i);
1279 //______________________________________________________________________
1280 void AliITSsimulationSDD::FillHistograms(){
1281 // fill 1D histograms from map
1285 for( Int_t i=0; i<fNofMaps; i++) {
1286 TH1F *hist =(TH1F *)fHis->UncheckedAt(i);
1287 Int_t nsamples = hist->GetNbinsX();
1288 for( Int_t j=0; j<nsamples; j++) {
1289 Double_t signal=fHitMap2->GetSignal(i,j);
1290 hist->Fill((Float_t)j,signal);
1294 //______________________________________________________________________
1295 void AliITSsimulationSDD::ResetHistograms(){
1296 // Reset histograms for this detector
1299 for (i=0;i<fNofMaps;i++ ) {
1300 if (fHis->At(i)) ((TH1F*)fHis->At(i))->Reset();
1303 //______________________________________________________________________
1304 TH1F *AliITSsimulationSDD::GetAnode(Int_t wing, Int_t anode) {
1305 // Fills a histogram from a give anode.
1307 if (!fHis) return 0;
1309 if(wing <=0 || wing > 2) {
1310 Warning("GetAnode","Wrong wing number: %d",wing);
1312 } // end if wing <=0 || wing >2
1313 if(anode <=0 || anode > fNofMaps/2) {
1314 Warning("GetAnode","Wrong anode number: %d",anode);
1316 } // end if ampde <=0 || andoe > fNofMaps/2
1318 Int_t index = (wing-1)*fNofMaps/2 + anode-1;
1319 return (TH1F*)(fHis->At(index));
1321 //______________________________________________________________________
1322 void AliITSsimulationSDD::WriteToFile(TFile *hfile) {
1323 // Writes the histograms to a file
1329 for(i=0; i<fNofMaps; i++) fHis->At(i)->Write(); //fAdcs[i]->Write();
1332 //______________________________________________________________________
1333 Float_t AliITSsimulationSDD::GetNoise() {
1334 // Returns the noise value
1335 //Bool_t do10to8=GetResp()->Do10to8();
1336 //noise will always be in the liniar part of the signal
1338 Int_t threshold = fT1[0];
1339 char opt1[20], opt2[20];
1340 AliITSCalibrationSDD* res = (AliITSCalibrationSDD*)GetCalibrationModel(fModule);
1342 res->GetParamOptions(opt1,opt2);
1344 Double_t noise,baseline;
1345 //GetBaseline(fModule);
1347 TCanvas *c2 = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("c2");
1348 if(c2) delete c2->GetPrimitive("noisehist");
1349 if(c2) delete c2->GetPrimitive("anode");
1350 else c2=new TCanvas("c2");
1352 c2->SetFillColor(0);
1354 TH1F *noisehist = new TH1F("noisehist","noise",100,0.,(float)2*threshold);
1355 TH1F *anode = new TH1F("anode","Anode Projection",fMaxNofSamples,0.,
1356 (float)fMaxNofSamples);
1358 for (i=0;i<fNofMaps;i++) {
1359 CompressionParam(i,decr,threshold);
1360 baseline = res->GetBaseline(i);
1361 noise = res->GetNoise(i);
1363 for (k=0;k<fMaxNofSamples;k++) {
1364 Float_t signal=(Float_t)fHitMap2->GetSignal(i,k);
1365 //if (signal <= (float)threshold) noisehist->Fill(signal-baseline);
1366 if (signal <= (float)(threshold+decr)) noisehist->Fill(signal);
1367 anode->Fill((float)k,signal);
1372 TF1 *gnoise = new TF1("gnoise","gaus",0.,threshold);
1373 noisehist->Fit("gnoise","RQ");
1376 Float_t mnoise = gnoise->GetParameter(1);
1377 cout << "mnoise : " << mnoise << endl;
1378 Float_t rnoise = gnoise->GetParameter(2);
1379 cout << "rnoise : " << rnoise << endl;
1383 //______________________________________________________________________
1384 void AliITSsimulationSDD::WriteSDigits(){
1385 // Fills the Summable digits Tree
1386 static AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
1388 for( Int_t i=0; i<fNofMaps; i++ ) {
1389 if( !fAnodeFire[i] ) continue;
1390 for( Int_t j=0; j<fMaxNofSamples; j++ ) {
1391 Double_t sig = fHitMap2->GetSignal( i, j );
1393 Int_t jdx = j*fScaleSize;
1394 Int_t index = fpList->GetHitIndex( i, j );
1395 AliITSpListItem pItemTmp2( fModule, index, 0. );
1396 // put the fScaleSize analog digits in only one
1397 for( Int_t ik=0; ik<fScaleSize; ik++ ) {
1398 AliITSpListItem *pItemTmp = fpList->GetpListItem(i,jdx+ik);
1399 if( pItemTmp == 0 ) continue;
1400 pItemTmp2.Add( pItemTmp );
1402 pItemTmp2.AddSignalAfterElect( fModule, index, sig );
1403 pItemTmp2.AddNoise(fModule,index,fHitNoiMap2->GetSignal(i,j));
1404 aliITS->AddSumDigit( pItemTmp2 );
1405 } // end if (sig > 0.2)
1410 //______________________________________________________________________
1411 void AliITSsimulationSDD::PrintStatus() const {
1412 // Print SDD simulation Parameters
1414 cout << "**************************************************" << endl;
1415 cout << " Silicon Drift Detector Simulation Parameters " << endl;
1416 cout << "**************************************************" << endl;
1417 cout << "Flag for Perpendicular tracks: " << (Int_t) fFlag << endl;
1418 cout << "Flag for noise checking: " << (Int_t) fCheckNoise << endl;
1419 cout << "Flag to switch off electronics: " << (Int_t) fDoFFT << endl;
1420 cout << "Number pf Anodes used: " << fNofMaps << endl;
1421 cout << "Number of Time Samples: " << fMaxNofSamples << endl;
1422 cout << "Scale size factor: " << fScaleSize << endl;
1423 cout << "**************************************************" << endl;