X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=PHOS%2FAliPHOSDigitizer.cxx;h=f19f3bfbd7272f9e25b41da3bb84ebb69746fe6e;hp=0ac1e94b6e1093041d959799730f26da8f136cb1;hb=f9223476c1758e9109c8a971ded5a47ee32739e7;hpb=349db5e753204c5a23bed574647a89eb383b02db diff --git a/PHOS/AliPHOSDigitizer.cxx b/PHOS/AliPHOSDigitizer.cxx index 0ac1e94b6e1..f19f3bfbd72 100644 --- a/PHOS/AliPHOSDigitizer.cxx +++ b/PHOS/AliPHOSDigitizer.cxx @@ -15,15 +15,80 @@ /* $Id$ */ +/* History of cvs commits: + * + * $Log: AliPHOSDigitizer.cxx,v $ + * Revision 1.104 2007/12/18 09:08:18 hristov + * Splitting of the QA maker into simulation and reconstruction dependent parts (Yves) + * + * Revision 1.103 2007/11/07 11:25:06 schutz + * Comment out the QA checking before starting digitization + * + * Revision 1.102 2007/10/19 18:04:29 schutz + * The standalone QA data maker is called from AliSimulation and AliReconstruction outside the event loop; i.e. re-reading the data. The QA data making in the event loop has been commented out. + * + * Revision 1.101 2007/10/14 21:08:10 schutz + * Introduced the checking of QA results from previous step before entering the event loop + * + * Revision 1.100 2007/10/10 09:05:10 schutz + * Changing name QualAss to QA + * + * Revision 1.99 2007/09/30 17:08:20 schutz + * Introducing the notion of QA data acquisition cycle (needed by online) + * + * Revision 1.98 2007/09/26 14:22:17 cvetan + * Important changes to the reconstructor classes. Complete elimination of the run-loaders, which are now steered only from AliReconstruction. Removal of the corresponding Reconstruct() and FillESD() methods. + * + * Revision 1.97 2007/08/07 14:12:03 kharlov + * Quality assurance added (Yves Schutz) + * + * Revision 1.96 2007/04/28 10:43:36 policheh + * Dead channels simulation: digit energy sets to 0. + * + * Revision 1.95 2007/04/10 07:20:52 kharlov + * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1 + * + * Revision 1.94 2007/02/01 10:34:47 hristov + * Removing warnings on Solaris x86 + * + * Revision 1.93 2006/10/17 13:17:01 kharlov + * Replace AliInfo by AliDebug + * + * Revision 1.92 2006/08/28 10:01:56 kharlov + * Effective C++ warnings fixed (Timur Pocheptsov) + * + * Revision 1.91 2006/04/29 20:25:30 hristov + * Decalibration is implemented (Yu.Kharlov) + * + * Revision 1.90 2006/04/22 10:30:17 hristov + * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov) + * + * Revision 1.89 2006/04/11 15:22:59 hristov + * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla) + * + * Revision 1.88 2006/03/13 14:05:43 kharlov + * Calibration objects for EMC and CPV + * + * Revision 1.87 2005/08/24 15:33:49 kharlov + * Calibration data for raw digits + * + * Revision 1.86 2005/07/12 20:07:35 hristov + * Changes needed to run simulation and reconstrruction in the same AliRoot session + * + * Revision 1.85 2005/05/28 14:19:04 schutz + * Compilation warnings fixed by T.P. + * + */ + //_________________________________________________________________________ //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute) ////////////////////////////////////////////////////////////////////////////// -// This TTask performs digitization of Summable digits (in the PHOS case it is just +// This class performs digitization of Summable digits (in the PHOS case it is just // the sum of contributions from all primary particles into a given cell). // In addition it performs mixing of summable digits from different events. -// The name of the TTask is also the title of the branch that will contain +// The name of the class is also the title of the branch that will contain // the created SDigits -// The title of the TTAsk is the name of the file that contains the hits from +// The title of the class is the name of the file that contains the hits from // which the SDigits are created // // For each event two branches are created in TreeD: @@ -35,7 +100,7 @@ // // Use case: // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ; -// root[1] d->ExecuteTask() +// root[1] d->Digitize() // Warning in : object already instantiated // //Digitizes SDigitis in all events found in file galice.root // @@ -46,7 +111,7 @@ // // Reads another set of sdigits from galice2.root // root[3] d1->MixWith("galice3.root") // // Reads another set of sdigits from galice3.root -// root[4] d->ExecuteTask("deb timing") +// root[4] d->Digitize("deb timing") // // Reads SDigits from files galice1.root, galice2.root .... // // mixes them and stores produced Digits in file galice1.root // // deb - prints number of produced digits @@ -55,714 +120,860 @@ // // --- ROOT system --- -#include "TFile.h" #include "TTree.h" #include "TSystem.h" -#include "TROOT.h" -#include "TFolder.h" -#include "TObjString.h" -#include "TGeometry.h" #include "TBenchmark.h" +#include "TRandom.h" +#include "TMath.h" // --- Standard library --- -#include // --- AliRoot header files --- - -#include "AliRun.h" -#include "AliHeader.h" -#include "AliRunDigitizer.h" +#include +#include "AliLog.h" +#include "AliDigitizationInput.h" #include "AliPHOSDigit.h" -#include "AliPHOS.h" -#include "AliPHOSGetter.h" #include "AliPHOSDigitizer.h" -#include "AliPHOSSDigitizer.h" #include "AliPHOSGeometry.h" #include "AliPHOSTick.h" +#include "AliPHOSSimParam.h" +#include "AliPHOSCalibData.h" +#include "AliRunLoader.h" +#include "AliPHOSLoader.h" +#include "AliPHOSPulseGenerator.h" ClassImp(AliPHOSDigitizer) //____________________________________________________________________________ - AliPHOSDigitizer::AliPHOSDigitizer() +AliPHOSDigitizer::AliPHOSDigitizer() : + AliDigitizer("",""), + fDefaultInit(kTRUE), + fDigitsInRun(0), + fInit(kFALSE), + fInput(0), + fInputFileNames(0x0), + fEventNames(0x0), + fEmcCrystals(0), + fEventFolderName(""), + fFirstEvent(0), + fLastEvent(0), + fcdb(0x0), + fEventCounter(0), + fPulse(0), + fADCValuesLG(0), + fADCValuesHG(0) { // ctor - - fPinNoise = 0.004 ; - fEMCDigitThreshold = 0.012 ; - fCPVNoise = 0.01; - fCPVDigitThreshold = 0.09 ; - fTimeResolution = 0.5e-9 ; - fTimeSignalLength = 1.0e-9 ; - fDigitsInRun = 0 ; - fADCchanelEmc = 0.0015; // width of one ADC channel in GeV - fADCpedestalEmc = 0.005 ; // - fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC - - fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais' - fADCpedestalCpv = 0.012 ; // - fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC - - fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude - fManager = 0 ; // We work in the standalong mode + InitParameters() ; + fDigInput = 0 ; // We work in the standalong mode } //____________________________________________________________________________ -AliPHOSDigitizer::AliPHOSDigitizer(const char *headerFile,const char * name) +AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName, + TString eventFolderName): + AliDigitizer("PHOSDigitizer", alirunFileName), + fDefaultInit(kFALSE), + fDigitsInRun(0), + fInit(kFALSE), + fInput(0), + fInputFileNames(0x0), + fEventNames(0x0), + fEmcCrystals(0), + fEventFolderName(eventFolderName), + fFirstEvent(0), + fLastEvent(0), + fcdb(0x0), + fEventCounter(0), + fPulse(0), + fADCValuesLG(0), + fADCValuesHG(0) { // ctor - SetName(name) ; - SetTitle(headerFile) ; - fManager = 0 ; // We work in the standalong mode + InitParameters() ; Init() ; - + fDefaultInit = kFALSE ; + fDigInput = 0 ; // We work in the standalone mode + fcdb = new AliPHOSCalibData(-1); } //____________________________________________________________________________ -AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * ard):AliDigitizer(ard) +AliPHOSDigitizer::AliPHOSDigitizer(AliDigitizationInput * rd) : + AliDigitizer(rd,"PHOSDigitizer"), + fDefaultInit(kFALSE), + fDigitsInRun(0), + fInit(kFALSE), + fInput(0), + fInputFileNames(0x0), + fEventNames(0x0), + fEmcCrystals(0), + fEventFolderName(fDigInput->GetInputFolderName(0)), + fFirstEvent(0), + fLastEvent(0), + fcdb (0x0), + fEventCounter(0), + fPulse(0), + fADCValuesLG(0), + fADCValuesHG(0) + { - // ctor - SetName(""); //Will call init in the digitizing - SetTitle("aliroot") ; + // ctor Init() is called by RunDigitizer + fDigInput = rd ; + SetTitle(static_cast(fDigInput->GetInputStream(0))->GetFileName(0)); + InitParameters() ; + fDefaultInit = kFALSE ; + fcdb = new AliPHOSCalibData(-1); } //____________________________________________________________________________ AliPHOSDigitizer::~AliPHOSDigitizer() { - // dtor + // dtor + delete [] fInputFileNames ; + delete [] fEventNames ; + delete fPulse; + delete [] fADCValuesLG; + delete [] fADCValuesHG; + + if(fcdb){ delete fcdb ; fcdb=0;} } //____________________________________________________________________________ -void AliPHOSDigitizer::Digitize(const Int_t event) +void AliPHOSDigitizer::Digitize(Int_t event) { // Makes the digitization of the collected summable digits. // It first creates the array of all PHOS modules - // filled with noise (different for EMC, CPV and PPSD) and + // filled with noise (different for EMC, and CPV) and // then adds contributions from SDigits. // This design avoids scanning over the list of digits to add // contribution to new SDigits only. - AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; - TClonesArray * digits = gime->Digits(GetName()) ; - + Bool_t toMakeNoise = kTRUE ; //Do not create noisy digits if merge with real data + + //First stream + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); + + Int_t readEvent = event ; + if (fDigInput) + readEvent = static_cast(fDigInput->GetInputStream(0))->GetCurrentEventNumber() ; + AliDebug(1,Form("Adding event %d from input stream 0 %s %s", + readEvent, GetTitle(), fEventFolderName.Data())) ; + rl->GetEvent(readEvent) ; + phosLoader->CleanSDigits() ; + phosLoader->LoadSDigits("READ") ; + + //Prepare Output + TClonesArray * digits = phosLoader->Digits() ; + if( !digits ) { + phosLoader->MakeDigitsArray() ; + digits = phosLoader->Digits() ; + } + digits->Clear() ; - const AliPHOSGeometry *geom = gime->PHOSGeometry() ; - + // + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; //Making digits with noise, first EMC - Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ(); + //Check which PHOS modules are present + Bool_t isPresent[5] ; + TString volpath ; + Int_t nmod=0 ; + for(Int_t i=0; i<5; i++){ + volpath = "/ALIC_1/PHOS_"; + volpath += i+1; + if (gGeoManager->CheckPath(volpath.Data())) { + isPresent[i]=1 ; + nmod++ ; + } + else{ + isPresent[i]=0 ; + } + } + + Int_t nEMC = nmod*geom->GetNPhi()*geom->GetNZ(); Int_t nCPV ; Int_t absID ; - TString name = geom->GetName() ; - nCPV = nEMC + geom->GetNumberOfCPVPadsZ()*geom->GetNumberOfCPVPadsPhi()* - geom->GetNModules() ; + //check if CPV exists + Bool_t isCPVpresent=0 ; + for(Int_t i=1; i<=5 && !isCPVpresent; i++){ + volpath = "/ALIC_1/PHOS_"; + volpath += i; + volpath += "/PCPV_1"; + if (gGeoManager->CheckPath(volpath.Data())) + isCPVpresent=1 ; + } + + if(isCPVpresent){ + nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * nmod ; + } + else{ + nCPV = nEMC ; + } digits->Expand(nCPV) ; - // get first the sdigitizer from the tasks list (must have same name as the digitizer) - const AliPHOSSDigitizer * sDigitizer = gime->SDigitizer(GetName()); - if ( !sDigitizer) { - cerr << "ERROR: AliPHOSDigitizer::Digitize -> SDigitizer with name " << GetName() << " not found " << endl ; - abort() ; - } - - // loop through the sdigits posted to the White Board and add them to the noise - TCollection * folderslist = gime->SDigitsFolder()->GetListOfFolders() ; - TIter next(folderslist) ; - TFolder * folder = 0 ; - TClonesArray * sdigits = 0 ; - Int_t input = 0 ; - TObjArray * sdigArray = new TObjArray(2) ; - while ( (folder = (TFolder*)next()) ) - if ( (sdigits = (TClonesArray*)folder->FindObject(GetName()) ) ) { - TString fileName(folder->GetName()) ; - fileName.ReplaceAll("_","/") ; - cout << "INFO: AliPHOSDigitizer::Digitize -> Adding SDigits " - << GetName() << " from " << fileName << endl ; - sdigArray->AddAt(sdigits, input) ; - input++ ; + //take all the inputs to add together and load the SDigits + TObjArray * sdigArray = new TObjArray(fInput) ; + sdigArray->AddAt(phosLoader->SDigits(), 0) ; + + for(Int_t i = 1 ; i < fInput ; i++){ + TString tempo(fEventNames[i]) ; + tempo += i ; + AliRunLoader* rl2 = AliRunLoader::GetRunLoader(tempo) ; + if(!rl2){ + rl2 = AliRunLoader::Open(fInputFileNames[i], tempo) ; + if (!rl2) { + Fatal("AliPHOSDigitizer", "Could not find the Run Loader for %s - %s",fInputFileNames[i].Data(), tempo.Data()) ; + return ; + } + rl2->LoadHeader(); + } + AliPHOSLoader * phosLoader2 = static_cast(rl2->GetLoader("PHOSLoader")); + + if(fDigInput){ + readEvent = static_cast(fDigInput->GetInputStream(i))->GetCurrentEventNumber() ; + } + TClonesArray * digs ; + if(AliPHOSSimParam::GetInstance()->IsStreamDigits(i)){ //This is Digits Stream + AliInfo(Form("Adding event %d from input stream %d %s %s", + readEvent, i, fInputFileNames[i].Data(), tempo.Data())) ; + rl2->GetEvent(readEvent) ; + phosLoader2->CleanDigits() ; + phosLoader2->LoadDigits("READ") ; + digs = phosLoader2->Digits() ; + toMakeNoise=0 ; //Do not add noise, it is already in stream } + else{ + AliInfo(Form("Adding event %d (SDigits) from input stream %d %s %s", + readEvent, i, fInputFileNames[i].Data(), tempo.Data())) ; + rl2->GetEvent(readEvent) ; + phosLoader2->CleanSDigits() ; + phosLoader2->LoadSDigits("READ") ; + digs = phosLoader2->SDigits() ; + } + sdigArray->AddAt(digs, i) ; + } - //Find the first crystall with signal + //Find the first crystal with signal Int_t nextSig = 200000 ; - Int_t i; - for(i=0; iAt(i) ; + TClonesArray * sdigits ; + for(Int_t i = 0 ; i < fInput ; i++){ + sdigits = static_cast(sdigArray->At(i)) ; if ( !sdigits->GetEntriesFast() ) continue ; - Int_t curNext = ((AliPHOSDigit *)sdigits->At(0))->GetId() ; - if(curNext < nextSig) - nextSig = curNext ; + Int_t curNext = static_cast(sdigits->At(0))->GetId() ; + if(curNext < nextSig) + nextSig = curNext ; } - - TArrayI index(input) ; + + TArrayI index(fInput) ; index.Reset() ; //Set all indexes to zero - + AliPHOSDigit * digit ; AliPHOSDigit * curSDigit ; - - TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ; - + +// TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ; + //Put Noise contribution - for(absID = 1; absID <= nEMC; absID++){ - Float_t noise = gRandom->Gaus(0., fPinNoise) ; - new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ; - //look if we have to add signal? - digit = (AliPHOSDigit *) digits->At(absID-1) ; - - if(absID==nextSig){ + Double_t apdNoise = 0. ; + if(toMakeNoise) + apdNoise = AliPHOSSimParam::GetInstance()->GetAPDNoise() ; + + Int_t emcpermod=geom->GetNPhi()*geom->GetNZ(); + Int_t idigit= 0; + for(Int_t imod=0; imod<5; imod++){ + if(!isPresent[imod]) + continue ; + Int_t firstAbsId=imod*emcpermod+1 ; + Int_t lastAbsId =(imod+1)*emcpermod ; + for(absID = firstAbsId ; absID <= lastAbsId ; absID++){ + Float_t noise = gRandom->Gaus(0.,apdNoise) ; + new((*digits)[idigit]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ; + //look if we have to add signal? + digit = static_cast(digits->At(idigit)) ; + idigit++ ; + + if(absID==nextSig){ //Add SDigits from all inputs - ticks->Clear() ; - Int_t contrib = 0 ; - Float_t a = digit->GetAmp() ; - Float_t b = TMath::Abs( a /fTimeSignalLength) ; - //Mark the beginnign of the signal - new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b); - //Mark the end of the ignal - new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b); - - //loop over inputs - for(i=0; iAt(i))->GetEntriesFast() > index[i] ) - curSDigit = (AliPHOSDigit*)((TClonesArray *)sdigArray->At(i))->At(index[i]) ; - else - curSDigit = 0 ; - //May be several digits will contribute from the same input - while(curSDigit && curSDigit->GetId() == absID){ - //Shift primary to separate primaries belonging different inputs - Int_t primaryoffset ; - if(fManager) - primaryoffset = fManager->GetMask(i) ; - else - primaryoffset = 10000000*i ; - curSDigit->ShiftPrimary(primaryoffset) ; - - a = curSDigit->GetAmp() ; - b = a /fTimeSignalLength ; - new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b); - new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b); - - *digit = *digit + *curSDigit ; //add energies - - index[i]++ ; - if(((TClonesArray *)sdigArray->At(i))->GetEntriesFast() > index[i] ) - curSDigit = (AliPHOSDigit*)((TClonesArray *)sdigArray->At(i))->At(index[i]) ; +// ticks->Clear() ; +// Int_t contrib = 0 ; + +//New Timing model is necessary +// Float_t a = digit->GetEnergy() ; +// Float_t b = TMath::Abs( a / fTimeSignalLength) ; +// //Mark the beginning of the signal +// new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b); +// //Mark the end of the signal +// new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b); + +// Calculate time as time of the largest digit + Float_t time = digit->GetTime() ; + Float_t eTime= digit->GetEnergy() ; + + //loop over inputs + for(Int_t i = 0 ; i < fInput ; i++){ + if( static_cast(sdigArray->At(i))->GetEntriesFast() > index[i] ){ + curSDigit = static_cast(static_cast(sdigArray->At(i))->At(index[i])) ; + if(AliPHOSSimParam::GetInstance()->IsStreamDigits(i)){ //This is Digits Stream + curSDigit->SetEnergy(Calibrate(curSDigit->GetEnergy(),curSDigit->GetId())) ; + curSDigit->SetTime(CalibrateT(curSDigit->GetTime(),curSDigit->GetId())) ; + } + } else curSDigit = 0 ; - } - } - - //calculate and set time - Float_t time = FrontEdgeTime(ticks) ; - digit->SetTime(time) ; - - //Find next signal module - nextSig = 200000 ; - for(i=0; iAt(i)) ; - Int_t curNext = nextSig ; - if(sdigits->GetEntriesFast() > index[i] ){ - curNext = ((AliPHOSDigit *) sdigits->At(index[i]))->GetId() ; + //May be several digits will contribute from the same input + while(curSDigit && curSDigit->GetId() == absID){ + //Shift primary to separate primaries belonging different inputs + Int_t primaryoffset ; + if(fDigInput) + primaryoffset = fDigInput->GetMask(i) ; + else + primaryoffset = 10000000*i ; + curSDigit->ShiftPrimary(primaryoffset) ; - } - if(curNext < nextSig) nextSig = curNext ; +//New Timing model is necessary +// a = curSDigit->GetEnergy() ; +// b = a /fTimeSignalLength ; +// new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b); +// new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b); + if(curSDigit->GetEnergy()>eTime){ + eTime=curSDigit->GetEnergy() ; + time=curSDigit->GetTime() ; + } + *digit += *curSDigit ; //add energies + + index[i]++ ; + if( static_cast(sdigArray->At(i))->GetEntriesFast() > index[i] ) + curSDigit = static_cast(static_cast(sdigArray->At(i))->At(index[i])) ; + else + curSDigit = 0 ; + } + } + + //calculate and set time +//New Timing model is necessary +// Float_t time = FrontEdgeTime(ticks) ; + digit->SetTime(time) ; + + //Find next signal module + nextSig = 200000 ; + for(Int_t i = 0 ; i < fInput ; i++){ + sdigits = static_cast(sdigArray->At(i)) ; + Int_t curNext = nextSig ; + if(sdigits->GetEntriesFast() > index[i] ){ + curNext = static_cast(sdigits->At(index[i]))->GetId() ; + } + if(curNext < nextSig) nextSig = curNext ; + } } } } - - ticks->Delete() ; - delete ticks ; - //Now CPV digits (different noise and no timing) - for(absID = nEMC+1; absID <= nCPV; absID++){ - Float_t noise = gRandom->Gaus(0., fCPVNoise) ; - new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ; - //look if we have to add signal? - if(absID==nextSig){ - digit = (AliPHOSDigit *) digits->At(absID-1) ; - //Add SDigits from all inputs - for(i=0; iAt(i))->GetEntriesFast() > index[i] ) - curSDigit = (AliPHOSDigit*)((TClonesArray *)sdigArray->At(i))->At(index[i]) ; - else - curSDigit = 0 ; - - //May be several digits will contribute from the same input - while(curSDigit && curSDigit->GetId() == absID){ - //Shift primary to separate primaries belonging different inputs - Int_t primaryoffset ; - if(fManager) - primaryoffset = fManager->GetMask(i) ; - else - primaryoffset = 10000000*i ; - curSDigit->ShiftPrimary(primaryoffset) ; - - //add energies - *digit = *digit + *curSDigit ; - index[i]++ ; - if(((TClonesArray *)sdigArray->At(i))->GetEntriesFast() > index[i] ) - curSDigit = (AliPHOSDigit*)((TClonesArray *)sdigArray->At(i))->At(index[i]) ; - else - curSDigit = 0 ; - } - } + //Apply non-linearity + if(AliPHOSSimParam::GetInstance()->IsCellNonlinearityOn()){ //Apply non-lineairyt on cell level + const Double_t aNL = AliPHOSSimParam::GetInstance()->GetCellNonLineairyA() ; + const Double_t bNL = AliPHOSSimParam::GetInstance()->GetCellNonLineairyB() ; + const Double_t cNL = AliPHOSSimParam::GetInstance()->GetCellNonLineairyC() ; + for(Int_t i = 0 ; i < nEMC ; i++){ + digit = static_cast( digits->At(i) ) ; + Double_t e= digit->GetEnergy() ; + // version(1) digit->SetEnergy(e*(1+a*TMath::Exp(-e/b))) ; + digit->SetEnergy(e*cNL*(1.+aNL*TMath::Exp(-e*e/2./bNL/bNL))) ; //Better agreement with data... + } + } - //Find next signal module - nextSig = 200000 ; - for(i=0; iAt(i) ; - Int_t curNext = nextSig ; - if(sdigits->GetEntriesFast() > index[i] ) - curNext = ((AliPHOSDigit *) sdigits->At(index[i]))->GetId() ; - if(curNext < nextSig) nextSig = curNext ; - } - - } + + //distretize energy if necessary + if(AliPHOSSimParam::GetInstance()->IsEDigitizationOn()){ + Float_t adcW=AliPHOSSimParam::GetInstance()->GetADCchannelW() ; + for(Int_t i = 0 ; i < nEMC ; i++){ + digit = static_cast( digits->At(i) ) ; + digit->SetEnergy(adcW*ceil(digit->GetEnergy()/adcW)) ; + } + } + + //Apply decalibration if necessary + for(Int_t i = 0 ; i < nEMC ; i++){ + digit = static_cast( digits->At(i) ) ; + Decalibrate(digit) ; } - delete sdigArray ; //We should not delete its contents +// ticks->Delete() ; +// delete ticks ; + //Now CPV digits (different noise and no timing) + Int_t cpvpermod = geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() ; + Int_t nEMCtotal=emcpermod*5 ; + Float_t cpvNoise = AliPHOSSimParam::GetInstance()->GetCPVNoise() ; + if(isCPVpresent){ //CPV is present in current geometry + for(Int_t imod=0; imod<5; imod++){ //module is present in current geometry + if(!isPresent[imod]) + continue ; + Int_t firstAbsId=nEMCtotal+imod*cpvpermod+1 ; + Int_t lastAbsId =nEMCtotal+(imod+1)*cpvpermod ; + for(absID = firstAbsId; absID <= lastAbsId; absID++){ + Float_t noise = gRandom->Gaus(0., cpvNoise) ; + new((*digits)[idigit]) AliPHOSDigit( -1,absID,noise, TimeOfNoise() ) ; + idigit++ ; + //look if we have to add signal? + if(absID==nextSig){ + digit = static_cast(digits->At(idigit-1)) ; + //Add SDigits from all inputs + for(Int_t i = 0 ; i < fInput ; i++){ + if( static_cast(sdigArray->At(i))->GetEntriesFast() > index[i] ) + curSDigit = static_cast( static_cast(sdigArray->At(i))->At(index[i])) ; + else + curSDigit = 0 ; + + //May be several digits will contribute from the same input + while(curSDigit && curSDigit->GetId() == absID){ + //Shift primary to separate primaries belonging different inputs + Int_t primaryoffset ; + if(fDigInput) + primaryoffset = fDigInput->GetMask(i) ; + else + primaryoffset = 10000000*i ; + curSDigit->ShiftPrimary(primaryoffset) ; + + //add energies + *digit += *curSDigit ; + index[i]++ ; + if( static_cast(sdigArray->At(i))->GetEntriesFast() > index[i] ) + curSDigit = static_cast( static_cast(sdigArray->At(i))->At(index[i]) ) ; + else + curSDigit = 0 ; + } + } + + //Find next signal module + nextSig = 200000 ; + for(Int_t i = 0 ; i < fInput ; i++){ + sdigits = static_cast(sdigArray->At(i)) ; + Int_t curNext = nextSig ; + if(sdigits->GetEntriesFast() > index[i] ) + curNext = static_cast( sdigits->At(index[i]) )->GetId() ; + if(curNext < nextSig) nextSig = curNext ; + } + + } + } + } + } + + delete sdigArray ; //We should not delete its contents + Int_t relId[4]; + + //set amplitudes in bad channels to zero + + for(Int_t i = 0 ; i GetEntriesFast(); i++){ + digit = static_cast( digits->At(i) ) ; + geom->AbsToRelNumbering(digit->GetId(),relId); + if(relId[1] == 0) // Emc + if(fcdb->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.); + } + //remove digits below thresholds - for(i = 0; i < nEMC ; i++){ - digit = (AliPHOSDigit*) digits->At(i) ; - if(sDigitizer->Calibrate( digit->GetAmp() ) < fEMCDigitThreshold) + Float_t emcThreshold = AliPHOSSimParam::GetInstance()->GetEmcDigitsThreshold() ; + for(Int_t i = 0 ; i < nEMC ; i++){ + digit = static_cast( digits->At(i) ) ; + + if(digit->GetEnergy() < emcThreshold){ digits->RemoveAt(i) ; - else - digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ; - } + continue ; + } + geom->AbsToRelNumbering(digit->GetId(),relId); - for(i = nEMC; i < nCPV ; i++) - if(sDigitizer->Calibrate(((AliPHOSDigit*)digits->At(i))->GetAmp()) < fCPVDigitThreshold) +// digit->SetEnergy(TMath::Ceil(digit->GetEnergy())-0.9999) ; + + Float_t tres = TimeResolution(digit->GetEnergy()) ; + digit->SetTime(gRandom->Gaus(digit->GetTime(), tres) ) ; + + fPulse->Reset(); + fPulse->SetAmplitude(digit->GetEnergy()/ + fcdb->GetADCchannelEmc(relId[0],relId[3],relId[2])); + fPulse->SetTZero(digit->GetTimeR()); + fPulse->MakeSamples(); + fPulse->GetSamples(fADCValuesHG, fADCValuesLG) ; + Int_t nSamples = fPulse->GetRawFormatTimeBins(); + digit->SetALTROSamplesHG(nSamples,fADCValuesHG); + digit->SetALTROSamplesLG(nSamples,fADCValuesLG); + } + + Float_t cpvDigitThreshold = AliPHOSSimParam::GetInstance()->GetCpvDigitsThreshold() ; + for(Int_t i = nEMC; i < nCPV ; i++){ + if( static_cast(digits->At(i))->GetEnergy() < cpvDigitThreshold ) digits->RemoveAt(i) ; + } digits->Compress() ; - Int_t ndigits = digits->GetEntriesFast() ; - digits->Expand(ndigits) ; - + //Set indexes in list of digits and make true digitization of the energy - for (i = 0 ; i < ndigits ; i++) { - AliPHOSDigit * digit = (AliPHOSDigit *) digits->At(i) ; + for (Int_t i = 0 ; i < ndigits ; i++) { + digit = static_cast( digits->At(i) ) ; digit->SetIndexInList(i) ; - Float_t energy = sDigitizer->Calibrate(digit->GetAmp()) ; - digit->SetAmp(DigitizeEnergy(energy,digit->GetId()) ) ; + if(digit->GetId() > fEmcCrystals){ //digitize CPV only + digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ; + } } } //____________________________________________________________________________ -Int_t AliPHOSDigitizer::DigitizeEnergy(Float_t energy, Int_t absId) +Float_t AliPHOSDigitizer::Calibrate(Float_t amp,Int_t absId){ + //Apply calibration + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absolute ID + Int_t relId[4]; + geom->AbsToRelNumbering(absId,relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + if(relId[1]==0){ //This Is EMC + Float_t calibration = fcdb->GetADCchannelEmc(module,column,row); + return amp*calibration ; + } + return 0 ; +} +//____________________________________________________________________________ +void AliPHOSDigitizer::Decalibrate(AliPHOSDigit *digit) { - Int_t chanel ; - if(absId <= fEmcCrystals){ //digitize as EMC - chanel = (Int_t) TMath::Ceil((energy - fADCpedestalEmc)/fADCchanelEmc) ; - if(chanel > fNADCemc ) chanel = fNADCemc ; + // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB + + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absolute ID + Int_t relId[4]; + geom->AbsToRelNumbering(digit->GetId(),relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + if(relId[1]==0){ //This Is EMC + Float_t decalib = fcdb->GetADCchannelEmcDecalib(module,column,row); // O(1) + Float_t calibration = fcdb->GetADCchannelEmc(module,column,row)*decalib; + Float_t energy = digit->GetEnergy()/calibration; + digit->SetEnergy(energy); //Now digit measures E in ADC counts + Float_t time = digit->GetTime() ; + time-=fcdb->GetTimeShiftEmc(module,column,row); + digit->SetTime(time) ; } - else{ //Digitize as CPV - chanel = (Int_t) TMath::Ceil((energy - fADCpedestalCpv)/fADCchanelCpv) ; - if(chanel > fNADCcpv ) chanel = fNADCcpv ; +} +//____________________________________________________________________________ +Float_t AliPHOSDigitizer::CalibrateT(Float_t time,Int_t absId){ + //Apply time calibration + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absolute ID + Int_t relId[4]; + geom->AbsToRelNumbering(absId,relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + time += fcdb->GetTimeShiftEmc(module,column,row); + return time ; +} +//____________________________________________________________________________ +Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId) +{ + // Returns digitized value of the CPV charge in a pad absId + + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + //Determine rel.position of the cell absId + Int_t relId[4]; + geom->AbsToRelNumbering(absId,relId); + Int_t module=relId[0]; + Int_t row =relId[2]; + Int_t column=relId[3]; + + Int_t channel = 0; + + if(absId > fEmcCrystals){ //digitize CPV only + + //reading calibration data for cell absId. + Float_t adcPedestalCpv = fcdb->GetADCpedestalCpv(module,column,row); + Float_t adcChanelCpv = fcdb->GetADCchannelCpv( module,column,row); + + channel = (Int_t) TMath::Ceil((charge - adcPedestalCpv)/adcChanelCpv) ; + Int_t nMax = AliPHOSSimParam::GetInstance()->GetNADCcpv() ; + if(channel > nMax ) channel = nMax ; } - return chanel ; + return channel ; } + //____________________________________________________________________________ -void AliPHOSDigitizer::Exec(Option_t *option) +void AliPHOSDigitizer::Digitize(Option_t *option) { - // Managing method + // Steering method to process digitization for events + // in the range from fFirstEvent to fLastEvent. + // This range is optionally set by SetEventRange(). + // if fLastEvent=-1, then process events until the end. + // by default fLastEvent = fFirstEvent (process only one event) + + if (!fInit) { // to prevent overwrite existing file + AliError(Form("Give a version name different from %s", + fEventFolderName.Data() )) ; + return ; + } - if(strcmp(GetName(), "") == 0 ) - Init() ; - if (strstr(option,"print")) { - Print(""); + Print(); return ; } - if(strstr(option,"tim")) - gBenchmark->Start("PHOSDigitizer"); - - AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; + if(strstr(option,"tim")) + gBenchmark->Start("PHOSDigitizer"); - Int_t nevents ; + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); - TTree * treeD ; + if (fLastEvent == -1) + fLastEvent = rl->GetNumberOfEvents() - 1 ; + else if (fDigInput) + fLastEvent = fFirstEvent ; + + Int_t nEvents = fLastEvent - fFirstEvent + 1; - if(fManager){ - treeD = fManager->GetTreeD() ; - nevents = 1 ; // Will process only one event - } - else { - gAlice->GetEvent(0) ; - nevents = (Int_t) gAlice->TreeE()->GetEntries() ; - treeD=gAlice->TreeD() ; - } - + Int_t ievent ; - //Check, if this branch already exits - if (treeD) { - TObjArray * lob = (TObjArray*)treeD->GetListOfBranches() ; - TIter next(lob) ; - TBranch * branch = 0 ; - Bool_t phosfound = kFALSE, digitizerfound = kFALSE ; - - while ( (branch = (TBranch*)next()) && (!phosfound || !digitizerfound) ) { - if ( (strcmp(branch->GetName(), "PHOS")==0) && - (strcmp(branch->GetTitle(), GetName())==0) ) - phosfound = kTRUE ; - - else if ( (strcmp(branch->GetName(), "AliPHOSDigitizer")==0) && - (strcmp(branch->GetTitle(), GetName())==0) ) - digitizerfound = kTRUE ; - } - - if ( phosfound ) { - cerr << "WARNING: AliPHOSDigitizer -> Digits branch with name " << GetName() - << " already exits" << endl ; - return ; - } - if ( digitizerfound ) { - cerr << "WARNING: AliPHOSDigitizer -> Digitizer branch with name " << GetName() - << " already exits" << endl ; - return ; - } - } + for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) { + fEventCounter++ ; - Int_t ievent ; - - for(ievent = 0; ievent < nevents; ievent++){ - - if(fManager){ - Int_t input ; - for(input = 0 ; input < fManager->GetNinputs(); input ++){ - TTree * treeS = fManager->GetInputTreeS(input) ; - if(!treeS){ - cerr << "AliPHOSDigitizer -> No Input " << endl ; - return ; - } - gime->ReadTreeS(treeS,input) ; - } - } - else - gime->Event(ievent,"S") ; - Digitize(ievent) ; //Add prepared SDigits to digits and add the noise - - WriteDigits(ievent) ; - + + WriteDigits() ; + if(strstr(option,"deb")) PrintDigits(option); //increment the total number of Digits per run - fDigitsInRun += gime->Digits()->GetEntriesFast() ; - } - + fDigitsInRun += phosLoader->Digits()->GetEntriesFast() ; + } + if(strstr(option,"tim")){ gBenchmark->Stop("PHOSDigitizer"); - cout << "AliPHOSDigitizer:" << endl ; - cout << " took " << gBenchmark->GetCpuTime("PHOSDigitizer") << " seconds for Digitizing " - << gBenchmark->GetCpuTime("PHOSDigitizer")/nevents << " seconds per event " << endl ; - cout << endl ; - } - + TString message ; + message = " took %f seconds for Digitizing %f seconds per event\n" ; + AliInfo(Form( message.Data(), + gBenchmark->GetCpuTime("PHOSDigitizer"), + gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents )); + } } - //____________________________________________________________________________ -Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) -{ // - ticks->Sort() ; //Sort in accordance with times of ticks - TIter it(ticks) ; - AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ; - Float_t time = ctick->CrossingTime(fTimeThreshold) ; - - AliPHOSTick * t ; - while((t=(AliPHOSTick*) it.Next())){ - if(t->GetTime() < time) //This tick starts before crossing - *ctick+=*t ; - else - return time ; - - time = ctick->CrossingTime(fTimeThreshold) ; - } - return time ; +Float_t AliPHOSDigitizer::TimeResolution(Float_t e){ + //calculate TOF resolution using beam-test resutls + Float_t a=AliPHOSSimParam::GetInstance()->GetTOFa() ; + Float_t b=AliPHOSSimParam::GetInstance()->GetTOFb() ; + return TMath::Sqrt(a*a+b*b/e) ; } + +////____________________________________________________________________________ +//Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const +//{ +// // Returns the shortest time among all time ticks +// +// ticks->Sort() ; //Sort in accordance with times of ticks +// TIter it(ticks) ; +// AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ; +// Float_t time = ctick->CrossingTime(fTimeThreshold) ; +// +// AliPHOSTick * t ; +// while((t=(AliPHOSTick*) it.Next())){ +// if(t->GetTime() < time) //This tick starts before crossing +// *ctick+=*t ; +// else +// return time ; +// +// time = ctick->CrossingTime(fTimeThreshold) ; +// } +// return time ; +//} + //____________________________________________________________________________ Bool_t AliPHOSDigitizer::Init() { - fPinNoise = 0.004 ; - fEMCDigitThreshold = 0.012 ; - fCPVNoise = 0.01; - fCPVDigitThreshold = 0.09 ; - fTimeResolution = 0.5e-9 ; - fTimeSignalLength = 1.0e-9 ; - fDigitsInRun = 0 ; - fADCchanelEmc = 0.0015; // width of one ADC channel in GeV - fADCpedestalEmc = 0.005 ; // - fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC - - fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais' - fADCpedestalCpv = 0.012 ; // - fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC - - fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude - - if(fManager) - SetTitle("aliroot") ; - else if (strcmp(GetTitle(),"")==0) - SetTitle("galice.root") ; - - if( strcmp(GetName(), "") == 0 ) - SetName("Default") ; - - AliPHOSGetter * gime = AliPHOSGetter::GetInstance(GetTitle(), GetName()) ; - if ( gime == 0 ) { - cerr << "ERROR: AliPHOSDigitizer::Init -> Could not obtain the Getter object !" << endl ; - return kFALSE; - } - - const AliPHOSGeometry * geom = gime->PHOSGeometry() ; - fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ; - - // Post Digits to the white board - gime->PostDigits(GetName() ) ; - - // Post Digitizer to the white board - gime->PostDigitizer(this) ; - - //Mark that we will use current header file - if(!fManager){ - gime->PostSDigits(GetName(),GetTitle()) ; - gime->PostSDigitizer(GetName(),GetTitle()) ; + // Makes all memory allocations + fInit = kTRUE ; + + AliPHOSGeometry *geom; + if (!(geom = AliPHOSGeometry::GetInstance())) + geom = AliPHOSGeometry::GetInstance("IHEP",""); +// const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ; + + fFirstEvent = 0 ; + fLastEvent = fFirstEvent ; + if (fDigInput) + fInput = fDigInput->GetNinputs() ; + else + fInput = 1 ; + + fInputFileNames = new TString[fInput] ; + fEventNames = new TString[fInput] ; + fInputFileNames[0] = GetTitle() ; + fEventNames[0] = fEventFolderName.Data() ; + Int_t index ; + for (index = 1 ; index < fInput ; index++) { + fInputFileNames[index] = static_cast(fDigInput->GetInputStream(index))->GetFileName(0); + TString tempo = fDigInput->GetInputFolderName(index) ; + fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fDigInput } - return kTRUE ; -} - -//__________________________________________________________________ -void AliPHOSDigitizer::MixWith(const char* headerFile) -{ - // Allows to produce digits by superimposing background and signal event. - // It is assumed, that headers file with SIGNAL events is opened in - // the constructor. - // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed - // Thus we avoid writing (changing) huge and expensive - // backgound files: all output will be writen into SIGNAL, i.e. - // opened in constructor file. - // - // One can open as many files to mix with as one needs. - // However only Sdigits with the same name (i.e. constructed with the same SDigitizer) - // can be mixed. - - if( strcmp(GetName(), "") == 0 ) - Init() ; - if(fManager){ - cout << "Can not use this method under AliRunDigitizer " << endl ; - return ; - } - - // check if the specified SDigits do not already exist on the White Board: - // //Folders/RunMC/Event/Data/PHOS/SDigits/headerFile/sdigitsname - - TString path = "Folders/RunMC/Event/Data/PHOS/SDigits" ; - path += headerFile ; - path += "/" ; - path += GetName() ; - if ( gROOT->FindObjectAny(path.Data()) ) { - cerr << "WARNING: AliPHOSDigitizer::MixWith -> Entry already exists, do not add" << endl ; - return; + //to prevent cleaning of this object while GetEvent is called + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + if(!rl){ + AliRunLoader::Open(GetTitle(), fEventFolderName) ; } - - AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; - gime->PostSDigits(GetName(),headerFile) ; - - // check if the requested file is already open or exist and if SDigits Branch exist - TFile * file = (TFile*)gROOT->FindObject(headerFile); - if ( !file ) { - file = new TFile(headerFile, "READ") ; - if (!file) { - cerr << "ERROR: AliPHOSDigitizer::MixWith -> File " << headerFile << " does not exist!" << endl ; - return ; - } - } - + return fInit ; } -//__________________________________________________________________ -void AliPHOSDigitizer::SetSplitFile(const TString splitFileName) const +//____________________________________________________________________________ +void AliPHOSDigitizer::InitParameters() { - // Diverts the Digits in a file separate from the hits file - - // I guess it is not going to work if we do merging - if (fManager) { - cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> Not yet available in case of merging activated " << endl ; - return ; - } + // Set initial parameters Digitizer - TDirectory * cwd = gDirectory ; - TFile * splitFile = gAlice->InitTreeFile("D",splitFileName.Data()); - splitFile->cd() ; - if ( !splitFile->Get("gAlice") ) - gAlice->Write(); - - TTree *treeE = gAlice->TreeE(); - if (!treeE) { - cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> No TreeE found "<Get("TreeE") ) { - AliHeader *header = new AliHeader(); - treeE->SetBranchAddress("Header", &header); - treeE->SetBranchStatus("*",1); - TTree *treeENew = treeE->CloneTree(); - treeENew->Write(); - } - - // copy AliceGeom - if ( !splitFile->Get("AliceGeom") ) { - TGeometry *AliceGeom = static_cast(cwd->Get("AliceGeom")); - if (!AliceGeom) { - cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> AliceGeom was not found in the input file "<Write(); - } - - cwd->cd() ; - gAlice->MakeTree("D",splitFile); - cout << "INFO: AliPHOSDigitizer::SetSPlitMode -> Digits will be stored in " << splitFileName.Data() << endl ; + fDigitsInRun = 0 ; + SetEventRange(0,-1) ; + fPulse = new AliPHOSPulseGenerator(); + fADCValuesLG = new Int_t[fPulse->GetRawFormatTimeBins()]; + fADCValuesHG = new Int_t[fPulse->GetRawFormatTimeBins()]; + } //__________________________________________________________________ -void AliPHOSDigitizer::Print(Option_t* option)const { +void AliPHOSDigitizer::Print(const Option_t *)const +{ // Print Digitizer's parameters - if( strcmp(GetName(), "") != 0 ){ - - cout << "------------------- "<< GetName() << " -------------" << endl ; - cout << "Digitizing sDigits from file(s): " <FindObjectAny("Folders/RunMC/Event/Data/PHOS/SDigits"))->GetListOfFolders() ; - TIter next(folderslist) ; - TFolder * folder = 0 ; + Int_t nStreams ; + if (fDigInput) + nStreams = GetNInputStreams() ; + else + nStreams = fInput ; - while ( (folder = (TFolder*)next()) ) { - if ( folder->FindObject(GetName()) ) - cout << "Adding SDigits " << GetName() << " from " << folder->GetName() << endl ; + Int_t index = 0 ; + for (index = 0 ; index < nStreams ; index++) { + TString tempo(fEventNames[index]) ; + tempo += index ; + printf ("Adding SDigits from %s \n", fInputFileNames[index].Data()) ; } - cout << endl ; - cout << "Writing digits to " << GetTitle() << endl ; - cout << endl ; - cout << "With following parameters: " << endl ; - cout << " Electronics noise in EMC (fPinNoise) = " << fPinNoise << endl ; - cout << " Threshold in EMC (fEMCDigitThreshold) = " << fEMCDigitThreshold << endl ; ; - cout << " Noise in CPV (fCPVNoise) = " << fCPVNoise << endl ; - cout << " Threshold in CPV (fCPVDigitThreshold) = " << fCPVDigitThreshold << endl ; - cout << "---------------------------------------------------" << endl ; + // printf("\nWith following parameters:\n") ; + // printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ; + // printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ; + // printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ; + // printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ; + printf(" ---------------------------------------------------\n") ; } else - cout << "AliPHOSDigitizer not initialized " << endl ; + AliInfo(Form("AliPHOSDigitizer not initialized" )) ; } //__________________________________________________________________ -void AliPHOSDigitizer::PrintDigits(Option_t * option){ + void AliPHOSDigitizer::PrintDigits(Option_t * option) +{ // Print a table of digits + - AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; - TClonesArray * digits = gime->Digits() ; + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); + TClonesArray * digits = phosLoader->Digits() ; + const AliPHOSGeometry *geom = AliPHOSGeometry::GetInstance() ; + + AliInfo(Form("%d", digits->GetEntriesFast())) ; + printf("\nevent %d", gAlice->GetEvNumber()) ; + printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ; - cout << "AliPHOSDigitiser: event " << gAlice->GetEvNumber() << endl ; - cout << " Number of entries in Digits list " << digits->GetEntriesFast() << endl ; - cout << endl ; - if(strstr(option,"all")||strstr(option,"EMC")){ - + + if(strstr(option,"all")||strstr(option,"EMC")){ //loop over digits AliPHOSDigit * digit; - cout << "EMC digits (with primaries): " << endl ; - cout << "Digit Id Amplitude Index Nprim Primaries list " << endl; - Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ; + printf("\nEMC digits (with primaries):\n") ; + printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ; + Int_t maxEmc = geom->GetNModules()*geom->GetNCristalsInModule() ; Int_t index ; for (index = 0 ; (index < digits->GetEntriesFast()) && - (((AliPHOSDigit * ) digits->At(index))->GetId() <= maxEmc) ; index++) { + (static_cast(digits->At(index))->GetId() <= maxEmc) ; index++) { digit = (AliPHOSDigit * ) digits->At(index) ; - if(digit->GetNprimary() == 0) continue; - cout << setw(6) << digit->GetId() << " " << setw(10) << digit->GetAmp() << " " - << setw(6) << digit->GetIndexInList() << " " - << setw(5) << digit->GetNprimary() <<" "; - + if(digit->GetNprimary() == 0) + continue; +// printf("%6d %8d %6.5e %4d %2d :", +// digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK + printf("%6d %.4f %6.5e %4d %2d :", + digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; Int_t iprimary; - for (iprimary=0; iprimaryGetNprimary(); iprimary++) - cout << setw(5) << digit->GetPrimary(iprimary+1) << " "; - cout << endl; - } - cout << endl; + for (iprimary=0; iprimaryGetNprimary(); iprimary++) { + printf("%d ",digit->GetPrimary(iprimary+1) ) ; + } + printf("\n") ; + } } - + if(strstr(option,"all")||strstr(option,"CPV")){ //loop over CPV digits AliPHOSDigit * digit; - cout << "CPV digits: " << endl ; - cout << "Digit Id Amplitude Index Nprim Primaries list " << endl; - Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ; + printf("\nCPV digits (with primaries):\n") ; + printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ; + Int_t maxEmc = geom->GetNModules()*geom->GetNCristalsInModule() ; Int_t index ; for (index = 0 ; index < digits->GetEntriesFast(); index++) { digit = (AliPHOSDigit * ) digits->At(index) ; if(digit->GetId() > maxEmc){ - cout << setw(6) << digit->GetId() << " " << setw(10) << digit->GetAmp() << " " - << setw(6) << digit->GetIndexInList() << " " - << setw(5) << digit->GetNprimary() <<" "; - + printf("%6d %8d %4d %2d :", + digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ; Int_t iprimary; - for (iprimary=0; iprimaryGetNprimary(); iprimary++) - cout << setw(5) << digit->GetPrimary(iprimary+1) << " "; - cout << endl; - } + for (iprimary=0; iprimaryGetNprimary(); iprimary++) { + printf("%d ",digit->GetPrimary(iprimary+1) ) ; + } + printf("\n") ; + } } } } //__________________________________________________________________ -void AliPHOSDigitizer::SetSDigitsBranch(const char* title) -{ - // we set title (comment) of the SDigits branch in the first! header file - if( strcmp(GetName(), "") == 0 ) - Init() ; - - AliPHOSGetter::GetInstance()->SDigits()->SetName(title) ; - -} -//__________________________________________________________________ -Float_t AliPHOSDigitizer::TimeOfNoise(void) +Float_t AliPHOSDigitizer::TimeOfNoise(void) const { // Calculates the time signal generated by noise - //to be rewritten, now returns just big number - return 1. ; - + //PH Info("TimeOfNoise", "Change me") ; + return gRandom->Rndm() * 1.28E-5; } -//____________________________________________________________________________ -void AliPHOSDigitizer::Reset() -{ - // sets current event number to the first simulated event - - if( strcmp(GetName(), "") == 0 ) - Init() ; - // Int_t inputs ; -// for(inputs = 0; inputs < fNinputs ;inputs++) -// fIevent->AddAt(-1, inputs ) ; +//__________________________________________________________________ +void AliPHOSDigitizer::Unload() +{ + + Int_t i ; + for(i = 1 ; i < fInput ; i++){ + TString tempo(fEventNames[i]) ; + tempo += i ; + AliRunLoader* rl = AliRunLoader::GetRunLoader(tempo) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); + phosLoader->UnloadSDigits() ; + } + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); + phosLoader->UnloadDigits() ; } //____________________________________________________________________________ -void AliPHOSDigitizer::WriteDigits(Int_t event) +void AliPHOSDigitizer::WriteDigits() { // Makes TreeD in the output file. @@ -773,32 +984,24 @@ void AliPHOSDigitizer::WriteDigits(Int_t event) // and branch "AliPHOSDigitizer", with the same title to keep all the parameters // and names of files, from which digits are made. - AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; - const TClonesArray * digits = gime->Digits(GetName()) ; - TTree * treeD ; - + AliRunLoader* rl = AliRunLoader::GetRunLoader(fEventFolderName) ; + AliPHOSLoader * phosLoader = static_cast(rl->GetLoader("PHOSLoader")); - if(fManager) - treeD = fManager->GetTreeD() ; - else { - if (!gAlice->TreeD() ) - gAlice->MakeTree("D"); - treeD = gAlice->TreeD(); - } - - + const TClonesArray * digits = phosLoader->Digits() ; + TTree * treeD = phosLoader->TreeD(); + if(!treeD){ + phosLoader->MakeTree("D"); + treeD = phosLoader->TreeD(); + } + // -- create Digits branch Int_t bufferSize = 32000 ; - TBranch * digitsBranch = treeD->Branch("PHOS",&digits,bufferSize); - digitsBranch->SetTitle(GetName()); - - // -- Create Digitizer branch - Int_t splitlevel = 0 ; - AliPHOSDigitizer * d = gime->Digitizer(GetName()) ; - TBranch * digitizerBranch = treeD->Branch("AliPHOSDigitizer", "AliPHOSDigitizer", &d,bufferSize,splitlevel); - digitizerBranch->SetTitle(GetName()); - + TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize); + digitsBranch->SetTitle(fEventFolderName); digitsBranch->Fill() ; - treeD->AutoSave() ; - + + phosLoader->WriteDigits("OVERWRITE"); + + Unload() ; + }