1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///_________________________________________________________________________
20 /// This class constructs Digits out of Hits
24 // --- Standard library ---
26 // --- ROOT system ---
30 #include <TGeoManager.h>
31 #include <TGeoPhysicalNode.h>
32 #include <AliGeomManager.h>
35 // --- AliRoot header files ---
36 #include "AliVZEROConst.h"
39 #include "AliVZEROhit.h"
40 #include "AliRunLoader.h"
41 #include "AliLoader.h"
42 #include "AliGRPObject.h"
43 #include "AliRunDigitizer.h"
44 #include "AliCDBManager.h"
45 #include "AliCDBStorage.h"
46 #include "AliCDBEntry.h"
47 #include "AliVZEROCalibData.h"
49 #include "AliVZEROdigit.h"
50 #include "AliVZERODigitizer.h"
52 ClassImp(AliVZERODigitizer)
54 AliVZERODigitizer::AliVZERODigitizer()
56 fCalibData(GetCalibData()),
57 fPhotoCathodeEfficiency(0.18),
59 fPMGain(TMath::Power((fPMVoltage / 112.5) ,7.04277)),
66 // default constructor
71 // fPhotoCathodeEfficiency = 0.18;
72 // fPMVoltage = 768.0;
73 // fPMGain = TMath::Power((fPMVoltage / 112.5) ,7.04277);
75 // fCalibData = GetCalibData();
78 //____________________________________________________________________________
79 AliVZERODigitizer::AliVZERODigitizer(AliRunDigitizer* manager)
80 :AliDigitizer(manager),
81 fCalibData(GetCalibData()),
82 fPhotoCathodeEfficiency(0.18),
84 fPMGain(TMath::Power((fPMVoltage / 112.5) ,7.04277)),
96 // fPhotoCathodeEfficiency = 0.18;
97 // fPMVoltage = 768.0;
98 // fPMGain = TMath::Power( (fPMVoltage / 112.5) ,7.04277 );
100 // fCalibData = GetCalibData();
104 //____________________________________________________________________________
105 AliVZERODigitizer::~AliVZERODigitizer()
116 //_____________________________________________________________________________
117 Bool_t AliVZERODigitizer::Init()
119 // Initialises the digitizer
121 // Initialises the Digit array
122 fDigits = new TClonesArray ("AliVZEROdigit", 1000);
124 // TGeoHMatrix *im = AliGeomManager::GetMatrix("VZERO/V0C");
131 //____________________________________________________________________________
132 void AliVZERODigitizer::Exec(Option_t* /*option*/)
134 // Creates digits from hits
136 Float_t map[80]; // 48 values on V0C + 32 on V0A
137 // Int_t pmNumber[80];
138 Float_t adc[64]; // 32 PMs on V0C + 32 PMs on V0A
139 Float_t time[80], time_ref[80], time2[64];
140 Float_t adc_gain[80];
141 Float_t adc_pedestal[64],adc_sigma[64];
143 Float_t pmGain_smeared[64];
146 // Smearing of the PM tubes intrinsic characteristics
148 for(Int_t ii=0; ii<64; ii++) {
149 pmGain_smeared[ii] = gRandom->Gaus(fPMGain, fPMGain/5); }
151 // Retrieval of ADC gain values and pedestal information from local CDB
152 // I use only the first 64th values of the calibration array in CDB
153 // as I have no beam burst structure - odd or even beam burst number
155 // Reminder : We have 16 scintillating cells mounted on 8 PMs
156 // on Ring 3 and Ring 4 in V0C - added to produce ADC outputs
159 for(Int_t i=0; i<16; i++) {
160 adc_gain[i] = fCalibData->GetGain(i);
161 cPM[i] = fPhotoCathodeEfficiency * pmGain_smeared[i];
164 for(Int_t j=16; j<48; j=j+2) {
166 adc_gain[j] = fCalibData->GetGain(i);
167 adc_gain[j+1] = fCalibData->GetGain(i);
168 cPM[j] = fPhotoCathodeEfficiency * pmGain_smeared[i];
169 cPM[j+1] = fPhotoCathodeEfficiency * pmGain_smeared[i];
172 for(Int_t i=48; i<80; i++){
173 adc_gain[i] = fCalibData->GetGain(i-16);
174 cPM[i] = fPhotoCathodeEfficiency * pmGain_smeared[i-16];
177 for(Int_t i=0; i<64; i++){
178 adc_pedestal[i] = fCalibData->GetPedestal(i);
179 adc_sigma[i] = fCalibData->GetSigma(i);
182 // for(Int_t i=0; i<64; i++) { printf(" i = %d pedestal = %f sigma = %f \n\n",
183 // i, adc_pedestal[i], adc_sigma[i] );}
185 AliRunLoader* outRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
187 Error("Exec", "Can not get output Run Loader");
191 AliLoader* outLoader = outRunLoader->GetLoader("VZEROLoader");
194 Error("Exec", "Can not get output VZERO Loader");
198 const char* mode = "update";
199 if(outRunLoader->GetEventNumber() == 0) mode = "recreate";
200 outLoader->LoadDigits(mode);
202 if (!outLoader->TreeD()) outLoader->MakeTree("D");
203 outLoader->MakeDigitsContainer();
204 TTree* treeD = outLoader->TreeD();
205 Int_t bufsize = 16000;
206 treeD->Branch("VZERODigit", &fDigits, bufsize);
208 for (Int_t iInput = 0; iInput < fManager->GetNinputs(); iInput++) {
209 AliRunLoader* runLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
210 AliLoader* loader = runLoader->GetLoader("VZEROLoader");
212 Error("Exec", "Can not get VZERO Loader for input %d", iInput);
216 if (!runLoader->GetAliRun()) runLoader->LoadgAlice();
218 AliVZERO* vzero = (AliVZERO*) runLoader->GetAliRun()->GetDetector("VZERO");
220 Error("Exec", "No VZERO detector for input %d", iInput);
225 TTree* treeH = loader->TreeH();
227 Error("Exec", "Cannot get TreeH for input %d", iInput);
231 for(Int_t i=0; i<80; i++) {map[i] = 0; time[i] = 0.0;}
233 TClonesArray* hits = vzero->Hits();
235 // Now makes Digits from hits
237 Int_t nTracks = (Int_t) treeH->GetEntries();
238 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
239 for (Int_t i=0; i<80; i++) {time_ref[i] = 999999.0;}
241 treeH->GetEvent(iTrack);
242 Int_t nHits = hits->GetEntriesFast();
243 for (Int_t iHit = 0; iHit < nHits; iHit++) {
244 AliVZEROhit* hit = (AliVZEROhit *)hits->UncheckedAt(iHit);
245 Int_t nPhot = hit->Nphot();
246 Int_t cell = hit->Cell();
247 map[cell] += Float_t(nPhot);
248 Float_t dt_scintillator = gRandom->Gaus(0,0.7);
249 Float_t t = dt_scintillator + 1e9*hit->Tof();
251 if(t < time_ref[cell]) time_ref[cell] = t;
252 time[cell] = TMath::Min(t,time_ref[cell]);
257 loader->UnloadHits();
261 // Now builds the scintillator cell response (80 cells i.e. 80 responses)
263 for (Int_t i=0; i<80; i++) {
264 Float_t q1 = Float_t ( map[i] )* cPM[i] * kQe;
265 Float_t noise = gRandom->Gaus(10.5,3.22);
266 Float_t pmResponse = q1/kC*TMath::Power(ktheta/kthau,1/(1-ktheta/kthau))
268 if(fCollisionMode >0) adc_gain[i] = adc_gain[i]/70.0; // reduce dynamics in Ion Collision Mode
269 map[i] = pmResponse * adc_gain[i];
270 Float_t MIP = 1.0/fCalibData->GetMIPperADC(GetPMNumber(i));
271 if(fCollisionMode >0) MIP=2.0;
272 // printf("cell = %d, ADC = %d, TDC = %f \n",i,map[i], time[i]*10.0 );
273 if(map[i] > (MIP/2.) )
274 {map[i] = gRandom->Gaus(map[i], (MIP/6.) );}
277 // Now transforms 80 cell responses into 64 photomultiplier responses
278 // Also adds the ADC pedestals taken out of the calibration data base
280 for (Int_t j=0; j<16; j++){
281 adc[j] = map [j] + gRandom->Gaus(adc_pedestal[j], adc_sigma[j]);
284 for (Int_t j=48; j<80; j++){
285 adc[j-16] = map [j] + gRandom->Gaus(adc_pedestal[j-16],adc_sigma[j-16]);
286 time2[j-16]= time[j]; }
288 for (Int_t j=0; j<16; j++){
289 adc[16+j] = map [16+2*j]+ map [16+2*j+1] + gRandom->Gaus(adc_pedestal[16+j], adc_sigma[16+j]);
290 Float_t min_time = TMath::Min(time [16+2*j],time [16+2*j+1]);
291 time2[16+j] = min_time;
292 if(min_time==0.0){time2[16+j]=TMath::Max(time[16+2*j],time[16+2*j+1]);}
296 // Now add digits to the digit Tree
298 for (Int_t i=0; i<64; i++) {
300 // printf(" Event, cell, adc, tof = %d %d %d %d\n",
301 // outRunLoader->GetEventNumber(),i, adc[i], Int_t((time2[i]*10.0) +0.5));
302 // multiply by 10 to have 100 ps per channel :
304 AddDigit(i, adc[i], (time2[i]*10.0) ) ;
308 outLoader->WriteDigits("OVERWRITE");
309 outLoader->UnloadDigits();
313 //____________________________________________________________________________
314 void AliVZERODigitizer::AddDigit(Int_t PMnumber, Float_t adc, Float_t time)
319 TClonesArray &ldigits = *fDigits;
321 if (((Int_t) gRandom->Uniform(2))<1) integrator = kFALSE;
322 else integrator = kTRUE;
324 new(ldigits[fNdigits++]) AliVZEROdigit(PMnumber,adc,time,0,kFALSE,kFALSE,integrator);
327 //____________________________________________________________________________
328 void AliVZERODigitizer::ResetDigit()
334 if (fDigits) fDigits->Delete();
337 //____________________________________________________________________________
338 void AliVZERODigitizer::GetCollisionMode()
340 // Retrieves the collision mode from GRP data
342 // Initialization of the GRP entry
344 Int_t run = AliCDBManager::Instance()->GetRun();
346 // printf("\n ++++++ Run Number retrieved as %d \n",run);
348 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data",run);
349 AliGRPObject* grpData = 0x0;
352 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
355 grpData = new AliGRPObject();
356 grpData->ReadValuesFromMap(m);
359 grpData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
362 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
365 if(!grpData) AliError("No GRP entry found in OCDB!");
367 // Retrieval of collision mode
369 TString beamType = grpData->GetBeamType();
370 if(beamType==AliGRPObject::GetInvalidString()){
371 AliError("GRP/GRP/Data entry: missing value for the beam type !");
372 AliError("\t VZERO cannot retrieve beam type\n");
376 if( (beamType.CompareTo("P-P") ==0) || (beamType.CompareTo("p-p") ==0) ){
379 else if( (beamType.CompareTo("Pb-Pb") ==0) || (beamType.CompareTo("A-A") ==0) ){
383 fBeamEnergy = grpData->GetBeamEnergy();
384 if(fBeamEnergy==AliGRPObject::GetInvalidFloat()) {
385 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
389 // printf("\n ++++++ Beam type and collision mode retrieved as %s %d @ %1.3f GeV ++++++\n\n",beamType.Data(), fCollisionMode, fBeamEnergy);
393 //____________________________________________________________________________
394 AliVZEROCalibData* AliVZERODigitizer::GetCalibData() const
397 AliCDBManager *man = AliCDBManager::Instance();
399 AliCDBEntry *entry=0;
401 entry = man->Get("VZERO/Calib/Data");
404 // AliWarning("Load of calibration data from default storage failed!");
405 // AliWarning("Calibration data will be loaded from local storage ($ALICE_ROOT)");
406 // Int_t runNumber = man->GetRun();
407 // entry = man->GetStorage("local://$ALICE_ROOT/OCDB")
408 // ->Get("VZERO/Calib/Data",runNumber);
412 // Retrieval of data in directory VZERO/Calib/Data:
415 AliVZEROCalibData *calibdata = 0;
417 if (entry) calibdata = (AliVZEROCalibData*) entry->GetObject();
418 if (!calibdata) AliFatal("No calibration data from calibration database !");
424 //____________________________________________________________________________
425 Int_t AliVZERODigitizer::GetPMNumber(Int_t cell) const
429 Int_t pmNumber[80] = { 0, 1, 2, 3, 4, 5, 6, 7,
430 8, 9, 10, 11, 12, 13, 14, 15,
431 16, 16, 17, 17, 18, 18, 19, 19, 20, 20, 21, 21, 22, 22, 23, 23,
432 24, 24, 25, 25, 26, 26, 27, 27, 28, 28, 29, 29, 30, 30, 31, 31,
433 32, 33, 34, 35, 36, 37, 38, 39,
434 40, 41, 42, 43, 44, 45, 46, 47,
435 48, 49, 50, 51, 52, 53, 54, 55,
436 56, 57, 58, 59, 60, 61, 62, 63 };
438 return pmNumber[cell];