66327acfa228b0048918a701516b461e179421ef
[u/mrichter/AliRoot.git] / START / AliSTARTDigitizer.cxx
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
18 #include <TTree.h> 
19 #include <TFile.h>
20 #include <TDirectory.h>
21 #include <TRandom.h>
22 #include <TArrayI.h>
23 #include <TError.h>
24 #include <TH1F.h>
25
26 #include "AliLog.h"
27 #include "AliSTARTDigitizer.h"
28 #include "AliSTART.h"
29 #include "AliSTARThit.h"
30 #include "AliSTARTdigit.h"
31 #include "AliRunDigitizer.h"
32 #include <AliDetector.h>
33 #include "AliRun.h"
34 #include <AliLoader.h>
35 #include <AliRunLoader.h>
36 #include <stdlib.h>
37 #include <Riostream.h>
38 #include <Riostream.h>
39
40 ClassImp(AliSTARTDigitizer)
41
42 //___________________________________________
43   AliSTARTDigitizer::AliSTARTDigitizer()  :AliDigitizer()
44 {
45 // Default ctor - don't use it
46   ;
47 }
48
49 //___________________________________________
50 AliSTARTDigitizer::AliSTARTDigitizer(AliRunDigitizer* manager) 
51   :AliDigitizer(manager),
52    fSTART(0),
53    fHits(0),
54    fdigits(0),
55    ftimeTDC(0),
56    fADC(0),
57    fEff(0)
58 {
59   //    cout<<"AliSTARTDigitizer::AliSTARTDigitizer"<<endl;
60 // ctor which should be used
61
62   AliDebug(1,"processed");
63
64   fSTART = 0;
65   fPhotons = 0;
66   fHits = 0;
67   fdigits = 0;
68
69   ftimeTDC = new TArrayI(24); 
70   fADC = new TArrayI(24); 
71
72   TFile* file = TFile::Open("$ALICE_ROOT/START/PMTefficiency.root");
73   fEff = (TH1F*) file->Get("hEff")->Clone();
74   fEff->SetDirectory(NULL);
75   file->Close();
76   delete file;
77 }
78
79 //------------------------------------------------------------------------
80 AliSTARTDigitizer::~AliSTARTDigitizer()
81 {
82 // Destructor
83
84   AliDebug(1,"START");
85
86   delete ftimeTDC;
87   delete fADC;
88   delete fEff;
89 }
90
91  //------------------------------------------------------------------------
92 Bool_t AliSTARTDigitizer::Init()
93 {
94 // Initialization
95   AliDebug(1," Init");
96  return kTRUE;
97 }
98  
99
100 //---------------------------------------------------------------------
101
102 void AliSTARTDigitizer::Exec(Option_t* /*option*/)
103 {
104
105   /*
106     Produde digits from hits
107         digits is TObject and includes
108         We are writing array if left & right  TDC
109         left & right  ADC (will need for slow simulation)
110         TOF first particle left & right
111         mean time and time difference (vertex position)
112         
113   */
114
115   AliRunLoader *inRL, *outRL;//in and out Run Loaders
116   AliLoader *pInStartLoader, *pOutStartLoader;// in and out STARTLoaders
117
118   outRL = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
119   pOutStartLoader = outRL->GetLoader("STARTLoader");
120
121   AliDebug(1,"start...");
122
123   //
124   // From hits to digits
125   //
126   Int_t hit, nhits;
127   Float_t meanTime;
128   Int_t countE[24];
129   Int_t volume,pmt,tr,sumRight;
130   Int_t  bestRightTDC,bestLeftTDC;
131   Float_t time[24]={24*0};
132   Float_t besttime[24]={24*0};
133   Float_t timeGaus[37]={24*0};
134   Float_t channelWidth=25.; //ps
135   
136   AliSTARThit  *startHit;
137   TBranch *brHits=0;
138
139   pOutStartLoader->LoadDigits("UPDATE");//probably it is necessary to load them before
140   fdigits= new AliSTARTdigit();
141   pOutStartLoader->GetDigitsDataLoader()->GetBaseLoader(0)->Post(fdigits);
142
143   Int_t nFiles=fManager->GetNinputs();
144   for (Int_t inputFile=0; inputFile<nFiles;  inputFile++) {
145     if (inputFile < nFiles-1) {
146       AliWarning(Form("ignoring input stream %d", inputFile));
147       continue;
148
149     }
150
151     Float_t besttimeright=9999.;
152     Float_t besttimeleft=9999.;
153     Float_t timeDiff;
154     sumRight=0;
155     for (Int_t i0=0; i0<24; i0++)
156       {
157         time[i0]=9999;  besttime[i0]=9999;      countE[i0]=0;
158       }
159
160     inRL = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
161     if (!inRL->GetAliRun()) inRL->LoadgAlice();
162     AliSTART *fSTART  = (AliSTART*) inRL->GetAliRun()->GetDetector("START");
163     pInStartLoader = inRL->GetLoader("STARTLoader");
164     pInStartLoader->LoadHits("READ");//probably it is necessary to load them before
165     TClonesArray *fHits = fSTART->Hits ();
166
167     TTree *th = pInStartLoader->TreeH();
168     brHits = th->GetBranch("START");
169     if (brHits) {
170       fSTART->SetHitsAddressBranch(brHits);
171     }else{
172        AliError("Branch START hit not found");
173       exit(111);
174     } 
175     Int_t ntracks    = (Int_t) th->GetEntries();
176
177     if (ntracks<=0) return;
178     // Start loop on tracks in the hits containers
179     for (Int_t track=0; track<ntracks;track++) {
180       brHits->GetEntry(track);
181       nhits = fHits->GetEntriesFast();
182        for (hit=0;hit<nhits;hit++) 
183         {
184           startHit   = (AliSTARThit*) fHits->UncheckedAt(hit);
185           if (!startHit) {
186             AliError("The unchecked hit doesn't exist");
187             break;
188           }
189           pmt=startHit->Pmt();
190           Int_t numpmt=pmt-1;
191           Float_t e=startHit->Etot();
192           volume = startHit->Volume();
193           if(RegisterPhotoE(e)) countE[numpmt]++;
194           besttime[numpmt] = startHit->Time();
195           if(besttime[numpmt]<time[numpmt])
196             {
197               time[numpmt]=besttime[numpmt];
198             }
199           
200         } //hits loop
201     } //track loop
202     
203     //best time right&left   
204     for (Int_t ipmt=0; ipmt<12; ipmt++)
205       {
206         timeGaus[ipmt]=gRandom->Gaus(time[ipmt],0.025);
207         if(timeGaus[ipmt]<besttimeleft) besttimeleft=timeGaus[ipmt]; //timeleft
208       }
209     for ( Int_t ipmt=12; ipmt<24; ipmt++)
210       {
211         timeGaus[ipmt]=gRandom->Gaus(time[ipmt],0.025);
212         if(timeGaus[ipmt]<besttimeright)  besttimeright=timeGaus[ipmt]; //timeright
213       }// besttime
214         
215     //folding with experimental time distribution
216       Float_t c = 29.9792; // mm/ns
217     Float_t koef=(350.-69.7)/c;
218     Float_t  besttimeleftR= besttimeleft;
219     besttimeleft=koef+besttimeleft;
220     bestLeftTDC=Int_t (besttimeleftR*1000/channelWidth);
221     bestRightTDC=Int_t (besttimeright*1000/channelWidth);
222     timeDiff=(c*(besttimeright-besttimeleftR)-(350.-69.7))/(2*c);
223     meanTime=(besttimeright+besttimeleftR)/2.;
224     Float_t ds=(c*(besttimeright-besttimeleftR)-(350.-69.7))/2;
225     AliDebug(2,Form(" timediff in ns %f  z= %f real point%f",timeDiff,timeDiff*c,ds));
226
227  
228     // Time to TDC signal
229     Int_t iTimeAv=Int_t (meanTime*1000/channelWidth); 
230     // time  channel numbres 
231     //       fill digits
232     fdigits->SetTimeBestLeft(bestLeftTDC);
233     fdigits->SetTimeBestRight(bestRightTDC);
234     fdigits->SetMeanTime(iTimeAv);
235     
236     //ADC features
237   
238     for (Int_t i=0; i<24; i++)
239       {
240
241         //  fill TDC
242         tr= Int_t (timeGaus[i]*1000/channelWidth); 
243         if(timeGaus[i]>100) tr=0;
244         ftimeTDC->AddAt(tr,i);
245         
246         //fill ADC
247         Int_t al= countE[i]; //1024 - channel width shoul be change
248         fADC->AddAt(al,i);
249       } //pmt loop
250
251     fdigits->SetTime(*ftimeTDC);
252     fdigits->SetADC(*fADC);
253      
254     pInStartLoader->UnloadHits();
255     
256   } //input streams loop
257      
258   pOutStartLoader->WriteDigits("OVERWRITE");
259   pOutStartLoader->UnloadDigits();
260 }
261
262
263 //------------------------------------------------------------------------
264 Bool_t AliSTARTDigitizer::RegisterPhotoE(Float_t e)
265 {
266
267   
268   //  Float_t hc=197.326960*1.e6; //mev*nm
269   Float_t hc=1.973*1.e-6; //gev*nm
270   Float_t lambda=hc/e;
271   Int_t bin=  fEff->GetXaxis()->FindBin(lambda);
272   Float_t eff=fEff->GetBinContent(bin);
273   Double_t  p = gRandom->Rndm();
274   if (p > eff)
275     return kFALSE;
276   
277   return kTRUE;
278 }
279 //----------------------------------------------------------------------------