1 /**************************************************************************
2 * Copyright(c) 1998-2000, 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 Revision 1.9 2002/10/14 14:57:32 hristov
19 Merging the VirtualMC branch to the main development branch (HEAD)
21 Revision 1.7.6.1 2002/07/24 10:07:52 alibrary
24 Revision 1.7 2001/11/02 15:37:26 hristov
25 Digitizer class created. Code cleaning and bug fixes (J.Chudoba)
27 Revision 1.5 2001/10/23 13:03:35 hristov
28 The access to several data members was changed from public to protected. The digitisation was adapted to the multi-event case (J.Chudoba)
30 Revision 1.4 2001/10/21 18:31:24 hristov
31 Several pointers were set to zero in the default constructors to avoid memory management problems
33 Revision 1.3 2001/05/16 14:57:20 alibrary
34 New files for folders and Stack
36 Revision 1.2 2001/03/14 18:16:08 jbarbosa
37 Corrected bug (more to correct).
38 File "points.dat" is no longer created.
40 Revision 1.1 2001/02/27 22:13:34 jbarbosa
41 Implementing merger class.
44 #include <Riostream.h>
47 #include <TObjArray.h>
49 #include <TDirectory.h>
50 #include <TParticle.h>
52 #include "AliRICHMerger.h"
53 #include "AliRICHChamber.h"
54 #include "AliHitMap.h"
55 #include "AliRICHHitMapA1.h"
57 #include "AliRICHHit.h"
58 #include "AliRICHSDigit.h"
59 #include "AliRICHDigit.h"
60 #include "AliRICHTransientDigit.h"
64 ClassImp(AliRICHMerger)
66 //___________________________________________
67 AliRICHMerger::AliRICHMerger()
69 // Default constructor
84 //------------------------------------------------------------------------
85 AliRICHMerger::~AliRICHMerger()
88 if (fTrH1) delete fTrH1;
89 if (fHitsBgr) delete fHitsBgr;
90 if (fSDigitsBgr) delete fSDigitsBgr;
91 if (fHitMap) delete fHitMap;
92 if (fList) delete fList;
93 if (fBgrFile) delete fBgrFile;
96 //------------------------------------------------------------------------
97 Bool_t AliRICHMerger::Exists(const AliRICHSDigit *padhit)
99 return (fHitMap[fNch]->TestHit(padhit->PadX(),padhit->PadY()));
102 //------------------------------------------------------------------------
103 void AliRICHMerger::Update(AliRICHSDigit *padhit)
105 AliRICHTransientDigit *pdigit =
106 static_cast<AliRICHTransientDigit*>(
107 fHitMap[fNch]->GetHit(padhit->PadX(),padhit->PadY()));
111 Int_t iqpad = Int_t(padhit->QPad()); // charge per pad
112 pdigit->AddSignal(iqpad);
113 pdigit->AddPhysicsSignal(iqpad);
115 // update list of tracks
125 pdigit->UpdateTrackList(track,charge);
128 //------------------------------------------------------------------------
129 void AliRICHMerger::CreateNew(AliRICHSDigit *padhit)
131 fList->AddAtAndExpand(
132 new AliRICHTransientDigit(fNch,fDigits),fCounter);
133 fHitMap[fNch]->SetHit(padhit->PadX(),padhit->PadY(),fCounter);
135 AliRICHTransientDigit* pdigit =
136 static_cast<AliRICHTransientDigit*>(fList->Last());
141 charge = padhit->QPad();
146 pdigit->AddToTrackList(track,charge);
151 //------------------------------------------------------------------------
152 void AliRICHMerger::Init()
156 if (fMerge && !fBgrFile) fBgrFile = InitBgr();
161 //------------------------------------------------------------------------
162 TFile* AliRICHMerger::InitBgr()
164 // Initialise background event
165 TFile *file = new TFile(fFnBgr);
166 // add error checking later
167 printf("\n AliRICHMerger has opened %s file with background event \n", fFnBgr);
168 fHitsBgr = new TClonesArray("AliRICHHit",1000);
169 fSDigitsBgr = new TClonesArray("AliRICHSDigit",1000);
173 //------------------------------------------------------------------------
174 void AliRICHMerger::Digitise(Int_t nev, Int_t flag)
177 // keep galice.root for signal and name differently the file for
178 // background when add! otherwise the track info for signal will be lost !
182 AliRICHChamber* iChamber;
183 AliSegmentation* segmentation;
185 fList = new TObjArray;
187 AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH");
188 fHitMap= new AliHitMap* [kNCH];
193 // Get Hits Tree header from file
194 if(fHitsBgr) fHitsBgr->Clear();
195 if(fSDigitsBgr) fSDigitsBgr->Clear();
196 if(fTrH1) delete fTrH1;
200 sprintf(treeName,"TreeH%d",fEvNrBgr);
201 fTrH1 = (TTree*)gDirectory->Get(treeName);
203 printf("ERROR: cannot find Hits Tree for event:%d\n",fEvNrBgr);
206 // Set branch addresses
209 sprintf(branchname,"%s",pRICH->GetName());
210 if (fTrH1 && fHitsBgr) {
211 branch = fTrH1->GetBranch(branchname);
212 if (branch) branch->SetAddress(&fHitsBgr);
214 if (fTrH1 && fSDigitsBgr) {
215 branch = fTrH1->GetBranch("RICHSDigits");
216 if (branch) branch->SetAddress(&fSDigitsBgr);
220 for (Int_t i =0; i<kNCH; i++) {
221 iChamber= &(pRICH->Chamber(i));
222 segmentation=iChamber->GetSegmentationModel(1);
223 fHitMap[i] = new AliRICHHitMapA1(segmentation, fList);
231 TTree *treeH = gAlice->TreeH();
232 Int_t ntracks =(Int_t) treeH->GetEntries();
233 for (fTrack=0; fTrack<ntracks; fTrack++) {
235 treeH->GetEvent(fTrack);
238 for(AliRICHHit* mHit=(AliRICHHit*)pRICH->FirstHit(-1);
240 mHit=(AliRICHHit*)pRICH->NextHit())
243 fNch = mHit->Chamber()-1; // chamber number
244 Int_t trackIndex = mHit->Track();
246 cerr<<"AliRICHMerger: chamber nr. fNch out of range: "<<fNch<<endl;
247 cerr<<" track: "<<fTrack<<endl;
250 iChamber = &(pRICH->Chamber(fNch));
254 // If flag is set, create digits only for some particles
256 TParticle *current = (TParticle*)gAlice->Particle(trackIndex);
257 if (current->GetPdgCode() >= 50000050)
259 TParticle *motherofcurrent = (TParticle*)gAlice->Particle(current->GetFirstMother());
260 particle = motherofcurrent->GetPdgCode();
264 particle = current->GetPdgCode();
268 //printf("Flag:%d\n",flag);
269 //printf("Track:%d\n",track);
270 //printf("Particle:%d\n",particle);
275 if(TMath::Abs(particle)==211 || TMath::Abs(particle)==111)
279 if(TMath::Abs(particle)==321 || TMath::Abs(particle)==130 || TMath::Abs(particle)==310
280 || TMath::Abs(particle)==311)
283 if (flag == 3 && TMath::Abs(particle)==2212)
286 if (flag == 4 && TMath::Abs(particle)==13)
289 if (flag == 5 && TMath::Abs(particle)==11)
292 if (flag == 6 && TMath::Abs(particle)==2112)
296 //printf ("Particle: %d, Flag: %d, Digitise: %d\n",particle,flag,digitise);
303 // Loop over pad hits
304 for (AliRICHSDigit* mPad=
305 (AliRICHSDigit*)pRICH->FirstPad(mHit,pRICH->SDigits());
307 mPad=(AliRICHSDigit*)pRICH->NextPad(pRICH->SDigits()))
309 Int_t iqpad = mPad->QPad(); // charge per pad
310 fDigits[0]=mPad->PadX();
311 fDigits[1]=mPad->PadY();
314 fDigits[4]=mPad->HitNumber();
316 // build the list of fired pads and update the info
322 } //end loop over clusters
323 }// track type condition
327 // open the file with background
331 ntracks = (Int_t)fTrH1->GetEntries();
335 for (fTrack = 0; fTrack < ntracks; fTrack++) {
337 if (fHitsBgr) fHitsBgr->Clear();
338 if (fSDigitsBgr) fSDigitsBgr->Clear();
340 fTrH1->GetEvent(fTrack);
344 for(Int_t i = 0; i < fHitsBgr->GetEntriesFast(); ++i)
346 mHit = (AliRICHHit*) (*fHitsBgr)[i];
347 fNch = mHit->Chamber()-1; // chamber number
348 iChamber = &(pRICH->Chamber(fNch));
351 // Loop over pad hits
352 for (AliRICHSDigit* mPad =
353 (AliRICHSDigit*)pRICH->FirstPad(mHit,fSDigitsBgr);
355 mPad = (AliRICHSDigit*)pRICH->NextPad(fSDigitsBgr))
357 fDigits[0] = mPad->PadX();
358 fDigits[1] = mPad->PadY();
359 fDigits[2] = mPad->QPad();
361 fDigits[4] = -1; // set hit number to -1 for bgr
363 // build the list of fired pads and update the info
369 } //end loop over clusters
373 TTree *treeK = gAlice->TreeK();
376 if (treeK) file = treeK->GetCurrentFile();
380 Int_t tracks[kMAXTRACKSPERRICHDIGIT];
381 Int_t charges[kMAXTRACKSPERRICHDIGIT];
382 Int_t nentries=fList->GetEntriesFast();
384 for (Int_t nent=0;nent<nentries;nent++) {
385 AliRICHTransientDigit *transDigit=(AliRICHTransientDigit*)fList->At(nent);
386 if (transDigit==0) continue;
387 Int_t ich=transDigit->GetChamber();
388 Int_t q=transDigit->Signal();
389 iChamber=&(pRICH->Chamber(ich));
390 AliRICHResponse * response=iChamber->GetResponseModel();
391 Int_t adcmax= (Int_t) response->MaxAdc();
394 // add white noise and do zero-suppression and signal truncation (new electronics,old electronics gaus 1.2,0.2)
395 //printf("Treshold: %d\n",iChamber->fTresh->GetHitIndex(transDigit->PadX(),transDigit->PadY()));
396 Int_t pedestal = iChamber->fTresh->GetHitIndex(transDigit->PadX(),transDigit->PadY());
398 //printf("Pedestal:%d\n",pedestal);
400 Float_t treshold = (pedestal + 4*2.2);
402 Float_t meanNoise = gRandom->Gaus(2.2, 0.3);
403 Float_t noise = gRandom->Gaus(0, meanNoise);
404 q+=(Int_t)(noise + pedestal);
406 // magic number to be parametrised !!!
413 if ( q >= adcmax) q=adcmax;
414 fDigits[0]=transDigit->PadX();
415 fDigits[1]=transDigit->PadY();
417 fDigits[3]=transDigit->Physics();
418 fDigits[4]=transDigit->Hit();
420 Int_t nptracks = transDigit->GetNTracks();
422 // this was changed to accomodate the real number of tracks
423 if (nptracks > kMAXTRACKSPERRICHDIGIT) {
424 printf("Attention - tracks > 10 %d\n",nptracks);
425 nptracks=kMAXTRACKSPERRICHDIGIT;
428 printf("Attention - tracks > 2 %d \n",nptracks);
429 //printf("cat,ich,ix,iy,q %d %d %d %d %d \n",
430 //icat,ich,digits[0],digits[1],q);
432 for (Int_t tr=0;tr<nptracks;tr++) {
433 tracks[tr]=transDigit->GetTrack(tr);
434 charges[tr]=transDigit->GetCharge(tr);
435 } //end loop over list of tracks for one pad
436 if (nptracks < kMAXTRACKSPERRICHDIGIT ) {
437 for (Int_t t=nptracks; t<kMAXTRACKSPERRICHDIGIT; t++) {
444 //fprintf(points,"%4d, %4d, %4d\n",digits[0],digits[1],digits[2]);
447 pRICH->AddDigits(ich,tracks,charges,fDigits);
449 gAlice->TreeD()->Fill();
452 for(Int_t ii=0;ii<kNCH;++ii) {
459 TClonesArray *richDigits;
460 for (Int_t k=0;k<kNCH;k++) {
461 richDigits = pRICH->DigitsAddress(k);
462 Int_t ndigit=richDigits->GetEntriesFast();
463 printf ("Chamber %d digits %d \n",k,ndigit);
465 pRICH->ResetDigits(); /// ??? should it be here???
466 gAlice->TreeD()->Write(0,TObject::kOverwrite);
468 // gAlice->TreeD()->Reset();
469 // delete fList; // deleted in dtor
470 // gObjectTable->Print();