Correct treatment of baseline and noise in SDD digits (P. Cerello)
[u/mrichter/AliRoot.git] / ITS / AliITSresponseSDD.cxx
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 /* $Id$ */
17
18 #include <Riostream.h>
19 #include <TRandom.h>
20
21 #include "AliITSresponseSDD.h"
22
23 //////////////////////////////////////////////////
24 //  Response class for set:ITS                      //
25 //  Specific subdetector implementation             //
26 //  for silicon drift detectors                     //
27 //                                                  //
28 //                                                  //
29 //////////////////////////////////////////////////////
30
31 const Int_t AliITSresponseSDD::fgkModules;   
32 const Int_t AliITSresponseSDD::fgkChips;  
33 const Int_t AliITSresponseSDD::fgkChannels; 
34 const Int_t AliITSresponseSDD::fgkMaxAdcDefault = 1024;
35 const Float_t AliITSresponseSDD::fgkDynamicRangeDefault = 132.;
36 const Float_t AliITSresponseSDD::fgkfChargeLossDefault = 0;
37 const Float_t AliITSresponseSDD::fgkDiffCoeffDefault = 3.23;
38 const Float_t AliITSresponseSDD::fgkDiffCoeff1Default = 30.;
39 const Float_t AliITSresponseSDD::fgkTemperatureDefault = 296.;
40 const TString AliITSresponseSDD::fgkParam1Default = "same";
41 const TString AliITSresponseSDD::fgkParam2Default = "same";
42 const Float_t AliITSresponseSDD::fgkNoiseDefault = 10.;
43 const Float_t AliITSresponseSDD::fgkBaselineDefault = 20.;
44 const TString AliITSresponseSDD::fgkOptionDefault = "1D";
45 const Float_t AliITSresponseSDD::fgkMinValDefault  = 4;
46 const Float_t AliITSresponseSDD::fgkDriftSpeedDefault = 7.3;
47 const Float_t AliITSresponseSDD::fgkNsigmasDefault = 3.;
48 const Int_t AliITSresponseSDD::fgkNcompsDefault = 121;
49 //______________________________________________________________________
50 ClassImp(AliITSresponseSDD)
51
52   AliITSresponseSDD::AliITSresponseSDD(){
53   // default constructor
54   fGaus = 0;
55   SetDeadChannels();
56   SetMaxAdc(fgkMaxAdcDefault);
57   SetDiffCoeff(fgkDiffCoeffDefault,fgkDiffCoeff1Default);
58   SetDriftSpeed(fgkDriftSpeedDefault);
59   SetNSigmaIntegration(fgkNsigmasDefault);
60   SetNLookUp(fgkNcompsDefault);
61   // SetClock();
62   SetNoiseParam(fgkNoiseDefault,fgkBaselineDefault);
63   SetNoiseAfterElectronics();
64   SetJitterError();
65   SetElectronics();
66   SetDynamicRange(fgkDynamicRangeDefault);
67   SetChargeLoss(fgkfChargeLossDefault);
68   SetThresholds(fgkMinValDefault,0.);
69   SetParamOptions(fgkParam1Default.Data(),fgkParam2Default.Data());
70   SetTemperature(fgkTemperatureDefault);
71   SetZeroSupp(fgkOptionDefault);
72   SetDataType();
73   SetFilenames();
74   SetOutputOption();
75   SetDo10to8();
76   // set the default zero suppression parameters
77   fCPar[0]=(Int_t) fBaseline;
78   fCPar[1]=(Int_t) fBaseline;
79   fCPar[2]=(Int_t)(2.*fNoiseAfterEl + 0.5);
80   fCPar[3]=(Int_t)(2.*fNoiseAfterEl + 0.5);
81   fCPar[4]=0;
82   fCPar[5]=0;
83   fCPar[6]=0;
84   fCPar[7]=0;
85 }
86 //______________________________________________________________________
87 AliITSresponseSDD::AliITSresponseSDD(const char *dataType){
88   // constructor
89   fGaus = 0;
90   SetDeadChannels();
91   SetMaxAdc(fgkMaxAdcDefault);
92   SetDiffCoeff(fgkDiffCoeffDefault,fgkDiffCoeff1Default);
93   SetDriftSpeed(fgkDriftSpeedDefault);
94   SetNSigmaIntegration(fgkNsigmasDefault);
95   SetNLookUp(fgkNcompsDefault);
96   // SetClock();
97   SetNoiseParam(fgkNoiseDefault,fgkBaselineDefault);
98   SetNoiseAfterElectronics();
99   SetJitterError();
100   SetElectronics();
101   SetDynamicRange(fgkDynamicRangeDefault);
102   SetChargeLoss(fgkfChargeLossDefault);
103   SetThresholds(fgkMinValDefault,0.);
104   SetParamOptions(fgkParam1Default.Data(),fgkParam2Default.Data());
105   SetTemperature(fgkTemperatureDefault);
106   SetZeroSupp(fgkOptionDefault);
107   SetDataType(dataType);
108   SetFilenames();
109   SetOutputOption();
110   SetDo10to8();
111   // set the default zero suppression parameters
112   fCPar[0]=(Int_t) fBaseline;
113   fCPar[1]=(Int_t) fBaseline;
114   fCPar[2]=(Int_t)(2.*fNoiseAfterEl + 0.5);
115   fCPar[3]=(Int_t)(2.*fNoiseAfterEl + 0.5);
116   fCPar[4]=0;
117   fCPar[5]=0;
118   fCPar[6]=0;
119   fCPar[7]=0;
120 }
121 //______________________________________________________________________
122 AliITSresponseSDD::AliITSresponseSDD(const AliITSresponseSDD &ob) : AliITSresponse(ob) {
123   // Copy constructor
124   // Copies are not allowed. The method is protected to avoid misuse.
125   Error("AliITSresponseSDD","Copy constructor not allowed\n");
126 }
127
128 //______________________________________________________________________
129 AliITSresponseSDD& AliITSresponseSDD::operator=(const AliITSresponseSDD& /* ob */){
130   // Assignment operator
131   // Assignment is not allowed. The method is protected to avoid misuse.
132   Error("= operator","Assignment operator not allowed\n");
133   return *this;
134 }
135
136 //______________________________________________________________________
137 AliITSresponseSDD::~AliITSresponseSDD() { 
138
139   if(fGaus) delete fGaus;
140 }
141
142 //______________________________________________________________________
143 Int_t AliITSresponseSDD::Convert8to10(Int_t signal) const {
144   // Undo the lossive 10 to 8 bit compression.
145   // code from Davide C. and Albert W.
146   if(Do10to8()){  // kTRUE if the compression is active
147     if (signal < 0 || signal > 255) {
148       Warning("Convert8to10","out of range signal=%d",signal);
149       return 0;
150     } // end if signal <0 || signal >255
151
152     if (signal < 128) return signal;
153     if (signal < 192) {
154       if (TMath::Odd(signal)) return (128+((signal-128)<<1));
155       else  return (128+((signal-128)<<1)+1);
156     } // end if signal < 192
157     if (signal < 224) {
158       if (TMath::Odd(signal)) return (256+((signal-192)<<3)+3);
159       else  return (256+((signal-192)<<3)+4);
160     } // end if signal < 224
161     if (TMath::Odd(signal)) return (512+((signal-224)<<4)+7);
162     return (512+((signal-224)<<4)+8);
163   }
164   else {  
165     return signal;
166   }
167 }
168
169 //______________________________________________________________________
170 void AliITSresponseSDD::SetCompressParam(Int_t  cp[8]){
171   // set compression param
172
173   Int_t i;
174   for (i=0; i<8; i++) {
175     fCPar[i]=cp[i];
176     //printf("\n CompressPar %d %d \n",i,fCPar[i]);    
177   } // end for i
178 }
179 //______________________________________________________________________
180 void AliITSresponseSDD::GiveCompressParam(Int_t  cp[8]) const {
181   // give compression param
182
183   Int_t i;
184   for (i=0; i<8; i++) {
185     cp[i]=fCPar[i];
186   } // end for i
187 }
188 //______________________________________________________________________
189 void AliITSresponseSDD::SetNLookUp(Int_t p1){
190   // Set number of sigmas over which cluster disintegration is performed
191   fNcomps=p1;
192   fGaus = new TArrayF(fNcomps+1);
193   for(Int_t i=0; i<=fNcomps; i++) {
194     Float_t x = -fNsigmas + (2.*i*fNsigmas)/(fNcomps-1);
195     (*fGaus)[i] = exp(-((x*x)/2));
196     //     cout << "fGaus[" << i << "]: " << fGaus->At(i) << endl;
197   }
198 }
199 //______________________________________________________________________
200 void AliITSresponseSDD::SetDeadChannels(Int_t nmod, Int_t nchip, Int_t nchan){
201   // Set fGain to zero to simulate a random distribution of 
202   // dead modules, dead chips and single dead channels
203
204   for( Int_t m=0; m<fgkModules; m++ ) 
205     for( Int_t n=0; n<fgkChips; n++ ) 
206       for( Int_t p=0; p<fgkChannels; p++ ) 
207         fGain[m][n][p] = 1.;
208                  
209   fDeadModules  = nmod;  
210   fDeadChips    = nchip;  
211   fDeadChannels = nchan; 
212     
213   // nothing to do
214   if( nmod == 0 && nchip == 0 && nchan == 0 ) return;
215
216   if( nmod < 0 || nmod > fgkModules ) 
217     { 
218       cout << "Wrong number of dead modules: " << nmod << endl; 
219       return; 
220     }
221   Int_t nmax = (fgkModules-nmod)*fgkChips; 
222   if( nchip < 0 || nchip > nmax ) 
223     { 
224       cout << "Wrong number of dead chips: " << nchip << endl; 
225       return; 
226     }
227   nmax = ((fgkModules - nmod)*fgkChips - nchip)*fgkChannels; 
228   if( nchan < 0 || nchan > nmax ) 
229     { 
230       cout << "Wrong number of dead channels: " << nchan << endl; 
231       return; 
232     }
233   
234   TRandom *gran = new TRandom();
235   
236   //  cout << "modules" << endl;
237   Int_t * mod = new Int_t [nmod];
238   Int_t i; //loop variable
239   for( i=0; i<nmod; i++ ) 
240     {
241       mod[i] = (Int_t) (1.+fgkModules*gran->Uniform());
242       cout << i+1 << ": Dead module nr: " << mod[i] << endl;
243       for(Int_t n=0; n<fgkChips; n++)
244         for(Int_t p=0; p<fgkChannels; p++)
245           fGain[mod[i]-1][n][p] = 0.;
246     }
247
248   //  cout << "chips" << endl;
249   Int_t * chip     = new Int_t[nchip];
250   Int_t * chipMod = new Int_t[nchip];
251   i = 0;
252   while( i < nchip ) 
253     {
254       Int_t module = (Int_t) (fgkModules*gran->Uniform() + 1.);
255       if( module <=0 || module > fgkModules ) 
256         cout << "Wrong module: " << module << endl;
257       Int_t flagMod = 0;
258       for( Int_t k=0; k<nmod; k++ ) 
259         if( module == mod[k] ) { flagMod = 1; break; }
260       if( flagMod == 1 ) continue;
261         
262       Int_t chi = (Int_t) (fgkChips*gran->Uniform() + 1.);
263       if( chi <=0 || chi > fgkChips ) cout << "Wrong chip: " << chi << endl;
264       i++;
265       chip[i-1] = chi; 
266       chipMod[i-1] = module;
267       for( Int_t m=0; m<fgkChannels; m++ ) 
268         fGain[module-1][chi-1][m] = 0.;
269       cout << i << ": Dead chip nr. " << chip[i-1] << " in module nr: " 
270            << chipMod[i-1] << endl;
271     }
272
273   //  cout << "channels" << endl;
274   Int_t * channel      = new Int_t[nchan];
275   Int_t * channelChip = new Int_t[nchan];
276   Int_t * channelMod  = new Int_t[nchan];
277   i = 0;
278   while( i < nchan ) 
279     {
280       Int_t k; //loop variable
281       Int_t module = (Int_t) (fgkModules*gran->Uniform() + 1.);
282       if( module <=0 || module > fgkModules ) 
283         cout << "Wrong module: " << module << endl;
284       Int_t flagMod = 0;
285       for( k=0; k<nmod; k++ ) 
286         if( module == mod[k] ) { flagMod = 1; break; }
287       if( flagMod == 1 ) continue;
288       Int_t chipp = (Int_t) (fgkChips*gran->Uniform() + 1.);
289       if( chipp <=0 || chipp > fgkChips ) cout << "Wrong chip: "<< chipp<<endl;
290       Int_t flagChip = 0;
291       for( k=0; k<nchip; k++) 
292         if( chipp == chip[k] && module == chipMod[k] ) { 
293           flagChip = 1; break; }
294       if( flagChip == 1 ) continue;
295       i++;
296       channel[i-1] = (Int_t) (fgkChannels*gran->Uniform() + 1.); 
297       if( channel[i-1] <=0 || channel[i-1] > fgkChannels ) 
298         cout << "Wrong channel: " << channel[i-1] << endl;
299       channelChip[i-1] = chipp;
300       channelMod[i-1] = module;
301       fGain[module-1][chipp-1][channel[i-1]-1] = 0.;
302       cout << i << ": Dead channel nr. " << channel[i-1] << " in chip nr. " 
303            << channelChip[i-1] << " in module nr: " << channelMod[i-1] 
304            << endl;
305     }
306     
307   delete [] mod;
308   delete [] chip;
309   delete [] chipMod;
310   delete [] channel;
311   delete [] channelMod;
312   delete [] channelChip;
313 }
314 //______________________________________________________________________
315 void AliITSresponseSDD::PrintGains(){
316   //
317
318   if( GetDeadModules() == 0 && 
319       GetDeadChips() == 0 && 
320       GetDeadChannels() == 0 )
321     return;  
322
323   // Print Electronics Gains
324   cout << "**************************************************" << endl; 
325   cout << "             Print Electronics Gains              " << endl;
326   cout << "**************************************************" << endl;
327
328   // Print SDD electronic gains
329   for(Int_t t=0; t<fgkModules;t++)
330     for(Int_t u=0; u<fgkChips;u++)
331       for(Int_t v=0; v<fgkChannels;v++)
332         {
333           if( fGain[t][u][v] != 1.0 )
334             cout << "Gain for Module: " << t+1 << ", Chip " << u+1 << 
335               ", Channel " << v+1 << " = " << fGain[t][u][v] << endl;
336         }
337 }
338 //______________________________________________________________________
339 void AliITSresponseSDD::Print(){
340   // Print SDD response Parameters
341
342   cout << "**************************************************" << endl;
343   cout << "   Silicon Drift Detector Response Parameters    " << endl;
344   cout << "**************************************************" << endl;
345   cout << "Diffusion Coefficients: "<< fDiffCoeff<< ", "<<fDiffCoeff1 << endl;
346
347   cout << "Hardware compression parameters: " << endl; 
348   for(Int_t i=0; i<8; i++) cout << "fCPar[" << i << "] = " << fCPar[i] <<endl;
349   cout << "Noise before electronics (arbitrary units): " << fNoise << endl;
350   cout << "Baseline (ADC units): " << fBaseline << endl;
351   cout << "Noise after electronics (ADC units): " << fNoiseAfterEl << endl;
352
353   cout << "Dynamic Range: " << fDynamicRange << endl;
354   cout << "Charge Loss: " << fChargeLoss << endl;
355   cout << "Temperature: " << Temperature() << " K " << endl;
356   cout << "Drift Speed: " << fDriftSpeed << endl;
357   cout << "Electronics (1=PASCAL, 2=OLA): " << fElectronics << endl;
358
359   cout << "N. of Sigma for signal integration: " << fNsigmas << endl;
360   cout << "N. of bins in lookup table: " << fNcomps << endl;
361
362   cout << "Max. ADC Value: " << fMaxAdc << endl;
363   cout << "Min. Value: " << fMinVal << endl;
364
365   PrintGains();
366
367 }
368
369
370