]>
Commit | Line | Data |
---|---|---|
19f796ed | 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 | ///_________________________________________________________________________ | |
19 | /// | |
20 | /// This class constructs Digits out of Hits | |
21 | /// | |
22 | /// | |
23 | ||
24 | // --- Standard library --- | |
25 | ||
26 | // --- ROOT system --- | |
27 | #include <TMath.h> | |
28 | #include <TTree.h> | |
29 | #include <TRandom.h> | |
30 | ||
31 | // --- AliRoot header files --- | |
32 | #include "AliRun.h" | |
33 | #include "AliACORDE.h" | |
34 | #include "AliACORDEhit.h" | |
35 | #include "AliRunLoader.h" | |
36 | #include "AliLoader.h" | |
37 | #include "AliRunDigitizer.h" | |
38 | #include "AliCDBManager.h" | |
39 | #include "AliCDBStorage.h" | |
40 | #include "AliCDBEntry.h" | |
41 | #include "AliACORDECalibData.h" | |
42 | #include "AliACORDEConstants.h" | |
43 | ||
44 | #include "AliACORDEdigit.h" | |
45 | #include "AliACORDEDigitizer.h" | |
46 | ||
47 | ClassImp(AliACORDEDigitizer) | |
48 | ||
49 | AliACORDEDigitizer::AliACORDEDigitizer() | |
50 | :AliDigitizer(), | |
51 | fCalibData(GetCalibData()), | |
52 | fNdigits(0), | |
53 | fDigits(0) | |
54 | ||
55 | { | |
56 | // default constructor | |
57 | } | |
58 | ||
59 | AliACORDEDigitizer::AliACORDEDigitizer(AliRunDigitizer* manager) | |
60 | :AliDigitizer(manager), | |
61 | fCalibData(GetCalibData()), | |
62 | fNdigits(0), | |
63 | fDigits(0) | |
64 | ||
65 | { | |
66 | // constructor | |
67 | ||
68 | } | |
69 | ||
70 | AliACORDEDigitizer::~AliACORDEDigitizer() | |
71 | { | |
72 | // destructor | |
73 | ||
74 | if (fDigits) { | |
75 | fDigits->Delete(); | |
76 | delete fDigits; | |
77 | fDigits=0; | |
78 | } | |
79 | } | |
80 | ||
81 | ||
82 | Bool_t AliACORDEDigitizer::Init() | |
83 | { | |
84 | // Initialises the digitizer | |
85 | ||
86 | // Initialises the Digit array | |
87 | fDigits = new TClonesArray ("AliACORDEdigit", 1000); | |
88 | ||
89 | return kTRUE; | |
90 | } | |
91 | ||
92 | void AliACORDEDigitizer::Exec(Option_t* /*option*/) | |
93 | { | |
94 | ||
95 | // Creates digits from hits | |
96 | ||
97 | // 1.- create and initialize temporal variables | |
98 | // 2.- get loaders to access hots and digits | |
99 | // 3.- loop over all input. | |
100 | // 3.1 for each input loop over all track | |
101 | // 3.2 for each track loop over all hits | |
102 | // 3.3 for each hit, check | |
103 | // if energy above threshold | |
104 | // if survives efficiency | |
105 | // then store the plastic, the time and track in temp arrays | |
106 | // if a hit already survived for this plastic take the hit | |
107 | // with the lowest time | |
108 | // 4.- loop over temporal array | |
109 | // if both plastic have a surviving hit and the time | |
110 | // difference is below the limit, add a new digit | |
111 | // | |
112 | ||
113 | ||
114 | // 1.- temporal variables | |
115 | Float_t emin = AliACORDEConstants::Instance()->HitEnergyThreshold(); | |
116 | Float_t td = AliACORDEConstants::Instance()->MaxHitTimeDifference(); | |
27172603 | 117 | Int_t modules[60]={0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; |
118 | Int_t moduls[60]={0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; | |
119 | Int_t mods; | |
19f796ed | 120 | Float_t PlasticTimes[2][60]; |
121 | Int_t PlasticTracks[2][60]; | |
122 | for (Int_t i=0;i<60;i++) { | |
123 | PlasticTimes[0][i]=-1.0; | |
124 | PlasticTimes[1][i]=-1.0; | |
125 | PlasticTracks[0][i]=-1; | |
126 | PlasticTracks[1][i]=-1; | |
127 | } | |
27172603 | 128 | |
19f796ed | 129 | // 2.- get loaders |
130 | AliRunLoader* outRunLoader = | |
131 | AliRunLoader::GetRunLoader(fManager->GetOutputFolderName()); | |
132 | if (!outRunLoader) { | |
133 | Error("Exec", "Can not get output Run Loader"); | |
134 | return;} | |
135 | ||
136 | AliLoader* outLoader = outRunLoader->GetLoader("ACORDELoader"); | |
137 | if (!outLoader) { | |
138 | Error("Exec", "Can not get output ACORDE Loader"); | |
139 | return;} | |
140 | ||
141 | outLoader->LoadDigits("update"); | |
142 | if (!outLoader->TreeD()) outLoader->MakeTree("D"); | |
143 | outLoader->MakeDigitsContainer(); | |
144 | TTree* treeD = outLoader->TreeD(); | |
145 | Int_t bufsize = 16000; | |
146 | treeD->Branch("ACORDEdigit", &fDigits, bufsize); | |
147 | ||
148 | // 3. loop over inputs | |
149 | for (Int_t iInput = 0; iInput < fManager->GetNinputs(); iInput++) { | |
150 | AliRunLoader* runLoader = | |
151 | AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput)); | |
152 | AliLoader* loader = runLoader->GetLoader("ACORDELoader"); | |
153 | if (!loader) { | |
154 | Error("Exec", "Can not get ACORDE Loader for input %d", iInput); | |
155 | continue;} | |
156 | ||
157 | if (!runLoader->GetAliRun()) runLoader->LoadgAlice(); | |
158 | ||
159 | AliACORDE* acorde = (AliACORDE*) runLoader->GetAliRun()->GetDetector("ACORDE"); | |
160 | if (!acorde) { | |
161 | Error("Exec", "No ACORDE detector for input %d", iInput); | |
162 | continue;} | |
163 | ||
164 | loader->LoadHits(); | |
165 | TTree* treeH = loader->TreeH(); | |
166 | if (!treeH) { | |
167 | Error("Exec", "Cannot get TreeH for input %d", iInput); | |
168 | continue; } | |
169 | ||
170 | TClonesArray* hits = acorde->Hits(); | |
171 | ||
172 | // 3.1 loop over all tracks | |
173 | Int_t nTracks = (Int_t) treeH->GetEntries(); | |
174 | for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) { | |
175 | acorde->ResetHits(); | |
176 | treeH->GetEvent(iTrack); | |
177 | Int_t nHits = hits->GetEntriesFast(); | |
178 | // 3.2 loop over hits | |
179 | for (Int_t iHit = 0; iHit < nHits; iHit++) { | |
180 | AliACORDEhit* hit = (AliACORDEhit *)hits->UncheckedAt(iHit); | |
181 | // 3.3 select hit | |
182 | // get hit info | |
183 | Float_t eloss_mev = hit->Eloss()*1000.0; | |
184 | Int_t module = hit->GetModule(); | |
27172603 | 185 | modules[module]=1; |
19f796ed | 186 | Int_t plastic = hit->GetPlastic(); |
187 | Float_t time_ns = hit->GetTime()*1e9; | |
188 | Float_t eff = TMath::Sqrt(fCalibData->GetEfficiency(module)); | |
19f796ed | 189 | // if enough energy and efficiency |
190 | if( eloss_mev > emin && gRandom->Uniform() < eff ) { | |
191 | // if first hit or earlier track | |
192 | if ((PlasticTimes[plastic-1][module-1] == -1.0) || | |
193 | (PlasticTimes[plastic-1][module-1] > time_ns) ) { | |
194 | PlasticTimes[plastic-1][module-1]= time_ns; | |
195 | PlasticTracks[plastic-1][module-1]= hit->GetTrack(); | |
196 | } | |
197 | } | |
198 | } // end of hit loop | |
199 | } // end of track loop | |
27172603 | 200 | for(Int_t i=0;i<60;i++){moduls[i]=modules[i];} |
19f796ed | 201 | |
202 | loader->UnloadHits(); | |
203 | ||
204 | } // end of input loop | |
205 | ||
206 | // 4.- loop over temporal arrays to add hits | |
19f796ed | 207 | Int_t tracks[3]={-1,-1,-1}; |
208 | for (Int_t i=0; i<60; i++) { | |
209 | // if both modules have a hit | |
19f796ed | 210 | // if time diff small enough |
211 | Float_t diff = TMath::Abs(PlasticTimes[0][i]-PlasticTimes[1][i]); | |
212 | if (diff < td) { | |
213 | tracks[0] = PlasticTracks[0][i]; | |
214 | if (PlasticTracks[0][i] != PlasticTracks[1][i]) | |
215 | tracks[1] = PlasticTracks[1][i]; | |
27172603 | 216 | if(moduls[i]==1) { |
217 | mods = i; | |
5b18f27b | 218 | // Float_t module_time = (PlasticTimes[0][i] > PlasticTimes[1][i] ? |
219 | // PlasticTimes[0][i] : PlasticTimes[1][i]); | |
220 | // AddDigit(tracks, mods, module_time); | |
221 | AddDigit(tracks, mods, 0); | |
27172603 | 222 | } |
19f796ed | 223 | } |
224 | } | |
19f796ed | 225 | treeD->Fill(); |
226 | outLoader->WriteDigits("OVERWRITE"); | |
227 | outLoader->UnloadDigits(); | |
228 | ResetDigit(); | |
229 | } | |
230 | ||
231 | //____________________________________________________________________________ | |
232 | ||
233 | void AliACORDEDigitizer::AddDigit(Int_t* track, Int_t module, Float_t time) | |
234 | { | |
235 | ||
236 | // Adds Digit | |
237 | ||
238 | TClonesArray &ldigits = *fDigits; | |
239 | new(ldigits[fNdigits++]) AliACORDEdigit(track,module,time); | |
240 | } | |
241 | ||
27172603 | 242 | void AliACORDEDigitizer::AddDigit(Int_t* modul,Float_t time) |
243 | { | |
244 | // MRC Adds Digit | |
245 | TClonesArray &ldigits = *fDigits; | |
246 | new(ldigits[fNdigits++]) AliACORDEdigit(modul,time); | |
247 | ||
19f796ed | 248 | |
27172603 | 249 | } |
19f796ed | 250 | void AliACORDEDigitizer::ResetDigit() |
251 | { | |
252 | // | |
253 | // Clears Digits | |
254 | // | |
255 | fNdigits = 0; | |
256 | if (fDigits) fDigits->Delete(); | |
257 | } | |
258 | ||
259 | ||
260 | AliACORDECalibData* AliACORDEDigitizer::GetCalibData() const | |
261 | ||
262 | { | |
263 | AliCDBManager *man = AliCDBManager::Instance(); | |
264 | ||
265 | AliCDBEntry *entry=0; | |
266 | ||
267 | entry = man->Get("ACORDE/Calib/Data"); | |
268 | ||
269 | if(!entry){ | |
270 | AliWarning("Load of calibration data from default storage failed!"); | |
271 | AliWarning("Calibration data will be loaded from local storage ($ALICE_ROOT)"); | |
272 | Int_t runNumber = man->GetRun(); | |
162637e4 | 273 | entry = man->GetStorage("local://$ALICE_ROOT/OCDB") |
19f796ed | 274 | ->Get("ACORDE/Calib/Data",runNumber); |
275 | ||
276 | } | |
277 | ||
278 | // Retrieval of data in directory ACORDE/Calib/Data: | |
279 | ||
280 | ||
281 | AliACORDECalibData *calibdata = 0; | |
282 | ||
283 | if (entry) calibdata = (AliACORDECalibData*) entry->GetObject(); | |
284 | if (!calibdata) AliError("No calibration data from calibration database !"); | |
285 | ||
286 | ||
287 | return calibdata; | |
288 | ||
289 | } | |
290 |