2 /**************************************************************************
3 * Copyright(c) 1998-2000, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
20 #include <TDirectory.h>
28 #include "AliSTARTDigitizer.h"
30 #include "AliSTARThit.h"
31 #include "AliSTARTdigit.h"
32 #include "AliRunDigitizer.h"
33 #include <AliDetector.h>
35 #include <AliLoader.h>
36 #include <AliRunLoader.h>
38 #include <Riostream.h>
39 #include <Riostream.h>
40 #include "AliSTARTParameters.h"
41 #include "AliCDBLocal.h"
42 #include "AliCDBStorage.h"
43 #include "AliCDBManager.h"
44 #include "AliCDBEntry.h"
46 ClassImp(AliSTARTDigitizer)
48 //___________________________________________
49 AliSTARTDigitizer::AliSTARTDigitizer() :AliDigitizer()
51 // Default ctor - don't use it
55 //___________________________________________
56 AliSTARTDigitizer::AliSTARTDigitizer(AliRunDigitizer* manager)
57 :AliDigitizer(manager),
66 // ctor which should be used
68 AliDebug(1,"processed");
74 ftimeCFD = new TArrayI(24);
75 fADC = new TArrayI(24);
76 ftimeLED = new TArrayI(24);
77 fADC0 = new TArrayI(24);
82 //------------------------------------------------------------------------
83 AliSTARTDigitizer::~AliSTARTDigitizer()
95 //------------------------------------------------------------------------
96 Bool_t AliSTARTDigitizer::Init()
104 //---------------------------------------------------------------------
106 void AliSTARTDigitizer::Exec(Option_t* /*option*/)
110 Produde digits from hits
111 digits is TObject and includes
112 We are writing array if left & right TDC
113 left & right ADC (will need for slow simulation)
114 TOF first particle left & right
115 mean time and time difference (vertex position)
120 AliRunLoader *outRL = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
121 AliLoader * pOutStartLoader = outRL->GetLoader("STARTLoader");
123 AliDebug(1,"start...");
126 // From hits to digits
130 Int_t volume, pmt, trCFD, trLED;
132 Int_t bestRightTDC, bestLeftTDC, qtCh;
133 Float_t time[24], besttime[24], timeGaus[24] ;
134 // Float_t channelWidth=25.; //ps
135 //Q->T-> coefficients !!!! should be asked!!!
136 Float_t gain[24],timeDelayCFD[24], timeDelayLED[24];
137 Int_t threshold =50; //photoelectrons
138 Float_t zdetA, zdetC;
139 TObjArray slewingLED;
140 TObjArray slewingRec;
141 AliSTARTParameters* param = AliSTARTParameters::Instance();
144 Int_t ph2Mip = param->GetPh2Mip();
145 Int_t channelWidth = param->GetChannelWidth() ;
146 Float_t delayVertex = param->GetTimeDelayTVD();
148 for (Int_t i=0; i<24; i++){
149 timeDelayCFD[i] = param->GetTimeDelayCFD(i);
150 timeDelayLED[i] = param->GetTimeDelayLED(i);
151 gain[i] = param->GetGain(i);
152 TGraph* gr = param ->GetSlew(i);
153 slewingLED.AddAtAndExpand(gr,i);
155 TGraph* gr1 = param ->GetSlewRec(i);
156 slewingRec.AddAtAndExpand(gr1,i);
158 TGraph* grEff = param ->GetPMTeff(i);
159 fEffPMT.AddAtAndExpand(grEff,i);
161 zdetC = param->GetZposition(0);
162 zdetA = param->GetZposition(1);
164 AliSTARThit *startHit;
167 Int_t nFiles=fManager->GetNinputs();
168 for (Int_t inputFile=0; inputFile<nFiles; inputFile++) {
169 if (inputFile < nFiles-1) {
170 AliWarning(Form("ignoring input stream %d", inputFile));
175 Float_t besttimeright=99999.;
176 Float_t besttimeleft=99999.;
177 Int_t pmtBestRight=9999;
178 Int_t pmtBestLeft=9999;
179 Int_t timeDiff=999, meanTime=0;
185 for (Int_t i0=0; i0<24; i0++)
187 time[i0]=besttime[i0]=timeGaus[i0]=99999; countE[i0]=0;
189 AliRunLoader * inRL = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
190 AliLoader * pInStartLoader = inRL->GetLoader("STARTLoader");
191 if (!inRL->GetAliRun()) inRL->LoadgAlice();
192 AliSTART *fSTART = (AliSTART*)inRL ->GetAliRun()->GetDetector("START");
195 pInStartLoader->LoadHits("READ");//probably it is necessary to load them before
196 TClonesArray *fHits = fSTART->Hits ();
197 TTree *th = pInStartLoader->TreeH();
198 brHits = th->GetBranch("START");
200 fSTART->SetHitsAddressBranch(brHits);
202 AliError("Branch START hit not found");
205 Int_t ntracks = (Int_t) th->GetEntries();
207 if (ntracks<=0) return;
208 // Start loop on tracks in the hits containers
209 for (Int_t track=0; track<ntracks;track++) {
210 brHits->GetEntry(track);
211 nhits = fHits->GetEntriesFast();
212 for (hit=0;hit<nhits;hit++)
214 startHit = (AliSTARThit*) fHits->UncheckedAt(hit);
216 AliError("The unchecked hit doesn't exist");
221 Double_t e=startHit->Etot();
222 volume = startHit->Volume();
224 if(e>0 && RegisterPhotoE(numpmt,e)) {
226 besttime[numpmt] = startHit->Time();
227 if(besttime[numpmt]<time[numpmt])
229 time[numpmt]=besttime[numpmt];
231 } //photoelectron accept
235 //spread time right&left by 25ps && besttime
236 Float_t c = 0.0299792; // cm/ps
238 Float_t koef=(zdetA-zdetC)/c; //correction position difference by cable
239 for (Int_t ipmt=0; ipmt<12; ipmt++){
240 if(countE[ipmt] > threshold) {
241 timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25)+koef;
242 if(timeGaus[ipmt]<besttimeleft){
243 besttimeleft=timeGaus[ipmt]; //timeleft
247 for ( Int_t ipmt=12; ipmt<24; ipmt++){
248 if(countE[ipmt] > threshold) {
249 timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25);
250 if(timeGaus[ipmt]<besttimeright) {
251 besttimeright=timeGaus[ipmt]; //timeright
256 //folding with alignmentz position distribution
258 bestLeftTDC=Int_t ((besttimeleft+1000*timeDelayCFD[pmtBestLeft])
261 bestRightTDC=Int_t ((besttimeright+1000*timeDelayCFD[pmtBestLeft])
263 timeDiff=Int_t (((besttimeleft-besttimeright)+1000*delayVertex)
265 meanTime=Int_t (((besttimeright+1000*timeDelayCFD[pmtBestLeft]+
266 besttimeleft+1000*timeDelayCFD[pmtBestLeft])/2.)
269 AliDebug(2,Form(" time in ns %f %f ",besttimeleft,besttimeright ));
270 for (Int_t i=0; i<24; i++)
272 Float_t al = countE[i];
273 if (al>threshold && timeGaus[i]<50000 ) {
276 // phe -> mV 0.3; 1MIP ->500phe -> ln (amp (mV)) = 5;
277 // max 200ns, HIJING mean 50000phe -> 15000mv -> ln = 15 (s zapasom)
279 qt= 50.*al*gain[i]/ph2Mip; // 50mv/Mip amp in mV
281 trCFD = Int_t (timeGaus[i] + 1000.*timeDelayCFD[i])/channelWidth;
282 trLED= Int_t (timeGaus[i] + 1000.*timeDelayLED[i]);
283 sl = ((TGraph*)slewingLED.At(i))->Eval(qt);
284 trLED = Int_t(( trLED + 1000*sl )/channelWidth);
285 qtCh=Int_t (1000.*TMath::Log(qt)) / channelWidth;
288 ftimeCFD->AddAt(Int_t (trCFD),i);
289 ftimeLED->AddAt(trLED,i);
290 sumMult += Int_t ((al/ph2Mip)*50) ;
294 if (sumMult > threshold){
295 fSumMult = Int_t (TMath::Log(sumMult));
297 fSTART->AddDigit(bestRightTDC,bestLeftTDC,meanTime,timeDiff,fSumMult,
298 ftimeCFD,fADC,ftimeLED,fADC0);
299 pOutStartLoader->UnloadHits();
300 } //input streams loop
303 pOutStartLoader->LoadDigits("UPDATE");
304 TTree *treeD = pOutStartLoader->TreeD();
306 pOutStartLoader->MakeTree("D");
307 treeD = pOutStartLoader->TreeD();
310 fSTART = (AliSTART*)outRL ->GetAliRun()->GetDetector("START");
311 // Make a branch in the tree
312 fSTART->MakeBranch("D");
315 pOutStartLoader->WriteDigits("OVERWRITE");
317 fSTART->ResetDigits();
318 pOutStartLoader->UnloadDigits();
323 //------------------------------------------------------------------------
324 Bool_t AliSTARTDigitizer::RegisterPhotoE(Int_t ipmt,Double_t energy)
328 // Float_t hc=197.326960*1.e6; //mev*nm
329 Double_t hc=1.973*1.e-6; //gev*nm
330 Float_t lambda=hc/energy;
331 Float_t eff = ((TGraph*) fEffPMT.At(ipmt))->Eval(lambda);
332 Double_t p = gRandom->Rndm();
340 //----------------------------------------------------------------------------