]>
Commit | Line | Data |
---|---|---|
451f5018 | 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: AliITSUDigitizer.cxx 52261 2011-10-23 15:46:57Z hristov $ */ | |
17 | /////////////////////////////////////////////////////////////////////////// | |
18 | //Piotr.Skowronski@cern.ch : // | |
19 | //Corrections applied in order to compile (only) // | |
20 | // with new I/O and folder structure // | |
21 | //To be implemented correctly by responsible // | |
22 | // // | |
23 | // Class used to steer // | |
24 | // the digitization for ITS // | |
25 | // // | |
26 | /////////////////////////////////////////////////////////////////////////// | |
27 | ||
28 | #include <stdlib.h> | |
29 | #include <TClonesArray.h> | |
30 | #include <TTree.h> | |
31 | #include <TBranch.h> | |
32 | ||
33 | #include "AliRun.h" | |
34 | #include "AliRunLoader.h" | |
35 | #include "AliLoader.h" | |
36 | #include "AliLog.h" | |
37 | #include "AliDigitizationInput.h" | |
38 | #include "AliITSUDigitizer.h" | |
39 | #include "AliITSUGeomTGeo.h" | |
40 | #include "AliITSUSimulation.h" | |
41 | #include "AliITSUSDigit.h" | |
42 | ||
43 | ClassImp(AliITSUDigitizer) | |
44 | ||
45 | //______________________________________________________________________ | |
46 | AliITSUDigitizer::AliITSUDigitizer() | |
47 | : fITS(0) | |
48 | ,fModActive(0) | |
49 | ,fInit(kFALSE) | |
50 | ,fRoif(-1) | |
51 | ,fRoiifile(0) | |
52 | ,fFlagFirstEv(kTRUE) | |
53 | { | |
54 | // Default constructor. | |
55 | } | |
56 | ||
57 | //______________________________________________________________________ | |
58 | AliITSUDigitizer::AliITSUDigitizer(AliDigitizationInput* digInp) | |
59 | :AliDigitizer(digInp) | |
60 | ,fITS(0) | |
61 | ,fModActive(0) | |
62 | ,fInit(kFALSE) | |
63 | ,fRoif(-1) | |
64 | ,fRoiifile(0) | |
65 | ,fFlagFirstEv(kTRUE) | |
66 | { | |
67 | // Standard constructor. | |
68 | } | |
69 | ||
70 | //______________________________________________________________________ | |
71 | AliITSUDigitizer::~AliITSUDigitizer() | |
72 | { | |
73 | // destructor. | |
74 | fITS = 0; // don't delete fITS. Done else where. | |
75 | delete[] fModActive; | |
76 | } | |
77 | ||
78 | //______________________________________________________________________ | |
79 | Bool_t AliITSUDigitizer::Init() | |
80 | { | |
81 | // Initialization. Set up region of interest, if switched on, and loads ITS and ITSgeom. | |
82 | // | |
83 | if (fInit) return kTRUE; | |
84 | // | |
85 | fInit = kTRUE; // Assume for now init will work. | |
86 | // | |
87 | if(!gAlice) { | |
88 | fITS = 0; | |
89 | fRoiifile = 0; | |
90 | fInit = kFALSE; | |
91 | AliFatal("gAlice not found"); | |
92 | } // end if | |
93 | // | |
94 | fITS = (AliITSU *)(gAlice->GetDetector("ITS")); | |
95 | if(!fITS){ | |
96 | fRoiifile = 0; | |
97 | fInit = kFALSE; | |
98 | AliFatal("ITS not found"); | |
99 | } | |
02d6eccc | 100 | if (!fITS->IsSimInitDone()) fITS->InitSimulation(); |
451f5018 | 101 | int nm = fITS->GetITSGeomTGeo()->GetNModules(); |
102 | fModActive = new Bool_t[nm]; | |
103 | for (Int_t i=nm;i--;) fModActive[i] = kTRUE; | |
104 | ||
105 | return fInit; | |
106 | } | |
107 | ||
108 | //______________________________________________________________________ | |
f9c7eb32 | 109 | void AliITSUDigitizer::Digitize(Option_t* /*opt*/) |
451f5018 | 110 | { |
111 | // Main digitization function. | |
112 | // | |
113 | if (!fInit) AliFatal("Init not successful, aborting."); | |
114 | // | |
115 | Int_t nfiles = GetDigInput()->GetNinputs(); | |
116 | Int_t event = GetDigInput()->GetOutputEventNr(); | |
117 | // | |
118 | TString loadname = Form("%sLoader",fITS->GetName()); | |
119 | // | |
120 | AliITSUGeomTGeo* geom = fITS->GetITSGeomTGeo(); | |
121 | Int_t nModules = geom->GetNModules(); | |
122 | Bool_t lmod; | |
123 | Int_t *fl = new Int_t[nfiles]; | |
124 | fl[0] = fRoiifile; | |
125 | int mask = 1; | |
126 | for (int id=0;id<nfiles;id++) if(id!=fRoiifile) fl[mask++] = id; | |
127 | // | |
128 | TClonesArray * sdig = new TClonesArray( "AliITSUSDigit",1000 ); | |
129 | // | |
130 | AliRunLoader *inRL = 0x0, *outRL = 0x0; | |
131 | AliLoader *ingime = 0x0, *outgime = 0x0; | |
132 | // | |
133 | outRL = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName()); | |
134 | if (!outRL) AliFatal("Can not get Output Run Loader"); | |
135 | // | |
136 | outRL->GetEvent(event); | |
137 | outgime = outRL->GetLoader(loadname); | |
138 | if ( outgime == 0x0) AliFatal("Can not get Output ITS Loader"); | |
139 | // | |
140 | outgime->LoadDigits("update"); | |
141 | if (outgime->TreeD() == 0x0) outgime->MakeTree("D"); | |
142 | // | |
143 | // Digitize | |
144 | fITS->MakeBranchInTreeD(outgime->TreeD()); | |
145 | if(fRoif!=0) AliDebug(1,"Region of Interest digitization selected"); | |
146 | else AliDebug(1,"No Region of Interest selected. Digitizing everything"); | |
147 | // | |
148 | for (int ifiles=0; ifiles<nfiles; ifiles++ ) { | |
149 | inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles])); | |
150 | ingime = inRL->GetLoader(loadname); | |
151 | if (ingime->TreeS() == 0x0) ingime->LoadSDigits(); | |
152 | } | |
153 | // | |
154 | for (int module=0; module<nModules; module++ ) { | |
155 | // | |
156 | if (!fRoif && !fModActive[module]) continue; | |
c92b1537 | 157 | int lr = geom->GetLayer(module); |
158 | AliITSUSimulation *sim = fITS->GetSimulationModel(lr); | |
159 | if (!sim) AliFatal(Form("The simulation model for layer %d is not available",lr)); | |
451f5018 | 160 | // |
161 | // Fill the module with the sum of SDigits | |
c92b1537 | 162 | sim->InitSimulationModule(fITS->GetModule(module), event, fITS->GetSegmentation(lr), fITS->GetResponseParam(lr)); |
451f5018 | 163 | // |
164 | for (int ifiles=0; ifiles<nfiles; ifiles++ ) { | |
165 | // | |
166 | if (!fRoif && !fModActive[module]) continue; | |
167 | inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles])); | |
168 | ingime = inRL->GetLoader(loadname); | |
169 | // | |
170 | TTree *treeS = ingime->TreeS(); | |
171 | fITS->SetTreeAddress(); | |
172 | // | |
173 | if( !treeS ) continue; | |
174 | TBranch *brchSDigits = treeS->GetBranch( fITS->GetName() ); | |
175 | if( brchSDigits ) brchSDigits->SetAddress( &sdig ); | |
176 | else { | |
177 | AliError(Form("branch ITS not found in TreeS, input file %d ", ifiles)); | |
178 | delete [] fl; | |
179 | return; | |
180 | } | |
181 | // | |
182 | sdig->Clear(); | |
183 | mask = GetDigInput()->GetMask(ifiles); | |
184 | brchSDigits->GetEvent( module ); | |
185 | lmod = sim->AddSDigitsToModule(sdig,mask); | |
186 | if(GetRegionOfInterest() && !ifiles) fModActive[module] = lmod; | |
187 | // | |
188 | } | |
189 | // Digitize current module sum(SDigits)->Digits | |
190 | sim->FinishSDigitiseModule(); | |
191 | // | |
192 | outgime->TreeD()->Fill(); // fills all branches - wasted disk space | |
193 | fITS->ResetDigits(); | |
194 | } // end for module | |
195 | // | |
196 | // fITS->WriteFOSignals(); | |
197 | outgime->TreeD()->AutoSave(); | |
198 | outgime->WriteDigits("OVERWRITE"); | |
199 | outgime->UnloadDigits(); | |
200 | for(int ifiles=0; ifiles<nfiles; ifiles++ ) { | |
201 | inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles])); | |
202 | ingime = inRL->GetLoader(loadname); | |
203 | ingime->UnloadSDigits(); | |
204 | } | |
205 | // | |
206 | delete[] fl; | |
207 | sdig->Clear(); | |
208 | delete sdig; | |
209 | for (Int_t i=nModules;i--;) fModActive[i] = kTRUE; | |
210 | // | |
211 | return; | |
212 | } | |
213 | ||
214 | //______________________________________________________________________ | |
215 | void AliITSUDigitizer::SetByRegionOfInterest(TTree *ts) | |
216 | { | |
217 | // Scans through the ITS branch of the SDigits tree, ts, for modules | |
218 | // which have SDigits in them. For these modules, a flag is set to | |
219 | // digitize only these modules. The value of fRoif determines how many | |
220 | // neighboring modules will also be turned on. fRoif=0 will turn on only | |
221 | // those modules with SDigits in them. fRoif=1 will turn on, in addition, | |
222 | // those modules that are +-1 module from the one with the SDigits. And | |
223 | // So on. This last feature is not supported yet. | |
224 | // Inputs: | |
225 | // TTree *ts The tree in which the existing SDigits will define the | |
226 | // region of interest. | |
227 | if (fRoif==0) return; | |
228 | if (ts==0) return; | |
229 | TBranch *brchSDigits = ts->GetBranch(fITS->GetName()); | |
230 | TClonesArray * sdig = new TClonesArray( "AliITSUSDigit",1000 ); | |
231 | // | |
232 | if( brchSDigits ) brchSDigits->SetAddress( &sdig ); | |
233 | else {AliError("Branch ITS not found in TreeS"); return;} | |
234 | // | |
235 | int nm = fITS->GetITSGeomTGeo()->GetNModules(); | |
236 | for (int m=0;m<nm;m++) { | |
237 | fModActive[m] = kFALSE; // Not active by default | |
238 | sdig->Clear(); | |
239 | brchSDigits->GetEvent(m); | |
240 | int ndig = sdig->GetEntries(); | |
241 | for(int i=0;i<ndig;i++) { | |
242 | // activate the necessary modules | |
243 | if ( ((AliITSUSDigit*)sdig->At(m))->GetSumSignal()>0.0 ) { // Must have non zero signal. | |
244 | fModActive[m] = kTRUE; | |
245 | break; | |
246 | } // end if | |
247 | } // end if. end for i. | |
248 | } // end for m | |
249 | AliDebug(1,"Digitization by Region of Interest selected"); | |
250 | sdig->Clear(); | |
251 | delete sdig; | |
252 | return; | |
253 | } |