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