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 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // ZDC event merging class //
21 ///////////////////////////////////////////////////////////////////////////////
27 #include <TDirectory.h>
31 // --- AliRoot header files
32 #include "AliZDCMerger.h"
34 #include "AliZDCHit.h"
35 #include "AliZDCMergedHit.h"
36 #include "AliZDCDigit.h"
37 #include "AliZDCFragment.h"
39 #include "AliDetector.h"
40 #include "AliHeader.h"
41 #include "AliGenEventHeader.h"
42 #include "AliGenHijingEventHeader.h"
44 ClassImp(AliZDCMerger)
46 //____________________________________________________________________________
47 AliZDCMerger::AliZDCMerger()
49 // Default constructor
70 //____________________________________________________________________________
71 AliZDCMerger::~AliZDCMerger()
74 if (fBgrFile) delete fBgrFile;
75 if (fTrHBgr) delete fTrHBgr;
76 if (fTrSDBgr) delete fTrSDBgr;
77 if (fHitsBgr) delete fHitsBgr;
78 if (fSpecnFile) delete fSpecnFile;
79 if (fSpecpFile) delete fSpecpFile;
82 //____________________________________________________________________________
83 void AliZDCMerger::InitMerging()
85 // Hits tree, impact parameter, num. of spectators n & p
86 // in background (full Hijing) event
89 Background(b, nspecn, nspecp);
91 // Production of nuclear fragments -> num. of FREE spectators n & p
92 Fragmentation(b, nspecn, nspecp, fFreeSpn, fFreeSpp);
94 // Extract from spectators distribution the signal events:
95 // NFreeSpectatorN spectator n & NFreeSpectatorP spectator p
99 //____________________________________________________________________________
100 void AliZDCMerger::Background(Float_t &fImpPar, Int_t &fSpecn, Int_t &fSpecp)
103 // --- Open the background file
104 if (fMerge && !fBgrFile) fBgrFile = OpenBgrFile();
106 // --- Read from the TreeE impact parameter (b),
107 // # of spectators n and p (fSpecn, fSpecp)
110 // // Get AliRun object from file or create it if not on file
111 // gAlice = (AliRun*)fBgrFile->Get("gAlice");
113 // gAlice = (AliRun*)fBgrFile->Get("gAlice");
114 // if (gAlice) printf("AliRun object found on file\n");
116 // printf("\n create new gAlice object");
117 // gAlice = new AliRun("gAlice","Alice test program");
121 AliHeader *header = gAlice->GetHeader();
122 AliGenEventHeader* mcHeader = header->GenEventHeader();
123 fImpPar = ((AliGenHijingEventHeader*) mcHeader)->ImpactParameter();
124 fSpecn = ((AliGenHijingEventHeader*) mcHeader)->Spectatorsn();
125 fSpecp = ((AliGenHijingEventHeader*) mcHeader)->Spectatorsp();
126 printf("\n HIJING simulation - b = %f fm, Nspecn = %d, Nspecp = %d\n",fImpPar,fSpecn,fSpecp);
129 //____________________________________________________________________________
130 TFile* AliZDCMerger::OpenBgrFile()
132 // Initialise background event
133 TFile *file = new TFile(fFnBgr,"UPDATE");
134 printf("\n AliZDCMerger --- Background event -> %s file opened \n", fFnBgr);
135 fHitsBgr = new TClonesArray("AliZDCHit",1000);
136 fMHits = new TClonesArray("AliZDCMergedHit",1000);
140 //____________________________________________________________________________
141 void AliZDCMerger::Fragmentation(Float_t fImpPar, Int_t fSpecn, Int_t fSpecp,
142 Int_t &fFreeSpn, Int_t &fFreeSpp)
144 Int_t j, zz[100], nn[100], nAlpha, Ztot, Ntot;
145 AliZDCFragment *frag = new AliZDCFragment(fImpPar);
146 for(j=0; j<=99; j++){
151 // Fragments generation
152 frag->GenerateIMF(zz, nAlpha);
157 frag->AttachNeutrons(zz, nn, Ztot, Ntot);
158 fFreeSpn = fSpecn-Ztot-2*nAlpha;
159 fFreeSpp = fSpecp-Ntot-2*nAlpha;
160 printf("\n Fragmentation -> FreeSpn = %d, FreeSpp = %d\n",fFreeSpn,fFreeSpp);
163 //____________________________________________________________________________
164 void AliZDCMerger::Mixing()
167 // printf("\n AliZDCMerger->Mixing\n");
169 // ### Background event Hits ###########################################
173 AliZDC *ZDC = (AliZDC *)gAlice->GetModule("ZDC");
174 // if(ZDC) printf("\n Ho trovato lo ZDC!\n");
177 if(fTrHBgr) delete fTrHBgr;
180 if(fTrSDBgr) delete fTrSDBgr;
183 // fNEvBgr = 0; // Let's suppose to have 1 full Hijing event per file
186 char treeBgrName[20];
187 sprintf(treeBgrName,"TreeH%d",fNEvBgr);
188 fTrHBgr = (TTree*)gDirectory->Get(treeBgrName);
190 printf("\n ERROR -> Can't find TreeH%d in background file\n",fNEvBgr);
197 sprintf(branchname,"%s",ZDC->GetName());
198 if(fTrHBgr && fHitsBgr){
199 // printf("\n fTrHBgr!=0 && fHitsBgr!=0\n");
200 branch = fTrHBgr->GetBranch(branchname);
201 if(branch) branch->SetAddress(&fHitsBgr);
204 Int_t ntracks = (Int_t) fTrHBgr->GetEntries();
205 // printf("\n --- ntracks = %d\n\n", ntracks);
207 Int_t itrack, nhits, ihit, j, sector[2];
209 AliZDCMergedHit *MHit;
211 TClonesArray &sdigits = *fMHits; // SDigits TCArray
215 for(itrack=0; itrack<ntracks; itrack++){
216 // printf(" itrack = %d", itrack);
217 // gAlice->ResetHits();
218 fTrHBgr->GetEvent(itrack);
220 // for(AliZDCHit* zdcHit=(AliZDCHit*)ZDC->FirstHit(-1);
222 // zdcHit = (AliZDCHit*)ZDC->NextHit()){
224 nhits = fHitsBgr->GetEntries();
225 // nhits = ZDCdp->GetNhits();
226 // printf(" nhits = %d \n", nhits);
227 for(ihit=0; ihit<nhits; ihit++){
228 zdcHit = (AliZDCHit*) fHitsBgr->UncheckedAt(ihit);
230 for(j=0; j<2; j++) sector[j] = zdcHit->GetVolume(j);
231 MHits[0] = zdcHit->GetPrimKinEn();
232 MHits[1] = zdcHit->GetXImpact();
233 MHits[2] = zdcHit->GetYImpact();
234 MHits[3] = zdcHit->GetSFlag();
235 MHits[4] = zdcHit->GetLightPMQ();
236 MHits[5] = zdcHit->GetLightPMC();
237 MHits[6] = zdcHit->GetEnergy();
238 MHit = new AliZDCMergedHit(sector, MHits);
240 new((*fMHits)[fNMhits]) AliZDCMergedHit(*MHit);
241 new (sdigits[fNMhits]) AliZDCMergedHit(*MHit);
247 // printf("\n fNMhits (after bckg) = %d, \n",fNMhits);
249 // AliZDCMergedHit *pippo7 = (AliZDCMergedHit*) fMHits->At(7);
253 // ### Signal event Hits ###########################################
259 // printf("\n fNMhits (after signal) = %d \n",fNMhits);
263 //____________________________________________________________________________
264 void AliZDCMerger::ExtractSignal(Int_t SpecType)
267 // printf("\n Entering in Extract Signal\n");
270 if(SpecType == 1){ // --- Signal for spectator neutrons
271 fFnSpecn = gSystem->ExpandPathName("$ALICE/$ALICE_LEVEL/ZDC/ZNsignalntu.root");
272 fSpecnFile = TFile::Open(fFnSpecn,"R");
274 // printf("\n --- ExtractSignal x n: file %s opened\n", fFnSpecn);
275 NumEvents = fFreeSpn;
277 else if(SpecType == 2){ // --- Signal for spectator protons
278 fFnSpecp = gSystem->ExpandPathName("$ALICE/$ALICE_LEVEL/ZDC/ZPsignalntu.root");
279 fSpecpFile = TFile::Open(fFnSpecp,"R");
281 // printf("\n --- ExtractSignal x p: file %s opened\n", fFnSpecp);
282 NumEvents = fFreeSpp;
284 // printf("\n # of free spectator = %d\n", NumEvents);
285 // printf("\n fNMhits (before adding signal) = %d\n",fNMhits);
287 TNtuple *ZDCSignal = (TNtuple*) gDirectory->Get("ZDCSignal");
288 Int_t nentries = (Int_t) ZDCSignal->GetEntries();
289 // printf("\n # entries = %d\n", nentries);
292 AliZDCMergedHit *MHit;
293 Float_t *entry, HitsSpec[7];
294 Int_t iev, rnd, Volume[2];
295 for(iev=0; iev<NumEvents; iev++){
296 rnd = (Int_t) (1000*gRandom->Rndm());
297 // printf("\n rnd = %d\n", rnd);
298 for(i=0; i<nentries; i++){
299 ZDCSignal->GetEvent(i);
300 entry = ZDCSignal->GetArgs();
302 // printf("\n entry[0] = %f\n", entry[0]);
303 for(k=0; k<2; k++) Volume[k] = (Int_t) entry[k+1];
305 HitsSpec[j] = entry[j+3];
307 MHit = new AliZDCMergedHit(Volume, HitsSpec);
308 new((*fMHits)[fNMhits++]) AliZDCMergedHit(*MHit);
311 else if(entry[0] > rnd) break;
316 // printf("\n fNMhits (after n signal) = %d\n",fNMhits);
319 else if(SpecType == 2){
320 // printf("\n fNMhits (after p signal) = %d\n",fNMhits);
326 //____________________________________________________________________________
327 void AliZDCMerger::Digitize(Int_t fNMhits, TClonesArray *fMHits)
330 // printf("\n AliZDCMerger->Digitize()");
332 AliZDC *ZDC = (AliZDC *)gAlice->GetModule("ZDC");
333 // if(ZDC) printf("\n Ho trovato lo ZDC!\n");
335 Int_t lightQ, lightC, sector[2], digit;
336 Int_t PMCZN = 0, PMCZP = 0, PMQZN[4], PMQZP[4], PMZEM1 = 0, PMZEM2 = 0;
343 AliZDCMergedHit *MHit;
345 // printf("\n fNMHits = %d\n", fNMhits);
347 for(imhit=0; imhit<fNMhits; imhit++){
349 MHit = (AliZDCMergedHit*) fMHits->UncheckedAt(imhit);
350 sector[0] = MHit->GetSector(0);
351 sector[1] = MHit->GetSector(1);
352 if((sector[1]!=1) && (sector[1]!=2) && (sector[1]!=3) && (sector[1]!=4)){
353 // printf("\n *** ERROR!!! sector[0] = %d, sector[1] = %d\n",
354 // sector[0], sector[1]);
357 lightQ = Int_t(MHit->GetLightPMQ());
358 lightC = Int_t(MHit->GetLightPMC());
359 // printf("\n imhit = %d -> DET. = %d, quad = %d,PMQ = %d, PMC = %d",
360 // imhit,sector[0], sector[1],lightQ, lightC);
362 if(sector[0] == 1){ //ZN
363 PMCZN = PMCZN + lightC;
364 PMQZN[sector[1]-1] = PMQZN[sector[1]-1] + lightQ;
366 else if(sector[0] == 2){ //ZP
367 PMCZP = PMCZP + lightC;
368 PMQZP[sector[1]-1] = PMQZP[sector[1]-1] + lightQ;
370 else if(sector[0] == 3){ //ZEM
371 if(sector[1] ==1) PMZEM1 = PMZEM1 + lightC;
372 else PMZEM2 = PMZEM2 + lightQ;
376 // ### Digits creation ###############################################
377 // Create digits for ZN
379 sector[0] = 1; // Detector = ZN
380 sector[1] = 0; // Common PM ADC
381 digit = Phe2ADCch(1, 0, PMCZN);
382 // printf("\n\n ZN ### PMCZN = %d ADCZN = %d",PMCZN, digit);
383 PedValue = AddPedestal();
385 // printf(" PedValue = %d",PedValue);
386 ZDC->AddDigit(sector, digit);
389 sector[1] = j+1; // Towers PM ADCs
390 digit = Phe2ADCch(1, j+1, PMQZN[j]);
391 // printf("\n PMQZN[%d] = %d phe ADCZN[%d] = %d ADCch",j,PMQZN[j],j,digit);
392 PedValue = AddPedestal();
394 // printf(" PedValue = %d",PedValue);
395 ZDC->AddDigit(sector, digit);
399 // Create digits for ZP
400 sector[0] = 2; // Detector = ZP
401 sector[1] = 0; // Common PM ADC
402 digit = Phe2ADCch(2, 0, PMCZP);
403 // printf("\n ZP --- PMCZP = %d phe ADCZP = %d ADCch",PMCZP,digit);
404 PedValue = AddPedestal();
406 printf(" PedValue = %d",PedValue);
407 ZDC->AddDigit(sector, digit);
409 sector[1] = j+1; // Towers PM ADCs
410 digit = Phe2ADCch(2, j+1, PMQZP[j]);
411 // printf("\n PMQZP[%d] = %d phe ADCZP[%d] = %d ADCch",j,PMQZP[j],j,digit);
412 PedValue = AddPedestal();
414 // printf(" PedValue = %d",PedValue);
415 ZDC->AddDigit(sector, digit);
419 // Create digits for ZEM
421 sector[1] = 1; // Detector = ZEM1
422 digit = Phe2ADCch(3, 1, PMZEM1);
423 // printf("\n ZEM *** PMZEM1 = %d phe ADCZEM1 = %d ADCch",PMZEM1,digit);
424 PedValue = AddPedestal();
426 // printf(" PedValue = %d\n",PedValue);
427 ZDC->AddDigit(sector, digit);
428 sector[1] = 2; // Detector = ZEM2
429 digit = Phe2ADCch(3, 2, PMZEM2);
430 // printf("\n ZEM *** PMZEM2 = %d phe ADCZEM2 = %d ADCch",PMZEM2,digit);
431 PedValue = AddPedestal();
433 // printf(" PedValue = %d\n",PedValue);
434 ZDC->AddDigit(sector, digit);
437 //_____________________________________________________________________________
438 Int_t AliZDCMerger::Phe2ADCch(Int_t Det, Int_t Quad, Int_t Light)
440 // Evaluation of the ADC channel corresponding to the light yield Light
442 if(gAlice->GetDebug() > 0){
443 // printf("\n Phe2ADCch -> Detector = %d, Quadrant = %d, Light = %d\n", Det, Quad, Light);
447 // Parameters for conversion of light yield in ADC channels
448 Float_t fPMGain[3][5]; // PM gain
449 Float_t fADCRes; // ADC conversion factor
454 fPMGain[i][j] = 100000.;
457 fADCRes = 0.00000064; // ADC Resolution: 250 fC/ADCch
459 ADCch = (Int_t) (Light*fPMGain[Det-1][Quad]*fADCRes);
464 //_____________________________________________________________________________
465 Int_t AliZDCMerger::AddPedestal()
467 // --- Pedestal value -> extracted from a gaussian distribution
468 // obtained from the beam test on the ZEM prototype (Aug. 2000)
471 Float_t PedMean = 50.;
472 Float_t PedWidth = 10.;
474 PedValue = (Int_t) gRandom->Gaus(PedMean,PedWidth);