#include <stdlib.h>
#include <TClonesArray.h>
#include <TTree.h>
#include <TBranch.h>
#include "AliRun.h"
#include "AliRunLoader.h"
#include "AliLoader.h"
#include "AliLog.h"
#include "AliDigitizationInput.h"
#include "AliITSDigitizer.h"
#include "AliITSgeom.h"
#include "AliITSgeomTGeo.h"
#include "AliITSsimulation.h"
ClassImp(AliITSDigitizer)
AliITSDigitizer::AliITSDigitizer() : AliDigitizer(),
fITS(0),
fModActive(0),
fInit(kFALSE),
fRoif(-1),
fRoiifile(0),
fFlagFirstEv(kTRUE){
fModActive = new Bool_t[AliITSgeomTGeo::GetNModules()];
for(Int_t i=0;i<AliITSgeomTGeo::GetNModules();i++) fModActive[i] = kTRUE;
}
AliITSDigitizer::AliITSDigitizer(AliDigitizationInput* digInp) : AliDigitizer(digInp),
fITS(0),
fModActive(0),
fInit(kFALSE),
fRoif(-1),
fRoiifile(0),
fFlagFirstEv(kTRUE){
fModActive = new Bool_t[AliITSgeomTGeo::GetNModules()];
for(Int_t i=0;i<AliITSgeomTGeo::GetNModules();i++) fModActive[i] = kTRUE;
}
AliITSDigitizer::~AliITSDigitizer(){
fITS = 0;
if(fModActive) delete[] fModActive;
}
Bool_t AliITSDigitizer::Init(){
fInit = kTRUE;
if(!gAlice) {
fITS = 0;
fRoiifile = 0;
fInit = kFALSE;
Warning("Init","gAlice not found");
return fInit;
}
fITS = (AliITS *)(gAlice->GetDetector("ITS"));
if(!fITS){
fRoiifile = 0;
fInit = kFALSE;
Warning("Init","ITS not found");
return fInit;
}
if(!fITS->GetITSgeom()){
fRoiifile = 0;
fInit = kFALSE;
Warning("Init","ITS geometry not found");
return fInit;
}
return fInit;
}
void AliITSDigitizer::Digitize(Option_t* opt){
char name[21] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
const char *all;
const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
strstr(opt,"SSD")};
if( !det[0] && !det[1] && !det[2] ) all = "All";
else all = 0;
Int_t nfiles = GetDigInput()->GetNinputs();
Int_t event = GetDigInput()->GetOutputEventNr();
AliITSsimulation *sim = 0;
if(fFlagFirstEv){
fITS->SetDefaults();
fITS->SetDefaultSimulation();
fFlagFirstEv=kFALSE;
}
if(!fInit){
Error("Exec","Init not successful, aborting.");
return;
}
snprintf(name,20,"%s",fITS->GetName());
Int_t size = fITS->GetITSgeom()->GetIndexMax();
Int_t module,id,ifiles,mask;
Bool_t lmod;
Int_t *fl = new Int_t[nfiles];
fl[0] = fRoiifile;
mask = 1;
for(id=0;id<nfiles;id++)
if(id!=fRoiifile)
{
fl[mask] = id;
mask++;
}
TClonesArray * sdig = new TClonesArray( "AliITSpListItem",1000 );
TString loadname(name);
loadname+="Loader";
AliRunLoader *inRL = 0x0, *outRL = 0x0;
AliLoader *ingime = 0x0, *outgime = 0x0;
outRL = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
if ( outRL == 0x0)
{
Error("Exec","Can not get Output Run Loader");
delete [] fl;
return;
}
outRL->GetEvent(event);
outgime = outRL->GetLoader(loadname);
if ( outgime == 0x0)
{
Error("Exec","Can not get Output ITS Loader");
delete [] fl;
return;
}
outgime->LoadDigits("update");
if (outgime->TreeD() == 0x0) outgime->MakeTree("D");
fITS->MakeBranchInTreeD(outgime->TreeD());
if(fRoif!=0) {
AliDebug(1,"Region of Interest digitization selected");
}
else {
AliDebug(1,"No Region of Interest selected. Digitizing everything");
}
for(ifiles=0; ifiles<nfiles; ifiles++ )
{
inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles]));
ingime = inRL->GetLoader(loadname);
if (ingime->TreeS() == 0x0) ingime->LoadSDigits();
}
for(module=0; module<size; module++ )
{
if(fRoif!=0) if(!fModActive[module]) continue;
id = fITS->GetITSgeom()->GetModuleType(module);
if(!all && !det[id]) continue;
sim = (AliITSsimulation*)fITS->GetSimulationModel(id);
if(!sim) {
Error( "Exec", "The simulation class was not instanciated!" );
exit(1);
}
sim->InitSimulationModule(module, event);
for(ifiles=0; ifiles<nfiles; ifiles++ )
{
if(fRoif!=0) if(!fModActive[module]) continue;
inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles]));
ingime = inRL->GetLoader(loadname);
TTree *treeS = ingime->TreeS();
fITS->SetTreeAddress();
if( !treeS ) continue;
TBranch *brchSDigits = treeS->GetBranch( name );
if( brchSDigits )
{
brchSDigits->SetAddress( &sdig );
} else {
Error( "Exec", "branch ITS not found in TreeS, input file %d ",
ifiles );
delete [] fl;
return;
}
sdig->Clear();
mask = GetDigInput()->GetMask(ifiles);
brchSDigits->GetEvent( module );
lmod = sim->AddSDigitsToModule(sdig,mask);
if(GetRegionOfInterest() && (ifiles==0))
{
fModActive[module] = lmod;
}
}
sim->FinishSDigitiseModule();
outgime->TreeD()->Fill();
fITS->ResetDigits();
}
fITS->WriteFOSignals();
outgime->TreeD()->AutoSave();
outgime->WriteDigits("OVERWRITE");
outgime->UnloadDigits();
for(ifiles=0; ifiles<nfiles; ifiles++ )
{
inRL = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(fl[ifiles]));
ingime = inRL->GetLoader(loadname);
ingime->UnloadSDigits();
}
delete[] fl;
sdig->Clear();
delete sdig;
for(Int_t i=0;i<fITS->GetITSgeom()->GetIndexMax();i++) fModActive[i] = kTRUE;
return;
}
void AliITSDigitizer::SetByRegionOfInterest(TTree *ts){
Int_t m,nm,i;
if(fRoif==0) return;
if(ts==0) return;
TBranch *brchSDigits = ts->GetBranch(fITS->GetName());
TClonesArray * sdig = new TClonesArray( "AliITSpListItem",1000 );
if( brchSDigits ) {
brchSDigits->SetAddress( &sdig );
} else {
Error( "SetByRegionOfInterest","branch ITS not found in TreeS");
return;
}
nm = fITS->GetITSgeom()->GetIndexMax();
for(m=0;m<nm;m++){
fModActive[m] = kFALSE;
sdig->Clear();
brchSDigits->GetEvent(m);
if(sdig->GetLast()>=0) for(i=0;i<sdig->GetLast();i++){
if(((AliITSpList*)sdig->At(m))->GetpListItem(i)->GetSignal()>0.0){
fModActive[m] = kTRUE;
break;
}
}
}
AliDebug(1,"Digitization by Region of Interest selected");
sdig->Clear();
delete sdig;
return;
}