Line data Source code
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 : // ZDC digitizer class //
21 : // //
22 : ///////////////////////////////////////////////////////////////////////////////
23 :
24 : // --- Standard libraries
25 : #include <stdlib.h>
26 : #include <stdio.h>
27 :
28 : // --- ROOT system
29 : #include <TTree.h>
30 : #include <TFile.h>
31 : #include <TNtuple.h>
32 : #include <TRandom.h>
33 :
34 : // --- AliRoot header files
35 : #include "AliRun.h"
36 : #include "AliHeader.h"
37 : #include "AliGenHijingEventHeader.h"
38 : #include "AliGenCocktailEventHeader.h"
39 : #include "AliDigitizationInput.h"
40 : #include "AliRunLoader.h"
41 : #include "AliLoader.h"
42 : #include "AliGRPObject.h"
43 : #include "AliCDBManager.h"
44 : #include "AliCDBEntry.h"
45 : #include "AliZDCSDigit.h"
46 : #include "AliZDCDigit.h"
47 : #include "AliZDCFragment.h"
48 : #include "AliZDCv3.h"
49 : #include "AliZDCDigitizer.h"
50 : #include "AliLog.h"
51 :
52 : class AliCDBStorage;
53 : class AliZDCPedestals;
54 :
55 12 : ClassImp(AliZDCDigitizer)
56 :
57 :
58 : //____________________________________________________________________________
59 0 : AliZDCDigitizer::AliZDCDigitizer() :
60 0 : fIsCalibration(0),
61 0 : fIsSignalInADCGate(kFALSE),
62 0 : fFracLostSignal(0.),
63 0 : fPedData(0),
64 0 : fSpectators2Track(kFALSE),
65 0 : fBeamEnergy(0.),
66 0 : fBeamType(""),
67 0 : fIspASystem(kFALSE),
68 0 : fIsRELDISgen(kFALSE),
69 0 : fSpectatorData(0x0),
70 0 : fSpectatorParam(1)
71 0 : {
72 : // Default constructor
73 0 : for(Int_t i=0; i<2; i++) fADCRes[i]=0.;
74 0 : }
75 :
76 : //____________________________________________________________________________
77 : AliZDCDigitizer::AliZDCDigitizer(AliDigitizationInput* digInput):
78 1 : AliDigitizer(digInput),
79 1 : fIsCalibration(0), //By default the simulation doesn't create calib. data
80 1 : fIsSignalInADCGate(kFALSE),
81 1 : fFracLostSignal(0.),
82 2 : fPedData(GetPedData()),
83 1 : fSpectators2Track(kFALSE),
84 1 : fBeamEnergy(0.),
85 1 : fBeamType(""),
86 1 : fIspASystem(kFALSE),
87 1 : fIsRELDISgen(kFALSE),
88 1 : fSpectatorData(NULL),
89 1 : fSpectatorParam(1)
90 5 : {
91 : // Get calibration data
92 1 : if(fIsCalibration!=0) printf("\n\t AliZDCDigitizer -> Creating calibration data (pedestals)\n");
93 12 : for(Int_t i=0; i<5; i++){
94 60 : for(Int_t j=0; j<5; j++)
95 25 : fPMGain[i][j] = 0.;
96 : }
97 2 : }
98 :
99 : //____________________________________________________________________________
100 : AliZDCDigitizer::~AliZDCDigitizer()
101 6 : {
102 : // Destructor
103 1 : if(fSpectatorData!=NULL){
104 0 : fSpectatorData->Close();
105 0 : delete fSpectatorData;
106 : }
107 3 : }
108 :
109 :
110 : //____________________________________________________________________________
111 : AliZDCDigitizer::AliZDCDigitizer(const AliZDCDigitizer &digitizer):
112 0 : AliDigitizer(),
113 0 : fIsCalibration(digitizer.fIsCalibration),
114 0 : fIsSignalInADCGate(digitizer.fIsSignalInADCGate),
115 0 : fFracLostSignal(digitizer.fFracLostSignal),
116 0 : fPedData(digitizer.fPedData),
117 0 : fSpectators2Track(digitizer.fSpectators2Track),
118 0 : fBeamEnergy(digitizer.fBeamEnergy),
119 0 : fBeamType(digitizer.fBeamType),
120 0 : fIspASystem(digitizer.fIspASystem),
121 0 : fIsRELDISgen(digitizer.fIsRELDISgen),
122 0 : fSpectatorData(digitizer.fSpectatorData),
123 0 : fSpectatorParam(digitizer.fSpectatorParam)
124 0 : {
125 : // Copy constructor
126 :
127 0 : for(Int_t i=0; i<5; i++){
128 0 : for(Int_t j=0; j<5; j++){
129 0 : fPMGain[i][j] = digitizer.fPMGain[i][j];
130 : }
131 : }
132 0 : for(Int_t i=0; i<2; i++) fADCRes[i] = digitizer.fADCRes[i];
133 :
134 :
135 0 : }
136 :
137 : //____________________________________________________________________________
138 : Bool_t AliZDCDigitizer::Init()
139 : {
140 : // Initialize the digitizer
141 :
142 :
143 : //printf(" **** Initializing AliZDCDigitizer fBeamEnergy = %1.0f GeV\n\n", fBeamEnergy);
144 3 : AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
145 1 : if(!entry) AliFatal("No calibration data loaded!");
146 : AliGRPObject* grpData = 0x0;
147 1 : if(entry){
148 3 : TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
149 1 : if(m){
150 : //m->Print();
151 0 : grpData = new AliGRPObject();
152 0 : grpData->ReadValuesFromMap(m);
153 0 : }
154 : else{
155 3 : grpData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
156 : }
157 1 : entry->SetOwner(0);
158 1 : AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
159 1 : }
160 1 : if(!grpData){
161 0 : AliError("No GRP entry found in OCDB! \n ");
162 0 : return kFALSE;
163 : }
164 :
165 2 : fBeamType = grpData->GetBeamType();
166 2 : if(fBeamType==AliGRPObject::GetInvalidString()){
167 0 : AliError("\t UNKNOWN beam type from GRP obj -> PMT gains not set in ZDC digitizer!!!\n");
168 0 : }
169 :
170 : // If beam energy is not set from Config.C (RELDIS / spectator generators)
171 1 : if(fBeamEnergy<0.01){
172 1 : fBeamEnergy = grpData->GetBeamEnergy();
173 1 : if(fBeamEnergy==AliGRPObject::GetInvalidFloat()){
174 0 : AliWarning("GRP/GRP/Data entry: missing value for the beam energy ! Using 0.");
175 0 : AliError("\t UNKNOWN beam type from GRP obj -> PMT gains not set in ZDC digitizer!!!\n");
176 0 : fBeamEnergy = 0.;
177 0 : }
178 : }
179 :
180 : /*if(fIspASystem){
181 : fBeamType = "p-A";
182 : AliInfo(" AliZDCDigitizer -> Manually setting beam type to p-A\n");
183 : }*/
184 :
185 : // Setting beam type for spectator generator and RELDIS generator
186 1 : if(((fBeamType.CompareTo("UNKNOWN")) == 0) || fIsRELDISgen){
187 1 : fBeamType = "A-A";
188 1 : AliInfo(" AliZDCDigitizer -> Manually setting beam type to A-A\n");
189 1 : }
190 1 : printf("\n\t AliZDCDigitizer -> beam type %s - beam energy = %f GeV\n", fBeamType.Data(), fBeamEnergy);
191 1 : if(fSpectators2Track) printf("\n\t AliZDCDigitizer -> spectator signal added at digit level\n\n");
192 :
193 1 : if(fBeamEnergy>0.1){
194 1 : ReadPMTGains();
195 : //CalculatePMTGains();
196 1 : }
197 : else{
198 0 : AliWarning("\n Beam energy is 0 -> ZDC PMT gains can't be set -> NO ZDC DIGITS!!!\n");
199 : }
200 :
201 : // ADC Caen V965
202 1 : fADCRes[0] = 0.0000008; // ADC Resolution high gain: 200 fC/adcCh
203 1 : fADCRes[1] = 0.0000064; // ADC Resolution low gain: 25 fC/adcCh
204 :
205 2 : if(fSpectatorParam==1) {printf("\t AliZDCDigitizer -> spectator kinematic from AliGenZDC generator\n\n");}
206 0 : else if(fSpectatorParam==2) {printf("\t AliZDCDigitizer -> spectator kinematic properties from HIJING\n\n");}
207 : //{printf(" AliZDCDigitizer -> spectator kinematic properties from HIJING generator\n\n");}
208 :
209 1 : return kTRUE;
210 1 : }
211 :
212 : //____________________________________________________________________________
213 : void AliZDCDigitizer::Digitize(Option_t* /*option*/)
214 : {
215 : // Execute digitization
216 :
217 : // ------------------------------------------------------------
218 : // !!! 2nd ZDC set added
219 : // *** 1st 3 arrays are digits from REAL (simulated) hits
220 : // *** last 2 are copied from simulated digits
221 : // --- pm[0][...] = light in ZN side C [C, Q1, Q2, Q3, Q4]
222 : // --- pm[1][...] = light in ZP side C [C, Q1, Q2, Q3, Q4]
223 : // --- pm[2][...] = light in ZEM [x, 1, 2, x, x]
224 : // --- pm[3][...] = light in ZN side A [C, Q1, Q2, Q3, Q4] ->NEW!
225 : // --- pm[4][...] = light in ZP side A [C, Q1, Q2, Q3, Q4] ->NEW!
226 : // ------------------------------------------------------------
227 8 : Float_t pm[5][5];
228 48 : for(Int_t iSector1=0; iSector1<5; iSector1++)
229 240 : for(Int_t iSector2=0; iSector2<5; iSector2++){
230 100 : pm[iSector1][iSector2] = 0;
231 : }
232 :
233 : // ------------------------------------------------------------
234 : // ### Out of time ADC added (22 channels)
235 : // --- same codification as for signal PTMs (see above)
236 : // ------------------------------------------------------------
237 : // Float_t pmoot[5][5];
238 : // for(Int_t iSector1=0; iSector1<5; iSector1++)
239 : // for(Int_t iSector2=0; iSector2<5; iSector2++){
240 : // pmoot[iSector1][iSector2] = 0;
241 : // }
242 :
243 : // impact parameter and number of spectators
244 : Float_t impPar = 0;
245 : Int_t specNTarg = 0, specPTarg = 0;
246 : Int_t specNProj = 0, specPProj = 0;
247 : Float_t signalTime0 = 0.;
248 :
249 : // loop over input streams
250 16 : for(Int_t iInput = 0; iInput<fDigInput->GetNinputs(); iInput++){
251 :
252 : // get run loader and ZDC loader
253 : AliRunLoader* runLoader =
254 4 : AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(iInput));
255 4 : AliLoader* loader = runLoader->GetLoader("ZDCLoader");
256 4 : if(!loader) continue;
257 :
258 : // load sdigits
259 4 : loader->LoadSDigits();
260 4 : TTree* treeS = loader->TreeS();
261 4 : if(!treeS) continue;
262 4 : AliZDCSDigit sdigit;
263 4 : AliZDCSDigit* psdigit = &sdigit;
264 4 : treeS->SetBranchAddress("ZDC", &psdigit);
265 :
266 : // loop over sdigits
267 18 : for(Int_t iSDigit=0; iSDigit<treeS->GetEntries(); iSDigit++){
268 2 : treeS->GetEntry(iSDigit);
269 : //
270 2 : if(!psdigit) continue;
271 4 : if((sdigit.GetSector(1) < 0) || (sdigit.GetSector(1) > 4)){
272 0 : AliError(Form("\nsector[0] = %d, sector[1] = %d\n",
273 : sdigit.GetSector(0), sdigit.GetSector(1)));
274 : continue;
275 : }
276 : // Checking if signal is inside ADC gate
277 4 : if(iSDigit==0) signalTime0 = sdigit.GetTrackTime();
278 : else{
279 : // Assuming a signal lenght of 20 ns, signal is in gate if
280 : // signal ENDS in signalTime0+50. -> sdigit.GetTrackTime()+20<=signalTime0+50.
281 0 : if(sdigit.GetTrackTime()<=signalTime0+30.) fIsSignalInADCGate = kTRUE;
282 0 : if(sdigit.GetTrackTime()>signalTime0+30.){
283 0 : fIsSignalInADCGate = kFALSE;
284 : // Vedi quaderno per spiegazione approx. usata nel calcolo della fraz. di segnale perso
285 0 : fFracLostSignal = (sdigit.GetTrackTime()-30)*(sdigit.GetTrackTime()-30)/280.;
286 0 : }
287 : }
288 :
289 2 : Float_t sdSignal = sdigit.GetLightPM();
290 2 : if(fIsSignalInADCGate == kFALSE){
291 10 : AliDebug(2,Form("\t Signal time %f -> fraction %f of ZDC signal for det.(%d, %d) out of ADC gate\n",
292 : sdigit.GetTrackTime(),fFracLostSignal,sdigit.GetSector(0),sdigit.GetSector(1)));
293 2 : sdSignal = (1-fFracLostSignal)*sdSignal;
294 2 : }
295 :
296 2 : pm[(sdigit.GetSector(0))-1][sdigit.GetSector(1)] += sdigit.GetLightPM();
297 : //Ch. debug
298 : /*printf("\t Detector %d, Tower %d -> pm[%d][%d] = %.0f \n",
299 : sdigit.GetSector(0), sdigit.GetSector(1),sdigit.GetSector(0)-1,
300 : sdigit.GetSector(1), pm[sdigit.GetSector(0)-1][sdigit.GetSector(1)]);
301 : */
302 :
303 2 : }
304 :
305 4 : loader->UnloadSDigits();
306 :
307 : // get the impact parameter and the number of spectators in case of hijing
308 8 : if(!runLoader->GetAliRun()) runLoader->LoadgAlice();
309 4 : AliHeader* header = runLoader->GetHeader();
310 4 : if(!header) continue;
311 4 : AliGenEventHeader* genHeader = header->GenEventHeader();
312 4 : if(!genHeader) continue;
313 : AliGenHijingEventHeader *hijingHeader = 0;
314 12 : if(genHeader->InheritsFrom(AliGenHijingEventHeader::Class())) hijingHeader = dynamic_cast <AliGenHijingEventHeader*> (genHeader);
315 12 : else if(genHeader->InheritsFrom(AliGenCocktailEventHeader::Class())){
316 4 : TList* listOfHeaders = ((AliGenCocktailEventHeader*) genHeader)->GetHeaders();
317 4 : if(listOfHeaders){
318 64 : for(Int_t iH = 0; iH < listOfHeaders->GetEntries(); ++iH) {
319 64 : AliGenEventHeader *currHeader = dynamic_cast <AliGenEventHeader *> (listOfHeaders->At(iH));
320 64 : if (currHeader && currHeader->InheritsFrom(AliGenHijingEventHeader::Class())) {
321 0 : hijingHeader = dynamic_cast <AliGenHijingEventHeader*> (currHeader);
322 0 : break;
323 : }
324 16 : }
325 4 : }
326 : else{
327 0 : printf(" No list of headers from generator \n");
328 : }
329 4 : }
330 4 : if(!hijingHeader){
331 4 : printf(" No HIJING header found in list of headers from generator\n");
332 : }
333 :
334 4 : if(hijingHeader && (!hijingHeader->GetSpectatorsInTheStack())){
335 0 : impPar = hijingHeader->ImpactParameter();
336 0 : Int_t freeSpecNProj=0, freeSpecPProj=0;
337 0 : Int_t freeSpecNTarg=0, freeSpecPTarg=0;
338 0 : if(!hijingHeader->GetFragmentationFromData()){ //No. of free spectators from our fragmentation model
339 0 : specNProj = hijingHeader->ProjSpectatorsn();
340 0 : specPProj = hijingHeader->ProjSpectatorsp();
341 0 : specNTarg = hijingHeader->TargSpectatorsn();
342 0 : specPTarg = hijingHeader->TargSpectatorsp();
343 0 : if(specNProj!=0 || specPProj!=0) Fragmentation(impPar, specNProj, specPProj, freeSpecNProj, freeSpecPProj);
344 0 : if(specNTarg!=0 || specPTarg!=0) Fragmentation(impPar, specNTarg, specPTarg, freeSpecNTarg, freeSpecPTarg);
345 : }
346 : else{ //No. of free spectators from data driven correction (written in HIJING ev. header)
347 0 : freeSpecNProj = hijingHeader->GetFreeProjSpecn();
348 0 : freeSpecPProj = hijingHeader->GetFreeProjSpecp();
349 0 : freeSpecNTarg = hijingHeader->GetFreeTargSpecn();
350 0 : freeSpecPTarg = hijingHeader->GetFreeTargSpecp();
351 : }
352 0 : printf("\t AliZDCDigitizer -> Adding spectator signal for PROJECTILE: %d free n and %d free p\n",freeSpecNProj,freeSpecPProj);
353 0 : printf("\t AliZDCDigitizer -> Adding spectator signal for TARGET: %d free n and %d free p\n",freeSpecNTarg,freeSpecPTarg);
354 : //
355 0 : if(fSpectatorParam==1){
356 0 : if(freeSpecNProj!=0) SpectatorSignal(1, freeSpecNProj, pm);
357 0 : if(freeSpecPProj!=0) SpectatorSignal(2, freeSpecPProj, pm);
358 0 : if(freeSpecNTarg!=0) SpectatorSignal(3, freeSpecNTarg, pm);
359 0 : if(freeSpecPTarg!=0) SpectatorSignal(4, freeSpecPTarg, pm);
360 : }
361 0 : else if(fSpectatorParam==2){
362 0 : SpectatorsFromHijing(1, freeSpecNProj, pm);
363 0 : SpectatorsFromHijing(2, freeSpecPProj, pm);
364 0 : SpectatorsFromHijing(3, freeSpecNTarg, pm);
365 0 : SpectatorsFromHijing(4, freeSpecPTarg, pm);
366 : }
367 0 : }
368 8 : }
369 :
370 :
371 : // get the output run loader and loader
372 : AliRunLoader* runLoader =
373 4 : AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
374 4 : AliLoader* loader = runLoader->GetLoader("ZDCLoader");
375 4 : if(!loader) {
376 0 : AliError("no ZDC loader found");
377 0 : return;
378 : }
379 :
380 : // create the output tree
381 : const char* mode = "update";
382 5 : if(runLoader->GetEventNumber() == 0) mode = "recreate";
383 4 : loader->LoadDigits(mode);
384 4 : loader->MakeTree("D");
385 4 : TTree* treeD = loader->TreeD();
386 4 : AliZDCDigit digit;
387 4 : AliZDCDigit* pdigit = &digit;
388 : const Int_t kBufferSize = 4000;
389 4 : treeD->Branch("ZDC", "AliZDCDigit", &pdigit, kBufferSize);
390 :
391 : // Create digits
392 4 : Int_t sector[2];
393 4 : Int_t digi[2], digioot[2];
394 48 : for(sector[0]=1; sector[0]<6; sector[0]++){
395 360 : for(sector[1]=0; sector[1]<5; sector[1]++){
396 256 : if((sector[0]==3) && ((sector[1]<1) || (sector[1]>2))) continue;
397 528 : for(Int_t res=0; res<2; res++){
398 352 : digi[res] = Phe2ADCch(sector[0], sector[1], pm[sector[0]-1][sector[1]], res)
399 352 : + Pedestal(sector[0], sector[1], res);
400 : }
401 176 : new(pdigit) AliZDCDigit(sector, digi);
402 100 : treeD->Fill();
403 :
404 : //Ch. debug
405 : //printf("\t DIGIT added -> det %d quad %d - digi[0,1] = [%d, %d]\n",
406 : // sector[0], sector[1], digi[0], digi[1]); // Chiara debugging!
407 :
408 : }
409 : } // Loop over detector
410 : // Adding in-time digits for 2 reference PTM signals (after signal ch.)
411 : // (for the moment the ref. signal is completely invented assuming a PMgain of 5*10^4!)
412 4 : Int_t sectorRef[2];
413 4 : sectorRef[1] = 5;
414 4 : Int_t sigRef[2];
415 : // Reference signal are set to 100 (high gain chain) and 800 (low gain chain)
416 8 : if(fIsCalibration==0) {sigRef[0]=100; sigRef[1]=800;}
417 0 : else {sigRef[0]=0; sigRef[1]=0;} // calibration -> simulation of pedestal values
418 : //
419 24 : for(Int_t iref=0; iref<2; iref++){
420 8 : sectorRef[0] = 3*iref+1;
421 48 : for(Int_t res=0; res<2; res++){
422 32 : sigRef[res] += Pedestal(sectorRef[0], sectorRef[1], res);
423 : }
424 16 : new(pdigit) AliZDCDigit(sectorRef, sigRef);
425 8 : treeD->Fill();
426 :
427 : //Ch. debug
428 : //printf("\t RefDigit added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n",
429 : // sectorRef[0], sectorRef[1], sigRef[0], sigRef[1]); // Chiara debugging!
430 : }
431 : //
432 : // --- Adding digits for out-of-time channels after signal digits
433 48 : for(sector[0]=1; sector[0]<6; sector[0]++){
434 360 : for(sector[1]=0; sector[1]<5; sector[1]++){
435 256 : if((sector[0]==3) && ((sector[1]<1) || (sector[1]>2))) continue;
436 528 : for(Int_t res=0; res<2; res++){
437 352 : digioot[res] = Pedestal(sector[0], sector[1], res); // out-of-time ADCs
438 : }
439 176 : new(pdigit) AliZDCDigit(sector, digioot);
440 100 : treeD->Fill();
441 :
442 : //Ch. debug
443 : //printf("\t DIGIToot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n",
444 : // sector[0], sector[1], digioot[0], digioot[1]); // Chiara debugging!
445 : }
446 : }
447 : // Adding out-of-time digits for 2 reference PTM signals (after out-of-time ch.)
448 4 : Int_t sigRefoot[2];
449 24 : for(Int_t iref=0; iref<2; iref++){
450 8 : sectorRef[0] = 3*iref+1;
451 48 : for(Int_t res=0; res<2; res++){
452 32 : sigRefoot[res] = Pedestal(sectorRef[0], sectorRef[1], res);
453 : }
454 16 : new(pdigit) AliZDCDigit(sectorRef, sigRefoot);
455 8 : treeD->Fill();
456 : //Ch. debug
457 : //printf("\t RefDigitoot added -> det = %d, quad = %d - digi[0,1] = [%d, %d]\n",
458 : // sectorRef[0], sectorRef[1], sigRefoot[0], sigRefoot[1]); // Chiara debugging!
459 :
460 : }
461 : //printf("\t AliZDCDigitizer -> TreeD has %d entries\n",(Int_t) treeD->GetEntries());
462 :
463 : // write the output tree
464 4 : loader->WriteDigits("OVERWRITE");
465 4 : loader->UnloadDigits();
466 8 : }
467 :
468 :
469 : //_____________________________________________________________________________
470 : void AliZDCDigitizer::ReadPMTGains()
471 : {
472 : // Read PMT gain from an external file
473 :
474 2 : char *fname = gSystem->ExpandPathName("$ALICE_ROOT/ZDC/PMTGainsdata.txt");
475 1 : FILE *fdata = fopen(fname,"r");
476 1 : if(fdata==NULL){
477 0 : AliWarning(" Can't open file $ALICE_ROOT/ZDC/PMTGainsdata.txt to read ZDC PMT Gains\n");
478 0 : AliWarning(" -> ZDC signal will be pedestal!!!!!!!!!!!!\n\n");
479 0 : return;
480 : }
481 : int read=1;
482 1 : Float_t data[5];
483 1 : Int_t beam[12], det[12];
484 1 : Float_t gain[12], aEne[12], bEne[12];
485 26 : for(int ir=0; ir<12; ir++){
486 144 : for(int ic=0; ic<5; ic++){
487 60 : read = fscanf(fdata,"%f ",&data[ic]);
488 60 : if(read==0) AliDebug(3, " Error in reading PMT gains from external file ");
489 : }
490 12 : beam[ir] = (int) (data[0]);
491 12 : det[ir] = (int) (data[1]);
492 12 : gain[ir] = data[2];
493 12 : aEne[ir] = data[3];
494 12 : bEne[ir] = data[4];
495 : }
496 1 : fclose(fdata);
497 :
498 1 : if(((fBeamType.CompareTo("A-A")) == 0)){
499 26 : for(int i=0; i<12; i++){
500 12 : if(beam[i]==1){
501 6 : Float_t scalGainFactor = fBeamEnergy/2760.;
502 11 : if(det[i]!=31 && det[i]!=32){
503 48 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = gain[i]/(aEne[i]*scalGainFactor);
504 4 : }
505 : else{
506 12 : for(int iq=1; iq<3; iq++) fPMGain[2][iq] = gain[i]/(aEne[i]*scalGainFactor);
507 : }
508 6 : }
509 : }
510 : //
511 1 : printf("\n AliZDCDigitizer::ReadPMTGains -> ZDC PMT gains for %s @ %1.0f+%1.0f A GeV: ZN(%1.0f), ZP(%1.0f), ZEM(%1.0f)\n",
512 1 : fBeamType.Data(),fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1]);
513 1 : }
514 0 : else if(((fBeamType.CompareTo("p-A")) == 0) || ((fBeamType.CompareTo("P-A")) == 0)){
515 0 : for(int i=0; i<12; i++){
516 0 : if(beam[i]==0 && fBeamEnergy!=0.){
517 0 : if(det[i]==1 || det[i]==2){ // Setting ZNC/ZPC gains as for p-p
518 0 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = gain[i]*(aEne[i]/fBeamEnergy+bEne[i]);
519 0 : }
520 : }
521 0 : if(beam[i]==1){
522 0 : Float_t scalGainFactor = fBeamEnergy/2760.;
523 : Float_t npartScalingFactor = 208./15.;
524 0 : if(det[i]==4 || det[i]==5){ // Setting ZNA/ZPA gains as for Pb-Pb
525 0 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = npartScalingFactor*gain[i]/(aEne[i]*scalGainFactor);
526 0 : }
527 0 : else if(det[i]==31 || det[i]==32){ // Setting ZEM gains as for Pb-Pb
528 0 : for(int iq=1; iq<3; iq++) fPMGain[2][iq] = npartScalingFactor*gain[i]/(aEne[i]*scalGainFactor);
529 0 : }
530 0 : }
531 : }
532 0 : printf("\n AliZDCDigitizer::ReadPMTGains -> ZDC PMT gains for %s: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
533 0 : fBeamType.Data(),fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
534 0 : }
535 0 : else if(((fBeamType.CompareTo("A-p")) == 0) || ((fBeamType.CompareTo("A-P")) == 0)){
536 0 : for(int i=0; i<12; i++){
537 0 : if(beam[i]==0 && fBeamEnergy!=0.){ // Setting ZNA/ZPA gains as for p-p
538 0 : if(det[i]==4 || det[i]==5){
539 0 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = gain[i]*(aEne[i]/fBeamEnergy+bEne[i]);
540 0 : }
541 : // Setting ZEM gains as for p-p
542 0 : else if(det[i] == 31) fPMGain[2][1] = gain[i]*(aEne[i]-fBeamEnergy*bEne[i]);
543 0 : else if(det[i] == 32) fPMGain[2][2] = gain[i]*(aEne[i]-fBeamEnergy*bEne[i]);
544 : }
545 0 : if(beam[i]==1){
546 0 : Float_t scalGainFactor = fBeamEnergy/2760.;
547 : Float_t npartScalingFactor = 208./15.;
548 0 : if(det[i]==1 || det[i]==2){ // Setting ZNC/ZPC gains as for Pb-Pb
549 0 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = npartScalingFactor*gain[i]/(aEne[i]*scalGainFactor);
550 0 : }
551 0 : }
552 : }
553 0 : printf("\n AliZDCDigitizer::ReadPMTGains -> ZDC PMT gains for %s: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
554 0 : fBeamType.Data(),fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
555 0 : }
556 : else{ // Setting pp gains as default value
557 0 : for(int i=0; i<12; i++){
558 0 : if(beam[i]==0 && fBeamEnergy!=0.){
559 0 : if(det[i]!=31 && det[i]!=32){
560 0 : for(Int_t j=0; j<5; j++) fPMGain[det[i]-1][j] = gain[i]*(aEne[i]/fBeamEnergy+bEne[i]);
561 0 : }
562 0 : else if(det[i] == 31) fPMGain[2][1] = gain[i]*(aEne[i]-fBeamEnergy*bEne[i]);
563 0 : else if(det[i] == 32) fPMGain[2][2] = gain[i]*(aEne[i]-fBeamEnergy*bEne[i]);
564 : }
565 : }
566 : //
567 0 : printf("\n AliZDCDigitizer::ReadPMTGains -> ZDC PMT gains for %s @ %1.0f+%1.0f GeV: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
568 0 : fBeamType.Data(),fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
569 : }
570 :
571 2 : }
572 :
573 : //_____________________________________________________________________________
574 : void AliZDCDigitizer::CalculatePMTGains()
575 : {
576 : // Calculate PMT gain according to beam type and beam energy
577 0 : if( ((fBeamType.CompareTo("P-P")) == 0) || ((fBeamType.CompareTo("p-p"))) ){
578 : // PTM gains rescaled to beam energy for p-p
579 : // New correction coefficients for PMT gains needed
580 : // to reproduce experimental spectra (from Grazia Jul 2010)
581 0 : if(fBeamEnergy != 0){
582 0 : for(Int_t j = 0; j < 5; j++){
583 0 : fPMGain[0][j] = 1.515831*(661.444/fBeamEnergy+0.000740671)*10000000;
584 0 : fPMGain[1][j] = 0.674234*(864.350/fBeamEnergy+0.00234375)*10000000;
585 0 : fPMGain[3][j] = 1.350938*(661.444/fBeamEnergy+0.000740671)*10000000;
586 0 : fPMGain[4][j] = 0.678597*(864.350/fBeamEnergy+0.00234375)*10000000;
587 : }
588 0 : fPMGain[2][1] = 0.869654*(1.32312-0.000101515*fBeamEnergy)*10000000;
589 0 : fPMGain[2][2] = 1.030883*(1.32312-0.000101515*fBeamEnergy)*10000000;
590 : //
591 0 : printf("\n AliZDCDigitizer::CalculatePMTGains -> ZDC PMT gains for p-p @ %1.0f+%1.0f GeV: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
592 0 : fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
593 :
594 0 : }
595 : }
596 0 : else if(((fBeamType.CompareTo("A-A")) == 0)){
597 : // PTM gains for Pb-Pb @ 2.7+2.7 A TeV ***************
598 : // rescaled for Pb-Pb @ 1.38+1.38 A TeV ***************
599 : // Values corrected after 2010 Pb-Pb data taking (7/2/2011 - Ch.)
600 : // Experimental data compared to EMD simulation for single nucleon peaks:
601 : // ZN gains must be divided by 4, ZP gains by 10!
602 0 : Float_t scalGainFactor = fBeamEnergy/2760.;
603 0 : for(Int_t j = 0; j < 5; j++){
604 0 : fPMGain[0][j] = 50000./(4*scalGainFactor); // ZNC
605 0 : fPMGain[1][j] = 100000./(5*scalGainFactor); // ZPC
606 0 : fPMGain[2][j] = 100000./scalGainFactor; // ZEM
607 0 : fPMGain[3][j] = 50000./(4*scalGainFactor); // ZNA
608 0 : fPMGain[4][j] = 100000./(5*scalGainFactor); // ZPA
609 : }
610 0 : printf("\n AliZDCDigitizer::CalculatePMTGains -> ZDC PMT gains for Pb-Pb @ %1.0f+%1.0f A GeV: ZN(%1.0f), ZP(%1.0f), ZEM(%1.0f)\n",
611 0 : fBeamEnergy, fBeamEnergy, fPMGain[0][0], fPMGain[1][0], fPMGain[2][1]);
612 0 : }
613 0 : else if(((fBeamType.CompareTo("p-A")) == 0) || ((fBeamType.CompareTo("P-A"))) ){
614 : // PTM gains for Pb-Pb @ 1.38+1.38 A TeV on side A
615 : // PTM gains rescaled to beam energy for p-p on side C
616 : // WARNING! Energies are set by hand for 2011 pA RUN!!!
617 0 : Float_t scalGainFactor = fBeamEnergy/2760.;
618 : Float_t npartScalingFactor = 208./15.;
619 :
620 0 : for(Int_t j = 0; j < 5; j++){
621 0 : fPMGain[0][j] = 1.515831*(661.444/fBeamEnergy+0.000740671)*10000000; //ZNC (p)
622 0 : fPMGain[1][j] = 0.674234*(864.350/fBeamEnergy+0.00234375)*10000000; //ZPC (p)
623 0 : if(j<2) fPMGain[2][j] = npartScalingFactor*100000./scalGainFactor; // ZEM (Pb)
624 : // Npart max scales from 400 in Pb-Pb to ~8 in pPb -> *40.
625 0 : fPMGain[3][j] = npartScalingFactor*50000/(4*scalGainFactor); // ZNA (Pb)
626 0 : fPMGain[4][j] = npartScalingFactor*100000/(5*scalGainFactor); // ZPA (Pb)
627 : }
628 0 : printf("\n AliZDCDigitizer::CalculatePMTGains -> ZDC PMT gains for p-Pb: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
629 0 : fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
630 0 : }
631 0 : else if(((fBeamType.CompareTo("A-p")) == 0) || ((fBeamType.CompareTo("A-P"))) ){
632 : // PTM gains for Pb-Pb @ 1.38+1.38 A TeV on side A
633 : // PTM gains rescaled to beam energy for p-p on side C
634 0 : Float_t scalGainFactor = fBeamEnergy/2760.;
635 : Float_t npartScalingFactor = 208./15.;
636 :
637 0 : for(Int_t j = 0; j < 5; j++){
638 0 : fPMGain[3][j] = 1.350938*(661.444/fBeamEnergy+0.000740671)*10000000; //ZNA (p)
639 0 : fPMGain[4][j] = 0.678597*(864.350/fBeamEnergy+0.00234375)*10000000; //ZPA (p)
640 : // Npart max scales from 400 in Pb-Pb to ~8 in pPb -> *40.
641 0 : fPMGain[1][j] = npartScalingFactor*50000/(4*scalGainFactor); // ZNC (Pb)
642 0 : fPMGain[2][j] = npartScalingFactor*100000/(5*scalGainFactor); // ZPC (Pb)
643 : }
644 0 : fPMGain[2][1] = 0.869654*(1.32312-0.000101515*fBeamEnergy)*10000000; // ZEM (pp)
645 0 : fPMGain[2][2] = 1.030883*(1.32312-0.000101515*fBeamEnergy)*10000000; // ZEM (pp)
646 0 : printf("\n AliZDCDigitizer::CalculatePMTGains -> ZDC PMT gains for Pb-p: ZNC(%1.0f), ZPC(%1.0f), ZEM(%1.0f), ZNA(%1.0f) ZPA(%1.0f)\n",
647 0 : fPMGain[0][0], fPMGain[1][0], fPMGain[2][1], fPMGain[3][0], fPMGain[4][0]);
648 0 : }
649 0 : }
650 :
651 : //_____________________________________________________________________________
652 : void AliZDCDigitizer::Fragmentation(Float_t impPar, Int_t specN, Int_t specP,
653 : Int_t &freeSpecN, Int_t &freeSpecP) const
654 : {
655 : // simulate fragmentation of spectators
656 :
657 0 : AliZDCFragment frag(impPar);
658 :
659 : // Fragments generation
660 0 : frag.GenerateIMF();
661 0 : Int_t nAlpha = frag.GetNalpha();
662 :
663 : // Attach neutrons
664 0 : frag.AttachNeutrons();
665 0 : Int_t ztot = frag.GetZtot();
666 0 : Int_t ntot = frag.GetNtot();
667 :
668 : // Removing fragments and alpha pcs
669 0 : freeSpecN = specN-ntot-2*nAlpha;
670 0 : freeSpecP = specP-ztot-2*nAlpha;
671 :
672 : // Removing deuterons
673 0 : Int_t ndeu = (Int_t) (freeSpecN*frag.DeuteronNumber());
674 0 : freeSpecN -= ndeu;
675 0 : freeSpecP -= ndeu;
676 : //
677 0 : if(freeSpecN<0) freeSpecN=0;
678 0 : if(freeSpecP<0) freeSpecP=0;
679 0 : AliDebug(2, Form("FreeSpn = %d, FreeSpp = %d", freeSpecN, freeSpecP));
680 0 : }
681 :
682 : //_____________________________________________________________________________
683 : void AliZDCDigitizer::SpectatorSignal(Int_t specType, Int_t numEvents, Float_t pm[5][5])
684 : {
685 : // add signal of the spectators
686 :
687 0 : if(!fSpectatorData) fSpectatorData = TFile::Open("$ALICE_ROOT/ZDC/SpectatorSignal.root");
688 0 : if(!fSpectatorData || !fSpectatorData->IsOpen()) {
689 0 : AliError((" No file $ALICE_ROOT/ZDC/SpectatorSignal.root found -> No spectators added!!!\n"));
690 0 : return;
691 : }
692 :
693 0 : TNtuple* zdcSignal=0x0;
694 :
695 0 : Float_t sqrtS = 2*fBeamEnergy;
696 : //
697 0 : if(TMath::Abs(sqrtS-5500) < 100.){
698 0 : AliInfo(" Extracting signal from SpectatorSignal/energy5500 ");
699 0 : fSpectatorData->cd("energy5500");
700 : //
701 0 : if(specType == 1) { // --- Signal for projectile spectator neutrons
702 0 : fSpectatorData->GetObject("energy5500/ZNCSignal;1",zdcSignal);
703 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorSignal.root file");
704 : }
705 0 : else if(specType == 2) { // --- Signal for projectile spectator protons
706 0 : fSpectatorData->GetObject("energy5500/ZPCSignal;1",zdcSignal);
707 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorSignal.root file");
708 : }
709 0 : else if(specType == 3) { // --- Signal for target spectator neutrons
710 0 : fSpectatorData->GetObject("energy5500/ZNASignal;1",zdcSignal);
711 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorSignal.root file");
712 : }
713 0 : else if(specType == 4) { // --- Signal for target spectator protons
714 0 : fSpectatorData->GetObject("energy5500/ZPASignal;1",zdcSignal);
715 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorSignal.root file");
716 : }
717 : }
718 0 : else if(TMath::Abs(sqrtS-5000) < 100.){
719 0 : AliInfo(" Extracting signal from SpectatorSignal/energy5020 ");
720 0 : fSpectatorData->cd("energy5020");
721 : //
722 0 : if(specType == 1) { // --- Signal for projectile spectator neutrons
723 0 : fSpectatorData->GetObject("energy5020/ZNCSignal;1",zdcSignal);
724 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorSignal.root file");
725 : }
726 0 : else if(specType == 2) { // --- Signal for projectile spectator protons
727 0 : fSpectatorData->GetObject("energy5020/ZPCSignal;1",zdcSignal);
728 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorSignal.root file");
729 : }
730 0 : else if(specType == 3) { // --- Signal for target spectator neutrons
731 0 : fSpectatorData->GetObject("energy5020/ZNASignal;1",zdcSignal);
732 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorSignal.root file");
733 : }
734 0 : else if(specType == 4) { // --- Signal for target spectator protons
735 0 : fSpectatorData->GetObject("energy5020/ZPASignal;1",zdcSignal);
736 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorSignal.root file");
737 : }
738 : }
739 0 : else if(TMath::Abs(sqrtS-2760) < 100.){
740 0 : AliInfo(" Extracting signal from SpectatorSignal/energy2760 ");
741 0 : fSpectatorData->cd("energy2760");
742 : //
743 0 : if(specType == 1) { // --- Signal for projectile spectator neutrons
744 0 : fSpectatorData->GetObject("energy2760/ZNCSignal;1",zdcSignal);
745 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorSignal.root file");
746 : }
747 0 : else if(specType == 2) { // --- Signal for projectile spectator protons
748 0 : fSpectatorData->GetObject("energy2760/ZPCSignal;1",zdcSignal);
749 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorSignal.root file");
750 : }
751 0 : else if(specType == 3) { // --- Signal for target spectator neutrons
752 0 : fSpectatorData->GetObject("energy2760/ZNASignal;1",zdcSignal);
753 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorSignal.root file");
754 : }
755 0 : else if(specType == 4) { // --- Signal for target spectator protons
756 0 : fSpectatorData->GetObject("energy2760/ZPASignal;1",zdcSignal);
757 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorSignal.root file");
758 : }
759 : }
760 :
761 0 : if(!zdcSignal){
762 0 : printf("\n No spectator signal available for ZDC digitization\n");
763 0 : return;
764 : }
765 :
766 0 : Int_t nentries = (Int_t) zdcSignal->GetEntries();
767 :
768 : Float_t *entry;
769 0 : Int_t rnd[125], volume[2]={0,0}, iev=0;
770 0 : for(int pl=0;pl<125;pl++) rnd[pl] = 0;
771 0 : if(numEvents > 125) {
772 0 : AliDebug(2,Form("numEvents (%d) is larger than 125", numEvents));
773 : numEvents = 125;
774 0 : }
775 0 : for(int pl=0;pl<numEvents;pl++){
776 0 : rnd[pl] = (Int_t) (9999*gRandom->Rndm());
777 0 : if(rnd[pl] >= 9999) rnd[pl] = 9998;
778 : //printf(" rnd[%d] = %d\n",pl,rnd[pl]);
779 : }
780 : // Sorting vector in ascending order with C function QSORT
781 0 : qsort((void*)rnd,numEvents,sizeof(Int_t),comp);
782 0 : do{
783 0 : for(int i=0; i<nentries; i++){
784 0 : zdcSignal->GetEvent(i);
785 0 : entry = zdcSignal->GetArgs();
786 0 : if(entry[0] == rnd[iev]){
787 0 : for(int k=0; k<2; k++) volume[k] = (Int_t) entry[k+1];
788 : //
789 0 : Float_t lightQ = entry[7];
790 0 : Float_t lightC = entry[8];
791 : //
792 0 : if(volume[0] != 3) { // ZN or ZP
793 0 : pm[volume[0]-1][0] += lightC;
794 0 : pm[volume[0]-1][volume[1]] += lightQ;
795 : //printf("\n pm[%d][0] = %.0f, pm[%d][%d] = %.0f\n",(volume[0]-1),pm[volume[0]-1][0],
796 : // (volume[0]-1),volume[1],pm[volume[0]-1][volume[1]]);
797 0 : }
798 : else {
799 0 : if(volume[1] == 1) pm[2][1] += lightC; // ZEM 1
800 0 : else pm[2][2] += lightQ; // ZEM 2
801 : //printf("\n pm[2][1] = %.0f, pm[2][2] = %.0f\n",pm[2][1],pm[2][2]);
802 : }
803 0 : }
804 0 : else if(entry[0] > rnd[iev]){
805 0 : iev++;
806 0 : continue;
807 : }
808 : }
809 0 : }while(iev<numEvents);
810 :
811 0 : }
812 :
813 : //_____________________________________________________________________________
814 : void AliZDCDigitizer::SpectatorsFromHijing(Int_t specType, Int_t numEvents, Float_t pm[5][5])
815 : {
816 : // Getting sigmal from spectators from Hijing+data driven correction
817 0 : if(!fSpectatorData) fSpectatorData = TFile::Open("$ALICE_ROOT/ZDC/SpectatorsFromData.root");
818 0 : if(!fSpectatorData || !fSpectatorData->IsOpen()) {
819 0 : AliError((" No file $ALICE_ROOT/ZDC/SpectatorsFromData.root found -> No spectators added!!!\n"));
820 0 : return;
821 : }
822 :
823 0 : TNtuple* zdcSignal=0x0;
824 :
825 0 : AliInfo(" Extracting signal from SpectatorsFromData.root");
826 : //
827 0 : if(specType == 1) { // --- Signal for projectile spectator neutrons
828 0 : fSpectatorData->GetObject("energy5020/ZNCSignal;1",zdcSignal);
829 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNCSignal from SpectatorsFromData.root file");
830 : }
831 0 : else if(specType == 2) { // --- Signal for projectile spectator protons
832 0 : fSpectatorData->GetObject("energy5020/ZPCSignal;1",zdcSignal);
833 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPCSignal from SpectatorsFromData.root file");
834 : }
835 0 : else if(specType == 3) { // --- Signal for target spectator neutrons
836 0 : fSpectatorData->GetObject("energy5020/ZNASignal;1",zdcSignal);
837 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZNASignal from SpectatorsFromData.root file");
838 : }
839 0 : else if(specType == 4) { // --- Signal for target spectator protons
840 0 : fSpectatorData->GetObject("energy5020/ZPASignal;1",zdcSignal);
841 0 : if(!zdcSignal) AliError(" PROBLEM!!! Can't retrieve ZPASignal from SpectatorsFromData.root file");
842 : }
843 :
844 0 : if(!zdcSignal){
845 0 : printf(" !!!!!!!!!!!!! No spectator signal available for ZDC digitization\n");
846 0 : return;
847 : }
848 :
849 0 : Int_t nentries = (Int_t) zdcSignal->GetEntries();
850 :
851 : Float_t *entry;
852 0 : Int_t rnd[125], volume[2]={0,0}, iev=0, base=nentries/8.;
853 0 : for(int pl=0;pl<125;pl++) rnd[pl] = 0;
854 0 : if(numEvents > 125) {
855 0 : AliDebug(2,Form("numEvents (%d) is larger than 125", numEvents));
856 : numEvents = 125;
857 0 : }
858 0 : for(int pl=0;pl<numEvents;pl++){
859 0 : rnd[pl] = (Int_t) (9999*gRandom->Rndm());
860 0 : if(rnd[pl] >= 9999) rnd[pl] = 9998;
861 : //printf(" rnd[%d] = %d\n",pl,rnd[pl]);
862 : }
863 : // Sorting vector in ascending order with C function QSORT
864 0 : qsort((void*)rnd,numEvents,sizeof(Int_t),comp);
865 0 : do{
866 0 : for(int i=0; i<nentries; i++){
867 0 : zdcSignal->GetEvent(i);
868 0 : entry = zdcSignal->GetArgs();
869 0 : if(entry[0] == rnd[iev]){
870 0 : for(int k=0; k<2; k++) volume[k] = (Int_t) entry[k+1];
871 : //
872 0 : Float_t lightQ = entry[7];
873 0 : Float_t lightC = entry[8];
874 : //
875 0 : if(volume[0] != 3) { // ZN or ZP
876 0 : pm[volume[0]-1][0] += lightC;
877 0 : pm[volume[0]-1][volume[1]] += lightQ;
878 : //printf("\n pm[%d][0] = %.0f, pm[%d][%d] = %.0f\n",(volume[0]-1),pm[volume[0]-1][0],
879 : // (volume[0]-1),volume[1],pm[volume[0]-1][volume[1]]);
880 0 : }
881 : else {
882 0 : if(volume[1] == 1) pm[2][1] += lightC; // ZEM 1
883 0 : else pm[2][2] += lightQ; // ZEM 2
884 : //printf("\n pm[2][1] = %.0f, pm[2][2] = %.0f\n",pm[2][1],pm[2][2]);
885 : }
886 0 : }
887 0 : else if(entry[0] > rnd[iev]){
888 0 : iev++;
889 0 : continue;
890 : }
891 : }
892 0 : }while(iev<numEvents);
893 :
894 0 : }
895 :
896 :
897 : //_____________________________________________________________________________
898 : Int_t AliZDCDigitizer::Phe2ADCch(Int_t Det, Int_t Quad, Float_t Light,
899 : Int_t Res) const
900 : {
901 : // Evaluation of the ADC channel corresponding to the light yield Light
902 352 : Int_t vADCch = (Int_t) (Light * fPMGain[Det-1][Quad] * fADCRes[Res]);
903 : // Ch. debug
904 : //printf("\t Phe2ADCch -> det %d quad %d - PMGain[%d][%d] %1.0f phe %1.2f ADC %d\n",
905 : // Det,Quad,Det-1,Quad,fPMGain[Det-1][Quad],Light,vADCch);
906 :
907 176 : return vADCch;
908 : }
909 :
910 : //_____________________________________________________________________________
911 : Int_t AliZDCDigitizer::Pedestal(Int_t Det, Int_t Quad, Int_t Res) const
912 : {
913 : // Returns a pedestal for detector det, PM quad, channel with res.
914 : //
915 : Float_t pedValue=0.;
916 : // Normal run
917 768 : if(fIsCalibration == 0){
918 : Int_t index=0, kNch=24;
919 384 : if(Quad!=5){
920 432 : if(Det==1) index = Quad+kNch*Res; // ZNC
921 352 : else if(Det==2) index = (Quad+5)+kNch*Res; // ZPC
922 224 : else if(Det==3) index = (Quad+9)+kNch*Res; // ZEM
923 240 : else if(Det==4) index = (Quad+12)+kNch*Res; // ZNA
924 160 : else if(Det==5) index = (Quad+17)+kNch*Res; // ZPA
925 : }
926 32 : else index = (Det-1)/3+22+kNch*Res; // Reference PMs
927 : //
928 384 : if(fPedData){
929 384 : Float_t meanPed = fPedData->GetMeanPed(index);
930 384 : Float_t pedWidth = fPedData->GetMeanPedWidth(index);
931 384 : pedValue = gRandom->Gaus(meanPed,pedWidth);
932 : //
933 : /*printf("\t AliZDCDigitizer::Pedestal -> det %d quad %d res %d - Ped[%d] = %d\n",
934 : Det, Quad, Res, index,(Int_t) pedValue); // Chiara debugging!
935 : */
936 384 : }
937 : else{
938 0 : printf(" AliZDCDigitizer::Pedestal -> No valid pedestal calibration object loaded!\n\n");
939 : }
940 384 : }
941 : // To create calibration object
942 : else{
943 0 : if(Res == 0) pedValue = gRandom->Gaus((35.+10.*gRandom->Rndm()),(0.5+0.2*gRandom->Rndm())); //High gain
944 0 : else pedValue = gRandom->Gaus((250.+100.*gRandom->Rndm()),(3.5+2.*gRandom->Rndm())); //Low gain
945 : }
946 :
947 384 : return (Int_t) pedValue;
948 : }
949 :
950 : //_____________________________________________________________________________
951 : AliCDBStorage* AliZDCDigitizer::SetStorage(const char *uri)
952 : {
953 :
954 : Bool_t deleteManager = kFALSE;
955 :
956 0 : AliCDBManager *manager = AliCDBManager::Instance();
957 0 : AliCDBStorage *defstorage = manager->GetDefaultStorage();
958 :
959 0 : if(!defstorage || !(defstorage->Contains("ZDC"))){
960 0 : AliWarning("No default storage set or default storage doesn't contain ZDC!");
961 0 : manager->SetDefaultStorage(uri);
962 : deleteManager = kTRUE;
963 0 : }
964 :
965 0 : AliCDBStorage *storage = manager->GetDefaultStorage();
966 :
967 0 : if(deleteManager){
968 0 : AliCDBManager::Instance()->UnsetDefaultStorage();
969 : defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
970 0 : }
971 :
972 0 : return storage;
973 : }
974 :
975 : //_____________________________________________________________________________
976 : AliZDCPedestals* AliZDCDigitizer::GetPedData() const
977 : {
978 :
979 : // Getting pedestal calibration object for ZDC set
980 : AliZDCPedestals *calibdata = 0x0;
981 3 : AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
982 1 : if(!entry) AliFatal("No calibration data loaded!");
983 : else{
984 3 : calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
985 1 : if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
986 : }
987 1 : return calibdata;
988 0 : }
989 :
|
