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 Revision 1.3 2002/03/25 10:48:55 nilsen
19 New ITS SDigit merging with region of interest cut. Update for changes in
20 AliDigitizer. Additional optimization should be done.
22 Revision 1.2 2002/03/15 17:26:40 nilsen
23 New SDigit version of ITS Digitizer.
25 Revision 1.1 2001/11/27 16:27:28 nilsen
26 Adding AliITSDigitizer class to do merging and digitization . Based on the
27 TTask method. AliITSDigitizer class added to the Makefile and ITSLinkDef.h
28 file. The following files required minor changes. AliITS, added functions
29 SetHitsAddressBranch, MakeBranchInTreeD and modified MakeBranchD.
30 AliITSsimulationSDD.cxx needed a Tree independent way of returning back to
31 the original Root Directory in function Compress1D. Now it uses gDirectory.
33 Revision 1.2 2002/03/01 E. Lopez
34 Digitization changed to start from SDigits instead of Hits.
35 The SDigits are reading as TClonesArray of AliITSpListItem
40 #include <TObjArray.h>
46 #include <AliRunDigitizer.h>
48 #include "AliITSDigitizer.h"
49 #include "AliITSpList.h"
50 #include "AliITSmodule.h"
51 #include "AliITSsimulation.h"
52 #include "AliITSDetType.h"
53 #include "AliITSgeom.h"
55 ClassImp(AliITSDigitizer)
57 //______________________________________________________________________
58 AliITSDigitizer::AliITSDigitizer() : AliDigitizer(){
59 // Default constructor. Assign fITS since it is never written out from
62 // Option_t * opt Not used
66 // A blank AliITSDigitizer class.
74 //______________________________________________________________________
75 AliITSDigitizer::AliITSDigitizer(AliRunDigitizer *mngr) : AliDigitizer(mngr){
76 // Standard constructor. Assign fITS since it is never written out from
79 // Option_t * opt Not used
83 // An AliItSDigitizer class.
91 //______________________________________________________________________
92 AliITSDigitizer::~AliITSDigitizer(){
93 // Default destructor.
95 // Option_t * opt Not used
101 fITS = 0; // don't delete fITS. Done else where.
102 if(fActive) delete[] fActive;
104 //______________________________________________________________________
105 Bool_t AliITSDigitizer::Init(){
106 // Initialization. Set up region of interest, if switched on, and
107 // loads ITS and ITSgeom.
115 fInit = kTRUE; // Assume for now init will work.
122 Warning("Init","gAlice not found");
125 fITS = (AliITS *)(gAlice->GetDetector("ITS"));
131 Warning("Init","ITS not found");
133 } else if(fITS->GetITSgeom()){
134 //cout << "fRoif,fRoiifile="<<fRoif<<" "<<fRoiifile<<endl;
135 fActive = new Bool_t[fITS->GetITSgeom()->GetIndexMax()];
141 Warning("Init","ITS geometry not found");
144 // fActive needs to be set to a default all kTRUE value
145 for(Int_t i=0;i<fITS->GetITSgeom()->GetIndexMax();i++) fActive[i] = kTRUE;
146 /* This will not work from Init. ts is aways returned as zero.
148 if(fRoif>=0 && fRoiifile>=0 && fRoiifile<GetManager()->GetNinputs()){
149 ts = GetManager()->GetInputTreeS(fRoiifile);
151 if(!gAlice) ts = gAlice->TreeS();
153 cout <<"The TTree TreeS needed to set by region not found."
154 " No region of interest cut will be applied."<< endl;
158 cout << "calling SetByReionOfInterest ts="<< ts <<endl;
159 SetByRegionOfInterest(ts);
164 //______________________________________________________________________
165 void AliITSDigitizer::Exec(Option_t* opt){
166 // Main digitization function.
168 // Option_t * opt list of sub detector to digitize. =0 all.
174 char name[20] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
176 const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
178 if( !det[0] && !det[1] && !det[2] ) all = "All";
180 AliITSsimulation *sim = 0;
181 AliITSDetType *iDetType = 0;
182 static Bool_t setDef = kTRUE;
185 Error("Exec","Init not succesfull, aborting.");
189 if( setDef ) fITS->SetDefaultSimulation();
191 sprintf(name,"%s",fITS->GetName());
193 Int_t nfiles = GetManager()->GetNinputs();
194 Int_t event = GetManager()->GetOutputEventNr();
195 Int_t size = fITS->GetITSgeom()->GetIndexMax();
196 Int_t module,id,ifiles,mask;
198 Int_t *fl = new Int_t[nfiles];
201 for(id=0;id<nfiles;id++) if(id!=fRoiifile){
202 // just in case fRoiifile!=0.
206 TClonesArray * sdig = new TClonesArray( "AliITSpListItem",1000 );
209 fITS->MakeBranchInTreeD(GetManager()->GetTreeD());
211 for(module=0; module<size; module++ ){
212 if(fActive && fRoif!=0) if(!fActive[module]) continue;
213 id = fITS->GetITSgeom()->GetModuleType(module);
214 if(!all && !det[id]) continue;
215 iDetType = fITS->DetType( id );
216 sim = (AliITSsimulation*)iDetType->GetSimulationModel();
218 Error( "Exec", "The simulation class was not instanciated!" );
222 // Fill the module with the sum of SDigits
223 sim->InitSimulationModule(module, event);
224 //cout << "Module=" << module;
225 for(ifiles=0; ifiles<nfiles; ifiles++ ){
226 if(fActive && fRoif!=0) if(!fActive[module]) continue;
227 //cout <<" fl[ifiles=" << ifiles << "]=" << fl[ifiles];
228 TTree *treeS = GetManager()->GetInputTreeS(fl[ifiles]);
229 if( !(treeS && fITS->GetSDigits()) ) continue;
230 TBranch *brchSDigits = treeS->GetBranch( name );
232 brchSDigits->SetAddress( &sdig );
234 Error( "Exec", "branch ITS not found in TreeS, input file %d ",
237 } // end if brchSDigits
239 mask = GetManager()->GetMask(ifiles);
240 // add summable digits to module
241 brchSDigits->GetEvent( module );
242 lmod = sim->AddSDigitsToModule(sdig,mask);
244 fActive[module] = lmod;
246 //cout << " fActive["<<module<<"]=";
247 //if(fActive[module]) cout << "kTRUE";
248 //else cout << "kFALSE";
250 //cout << " end ifiles loop" << endl;
251 // Digitize current module sum(SDigits)->Digits
252 sim->FinishSDigitiseModule();
254 // fills all branches - wasted disk space
255 GetManager()->GetTreeD()->Fill();
258 //cout << "end modules loop"<<endl;
260 GetManager()->GetTreeD()->AutoSave();
265 for(Int_t i=0;i<fITS->GetITSgeom()->GetIndexMax();i++) fActive[i] = kTRUE;
268 //______________________________________________________________________
269 void AliITSDigitizer::SetByRegionOfInterest(TTree *ts){
270 // Scans through the ITS branch of the SDigits tree, ts, for modules
271 // which have SDigits in them. For these modules, a flag is set to
272 // digitize only these modules. The value of fRoif determines how many
273 // neighboring modules will also be turned on. fRoif=0 will turn on only
274 // those modules with SDigits in them. fRoif=1 will turn on, in addition,
275 // those modules that are +-1 module from the one with the SDigits. And
276 // So on. This last feature is not supported yet.
278 // TTree *ts The tree in which the existing SDigits will define the
279 // region of interest.
288 TBranch *brchSDigits = ts->GetBranch(fITS->GetName());
289 TClonesArray * sdig = new TClonesArray( "AliITSpListItem",1000 );
290 //cout << "Region of Interest ts="<<ts<<" brchSDigits="<<brchSDigits<<" sdig="<<sdig<<endl;
293 brchSDigits->SetAddress( &sdig );
295 Error( "SetByRegionOfInterest","branch ITS not found in TreeS");
297 } // end if brchSDigits
299 nm = fITS->GetITSgeom()->GetIndexMax();
301 //cout << " fActive["<<m<<"]=";
302 fActive[m] = kFALSE; // Not active by default
304 brchSDigits->GetEvent(m);
305 if(sdig->GetLast()>=0) for(i=0;i<sdig->GetLast();i++){
306 // activate the necessary modules
307 if(((AliITSpList*)sdig->At(m))->GetpListItem(i)->GetSignal()>0.0){ // Must have non zero signal.
311 } // end if. end for i.
312 //cout << fActive[m];