"Setting correct gains for Pb-p in simulation"
[u/mrichter/AliRoot.git] / T0 / AliT0Digitizer.cxx
CommitLineData
dc7ca31d 1
2/**************************************************************************
3 * Copyright(c) 1998-2000, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
9edefa04 17/* $Id$ */
dc7ca31d 18
94249139 19/******************************************************************
20 * Produde digits from hits
21 * digits is TObject and includes
22 * We are writing array if C & A TDC
23 * C & A ADC (will need for slow simulation)
24 * TOF first particle C & A
25 * mean time and time difference (vertex position)
26 *
27 * Alla.Maevskaya@cern.ch
28 ****************************************************************/
29
30
dc7ca31d 31#include <TArrayI.h>
9edefa04 32#include <TFile.h>
dc7ca31d 33#include <TGraph.h>
9edefa04 34#include <TH1F.h>
35#include <TMath.h>
36#include <TRandom.h>
37#include <TTree.h>
dc7ca31d 38
39#include "AliLog.h"
40#include "AliT0Digitizer.h"
41#include "AliT0.h"
42#include "AliT0hit.h"
43#include "AliT0digit.h"
f21fc003 44#include "AliDigitizationInput.h"
dc7ca31d 45#include "AliRun.h"
46#include <AliLoader.h>
47#include <AliRunLoader.h>
48#include <stdlib.h>
dc7ca31d 49#include "AliT0Parameters.h"
dc7ca31d 50
51ClassImp(AliT0Digitizer)
52
53//___________________________________________
c41ceaac 54 AliT0Digitizer::AliT0Digitizer() :AliDigitizer(),
55 fT0(0),
56 fHits(0),
57 fdigits(0),
58 ftimeCFD(new TArrayI(24)),
59 ftimeLED (new TArrayI(24)),
60 fADC(new TArrayI(24)),
61 fADC0 (new TArrayI(24)),
62 fSumMult(0),
29d3e0eb 63 fAmpLED(0),
d0bcd1fb 64 fAmpQTC(0),
29d3e0eb 65 fParam(0)
dc7ca31d 66{
67// Default ctor - don't use it
68 ;
69}
70
71//___________________________________________
f21fc003 72AliT0Digitizer::AliT0Digitizer(AliDigitizationInput* digInput)
73 :AliDigitizer(digInput),
dc7ca31d 74 fT0(0),
75 fHits(0),
76 fdigits(0),
c41ceaac 77 ftimeCFD(new TArrayI(24)),
78 ftimeLED (new TArrayI(24)),
79 fADC(new TArrayI(24)),
80 fADC0 (new TArrayI(24)),
81 fSumMult(0),
29d3e0eb 82 fAmpLED(0),
d0bcd1fb 83 fAmpQTC(0),
29d3e0eb 84 fParam(0)
dc7ca31d 85{
86// ctor which should be used
87
88 AliDebug(1,"processed");
29d3e0eb 89 fParam = AliT0Parameters::Instance();
90 fParam->Init();
dc7ca31d 91
29d3e0eb 92 for (Int_t i=0; i<24; i++){
d0bcd1fb 93 TGraph* gr = fParam ->GetAmpLED(i);
4cde7903 94 if(gr) {
95 Int_t np = gr->GetN();
96 Double_t *x = gr->GetX();
97 Double_t *y = gr->GetY();
98
99 Double_t *x1 = new Double_t[np];
100 Double_t *y1 = new Double_t[np];
101 for (Int_t ii=0; ii<np; ii++) {
102 y1[ii]=y[np-ii-1];
103 x1[ii]=x[np-ii-1];
104 }
5c25947b 105 TGraph *grInverse = new TGraph(np,y1,x1);
106 fAmpLED.AddAtAndExpand(grInverse,i);
d097d860 107 if (x1) delete [] x1;
108 if (y1) delete [] y1;
4cde7903 109 }
29d3e0eb 110 }
d0bcd1fb 111 for (Int_t i=0; i<24; i++){
112 TGraph* grq = fParam ->GetQTC(i);
4cde7903 113 if(grq){
114 Int_t npq = grq->GetN();
115 Double_t *xq = grq->GetX();
116 Double_t *yq = grq->GetY();
117 Double_t *x1q = new Double_t[npq];
118 Double_t *y1q = new Double_t[npq];
119 for (Int_t ii=1; ii<npq; ii++) {
120 y1q[ii]=yq[ii-1];
121 x1q[ii]=xq[ii-1];
122 }
123 TGraph *grInverseQTC = new TGraph(npq,y1q,x1q);
124 fAmpQTC.AddAtAndExpand(grInverseQTC,i);
d097d860 125 if (x1q) delete [] x1q;
126 if (y1q) delete [] y1q;
5c25947b 127
d0bcd1fb 128 }
d0bcd1fb 129 }
dc7ca31d 130}
12e9daf9 131
4cde7903 132
dc7ca31d 133//------------------------------------------------------------------------
134AliT0Digitizer::~AliT0Digitizer()
135{
136// Destructor
137
138 AliDebug(1,"T0");
139
140 delete ftimeCFD;
141 delete fADC;
142 delete ftimeLED;
143 delete fADC0;
29d3e0eb 144 // delete fAmpLED;
dc7ca31d 145}
146
29d3e0eb 147//------------------------------------------------------------------------
dc7ca31d 148Bool_t AliT0Digitizer::Init()
149{
150// Initialization
151 AliDebug(1," Init");
152 return kTRUE;
153}
154
dc7ca31d 155//---------------------------------------------------------------------
f21fc003 156void AliT0Digitizer::Digitize(Option_t* /*option*/)
dc7ca31d 157{
158
159 /*
160 Produde digits from hits
161 digits is TObject and includes
c41ceaac 162 We are writing array if C & A TDC
163 C & A ADC (will need for slow simulation)
164 TOF first particle C & A
dc7ca31d 165 mean time and time difference (vertex position)
166
167 */
168
169 //output loader
f21fc003 170 AliRunLoader *outRL = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
dc7ca31d 171 AliLoader * pOutStartLoader = outRL->GetLoader("T0Loader");
172
173 AliDebug(1,"start...");
174 //input loader
175 //
176 // From hits to digits
177 //
74adb36a 178 Int_t hit, nhits;
dc7ca31d 179 Int_t countE[24];
180 Int_t volume, pmt, trCFD, trLED;
4cde7903 181 Float_t sl=0, qt=0;
af3bb1cc 182 Int_t bestATDC=0;
183 Int_t bestCTDC=0;
d0bcd1fb 184 Float_t qtCh=0;
dc7ca31d 185 Float_t time[24], besttime[24], timeGaus[24] ;
74adb36a 186 //Q->T-> coefficients !!!! should be asked!!!
187 Float_t timeDelayCFD[24];
dc7ca31d 188 Int_t threshold =50; //photoelectrons
189 Float_t zdetA, zdetC;
74adb36a 190 Int_t sumMultCoeff = 100;
740f0839 191 Int_t refpoint=0;
446d6ec4 192
29d3e0eb 193 Int_t ph2Mip = fParam->GetPh2Mip();
8955c6b4 194 Float_t channelWidth = fParam->GetChannelWidth() ;
29d3e0eb 195 Float_t delayVertex = fParam->GetTimeDelayTVD();
5325480c 196
d0bcd1fb 197
29d3e0eb 198 zdetC = TMath::Abs(fParam->GetZPosition("T0/C/PMT1"));
199 zdetA = TMath::Abs(fParam->GetZPosition("T0/A/PMT15"));
74adb36a 200
d0bcd1fb 201 // printf(" !!!!!Z det A = %f C = % f",zdetA,zdetC);
74adb36a 202 AliT0hit *startHit;
203 TBranch *brHits=0;
204
f21fc003 205 Int_t nFiles=fDigInput->GetNinputs();
dc7ca31d 206 for (Int_t inputFile=0; inputFile<nFiles; inputFile++) {
207 if (inputFile < nFiles-1) {
74adb36a 208 AliWarning(Form("ignoring input stream %d", inputFile));
209 continue;
210
dc7ca31d 211 }
4cde7903 212 Float_t slew = 0;
c41ceaac 213 Float_t besttimeC=99999.;
214 Float_t besttimeA=99999.;
af3bb1cc 215 Int_t pmtBestC=99999;
216 Int_t pmtBestA=99999;
217 Int_t timeDiff=99999, meanTime=99999;
dc7ca31d 218 Int_t sumMult =0; fSumMult=0;
c41ceaac 219 bestATDC = 99999; bestCTDC = 99999;
dc7ca31d 220
c41ceaac 221
dc7ca31d 222 ftimeCFD -> Reset();
223 fADC -> Reset();
224 fADC0 -> Reset();
225 ftimeLED ->Reset();
226 for (Int_t i0=0; i0<24; i0++)
227 {
af3bb1cc 228 time[i0]=besttime[i0]=timeGaus[i0]=999999; countE[i0]=0;
dc7ca31d 229 }
f21fc003 230 AliRunLoader * inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(inputFile));
dc7ca31d 231 AliLoader * pInStartLoader = inRL->GetLoader("T0Loader");
232 if (!inRL->GetAliRun()) inRL->LoadgAlice();
c2337900 233 fT0 = (AliT0*)inRL ->GetAliRun()->GetDetector("T0");
dc7ca31d 234
235 //read Hits
236 pInStartLoader->LoadHits("READ");//probably it is necessary to load them before
955401cc 237 fHits = fT0->Hits ();
dc7ca31d 238 TTree *th = pInStartLoader->TreeH();
239 brHits = th->GetBranch("T0");
240 if (brHits) {
241 fT0->SetHitsAddressBranch(brHits);
242 }else{
af3bb1cc 243 AliWarning("Branch T0 hit not found for this event");
244 // fT0->AddDigit(bestATDC,bestCTDC,meanTime,timeDiff,fSumMult, refpoint,
245 // ftimeCFD,fADC0,ftimeLED,fADC);
246 continue;
dc7ca31d 247 }
248 Int_t ntracks = (Int_t) th->GetEntries();
249
250 if (ntracks<=0) return;
251 // Start loop on tracks in the hits containers
252 for (Int_t track=0; track<ntracks;track++) {
253 brHits->GetEntry(track);
254 nhits = fHits->GetEntriesFast();
255 for (hit=0;hit<nhits;hit++)
256 {
257 startHit = (AliT0hit*) fHits->UncheckedAt(hit);
258 if (!startHit) {
259 AliError("The unchecked hit doesn't exist");
260 break;
261 }
262 pmt=startHit->Pmt();
263 Int_t numpmt=pmt-1;
264 Double_t e=startHit->Etot();
265 volume = startHit->Volume();
266
dc7ca31d 267 if(e>0 ) {
268 countE[numpmt]++;
269 besttime[numpmt] = startHit->Time();
270 if(besttime[numpmt]<time[numpmt])
271 {
272 time[numpmt]=besttime[numpmt];
273 }
274 } //photoelectron accept
275 } //hits loop
276 } //track loop
277
c41ceaac 278 //spread time A&C by 25ps && besttime
dc7ca31d 279 Float_t c = 0.0299792; // cm/ps
280
281 Float_t koef=(zdetA-zdetC)/c; //correction position difference by cable
282 for (Int_t ipmt=0; ipmt<12; ipmt++){
283 if(countE[ipmt] > threshold) {
284 timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25)+koef;
c41ceaac 285 if(timeGaus[ipmt]<besttimeC){
286 besttimeC=timeGaus[ipmt]; //timeC
287 pmtBestC=ipmt;}
9a5cba6d 288 }
dc7ca31d 289 }
9a5cba6d 290 for ( Int_t ipmt=12; ipmt<24; ipmt++){
dc7ca31d 291 if(countE[ipmt] > threshold) {
292 timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25);
c41ceaac 293 if(timeGaus[ipmt]<besttimeA) {
294 besttimeA=timeGaus[ipmt]; //timeA
295 pmtBestA=ipmt;}
dc7ca31d 296 }
297 }
dc7ca31d 298
6bfb23c3 299 timeDelayCFD[0] = fParam->GetTimeDelayCFD(0);
446d6ec4 300
dc7ca31d 301 for (Int_t i=0; i<24; i++)
302 {
303 Float_t al = countE[i];
304 if (al>threshold && timeGaus[i]<50000 ) {
305 //fill ADC
306 // QTC procedure:
307 // phe -> mV 0.3; 1MIP ->500phe -> ln (amp (mV)) = 5;
308 // max 200ns, HIJING mean 50000phe -> 15000mv -> ln = 15 (s zapasom)
309 // channel 25ps
d0bcd1fb 310 qt= al/ph2Mip; // 50mv/Mip amp in mV
71cb50f8 311 // before will we have calibration for high multiplicity
205262fe 312 // if (qt > 115.) qt =115.; //must be fix!!!
dc7ca31d 313 // fill TDC
29d3e0eb 314 timeDelayCFD[i] = fParam->GetTimeDelayCFD(i);
12e9daf9 315 trCFD = Int_t (timeGaus[i]/channelWidth + timeDelayCFD[i]);
d0bcd1fb 316
29d3e0eb 317 TGraph* gr = ((TGraph*)fAmpLED.At(i));
4cde7903 318 if(gr) sl = gr->Eval(qt);
29d3e0eb 319
d0bcd1fb 320 TGraph* gr1 = ((TGraph*)fAmpQTC.At(i));
4cde7903 321 if(gr1) qtCh = gr1->Eval(qt);
dc7ca31d 322 fADC0->AddAt(0,i);
d52fab81 323 if(qtCh)
40c34546 324 fADC->AddAt(Int_t(qtCh),i);
dc7ca31d 325 // sumMult += Int_t ((al*gain[i]/ph2Mip)*50) ;
d0bcd1fb 326 sumMult += Int_t (qtCh/sumMultCoeff) ;
dc7ca31d 327
74adb36a 328 // put slewing
71cb50f8 329 TGraph *fu=(TGraph*) fParam ->GetWalk(i);
4cde7903 330 if(fu) slew=fu->Eval(Float_t(qtCh));
331
446d6ec4 332 // trCFD=trCFD-Int_t(fMaxValue[i]-slew);
40c34546 333 trCFD = trCFD + Int_t(slew); //for the same channel as cosmic
74adb36a 334 ftimeCFD->AddAt(Int_t (trCFD),i);
d0bcd1fb 335 trLED = Int_t(trCFD + sl );
336 ftimeLED->AddAt(trLED,i);
4cde7903 337 AliDebug(1,Form(" pmt %i : delay %f time in ns %f time in channels %i LEd %i ", i, timeDelayCFD[i], timeGaus[i],trCFD, trLED ));
338 AliDebug(1,Form(" qt in MIP %f led-cfd in %f qt in channels %f ",qt, sl, qtCh));
74adb36a 339
dc7ca31d 340 }
341 } //pmt loop
342
9a5cba6d 343 //folding with alignmentz position distribution
344 if( besttimeC > 10000. && besttimeC <15000)
345 bestCTDC=Int_t ((besttimeC+timeDelayCFD[pmtBestC])
346 /channelWidth);
347
348 if( besttimeA > 10000. && besttimeA <15000)
349 bestATDC=Int_t ((besttimeA+timeDelayCFD[pmtBestA])
350 /channelWidth);
351
352 if (bestATDC < 99999 && bestCTDC < 99999)
353 {
354 timeDiff=Int_t (((besttimeC-besttimeA)+1000*delayVertex)
355 /channelWidth);
446d6ec4 356 meanTime=Int_t (((besttimeC+besttimeA)/2. )/channelWidth);
9a5cba6d 357 }
9a5cba6d 358
dc7ca31d 359 if (sumMult > threshold){
360 fSumMult = Int_t (1000.* TMath::Log(Double_t(sumMult) / Double_t(sumMultCoeff))
361 /channelWidth);
b5a9f753 362 AliDebug(10,Form("summult mv %i mult in chammens %i in ps %f ",
dc7ca31d 363 sumMult, fSumMult, fSumMult*channelWidth));
364 }
dc7ca31d 365
740f0839 366 fT0->AddDigit(bestATDC,bestCTDC,meanTime,timeDiff,fSumMult, refpoint,
c41ceaac 367 ftimeCFD,fADC0,ftimeLED,fADC);
dc7ca31d 368
740f0839 369
370 AliDebug(10,Form(" Digits wrote refpoint %i bestATDC %i bestCTDC %i meanTime %i timeDiff %i fSumMult %i ",refpoint ,bestATDC,bestCTDC,meanTime,timeDiff,fSumMult ));
dc7ca31d 371 pOutStartLoader->UnloadHits();
372 } //input streams loop
373
374 //load digits
375 pOutStartLoader->LoadDigits("UPDATE");
376 TTree *treeD = pOutStartLoader->TreeD();
377 if (treeD == 0x0) {
378 pOutStartLoader->MakeTree("D");
379 treeD = pOutStartLoader->TreeD();
380 }
381 treeD->Reset();
382 fT0 = (AliT0*)outRL ->GetAliRun()->GetDetector("T0");
383 // Make a branch in the tree
384 fT0->MakeBranch("D");
385 treeD->Fill();
386
387 pOutStartLoader->WriteDigits("OVERWRITE");
388
389 fT0->ResetDigits();
390 pOutStartLoader->UnloadDigits();
391
392}