New ITS code for new structure and simulations.
[u/mrichter/AliRoot.git] / ITS / AliITSsimulationSDD.cxx
CommitLineData
b0f5e3fc 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16
17#include <iostream.h>
18#include <stdlib.h>
19#include <stdio.h>
20#include <string.h>
21
22
23
24#include "AliRun.h"
25#include "AliITSetfSDD.h"
26#include "AliITSsimulationSDD.h"
27#include "AliITSHuffman.h"
28
29const Int_t kMaxNeighbours = 4;
30
31ClassImp(AliITSsimulationSDD)
32////////////////////////////////////////////////////////////////////////
33// Version: 0
34// Written by Piergiorgio Cerello
35// November 23 1999
36//
37// AliITSsimulationSDD is the simulation of SDDs.
38 //
39//Begin_Html
40/*
41<img src="picts/ITS/AliITShit_Class_Diagram.gif">
42</pre>
43<br clear=left>
44<font size=+2 color=red>
45<p>This show the relasionships between the ITS hit class and the rest of Aliroot.
46</font>
47<pre>
48*/
49//End_Html
50//_____________________________________________________________________________
51
52Int_t power(Int_t b, Int_t e) {
53 // copute b to the e power, where bothe b and e are Int_ts.
54 Int_t power = 1,i;
55 for(i=0; i<e; i++) power *= b;
56 return power;
57}
58
59//_____________________________________________
60
61void FastFourierTransform(AliITSetfSDD *alisddetf,Double_t *real,
62 Double_t *imag,Int_t direction) {
63 // Do a Fast Fourier Transform
64
65 Int_t samples = alisddetf->GetSamples();
66 Int_t l = (Int_t) ((log((Float_t) samples)/log(2.))+0.5);
67 Int_t m1 = samples;
68 Int_t m = samples/2;
69 Int_t m2 = samples/m1;
70 Int_t i,j,k;
71 for(i=1; i<=l; i++) {
72 for(j=0; j<samples; j += m1) {
73 Int_t p = 0;
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);
82 *(real+k) += xr;
83 *(imag+k) += xi;
84 p += m2;
85 }
86 }
87 m1 = m;
88 m /= 2;
89 m2 += m2;
90 }
91
92 for(j=0; j<samples; j++) {
93 Int_t j1 = j;
94 Int_t p = 0;
95 Int_t i1;
96 for(i1=1; i1<=l; i1++) {
97 Int_t j2 = j1;
98 j1 /= 2;
99 p = p + p + j2 - j1 - j1;
100 }
101 if(p >= j) {
102 Double_t xr = *(real+j);
103 Double_t xi = *(imag+j);
104 *(real+j) = *(real+p);
105 *(imag+j) = *(imag+p);
106 *(real+p) = xr;
107 *(imag+p) = xi;
108 }
109 }
110 if(direction == -1) {
111 for(i=0; i<samples; i++) {
112 *(real+i) /= samples;
113 *(imag+i) /= samples;
114 }
115 }
116 return;
117}
118//_____________________________________________________________________________
119
120AliITSsimulationSDD::AliITSsimulationSDD(){
121 // Default constructor
122
123 fResponse = 0;
124 fSegmentation = 0;
125 fHis = 0;
126 fD.Set(0);
127 fT1.Set(0);
128 fT2.Set(0);
129 fTol.Set(0);
130}
131//_____________________________________________________________________________
132AliITSsimulationSDD::AliITSsimulationSDD(AliITSsimulationSDD &source){
133 // Copy constructor to satify Coding roules only.
134 if(this==&source) return;
135 printf("Not allowed to make a copy of AliITSsimulationSDD "
136 "Using default creater instead\n");
137 AliITSsimulationSDD();
138}
139//_____________________________________________________________________________
140AliITSsimulationSDD& AliITSsimulationSDD::operator=(AliITSsimulationSDD
141 &source){
142 // Copy constructor to satify Coding roules only.
143 if(this==&source) return *this;
144 printf("Not allowed to make a = with AliITSsimulationSDD "
145 "Using default creater instead\n");
146 return *this ;
147}
148//_____________________________________________________________________________
149
150AliITSsimulationSDD::AliITSsimulationSDD(AliITSsegmentation *seg,
151 AliITSresponse *resp) {
152 // Constructor
153 fResponse = resp;
154 fSegmentation = seg;
155
156 fHitMap2 = new AliITSMapA2(fSegmentation);
157 fHitMap1 = new AliITSMapA1(fSegmentation);
158
159 //
160 fNofMaps=fSegmentation->Npz();
161 fMaxNofSamples=fSegmentation->Npx();
162
163 Float_t sddLength = fSegmentation->Dx();
164 Float_t sddWidth = fSegmentation->Dz();
165
166 Int_t dummy=0;
167 Float_t anodePitch = fSegmentation->Dpz(dummy);
168 Double_t timeStep = (Double_t)fSegmentation->Dpx(dummy);
169 Float_t driftSpeed=fResponse->DriftSpeed();
170
171 if(anodePitch*(fNofMaps/2) > sddWidth) {
172 Warning("AliITSsimulationSDD",
173 "Too many anodes %d or too big pitch %f \n",fNofMaps/2,anodePitch);
174 }
175
176 if(timeStep*fMaxNofSamples < sddLength/driftSpeed) {
177 Error("AliITSsimulationSDD",
178 "Time Interval > Allowed Time Interval: exit\n");
179 return;
180 }
181
182 fElectronics = new AliITSetfSDD(timeStep);
183
184 Option_t *opt1, *opt2;
185 fResponse->ParamOptions(opt1,opt2);
186 fParam=opt2;
187 char *same = strstr(opt1,"same");
188 if (same) {
189 fNoise.Set(0);
190 fBaseline.Set(0);
191 } else {
192 fNoise.Set(fNofMaps);
193 fBaseline.Set(fNofMaps);
194 }
195
196 //
197 Option_t *opt=fResponse->ZeroSuppOption();
198 if (strstr(fParam,"file") ) {
199 fD.Set(fNofMaps);
200 fT1.Set(fNofMaps);
201 if (strstr(opt,"2D")) {
202 fT2.Set(fNofMaps);
203 fTol.Set(0);
204 Init2D(); // desactivate if param change module by module
205 } else if(strstr(opt,"1D")) {
206 fT2.Set(2);
207 fTol.Set(2);
208 Init1D(); // desactivate if param change module by module
209 }
210 } else {
211 fD.Set(2);
212 fTol.Set(2);
213 fT1.Set(2);
214 fT2.Set(2);
215 SetCompressParam();
216 }
217
218
219 Bool_t write=fResponse->OutputOption();
220 if(write && strstr(opt,"2D")) MakeTreeB();
221
222 // call here if baseline does not change by module
223 // ReadBaseline();
224
225 fITS = (AliITS*)gAlice->GetModule("ITS");
226 Int_t size=fNofMaps*fMaxNofSamples;
227 fStream = new AliITSInStream(size);
228
229}
230
231
232//_____________________________________________________________________________
233
234AliITSsimulationSDD::~AliITSsimulationSDD() {
235 // destructor
236
237 delete fHitMap1;
238 delete fHitMap2;
239 delete fStream;
240
241 fD.Set(0);
242 fT1.Set(0);
243 fT2.Set(0);
244 fTol.Set(0);
245 fNoise.Set(0);
246 fBaseline.Set(0);
247
248 if (fHis) {
249 fHis->Delete();
250 delete fHis;
251 }
252}
253//_____________________________________________________________________________
254
255void AliITSsimulationSDD::DigitiseModule(AliITSmodule *mod,Int_t md,Int_t ev){
256 // create maps to build the lists of tracks
257 // for each digit
258
259 fModule=md;
260 fEvent=ev;
261
262 TObjArray *fHits = mod->GetHits();
263 Int_t nhits = fHits->GetEntriesFast();
264 if (!nhits) return;
265
266
267 TObjArray *list=new TObjArray;
268 static TClonesArray *padr=0;
269 if(!padr) padr=new TClonesArray("TVector",1000);
270 Int_t arg[5] = {0,0,0,0,0};
271 fHitMap1->SetArray(list);
272
273
274 Int_t NofAnodes=fNofMaps/2;
275
276 Float_t sddLength = fSegmentation->Dx();
277 Float_t sddWidth = fSegmentation->Dz();
278
279 Int_t dummy=0;
280 Float_t anodePitch = fSegmentation->Dpz(dummy);
281 Float_t timeStep = fSegmentation->Dpx(dummy);
282
283 Float_t driftSpeed=fResponse->DriftSpeed();
284
285 // Piergiorgio's part (apart for few variables which I made float
286 // when i thought that can be done
287
288 // Fill detector maps with GEANT hits
289 // loop over hits in the module
290
291 const Float_t kconv=1000000.; // GeV->KeV
292 Int_t ii;
293 for(ii=0; ii<nhits; ii++) {
294 AliITShit *hit = (AliITShit*) fHits->At(ii);
295 Int_t hitDetector = hit->GetDetector();
296 Float_t xL[3];
297 hit->GetPositionL(xL[0],xL[1],xL[2]);
298 // cout << "hit local coordinates: " << xL[0] << "," << xL[1] << "," << xL[2] << endl;
299 // Deposited energy in keV
300 Float_t avpath = 0.;
301 Float_t avanod = 0.;
302 Float_t depEnergy = kconv*hit->GetIonization();
303 AliITShit *hit1 = 0;
304 if(depEnergy == 0.) {
305 ii++;
306 Float_t xL1[3];
307 hit1 = (AliITShit*) fHits->At(ii);
308 hit1->GetPositionL(xL1[0],xL1[1],xL1[2]);
309 //cout << "hit1 local coordinates: " << xL1[0] << "," << xL1[1] << "," << xL1[2] << endl;
310 //cout << "radius1: " << TMath::Sqrt(xL1[0]*xL1[0]+xL1[1]*xL1[1]) << ", azimuth: " << TMath::ATan2(xL1[0],xL1[1]) << endl;
311 avpath = xL1[0];
312 avanod = xL1[2];
313 depEnergy = kconv*hit1->GetIonization();
314 }
315 Float_t avDrft = xL[0]+avpath;
316 Float_t avAnode = xL[2]+avanod;
317
318 if(avpath != 0.) avDrft /= 2.;
319 if(avanod != 0.) avAnode /= 2.;
320
321 Float_t driftPath = 10000.*avDrft;
322 //printf("sddLength %f avDrft driftPath %f %f\n",sddLength,avDrft, driftPath);
323 Int_t iWing = 2;
324 if(driftPath < 0) {
325 iWing = 1;
326 driftPath = -driftPath;
327 }
328 driftPath = sddLength-driftPath;
329 Int_t detector = 2*(hitDetector-1) + iWing;
330 if(driftPath < 0) {
331 cout << "Warning: negative drift path " << driftPath << endl;
332 continue;
333 }
334
335 // Drift Time
336 Float_t driftTime = driftPath/driftSpeed;
337 Int_t timeSample = (Int_t) (driftTime/timeStep + 1);
338 if(timeSample > fMaxNofSamples) {
339 cout << "Warning: Wrong Time Sample: " << timeSample << endl;
340 continue;
341 }
342
343 // Anode
344 Float_t xAnode = 10000.*(avAnode)/anodePitch + NofAnodes/2; // +1?
345 // Int_t iAnode = 0.5+xAnode; // xAnode?
346 if(xAnode*anodePitch > sddWidth || xAnode*anodePitch < 0.)
347 { cout << "Warning: Z = " << xAnode*anodePitch << endl; }
348 Int_t iAnode = (Int_t) (1.+xAnode); // xAnode?
349 // cout << "iAnode " << iAnode << endl;
350 if(iAnode < 0 || iAnode > NofAnodes) {
351 cout << "Warning: Wrong iAnode: " << iAnode << endl;
352 continue;
353 }
354
355
356 // work with the idtrack=entry number in the TreeH
357 // Int_t idtrack=mod->GetHitTrackIndex(ii);
358 // or store straight away the particle position in the array
359 // of particles :
360 Int_t idtrack = hit->GetTrack();
361
362 // Signal 2d Shape
363 Double_t qRef = (Double_t)fResponse->Qref();
364 Double_t diffCoeff = (Double_t)fResponse->DiffCoeff();
365
366 Double_t gamma = 1. + 0.155*depEnergy/qRef;
367 // Squared Sigma along the anodes
368 Double_t sigma2A = 2.*diffCoeff*driftTime*gamma;
369 Double_t sigmaT = TMath::Sqrt(sigma2A)/driftSpeed;
370
371 // Peak amplitude in nanoAmpere
372 Double_t eVpairs = 3.6;
373 Double_t amplitude = 160.*depEnergy/(timeStep*eVpairs*2.*acos(-1.)*sigmaT*TMath::Sqrt(sigma2A));
374
375 // Spread the charge
376 // Pixel index
377 Int_t ja = iAnode;
378 Int_t jt = timeSample;
379 // Sub-pixel index
380 Int_t nsplit = 8;
381 nsplit = (nsplit+1)/2*2;
382 // Sub-pixel size
383 Double_t aStep = anodePitch/nsplit;
384 Double_t tStep = timeStep/nsplit;
385 // Define SDD window corresponding to the hit
386 Int_t anodeWindow = (Int_t) (4.*TMath::Sqrt(sigma2A)/anodePitch + 1);
387 Int_t timeWindow = (Int_t) (4.*sigmaT/timeStep + 1);
388 Int_t jamin = (ja - anodeWindow/2 - 1)*nsplit + 1;
389 Int_t jamax = (ja + anodeWindow/2)*nsplit;
390 if(jamin <= 0) jamin = 1;
391 if(jamax > NofAnodes*nsplit) jamax = NofAnodes*nsplit;
392 Int_t jtmin = (jt - timeWindow/2 - 1)*nsplit + 1;
393 Int_t jtmax = (jt + timeWindow/2)*nsplit;
394 if(jtmin <= 0) jtmin = 1;
395 if(jtmax > fMaxNofSamples*nsplit) jtmax = fMaxNofSamples*nsplit;
396 Double_t rlAnode = log(aStep*amplitude);
397 // Spread the charge in the anode-time window
398 Int_t ka;
399 for(ka=jamin; ka <=jamax; ka++) {
400 Int_t ia = (ka-1)/nsplit + 1;
401 if(ia <= 0) { cout << "Warning: ia < 1: " << endl; continue; }
402 if(ia > NofAnodes) ia = NofAnodes;
403 Double_t aExpo = aStep*(ka)-xAnode*anodePitch;
404 Double_t anodeAmplitude = rlAnode - 0.5*aExpo*aExpo/sigma2A;
405 // Protect against overflows
406 if(anodeAmplitude > -87.3)
407 anodeAmplitude = exp(anodeAmplitude);
408 else
409 anodeAmplitude = 0;
410 if(anodeAmplitude) {
411 Double_t rlTime = log(tStep*anodeAmplitude);
412 Int_t kt;
413 for(kt=jtmin; kt<=jtmax; kt++) {
414 Int_t it = (kt-1)/nsplit+1;
415 if(it<=0) { cout << "Warning: it < 1: " << endl; continue; }
416 if(it>fMaxNofSamples) it = fMaxNofSamples;
417 Double_t tExpo = (tStep*(kt)-driftTime)/sigmaT;
418 Double_t timeAmplitude = rlTime - 0.5*tExpo*tExpo;
419 // Protect against overflows
420 if(timeAmplitude > -87.3)
421 timeAmplitude = exp(timeAmplitude);
422 else
423 timeAmplitude = 0;
424
425 Int_t index = ((detector+1)%2)*NofAnodes+ia-1;
426 // build the list of digits for this module
427 arg[0]=index;
428 arg[1]=it;
429 arg[2]=idtrack;
430 ListOfFiredCells(arg,timeAmplitude,list,padr);
431 } // loop over time in window
432 } // end if anodeAmplitude
433 } // loop over anodes in window
434 } // end loop over hits
435
436 Int_t nentries=list->GetEntriesFast();
437 // introduce the electronics effects and do zero-suppression if required
438 if (nentries) {
439 ChargeToSignal();
440
441 Option_t *opt=fResponse->ZeroSuppOption();
442 ZeroSuppression(opt);
443 }
444
445 // clean memory
446 list->Delete();
447 delete list;
448
449 padr->Delete();
450
451 fHitMap1->ClearMap();
452 fHitMap2->ClearMap();
453
454 //gObjectTable->Print();
455}
456
457
458//____________________________________________
459
460void AliITSsimulationSDD::ListOfFiredCells(Int_t *arg,Double_t timeAmplitude,
461 TObjArray *list,TClonesArray *padr){
462 // Returns the list of "fired" cells.
463
464 Int_t index=arg[0];
465 Int_t it=arg[1];
466 Int_t idtrack=arg[2];
467 Int_t counter=arg[3];
468 Int_t countadr=arg[4];
469
470 Int_t digits[3];
471
472 digits[0]=index;
473 digits[1]=it-1;
474 digits[2]=(Int_t)timeAmplitude;
475 Float_t phys;
476 if (idtrack >= 0) phys=(Float_t)timeAmplitude;
477 else phys=0;
478
479 Double_t charge=timeAmplitude;
480 AliITSTransientDigit* pdigit;
481 // build the list of fired cells and update the info
482 if (!fHitMap1->TestHit(index, it-1)) {
483
484 new((*padr)[countadr++]) TVector(2);
485 TVector &trinfo=*((TVector*) (*padr)[countadr-1]);
486 trinfo(0)=(Float_t)idtrack;
487 trinfo(1)=(Float_t)timeAmplitude;
488
489 list->AddAtAndExpand(
490 new AliITSTransientDigit(phys,digits),counter);
491
492 fHitMap1->SetHit(index, it-1, counter);
493 fHitMap2->SetHit(index, it-1, charge);
494 counter++;
495
496 pdigit=(AliITSTransientDigit*)list->
497 At(list->GetLast());
498 // list of tracks
499 TObjArray *trlist=(TObjArray*)pdigit->TrackList();
500 trlist->Add(&trinfo);
501
502 } else {
503 pdigit=
504 (AliITSTransientDigit*) fHitMap1->GetHit(index, it-1);
505 charge += fHitMap2->GetSignal(index,it-1);
506 fHitMap2->SetHit(index, it-1, charge);
507 // update charge
508 (*pdigit).fSignal=(Int_t)charge;
509 (*pdigit).fPhysics+=phys;
510 // update list of tracks
511 TObjArray* trlist=(TObjArray*)pdigit->TrackList();
512 Int_t lastentry=trlist->GetLast();
513 TVector *ptrkp=(TVector*)trlist->At(lastentry);
514 TVector &trinfo=*ptrkp;
515 Int_t lasttrack=Int_t(trinfo(0));
516 Float_t lastcharge=(trinfo(1));
517
518 if (lasttrack==idtrack ) {
519 lastcharge+=(Float_t)timeAmplitude;
520 trlist->RemoveAt(lastentry);
521 trinfo(0)=lasttrack;
522 trinfo(1)=lastcharge;
523 trlist->AddAt(&trinfo,lastentry);
524 } else {
525
526 new((*padr)[countadr++]) TVector(2);
527 TVector &trinfo=*((TVector*) (*padr)[countadr-1]);
528 trinfo(0)=(Float_t)idtrack;
529 trinfo(1)=(Float_t)timeAmplitude;
530
531 trlist->Add(&trinfo);
532 }
533
534#ifdef print
535 // check the track list - debugging
536 Int_t trk[50];
537 Float_t chtrk[50];
538 Int_t nptracks=trlist->GetEntriesFast();
539 if (nptracks > 2) {
540 Int_t tr;
541 for(tr=0;tr<nptracks;tr++) {
542 TVector *pptrkp=(TVector*)trlist->At(tr);
543 TVector &pptrk=*pptrkp;
544 trk[tr]=Int_t(pptrk(0));
545 chtrk[tr]=(pptrk(1));
546 printf("nptracks %d \n",nptracks);
547 // set printings
548 }
549 } // end if nptracks
550#endif
551 } // end if pdigit
552
553 arg[3]=counter;
554 arg[4]=countadr;
555
556
557}
558
559
560//____________________________________________
561
562void AliITSsimulationSDD::AddDigit(Int_t i, Int_t j, Int_t signal){
563 // Adds a Digit.
564 // tag with -1 signals coming from background tracks
565 // tag with -2 signals coming from pure electronic noise
566
567 Int_t digits[3], tracks[3];
568 Float_t phys, charges[3];
569
570 Int_t trk[20];
571 Float_t chtrk[20];
572
573 signal=Convert8to10(signal); // set a flag in case non-ZS are 10-bit
574 AliITSTransientDigit *obj = (AliITSTransientDigit*)fHitMap1->GetHit(i,j);
575 digits[0]=i;
576 digits[1]=j;
577 digits[2]=signal;
578 // printf("module anode, time, signal %d %d %d %d\n",fModule,i,j,signal);
579 if (!obj) {
580 phys=0;
581 Int_t k;
582 for(k=0;k<3;k++) {
583 tracks[k]=-2;
584 charges[k]=0;
585 }
586 fITS->AddDigit(1,phys,digits,tracks,charges);
587 } else {
588 phys=obj->fPhysics;
589 //printf("AddDigit - test: fCoord1 fCoord2 fSignal %d %d %d i j signal %d %d %d \n",obj->fCoord1,obj->fCoord2,obj->fSignal,i,j,signal);
590
591 TObjArray* trlist=(TObjArray*)obj->TrackList();
592 Int_t nptracks=trlist->GetEntriesFast();
593
594 if (nptracks > 20) {
595 cout<<"Attention - nptracks > 20 "<<nptracks<<endl;
596 nptracks=20;
597 }
598 Int_t tr;
599 for(tr=0;tr<nptracks;tr++) {
600 TVector &pp =*((TVector*)trlist->At(tr));
601 trk[tr]=Int_t(pp(0));
602 chtrk[tr]=(pp(1));
603 }
604 if (nptracks > 1) {
605 //printf("AddDigit: nptracks %d\n",nptracks);
606 SortTracks(trk,chtrk,nptracks);
607 }
608 Int_t i;
609 if (nptracks < 3 ) {
610 for(i=0; i<nptracks; i++) {
611 tracks[i]=trk[i];
612 charges[i]=chtrk[i];
613 }
614 for(i=nptracks; i<3; i++) {
615 tracks[i]=0;
616 charges[i]=0;
617 }
618 } else {
619 for(i=0; i<3; i++) {
620 tracks[i]=trk[i];
621 charges[i]=chtrk[i];
622 }
623 }
624
625 fITS->AddDigit(1,phys,digits,tracks,charges);
626
627 }
628
629}
630
631//____________________________________________
632
633void AliITSsimulationSDD::SortTracks(Int_t *tracks,Float_t *charges,Int_t ntr){
634 //
635 // Sort the list of tracks contributing to a given digit
636 // Only the 3 most significant tracks are acctually sorted
637 //
638
639 //
640 // Loop over signals, only 3 times
641 //
642
643
644 Float_t qmax;
645 Int_t jmax;
646 Int_t idx[3] = {-3,-3,-3};
647 Float_t jch[3] = {-3,-3,-3};
648 Int_t jtr[3] = {-3,-3,-3};
649 Int_t i,j,imax;
650
651 if (ntr<3) imax=ntr;
652 else imax=3;
653 for(i=0;i<imax;i++){
654 qmax=0;
655 jmax=0;
656
657 for(j=0;j<ntr;j++){
658
659 if((i == 1 && j == idx[i-1] )
660 ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue;
661
662 if(charges[j] > qmax) {
663 qmax = charges[j];
664 jmax=j;
665 }
666 }
667
668 if(qmax > 0) {
669 idx[i]=jmax;
670 jch[i]=charges[jmax];
671 jtr[i]=tracks[jmax];
672 }
673
674 }
675
676 for(i=0;i<3;i++){
677 if (jtr[i] == -3) {
678 charges[i]=0;
679 tracks[i]=0;
680 } else {
681 charges[i]=jch[i];
682 tracks[i]=jtr[i];
683 }
684 }
685
686}
687//____________________________________________
688void AliITSsimulationSDD::ChargeToSignal() {
689 // add baseline, noise, electronics and ADC saturation effects
690
691 Double_t InZR[fMaxNofSamples];
692 Double_t InZI[fMaxNofSamples];
693 Double_t OutZR[fMaxNofSamples];
694 Double_t OutZI[fMaxNofSamples];
695
696
697 Float_t maxadc = fResponse->MaxAdc();
698 Float_t TopValue = fResponse->MagicValue();
699 Float_t norm = maxadc/TopValue;
700
701
702 Option_t *opt1, *opt2;
703 fResponse->ParamOptions(opt1,opt2);
704 char *read = strstr(opt1,"file");
705
706 Float_t baseline, noise;
707
708 if (read) {
709 static Bool_t readfile=kTRUE;
710 //read baseline and noise from file
711 if (readfile) ReadBaseline();
712 readfile=kFALSE;
713 } else fResponse->GetNoiseParam(noise,baseline);
714
715 Float_t contrib=0;
716 Bool_t first=kTRUE;
717
718 TRandom *random = new TRandom();
719 Int_t i,k;
720 for(i=0;i<=fNofMaps;i++) {
721 if (read) GetAnodeBaseline(i,baseline,noise);
722 if (!first) FastFourierTransform(fElectronics,&InZR[0],&InZI[0],1);
723 for(k=0;k<fMaxNofSamples;k++) {
724 if (!first) {
725 // analog to digital ?
726 Double_t signal = OutZR[k]*norm;
727 if (signal > maxadc) signal = maxadc;
728 // back to analog: ?
729 signal /=norm;
730 //printf("ChargeToSignal: signal %f\n",signal);
731 fHitMap2->SetHit(i-1,k,signal);
732
733 Double_t rw = fElectronics->GetTraFunReal(k);
734 Double_t iw = fElectronics->GetTraFunImag(k);
735 OutZR[k] = InZR[k]*rw - InZI[k]*iw;
736 OutZI[k] = InZR[k]*iw + InZI[k]*rw;
737 if(i+1 < fNofMaps) InZR[k] = fHitMap2->GetSignal(i+1,k);
738 }
739
740 if (first) {
741 InZR[k] = fHitMap2->GetSignal(i,k);
742 }
743 InZI[k] = 0.;
744 // add baseline and noise
745 contrib = baseline + noise*random->Gaus();
746 InZR[k] += contrib;
747
748 } // loop over time
749
750 if (first) {
751 FastFourierTransform(fElectronics,&InZR[0],&InZI[0],1);
752 for(k=0; k<fMaxNofSamples; k++) {
753 Double_t rw = fElectronics->GetTraFunReal(k);
754 Double_t iw = fElectronics->GetTraFunImag(k);
755 OutZR[k] = InZR[k]*rw - InZI[k]*iw;
756 OutZI[k] = InZR[k]*iw + InZI[k]*rw;
757 InZR[k] = fHitMap2->GetSignal(i+1,k);
758 InZI[k] = 0.;
759 // add baseline and noise
760 contrib = baseline + noise*random->Gaus();
761 InZR[k] += contrib;
762 }
763 }
764 FastFourierTransform(fElectronics,&OutZR[0],&OutZI[0],-1);
765 first = kFALSE;
766 } // loop over anodes
767
768}
769
770//____________________________________________
771void AliITSsimulationSDD::GetAnodeBaseline(Int_t i,Float_t &baseline,
772 Float_t &noise){
773 // Returns the Baseline for a particular anode.
774 baseline=fBaseline[i];
775 noise=fNoise[i];
776
777}
778
779//____________________________________________
780void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl,
781 Int_t &th){
782 // Returns the compression alogirthm parameters
783 Int_t size = fD.GetSize();
784 if (size > 2 ) {
785 db=fD[i]; tl=fT1[i]; th=fT2[i];
786 } else {
787 if (size <= 2 && i>=fNofMaps/2) {
788 db=fD[1]; tl=fT1[1]; th=fT2[1];
789 } else {
790 db=fD[0]; tl=fT1[0]; th=fT2[0];
791 }
792 }
793}
794//____________________________________________
795void AliITSsimulationSDD::CompressionParam(Int_t i,Int_t &db,Int_t &tl){
796 // returns the compression alogirthm parameters
797 Int_t size = fD.GetSize();
798 if (size > 2 ) {
799 db=fD[i]; tl=fT1[i];
800 } else {
801 if (size <= 2 && i>=fNofMaps/2) {
802 db=fD[1]; tl=fT1[1];
803 } else {
804 db=fD[0]; tl=fT1[0];
805 }
806 }
807
808}
809//____________________________________________
810void AliITSsimulationSDD::SetCompressParam(){
811 // Sets the compression alogirthm parameters
812 Int_t cp[8],i;
813
814 fResponse->GiveCompressParam(cp);
815 for(i=0; i<2; i++) {
816 fD[i] =cp[i];
817 fT1[i] =cp[i+2];
818 fT2[i] =cp[i+4];
819 fTol[i]=cp[i+6];
820 /*
821 printf("\n i, fD, fT1, fT2, fTol %d %d %d %d %d\n",
822 i,fD[i],fT1[i],fT2[i],fTol[i]);
823 */
824 }
825}
826
827//____________________________________________
828void AliITSsimulationSDD::ReadBaseline(){
829 // read baseline and noise from file - either a .root file and in this
830 // case data should be organised in a tree with one entry for each
831 // module => reading should be done accordingly
832 // or a classic file and do smth. like this:
833 //
834 // Read baselines and noise for SDD
835 //
836
837
838 Int_t na,pos;
839 Float_t bl,n;
840 const char *kinput, *base,*kparam;
841 char *filtmp;
842
843 fResponse->Filenames(kinput,base,kparam);
844 fFileName=base;
845//
846 filtmp = gSystem->ExpandPathName(fFileName.Data());
847 FILE *bline = fopen(filtmp,"r");
848 printf("filtmp %s\n",filtmp);
849 na = 0;
850
851 if(bline) {
852 while(fscanf(bline,"%d %f %f",&pos, &bl, &n) != EOF) {
853 //printf("na, pos, bl, n %d %d %f %f\n",na, pos, bl, n);
854 if (pos != na+1) {
855 Error("ReadBaseline","Anode number not in increasing order!",
856 filtmp);
857 exit(1);
858 }
859 fBaseline[na]=bl;
860 fNoise[na]=n;
861 na++;
862 }
863 } else {
864 Error("ReadBaseline"," THE BASELINE FILE %s DOES NOT EXIST !",
865 filtmp);
866 exit(1);
867 } // end if(bline)
868
869 fclose(bline);
870 delete [] filtmp;
871
872
873}
874
875//____________________________________________
876Int_t AliITSsimulationSDD::Convert10to8(Int_t signal) {
877 // To the 10 to 8 bit lossive compression.
878 // code from Davide C. and Albert W.
879
880 if (signal < 128) return signal;
881 if (signal < 256) return (128+((signal-128)>>1));
882 if (signal < 512) return (192+((signal-256)>>3));
883 if (signal < 1024) return (224+((signal-512)>>4));
884 return 0;
885
886}
887
888//____________________________________________
889Int_t AliITSsimulationSDD::Convert8to10(Int_t signal) {
890 // Undo the lossive 10 to 8 bit compression.
891 // code from Davide C. and Albert W.
892 if (signal < 0 || signal > 255) {
893 printf("<Convert8to10> out of range %d \n",signal);
894 return 0;
895 }
896
897 if (signal < 128) return signal;
898 if (signal < 192) {
899 if (TMath::Odd(signal)) return (128+((signal-128)<<1));
900 else return (128+((signal-128)<<1)+1);
901 }
902 if (signal < 224) {
903 if (TMath::Odd(signal)) return (256+((signal-192)<<3)+3);
904 else return (256+((signal-192)<<3)+4);
905 }
906 if (TMath::Odd(signal)) return (512+((signal-224)<<4)+7);
907 else return (512+((signal-224)<<4)+7);
908 return 0;
909
910}
911
912//____________________________________________
913AliITSMap* AliITSsimulationSDD::HitMap(Int_t i){
914 //Return the correct map.
915 return ((i==0)? fHitMap1 : fHitMap2);
916}
917
918
919//____________________________________________
920void AliITSsimulationSDD::ZeroSuppression(Option_t *option) {
921 // perform the zero suppresion
922 if (strstr(option,"2D")) {
923 //Init2D(); // activate if param change module by module
924 Compress2D();
925 } else if (strstr(option,"1D")) {
926 //Init1D(); // activate if param change module by module
927 Compress1D();
928 } else StoreAllDigits();
929
930}
931
932//____________________________________________
933void AliITSsimulationSDD::Init2D(){
934 // read in and prepare arrays: fD, fT1, fT2
935 // savemu[nanodes], savesigma[nanodes]
936 // read baseline and noise from file - either a .root file and in this
937 // case data should be organised in a tree with one entry for each
938 // module => reading should be done accordingly
939 // or a classic file and do smth. like this ( code from Davide C. and
940 // Albert W.) :
941 //
942 // Read 2D zero-suppression parameters for SDD
943 //
944
945 if (!strstr(fParam,"file")) return;
946
947 Int_t na,pos,tempTh;
948 Float_t mu,sigma;
949 Float_t savemu[fNofMaps], savesigma[fNofMaps];
950 const char *kinput,*kbasel,*kpar;
951 char *filtmp;
952
953
954 Int_t minval = fResponse->MinVal();
955
956 fResponse->Filenames(kinput,kbasel,kpar);
957 fFileName=kpar;
958
959//
960 filtmp = gSystem->ExpandPathName(fFileName.Data());
961 FILE *param = fopen(filtmp,"r");
962 na = 0;
963
964 if(param) {
965 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
966 if (pos != na+1) {
967 Error("Init2D ","Anode number not in increasing order!",
968 filtmp);
969 exit(1);
970 }
971 savemu[na]=mu;
972 savesigma[na]=sigma;
973 if ((2.*sigma) < mu) {
974 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
975 mu = 2.0 * sigma;
976 } else fD[na] = 0;
977 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
978 if (tempTh < 0) tempTh=0;
979 fT1[na] = tempTh;
980 tempTh = (Int_t)floor(mu+3.0*sigma+0.5) - minval;
981 if (tempTh < 0) tempTh=0;
982 fT2[na] = tempTh;
983 na++;
984 } // end while
985
986 } else {
987 Error("Init2D "," THE FILE %s DOES NOT EXIST !",
988 filtmp);
989 exit(1);
990 } // end if(param)
991
992 fclose(param);
993 delete [] filtmp;
994
995
996
997}
998//____________________________________________
999void AliITSsimulationSDD::Compress2D(){
1000 //
1001 // simple ITS cluster finder -- online zero-suppression conditions
1002 //
1003 //
1004
1005 //printf("Compress2D!\n");
1006
1007 Int_t db,tl,th;
1008 Int_t minval = fResponse->MinVal();
1009 Bool_t write=fResponse->OutputOption();
1010
1011 Int_t nz, nl, nh, low, i, j;
1012
1013 for(i=0; i<fNofMaps; i++) {
1014 CompressionParam(i,db,tl,th);
1015 nz=0;
1016 nl=0;
1017 nh=0;
1018 low=0;
1019 for(j=0; j<fMaxNofSamples; j++) {
1020 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1021 signal -= db; // if baseline eq. is done here
1022 if (signal <= 0) {nz++; continue;}
1023 if ((signal - tl) < minval) low++;
1024 if ((signal - th) >= minval) {
1025 nh++;
1026 Bool_t cond=kTRUE;
1027 //printf("Compress2D : i j %d %d signal %d\n",i,j,signal);
1028 FindCluster(i,j,signal,minval,cond);
1029 } else if ((signal - tl) >= minval) nl++;
1030 } // loop time samples
1031 if (write) TreeB()->Fill(nz,nl,nh,low,i+1);
1032 //if (nz != 256 && low != 256) printf("i, nz, nl, nh low %d %d %d %d %d\n",i,nz,nl,nh,low);
1033 } // loop anodes
1034
1035 char hname[30];
1036 if (write) {
1037 sprintf(hname,"TNtuple%d_%d",fModule,fEvent);
1038 TreeB()->Write(hname);
1039 // reset tree
1040 TreeB()->Reset();
1041 }
1042
1043}
1044
1045//_____________________________________________________________________________
1046void AliITSsimulationSDD::FindCluster(Int_t i,Int_t j,Int_t signal,
1047 Int_t minval,Bool_t cond){
1048//
1049// Find clusters according to the online 2D zero-suppression algorithm
1050//
1051
1052 Bool_t high=kFALSE;
1053
1054 fHitMap2->FlagHit(i,j);
1055//
1056// check the online zero-suppression conditions
1057//
1058 Int_t nn;
1059 Int_t dbx,tlx,thx;
1060 Int_t Xlist[kMaxNeighbours], Ylist[kMaxNeighbours];
1061 fSegmentation->Neighbours(i,j,&nn,Xlist,Ylist);
1062 Int_t in,ix,iy;
1063 for(in=0; in<nn; in++) {
1064 ix=Xlist[in];
1065 iy=Ylist[in];
1066 if (fHitMap2->TestHit(ix,iy)==kUnused) {
1067 CompressionParam(ix,dbx,tlx,thx);
1068 Int_t qn = (Int_t)(fHitMap2->GetSignal(ix,iy));
1069 qn -= dbx; // if baseline eq. is done here
1070 if ((qn-tlx) < minval) {
1071 fHitMap2->FlagHit(ix,iy);
1072 continue;
1073 } else {
1074 if ((qn - thx) >= minval) high=kTRUE;
1075 if (cond) {
1076 signal = Convert10to8(signal);
1077 //printf("FindCl -cond : i j %d %d signal %d\n",i,j,signal);
1078 AddDigit(i,j,signal);
1079 }
1080 Int_t qns = Convert10to8(qn);
1081 //printf("FindCl : i j %d %d qns %d\n",i,j,qns);
1082 if (!high) AddDigit(ix,iy,qns);
1083 cond=kFALSE;
1084 if(!high) fHitMap2->FlagHit(ix,iy);
1085 }
1086 } // TestHit
1087 } // loop over neighbours
1088
1089}
1090
1091//____________________________________________
1092void AliITSsimulationSDD::Init1D(){
1093 // this is just a copy-paste of input taken from 2D algo
1094 // Torino people should give input
1095 //
1096 // Read 1D zero-suppression parameters for SDD
1097 //
1098
1099 if (!strstr(fParam,"file")) return;
1100
1101 Int_t na,pos,tempTh;
1102 Float_t mu,sigma;
1103 Float_t savemu[fNofMaps], savesigma[fNofMaps];
1104 const char *kinput,*kbasel,*kpar;
1105 char *filtmp;
1106
1107
1108 Int_t minval = fResponse->MinVal();
1109 fResponse->Filenames(kinput,kbasel,kpar);
1110 fFileName=kpar;
1111
1112// set first the disable and tol param
1113 SetCompressParam();
1114//
1115 filtmp = gSystem->ExpandPathName(fFileName.Data());
1116 FILE *param = fopen(filtmp,"r");
1117 na = 0;
1118
1119 if (param) {
1120 fscanf(param,"%d %d %d %d ", &fT2[0], &fT2[1], &fTol[0], &fTol[1]);
1121 while(fscanf(param,"%d %f %f",&pos, &mu, &sigma) != EOF) {
1122 if (pos != na+1) {
1123 Error("Init1D ","Anode number not in increasing order!",
1124 filtmp);
1125 exit(1);
1126 }
1127 savemu[na]=mu;
1128 savesigma[na]=sigma;
1129 if ((2.*sigma) < mu) {
1130 fD[na] = (Int_t)floor(mu - 2.0*sigma + 0.5);
1131 mu = 2.0 * sigma;
1132 } else fD[na] = 0;
1133 tempTh = (Int_t)floor(mu+2.25*sigma+0.5) - minval;
1134 if (tempTh < 0) tempTh=0;
1135 fT1[na] = tempTh;
1136 na++;
1137 } // end while
1138 } else {
1139 Error("Init1D "," THE FILE %s DOES NOT EXIST !",
1140 filtmp);
1141 exit(1);
1142 } // end if(param)
1143
1144 fclose(param);
1145 delete [] filtmp;
1146
1147
1148}
1149
1150//____________________________________________
1151void AliITSsimulationSDD::Compress1D(){
1152 // 1D zero-suppression algorithm (from Gianluca A.)
1153
1154 Int_t dis,tol,thres,decr,diff;
1155 //char *dfile=strstr(fParam,"file");
1156
1157 UChar_t *str=fStream->Stream();
1158 Int_t counter=0;
1159
1160 Int_t last=0,k,i,j;
1161 for(k=1; k<=2; k++) {
1162 tol = Tolerance(k-1);
1163 dis = Disable(k-1);
1164 for(i=0; i<fNofMaps/2; i++) {
1165 Bool_t firstSignal=kTRUE;
1166 CompressionParam(k*i,decr,thres);
1167 for(j=0; j<fMaxNofSamples; j++) {
1168 Int_t signal=(Int_t)(fHitMap2->GetSignal(k*i,j));
1169 signal -= decr; // if baseline eq.
1170 signal = Convert10to8(signal);
1171 if (signal < thres) {
1172 signal=0;
1173 diff=128;
1174 last=0;
1175 // write diff in the buffer for HuffT
1176 str[counter]=(UChar_t)diff;
1177 counter++;
1178 continue;
1179 }
1180 diff=signal-last;
1181 if (diff > 127) diff=127;
1182 if (diff < -128) diff=-128;
1183
1184 if (signal < dis) {
1185 if (tol==1 && (diff >= -2 && diff <= 1)) diff=0;
1186 if (tol==2 && (diff >= -4 && diff <= 3)) diff=0;
1187 if (tol==3 && (diff >= -16 && diff <= 15)) diff=0;
1188 AddDigit(k*i,j,last+diff);
1189 } else {
1190 AddDigit(k*i,j,signal);
1191 }
1192
1193 diff += 128;
1194 // write diff in the buffer used to compute Huffman tables
1195 if (firstSignal) str[counter]=(UChar_t)signal;
1196 else str[counter]=(UChar_t)diff;
1197 counter++;
1198
1199 last=signal;
1200 firstSignal=kFALSE;
1201 } // loop time samples
1202 } // loop anodes one half of detector
1203 }
1204
1205 // check
1206 fStream->CheckCount(counter);
1207
1208 // open file and write out the stream of diff's
1209
1210 static Bool_t open=kTRUE;
1211 static TFile *OutFile;
1212 Bool_t write = fResponse->OutputOption();
1213
1214 if (write ) {
1215 if(open) {
1216 SetFileName("stream.root");
1217 cout<<"filename "<<fFileName<<endl;
1218 OutFile=new TFile(fFileName,"recreate");
1219 cout<<"I have opened "<<fFileName<<" file "<<endl;
1220 }
1221 open=kFALSE;
1222 OutFile->cd();
1223 fStream->Write();
1224 } // endif write
1225
1226 fStream->ClearStream();
1227
1228 // back to galice.root file
1229
1230 TTree *fAli=gAlice->TreeK();
1231 TFile *file = 0;
1232
1233 if (fAli) file =fAli->GetCurrentFile();
1234 file->cd();
1235
1236
1237}
1238//____________________________________________
1239void AliITSsimulationSDD::StoreAllDigits(){
1240 // if non-zero-suppressed data
1241
1242 Int_t digits[3],i,j;
1243
1244 for(i=0; i<fNofMaps; i++) {
1245 for(j=0; j<fMaxNofSamples; j++) {
1246 Int_t signal=(Int_t)(fHitMap2->GetSignal(i,j));
1247 signal = Convert10to8(signal);
1248 signal = Convert8to10(signal); // ?
1249 digits[0]=i;
1250 digits[1]=j;
1251 digits[2]=signal;
1252 fITS->AddRealDigit(1,digits);
1253 }
1254 }
1255}
1256//____________________________________________
1257
1258void AliITSsimulationSDD::CreateHistograms(){
1259 // Creates histograms of maps for debugging
1260
1261 Int_t i;
1262 for(i=0;i<fNofMaps;i++) {
1263 TString *sddName = new TString("sdd_");
1264 Char_t candNum[4];
1265 sprintf(candNum,"%d",i+1);
1266 sddName->Append(candNum);
1267 (*fHis)[i] = new TH1F(sddName->Data(),"SDD maps",
1268 fMaxNofSamples,0.,(Float_t) fMaxNofSamples);
1269 delete sddName;
1270 }
1271
1272}
1273//____________________________________________
1274
1275void AliITSsimulationSDD::ResetHistograms(){
1276 //
1277 // Reset histograms for this detector
1278 //
1279 Int_t i;
1280 for(i=0;i<fNofMaps;i++ ) {
1281 if ((*fHis)[i]) ((TH1F*)(*fHis)[i])->Reset();
1282 }
1283
1284}
1285
1286
1287//____________________________________________
1288
1289TH1F *AliITSsimulationSDD::GetAnode(Int_t wing, Int_t anode) {
1290 // Fills a histogram from a give anode.
1291 if (!fHis) return 0;
1292
1293 if(wing <=0 || wing > 2) {
1294 cout << "Wrong wing number: " << wing << endl;
1295 return NULL;
1296 }
1297 if(anode <=0 || anode > fNofMaps/2) {
1298 cout << "Wrong anode number: " << anode << endl;
1299 return NULL;
1300 }
1301
1302 Int_t index = (wing-1)*fNofMaps/2 + anode-1;
1303 return (TH1F*)((*fHis)[index]);
1304}
1305
1306//____________________________________________
1307
1308void AliITSsimulationSDD::WriteToFile(TFile *hfile) {
1309 // Writes the histograms to a file
1310 if (!fHis) return;
1311
1312 hfile->cd();
1313 Int_t i;
1314 for(i=0; i<fNofMaps; i++) (*fHis)[i]->Write(); //fAdcs[i]->Write();
1315 return;
1316}
1317//____________________________________________
1318Float_t AliITSsimulationSDD::GetNoise(Float_t threshold) {
1319 // Returns the noise value
1320 if (!fHis) return 0.;
1321
1322 TH1F *noisehist = new TH1F("noisehist","noise",100,0.,threshold);
1323 Int_t i,k;
1324 for(i=0;i<fNofMaps;i++) {
1325 Int_t nOfBinsA = ((TH1F*)(*fHis)[i])->GetNbinsX();
1326 for(k=0;k<nOfBinsA;k++) {
1327 Float_t content = ((TH1F*)(*fHis)[i])->GetBinContent(k+1);
1328 if (content < threshold) noisehist->Fill(content);
1329 }
1330 }
1331 TF1 *gnoise = new TF1("gnoise","gaus",0.,threshold);
1332 noisehist->Fit("gnoise","RQ");
1333 noisehist->Draw();
1334 Float_t mnoise = gnoise->GetParameter(1);
1335 cout << "mnoise : " << mnoise << endl;
1336 Float_t rnoise = gnoise->GetParameter(2);
1337 cout << "rnoise : " << rnoise << endl;
1338 delete noisehist;
1339 return rnoise;
1340}
1341void AliITSsimulationSDD::Streamer(TBuffer &R__b)
1342{
1343 // Stream an object of class AliITSsimulationSDD.
1344
1345 if (R__b.IsReading()) {
1346 Version_t R__v = R__b.ReadVersion(); if (R__v) { }
1347 AliITSsimulation::Streamer(R__b);
1348 R__b >> fITS;
1349 R__b >> fHitMap1;
1350 R__b >> fHitMap2;
1351 R__b >> fStream;
1352 R__b >> fElectronics;
1353 fD.Streamer(R__b);
1354 fT1.Streamer(R__b);
1355 fT2.Streamer(R__b);
1356 fTol.Streamer(R__b);
1357 fBaseline.Streamer(R__b);
1358 fNoise.Streamer(R__b);
1359 R__b >> fTreeB;
1360 //R__b.ReadArray(fParam); // Not to be printed out?
1361 fFileName.Streamer(R__b);
1362 R__b >> fNofMaps;
1363 R__b >> fMaxNofSamples;
1364 R__b >> fModule;
1365 R__b >> fEvent;
1366 R__b >> fHis;
1367 } else {
1368 R__b.WriteVersion(AliITSsimulationSDD::IsA());
1369 AliITSsimulation::Streamer(R__b);
1370 R__b << fITS;
1371 R__b << fHitMap1;
1372 R__b << fHitMap2;
1373 R__b << fStream;
1374 R__b << fElectronics;
1375 fD.Streamer(R__b);
1376 fT1.Streamer(R__b);
1377 fT2.Streamer(R__b);
1378 fTol.Streamer(R__b);
1379 fBaseline.Streamer(R__b);
1380 fNoise.Streamer(R__b);
1381 R__b << fTreeB;
1382 //R__b.WriteArray(fParam, __COUNTER__); // Not to be printed out?
1383 fFileName.Streamer(R__b);
1384 R__b << fNofMaps;
1385 R__b << fMaxNofSamples;
1386 R__b << fModule;
1387 R__b << fEvent;
1388 R__b << fHis;
1389 }
1390}