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 //_________________________________________________________________________
19 // This is a TTask that constructs SDigits out of Hits
20 // A Summable Digits is the sum of all hits in a cell
21 // A threshold is applied
23 //-- Author: Alla Maevskaia(INR)
24 //////////////////////////////////////////////////////////////////////////////
26 // --- Standard library ---
27 #include <Riostream.h>
30 // --- ROOT system ---
37 #include <TVirtualMC.h>
39 // --- AliRoot header files ---
40 #include "AliConfig.h"
41 #include "AliDetector.h"
43 #include "AliFMDSDigitizer.h"
44 #include "AliFMDdigit.h"
45 #include "AliFMDhit.h"
47 #include "AliLoader.h"
49 #include "AliRunLoader.h"
52 ClassImp(AliFMDSDigitizer)
54 //____________________________________________________________________________
55 AliFMDSDigitizer::AliFMDSDigitizer():TTask("AliFMDSDigitizer","")
64 //____________________________________________________________________________
65 AliFMDSDigitizer::AliFMDSDigitizer(const char* HeaderFile,char *SdigitsFile ):TTask("AliFMDSDigitizer","")
67 fNevents = 0 ; // Number of events to digitize, 0 means all evens in current file
68 // add Task to //root/Tasks folder
69 fRunLoader = AliRunLoader::Open(HeaderFile);//Load event in default folder
70 if (fRunLoader == 0x0)
72 Fatal("AliFMDSDigitizer","Can not open session. Header File is %s ",HeaderFile);
73 return;//never reached
75 AliLoader* gime = fRunLoader->GetLoader("FMDLoader");
78 Fatal("AliFMDSDigitizer","Can not find FMD (loader) in specified event");
79 return;//never reached
81 //add Task to //root/Tasks folder
82 gime->PostSDigitizer(this);
85 //____________________________________________________________________________
86 AliFMDSDigitizer::~AliFMDSDigitizer()
89 AliLoader* gime = fRunLoader->GetLoader("FMDLoader");
90 gime->CleanSDigitizer();
93 //---------------------------------------------------------------------
94 void AliFMDSDigitizer::SetRingsSi1(Int_t ringsSi1)
98 void AliFMDSDigitizer::SetSectorsSi1(Int_t sectorsSi1)
100 fSectorsSi1=sectorsSi1;
102 void AliFMDSDigitizer::SetRingsSi2(Int_t ringsSi2)
106 void AliFMDSDigitizer::SetSectorsSi2(Int_t sectorsSi2)
108 fSectorsSi2=sectorsSi2;
111 //____________________________________________________________________________
112 void AliFMDSDigitizer::Exec(Option_t *option)
114 Int_t NumberOfRings[5]=
115 {fRingsSi1,fRingsSi2,fRingsSi1,fRingsSi2,fRingsSi1};
116 Int_t NumberOfSectors[5]=
117 {fSectorsSi1,fSectorsSi2,fSectorsSi1,fSectorsSi2,fSectorsSi1};
121 Error("Exec","Run Loader loader is NULL - Session not opened");
124 AliLoader* gime = fRunLoader->GetLoader("FMDLoader");
127 Fatal("AliFMDReconstruction","Can not find FMD (loader) in specified event");
128 return;//never reached
131 if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
132 if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
133 if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics("READ");
137 retval = gime->LoadHits("READ");
140 Error("Exec","Error occured while loading hits. Exiting.");
145 // Initialise Hit array
146 fHits = new TClonesArray ("AliFMDhit", 1000);
147 fSDigits = new TClonesArray ("AliFMDdigit", 1000);
149 AliFMD *FMD = (AliFMD *) gAlice->GetDetector ("FMD");
152 fNevents = (Int_t) fRunLoader->TreeE ()->GetEntries ();
154 for (Int_t ievent = 0; ievent < fNevents; ievent++)
156 fRunLoader->GetEvent (ievent);
158 TTree* TH = gime->TreeH();
161 Error("Exec","Can not get TreeH");
164 if (gime->TreeS () == 0) gime->MakeTree("S");
166 //Make branch for digits
167 FMD->MakeBranch ("S");
168 //Now made SDigits from hits, for PHOS it is the same
169 Int_t volume, sector, ring, charge;
171 Float_t de[10][50][150];
172 Int_t ivol, isec, iring;
176 TClonesArray *FMDhits = FMD->Hits ();
178 // Event ------------------------- LOOP
180 for (ivol = 0; ivol < 10; ivol++)
181 for (isec = 0; isec < 50; isec++)
182 for (iring = 0; iring < 150; iring++)
183 de[ivol][isec][iring] = 0;
187 FMDhits = FMD->Hits ();
190 Stat_t ntracks = TH->GetEntries ();
191 for (Int_t track = 0; track < ntracks; track++)
193 gAlice->ResetHits ();
194 nbytes += TH->GetEvent(track);
195 particle = fRunLoader->Stack()->Particle (track);
196 Int_t nhits = FMDhits->GetEntriesFast ();
198 for (hit = 0; hit < nhits; hit++)
200 fmdHit = (AliFMDhit *) FMDhits->UncheckedAt (hit);
202 volume = fmdHit->Volume ();
203 sector = fmdHit->NumberOfSector ();
204 ring = fmdHit->NumberOfRing ();
206 de[volume][sector][ring] += e;
213 Float_t I = 1.664 * 0.04 * 2.33 / 22400; // = 0.69e-6;
214 for (ivol = 1; ivol < 6; ivol++)
216 for (isec = 1; isec <= NumberOfSectors[ivol-1]; isec++)
218 for (iring = 1; iring <= NumberOfRings[ivol-1]; iring++)
223 charge = Int_t (de[ivol][isec][iring] / I);
226 //dinamic diapason from MIP(0.155MeV) to 30MIP(4.65MeV)
228 Float_t channelWidth = (22400 * 30) / 1024;
230 digit[4] = Int_t (digit[3] / channelWidth);
231 FMD->AddSDigit(digit);
237 gime->TreeS()->Reset();
238 gime->TreeS()->Fill();
239 gime->WriteSDigits("OVERWRITE");
245 //__________________________________________________________________
246 void AliFMDSDigitizer::SetSDigitsFile(char * file ){
247 if(!fSDigitsFile.IsNull())
248 cout << "Changing SDigits file from " <<(char *)fSDigitsFile.Data() << " to " << file << endl ;
251 //__________________________________________________________________
252 void AliFMDSDigitizer::Print(Option_t* option)const
254 cout << "------------------- "<< GetName() << " -------------" << endl ;
255 if(fSDigitsFile.IsNull())
256 cout << " Writing SDigitis to file galice.root "<< endl ;
258 cout << " Writing SDigitis to file " << (char*) fSDigitsFile.Data() << endl ;