#include "AliEMCALDigitizer.h"
#include "AliEMCAL.h"
#include "AliEMCALGeometry.h"
-//JLK
-//#include "AliEMCALHistoUtilities.h"
#include "AliEMCALRecParam.h"
#include "AliEMCALReconstructor.h"
#include "AliCDBManager.h"
//____________________________________________________________________________
AliEMCALClusterizerv1::AliEMCALClusterizerv1()
: AliEMCALClusterizer(),
- //JLK
- //fHists(0),fPointE(0),fPointL1(0),fPointL2(0),
- //fPointDis(0),fPointMult(0),fDigitAmp(0),fMaxE(0),
- //fMaxL1(0),fMaxL2(0),fMaxDis(0),
fGeom(0),
fDefaultInit(kFALSE),
fToUnfold(kFALSE),
if(!gMinuit)
gMinuit = new TMinuit(100) ;
- //JLK
- //fHists = BookHists();
}
//____________________________________________________________________________
while ( (digit = dynamic_cast<AliEMCALDigit *>(nextdigitC())) ) { // clean up digits
e = Calibrate(digit->GetAmp(), digit->GetId());
- //JLK
- //AliEMCALHistoUtilities::FillH1(fHists, 10, digit->GetAmp());
- //AliEMCALHistoUtilities::FillH1(fHists, 11, e);
if ( e < fMinECut || digit->GetTimeR() > fTimeCut )
digitsC->Remove(digit);
else
if(fGeom->CheckAbsCellId(digit->GetId()) && (Calibrate(digit->GetAmp(), digit->GetId()) > fECAClusteringThreshold ) ){
// start a new Tower RecPoint
if(fNumberOfECAClusters >= fRecPoints->GetSize()) fRecPoints->Expand(2*fNumberOfECAClusters+1) ;
+
AliEMCALRecPoint *recPoint = new AliEMCALRecPoint("") ;
fRecPoints->AddAt(recPoint, fNumberOfECAClusters) ;
recPoint = dynamic_cast<AliEMCALRecPoint *>(fRecPoints->At(fNumberOfECAClusters)) ;
} // while digit
delete digitsC ;
-
+
AliDebug(1,Form("total no of clusters %d from %d digits",fNumberOfECAClusters,fDigitsArr->GetEntriesFast()));
}
printf("Index Ene(GeV) Multi Module GX GY GZ lX lY lZ Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
}
Int_t index =0;
- Float_t maxE=0;
- Float_t maxL1=0;
- Float_t maxL2=0;
- Float_t maxDis=0;
-
- //JLK
- //AliEMCALHistoUtilities::FillH1(fHists, 12, double(fRecPoints->GetEntries()));
for (index = 0 ; index < fRecPoints->GetEntries() ; index++) {
AliEMCALRecPoint * rp = dynamic_cast<AliEMCALRecPoint * >(fRecPoints->At(index)) ;
rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(),
globalpos.X(), globalpos.Y(), globalpos.Z(), localpos.X(), localpos.Y(), localpos.Z(),
rp->GetDispersion(), lambda[0], lambda[1], nprimaries) ;
- /////////////
- if(rp->GetEnergy()>maxE){
- maxE=rp->GetEnergy();
- maxL1=lambda[0];
- maxL2=lambda[1];
- maxDis=rp->GetDispersion();
- }
- //JLK
- //fPointE->Fill(rp->GetEnergy());
- //fPointL1->Fill(lambda[0]);
- //fPointL2->Fill(lambda[1]);
- //fPointDis->Fill(rp->GetDispersion());
- //fPointMult->Fill(rp->GetMultiplicity());
- /////////////
if(strstr(option,"deb")){
for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
printf("%d ", primaries[iprimary] ) ;
}
}
- //JLK
- // fMaxE->Fill(maxE);
- // fMaxL1->Fill(maxL1);
- // fMaxL2->Fill(maxL2);
- // fMaxDis->Fill(maxDis);
-
if(strstr(option,"deb"))
printf("\n-----------------------------------------------------------------------\n");
}
}
-/*
-TList* AliEMCALClusterizerv1::BookHists()
-{
- //set up histograms for monitoring clusterizer performance
-
- gROOT->cd();
-
- fPointE = new TH1F("00_pointE","point energy", 2000, 0.0, 150.);
- fPointL1 = new TH1F("01_pointL1","point L1", 1000, 0.0, 3.);
- fPointL2 = new TH1F("02_pointL2","point L2", 1000, 0.0, 3.);
- fPointDis = new TH1F("03_pointDisp","point dispersion", 1000, 0.0, 10.);
- fPointMult = new TH1F("04_pointMult","#cell in point(cluster)", 101, -0.5, 100.5);
- fDigitAmp = new TH1F("05_digitAmp","Digit Amplitude", 2000, 0.0, 5000.);
- fMaxE = new TH1F("06_maxE","Max point energy", 2000, 0.0, 150.);
- fMaxL1 = new TH1F("07_maxL1","Largest (first) of eigenvalue of covariance matrix", 1000, 0.0, 3.);
- fMaxL2 = new TH1F("08_maxL2","Smalest (second) of eigenvalue of covariace matrix", 1000, 0.0, 3.);
- fMaxDis = new TH1F("09_maxDis","Point dispersion", 1000, 0.0, 10.); // 9
- //
- new TH1F("10_adcOfDigits","adc of digits(threshold control)", 1001, -0.5, 1000.5); // 10
- new TH1F("11_energyOfDigits","energy of digits(threshold control)", 1000, 0.0, 1.); // 11
- new TH1F("12_numberOfPoints","number of points(clusters)", 101, -0.5, 100.5); // 12
-
- return AliEMCALHistoUtilities::MoveHistsToList("EmcalClusterizerv1ControlHists", kFALSE);
-}
-
-void AliEMCALClusterizerv1::SaveHists(const char *fn)
-{
- AliEMCALHistoUtilities::SaveListOfHists(fHists, fn, kTRUE);
-}
-*/
-
//___________________________________________________________________
void AliEMCALClusterizerv1::PrintRecoInfo()
{
printf(" AliEMCALClusterizerv1::PrintRecoInfo() : version %s \n", Version() );
- //JLK
- //TH1F *h = (TH1F*)fHists->At(12);
- //if(h) {
- // printf(" ## Multiplicity of RecPoints ## \n");
- // for(int i=1; i<=h->GetNbinsX(); i++) {
- // int nbin = int((*h)[i]);
- // int mult = int(h->GetBinCenter(i));
- // if(nbin > 0) printf(" %i : %5.5i %6.3f %% \n", mult, nbin, 100.*nbin/h->GetEntries());
- // }
- // }
-
-}
-/*
-//___________________________________________________________________
-void AliEMCALClusterizerv1::DrawLambdasHists()
-{
- if(fMaxL1) {
- fMaxL1->Draw();
- if(fMaxL2) fMaxL2->Draw("same");
- if(fMaxDis) {
- fMaxDis->Draw("same");
- }
- }
}
-*/
static void UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Double_t & fret, Double_t * x, Int_t iflag) ;
// Chi^2 of the fit. Should be static to be passes to MINUIT
virtual const char * Version() const { return "clu-v1" ; }
-
- //JLK
- //TList* BookHists();
- //void SaveHists(const char *fn="reco.root"); //*MENU*
- //void DrawLambdasHists(); //*MENU*
void PrintRecoInfo(); //*MENU*
virtual void MakeClusters();
- //JLK
-/////////////////////
- // TList *fHists; //!
- // TH1F* fPointE; //histogram of point energy
- // TH1F* fPointL1; //histogram of point L1
- // TH1F* fPointL2; //histogram of point L2
- // TH1F* fPointDis; //histogram of point dispersion
- // TH1F* fPointMult; //histogram of point multiplicity
- // TH1F* fDigitAmp; //histogram of digit ADC Amplitude
- // TH1F* fMaxE; //histogram of maximum point energy
- // TH1F* fMaxL1; //histogram of largest (first) of eigenvalue of covariance matrix
- // TH1F* fMaxL2; //histogram of smalest (second) of eigenvalue of covariace matrix
- // TH1F* fMaxDis; //histogram of point dispersion
-///////////////////////
-
-
private:
AliEMCALClusterizerv1(const AliEMCALClusterizerv1 &); //copy ctor
AliEMCALClusterizerv1 & operator = (const AliEMCALClusterizerv1 &);
Float_t fTimeCut ; // Maximum time difference between the digits in ont EMC cluster
Float_t fMinECut; // Minimum energy for a digit to be a member of a cluster
- //JLK
- //ClassDef(AliEMCALClusterizerv1,6) // Clusterizer implementation version 1
ClassDef(AliEMCALClusterizerv1,7) // Clusterizer implementation version 1
};
//_________________________________________________________________________
//
//////////////////////////////////////////////////////////////////////////////
-// Class performs digitization of Summable digits
+// Class performs digitization of Summable digits from simulated data
//
// In addition it performs mixing of summable digits from different events.
//
#include <TTree.h>
#include <TSystem.h>
#include <TBenchmark.h>
-#include <TList.h>
-#include <TH1.h>
#include <TBrowser.h>
#include <TObjectTable.h>
#include <TRandom.h>
#include "AliEMCALSDigitizer.h"
#include "AliEMCALGeometry.h"
#include "AliEMCALTick.h"
-//JLK
-//#include "AliEMCALHistoUtilities.h"
ClassImp(AliEMCALDigitizer)
fEventFolderName(""),
fFirstEvent(0),
fLastEvent(0),
- //JLK
- //fControlHists(0),
- //fHists(0),
fCalibData(0x0)
{
// ctor
InitParameters() ;
- fManager = 0 ; // We work in the standalong mode
+ fManager = 0 ; // We work in the standalone mode
}
//____________________________________________________________________________
fEventFolderName(eventFolderName),
fFirstEvent(0),
fLastEvent(0),
- //JLK
- //fControlHists(0),
- //fHists(0),
fCalibData(0x0)
{
// ctor
InitParameters() ;
Init() ;
- fManager = 0 ; // We work in the standalong mode
+ fManager = 0 ; // We work in the standalone mode
}
//____________________________________________________________________________
fEventFolderName(d.fEventFolderName),
fFirstEvent(d.fFirstEvent),
fLastEvent(d.fLastEvent),
- //JLK
- //fControlHists(d.fControlHists),
- //fHists(d.fHists),
fCalibData(d.fCalibData)
{
// copyy ctor
fEventFolderName(0),
fFirstEvent(0),
fLastEvent(0),
- //JLK
- //fControlHists(0),
- //fHists(0),
fCalibData(0x0)
{
// ctor Init() is called by RunDigitizer
delete [] fInputFileNames ;
delete [] fEventNames ;
- //JLK
- //if(fHists) delete fHists;
}
//____________________________________________________________________________
// Makes the digitization of the collected summable digits
// for this it first creates the array of all EMCAL modules
- // filled with noise (different for EMC, CPV and PPSD) and
- // after that adds contributions from SDigits. This design
- // helps to avoid scanning over the list of digits to add
- // contribution of any new SDigit.
+ // filled with noise and after that adds contributions from
+ // SDigits. This design helps to avoid scanning over the
+ // list of digits to add contribution of any new SDigit.
+ //
+ // JLK 26-Jun-2008
+ // Note that SDigit energy info is stored as an amplitude, so we
+ // must call the Calibrate() method of the SDigitizer to convert it
+ // back to an energy in GeV before adding it to the Digit
+ //
static int nEMC=0; //max number of digits possible
AliRunLoader *rl = AliRunLoader::GetRunLoader();
rl->GetEvent(readEvent);
TClonesArray * digits = emcalLoader->Digits() ;
- digits->Delete() ; //JLK why is this created then deleted?
+ digits->Delete() ; //correct way to clear array when memory is
+ //allocated by objects stored in array
// Load Geometry
AliEMCALGeometry *geom = 0;
//Put Noise contribution
for(absID = 0; absID < nEMC; absID++){ // Nov 30, 2006 by PAI; was from 1 to nEMC
- Float_t amp = 0 ;
+ Float_t energy = 0 ;
// amplitude set to zero, noise will be added later
new((*digits)[absID]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ); // absID-1->absID
//look if we have to add signal?
// Calculate time as time of the largest digit
Float_t time = digit->GetTime() ;
- Float_t eTime= digit->GetAmp() ;
+ Float_t aTime= digit->GetAmp() ;
// loop over input
for(i = 0; i< fInput ; i++){ //loop over (possible) merge sources
//b = a /fTimeSignalLength ;
//new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime(),0, b);
//new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
- if(curSDigit->GetAmp()>eTime) {
- eTime = curSDigit->GetAmp();
+ if(curSDigit->GetAmp()>aTime) {
+ aTime = curSDigit->GetAmp();
time = curSDigit->GetTime();
}
- *digit = *digit + *curSDigit ; //add energies
+ *digit = *digit + *curSDigit ; //adds amplitudes of each digit
index[i]++ ;
if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
curSDigit = 0 ;
}
}
+ //Here we convert the summed amplitude to an energy in GeV
+ energy = sDigitizer->Calibrate(digit->GetAmp()) ; // GeV
// add fluctuations for photo-electron creation
- amp = sDigitizer->Calibrate(digit->GetAmp()) ; // GeV
- amp *= static_cast<Float_t>(gRandom->Poisson(fMeanPhotonElectron)) / static_cast<Float_t>(fMeanPhotonElectron) ;
+ energy *= static_cast<Float_t>(gRandom->Poisson(fMeanPhotonElectron)) / static_cast<Float_t>(fMeanPhotonElectron) ;
//calculate and set time
//New timing model needed - JLK 28-April-2008
}
}
// add the noise now
- amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
- digit->SetAmp(sDigitizer->Digitize(amp)) ;
- AliDebug(10,Form(" absID %5i amp %f nextSig %5i\n",
- absID, amp, nextSig));
- } // for(absID = 1; absID <= nEMC; absID++)
+ energy += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
+ // JLK 26-June-2008
+ //Now digitize the energy according to the sDigitizer method,
+ //which merely converts the energy to an integer. Later we will
+ //check that the stored value matches our allowed dynamic ranges
+ digit->SetAmp(sDigitizer->Digitize(energy)) ;
+ AliDebug(10,Form(" absID %5i energy %f nextSig %5i\n",
+ absID, energy, nextSig));
+ } // for(absID = 0; absID < nEMC; absID++)
//ticks->Delete() ;
//delete ticks ;
+ //JLK is it better to call Clear() here?
delete sdigArray ; //We should not delete its contents
//remove digits below thresholds
for(i = 0 ; i < nEMC ; i++){
digit = dynamic_cast<AliEMCALDigit*>( digits->At(i) ) ;
- Float_t threshold = fDigitThreshold ;
+ Float_t threshold = fDigitThreshold ; //this is in GeV
+ //need to calibrate digit amplitude to energy in GeV for comparison
if(sDigitizer->Calibrate( digit->GetAmp() ) < threshold)
digits->RemoveAt(i) ;
else
Int_t ndigits = digits->GetEntriesFast() ;
- //JLK
- //Set indexes in list of digits and fill hists.
- //AliEMCALHistoUtilities::FillH1(fHists, 0, Double_t(ndigits));
- Float_t energy=0., esum=0.;
+ //JLK 26-June-2008
+ //After we have done the summing and digitizing to create the
+ //digits, now we want to calibrate the resulting amplitude to match
+ //the dynamic range of our real data.
+ Float_t energy=0;
for (i = 0 ; i < ndigits ; i++) {
digit = dynamic_cast<AliEMCALDigit *>( digits->At(i) ) ;
digit->SetIndexInList(i) ;
energy = sDigitizer->Calibrate(digit->GetAmp()) ;
- esum += energy;
- digit->SetAmp(DigitizeEnergy(energy, digit->GetId()) ) ; // for what ??
- //JLK
- //AliEMCALHistoUtilities::FillH1(fHists, 2, double(digit->GetAmp()));
- //AliEMCALHistoUtilities::FillH1(fHists, 3, double(energy));
- //AliEMCALHistoUtilities::FillH1(fHists, 4, double(digit->GetId()));
- // if(digit->GetId() == nEMC) {
- // printf(" i %i \n", i );
- // digit->Dump();
- // assert(0);
- //}
+ digit->SetAmp(DigitizeEnergy(energy, digit->GetId()) ) ;
}
- //JLK
- //AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
+
}
// //_____________________________________________________________________
Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t AbsId)
{
+ // JLK 26-June-2008
// Returns digitized value of the energy in a cell absId
-
+ // using the calibration constants stored in the OCDB
+ // or default values if no CalibData object is found.
+ // This effectively converts everything to match the dynamic range
+ // of the real data we will collect
+ //
// Load Geometry
const AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
//to prevent cleaning of this object while GetEvent is called
emcalLoader->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
- //PH Print();
//Calibration instance
fCalibData = emcalLoader->CalibData();
return fInit ;
fNADCEC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC - 65536
fTimeThreshold = 0.001*10000000 ; // Means 1 MeV in terms of SDigits amplitude ??
- //JLK
- // hists. for control; no hists on default
- //fControlHists = 0;
- //fHists = 0;
+
}
//__________________________________________________________________
//__________________________________________________________________
void AliEMCALDigitizer::Browse(TBrowser* b)
{
- //JLK
- //if(fHists) b->Add(fHists);
TTask::Browse(b);
}
-
-/*
-//__________________________________________________________________
-TList *AliEMCALDigitizer::BookControlHists(int var)
-{
- // 22-nov-04
- // histograms for monitoring digitizer performance
-
- Info("BookControlHists"," started ");
- gROOT->cd();
- const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
- if(var>=1){
- new TH1F("hDigiN", "#EMCAL digits with fAmp > fDigitThreshold",
- fNADCEC+1, -0.5, Double_t(fNADCEC));
- new TH1F("HDigiSumEnergy","Sum.EMCAL energy from digi", 1000, 0.0, 200.);
- new TH1F("hDigiAmp", "EMCAL digital amplitude", fNADCEC+1, -0.5, Double_t(fNADCEC));
- new TH1F("hDigiEnergy","EMCAL cell energy", 2000, 0.0, 200.);
- new TH1F("hDigiAbsId","EMCAL absId cells with fAmp > fDigitThreshold ",
- geom->GetNCells(), 0.5, Double_t(geom->GetNCells())+0.5);
- }
-
- fHists = AliEMCALHistoUtilities::MoveHistsToList("EmcalDigiControlHists", kFALSE);
- fHists = 0; //huh? JLK 03-Mar-2006
- return fHists;
-}
-
-//__________________________________________________________________
-void AliEMCALDigitizer::SaveHists(const char* name, Bool_t kSingleKey, const char* opt)
-{
- AliEMCALHistoUtilities::SaveListOfHists(fHists, name, kSingleKey, opt);
-}
-*/
/* $Id$ */
//_________________________________________________________________________
-// Task Class for making Digits in EMCAL
+// Task Class for making Digits in EMCAL from SIMULATED DATA
//
//*-- Author: Sahal Yacoob (LBL)
// based on : AliPHOSDigit
// --- ROOT system ---
#include "TObjString.h"
-class TArrayI ;
class TClonesArray ;
-class TList;
class TBrowser;
// --- Standard library ---
#include "AliDigitizer.h"
#include "AliConfig.h"
#include "AliEMCALCalibData.h"
+
class AliEMCALSDigitizer ;
class AliRunDigitizer ;
}
virtual void Browse(TBrowser* b);
- //JLK
- // hists
- //void SetControlHists(Int_t var=0) {fControlHists=var;}
- //Int_t GetControlHist() const {return fControlHists;}
- //TList *GetListOfHists() {return fHists;}
- //TList* BookControlHists(int var=0);
- //void SaveHists(const char* name="RF/TRD1/Digitizations/DigiVar?",
- //Bool_t kSingleKey=kTRUE, const char* opt="RECREATE"); // *MENU*
private:
TString fEventFolderName; // skowron: name of EFN to read data from in stand alone mode
Int_t fFirstEvent; // first event to process
Int_t fLastEvent; // last event to process
- //JLK
- // Control hists
- //Int_t fControlHists; //!
- //TList *fHists; //!
+
AliEMCALCalibData * fCalibData; //Calibration data pointer
- //JLK
- //ClassDef(AliEMCALDigitizer,6) // description
ClassDef(AliEMCALDigitizer,7) // description
};
// A Summable Digits is the sum of all hits originating
// from one in one tower of the EMCAL
// A threshold for assignment of the primary to SDigit is applied
+//
+// JLK 26-Jun-2008 Added explanation:
+// SDigits need to hold the energy sum of the hits, but AliEMCALDigit
+// can (should) only store amplitude. Therefore, the SDigit energy is
+// "digitized" before being stored and must be "calibrated" back to an
+// energy before SDigits are summed to form true Digits
+//
// SDigits are written to TreeS, branch "EMCAL"
// AliEMCALSDigitizer with all current parameters is written
// to TreeS branch "AliEMCALSDigitizer".
// --- ROOT system ---
#include <TBenchmark.h>
-#include <TH1.h>
#include <TBrowser.h>
#include <Riostream.h>
#include <TMath.h>
#include "AliEMCALHit.h"
#include "AliEMCALSDigitizer.h"
#include "AliEMCALGeometry.h"
-//JLK
-//#include "AliEMCALHistoUtilities.h"
ClassImp(AliEMCALSDigitizer)
fFirstEvent(0),
fLastEvent(0),
fSampling(0.)
- //JLK
- //fControlHists(0),
- //fHists(0)
{
// ctor
InitParameters();
fFirstEvent(0),
fLastEvent(0),
fSampling(0.)
- //JLK
- //fControlHists(1),
- //fHists(0)
{
// ctor
Init();
InitParameters() ;
- //JLK
- //if(fControlHists) BookControlHists(1);
}
fFirstEvent(sd.fFirstEvent),
fLastEvent(sd.fLastEvent),
fSampling(sd.fSampling)
- //JLK
- //fControlHists(sd.fControlHists),
- //fHists(sd.fHists)
{
//cpy ctor
}
if (geom->GetSampling() == 0.) {
Fatal("InitParameters", "Sampling factor not set !") ;
}
+
+ //
+ //JLK 26-Jun-2008 THIS SHOULD HAVE BEEN EXPLAINED AGES AGO:
+ //
+ //In order to be able to store SDigit Energy info into
+ //AliEMCALDigit, we need to convert it temporarily to an ADC amplitude
+ //and later when summing SDigits to form digits, convert it back to
+ //energy. These fA and fB parameters accomplish this through the
+ //Digitize() and Calibrate() methods
+ //
// Initializes parameters
fA = 0;
fB = 1.e+6; // Changed 24 Apr 2007. Dynamic range now 2 TeV
AliDebug(2,Form(" Sampling = %f\n", fSampling));
AliDebug(2,Form("---------------------------------------------------\n"));
- // Print("");
}
AliDebug(2,Form(" Sampling = %f\n", fSampling));
AliDebug(2,Form("---------------------------------------------------\n"));
- // Print();
return ;
}
nSdigits = sdigits->GetEntriesFast() ;
fSDigitsInRun += nSdigits ;
-
- //JLK
- //Double_t e=0.,esum=0.;
- //AliEMCALHistoUtilities::FillH1(fHists, 0, double(sdigits->GetEntriesFast()));
+
for (i = 0 ; i < sdigits->GetEntriesFast() ; i++) {
AliEMCALDigit * sdigit = dynamic_cast<AliEMCALDigit *>(sdigits->At(i)) ;
sdigit->SetIndexInList(i) ;
-
- //JLK
- //AliEMCALHistoUtilities::FillH1(fHists, 2, double(sdigit->GetAmp()));
- //e = double(Calibrate(sdigit->GetAmp()));
- //esum += e;
- //AliEMCALHistoUtilities::FillH1(fHists, 3, e);
- //AliEMCALHistoUtilities::FillH1(fHists, 4, double(sdigit->GetId()));
}
- //JLK
- //if(esum>0.) AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
// Now write SDigits
//__________________________________________________________________
Int_t AliEMCALSDigitizer::Digitize(Float_t energy)const {
// Digitize the energy
- Double_t aSignal = fA + energy*fB;
- if (TMath::Abs(aSignal)>2147483647.0) {
- //PH 2147483647 is the max. integer
- //PH This apparently is a problem which needs investigation
- AliWarning(Form("Too big or too small energy %f",aSignal));
- aSignal = TMath::Sign((Double_t)2147483647,aSignal);
- }
- return (Int_t ) aSignal;
+ //
+ //JLK 26-Jun-2008 EXPLANATION LONG OVERDUE:
+ //
+ //We have to digitize the SDigit energy so that it can be stored in
+ //AliEMCALDigit, which has only an ADC amplitude member and
+ //(rightly) no energy member. This method converts the energy to an
+ //integer which can be re-converted back to an energy with the
+ //Calibrate(energy) method when it is time to create Digits from SDigits
+ //
+ Double_t aSignal = fA + energy*fB;
+ if (TMath::Abs(aSignal)>2147483647.0) {
+ //PH 2147483647 is the max. integer
+ //PH This apparently is a problem which needs investigation
+ AliWarning(Form("Too big or too small energy %f",aSignal));
+ aSignal = TMath::Sign((Double_t)2147483647,aSignal);
}
-
+ return (Int_t ) aSignal;
+}
//__________________________________________________________________
+Float_t AliEMCALSDigitizer::Calibrate(Int_t amp)const {
+ //
+ // Convert the amplitude back to energy in GeV
+ //
+ //JLK 26-Jun-2008 EXPLANATION LONG OVERDUE:
+ //
+ //We have to digitize the SDigit energy with the method Digitize()
+ //so that it can be stored in AliEMCALDigit, which has only an ADC
+ //amplitude member and (rightly) no energy member. This method is
+ //just the reverse of Digitize(): it converts the stored amplitude
+ //back to an energy value in GeV so that the SDigit energies can be
+ //summed before adding noise and creating digits out of them
+ //
+ return (Float_t)(amp - fA)/fB;
+
+}
+
+//__________________________________________________________________
void AliEMCALSDigitizer::Print1(Option_t * option)
{
Print();
}
}
delete tempo ;
- printf("\n** Sum %i : %10.3f GeV/c **\n ", isum, double(isum)*1.e-6);
+ printf("\n** Sum %i : %10.3f GeV/c **\n ", isum, Calibrate(isum));
} else printf("\n");
}
//____________________________________________________________________________
void AliEMCALSDigitizer::Browse(TBrowser* b)
{
- //JLK
- //if(fHists) b->Add(fHists);
TTask::Browse(b);
}
-
-/*
-//____________________________________________________________________________
-TList *AliEMCALSDigitizer::BookControlHists(int var)
-{
- //book histograms for monitoring sdigitization
- // 22-nov-04
- gROOT->cd();
- const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance() ;
- if(var>=1){
- AliDebug(1, " BookControlHists() in action ");
- new TH1F("HSDigiN", "#EMCAL sdigits ", 1001, -0.5, 1000.5);
- new TH1F("HSDigiSumEnergy","Sum.EMCAL energy", 1000, 0.0, 100.);
- new TH1F("HSDigiAmp", "EMCAL sdigits amplitude", 1000, 0., 2.e+9);
- new TH1F("HSDigiEnergy","EMCAL cell energy", 1000, 0.0, 100.);
- new TH1F("HSDigiAbsId","EMCAL absID for sdigits",
- geom->GetNCells(), 0.5, Double_t(geom->GetNCells())+0.5);
- }
-
- fHists = AliEMCALHistoUtilities::MoveHistsToList("EmcalSDigiControlHists", kFALSE);
- // fHists = 0; ??
-
- return fHists;
-}
-
-//____________________________________________________________________________
-void AliEMCALSDigitizer::SaveHists(const char* name, Bool_t kSingleKey, const char* opt)
-{
- AliEMCALHistoUtilities::SaveListOfHists(fHists, name, kSingleKey, opt);
-}
-*/
// --- ROOT system ---
#include "TTask.h"
-class TFile ;
-class TList;
class TBrowser;
-//class TBrowser;
// --- Standard library ---
AliEMCALSDigitizer(const AliEMCALSDigitizer & sd) ;
virtual ~AliEMCALSDigitizer(); // dtor
- Float_t Calibrate(Int_t amp)const {return (amp - fA)/fB ; }
- Int_t Digitize(Float_t energy)const;
+ Int_t Digitize(Float_t energy)const; //convert energy in GeV to int amplitude
+ Float_t Calibrate(Int_t amp)const; //opposite of Digitize()
+
virtual void Exec(Option_t *option);
Int_t GetSDigitsInRun() const {return fSDigitsInRun ;}
virtual void Print(Option_t *option="") const;
const AliEMCALSDigitizer & operator = (const AliEMCALSDigitizer & /*sd*/) {return *this ;}
virtual void Browse(TBrowser* b);
- //JLK
- // hists
- //void SetControlHists(Int_t var=0) {fControlHists=var;}
- //Int_t GetControlHist() const {return fControlHists;}
- //TList *GetListOfHists() {return fHists;}
- //TList* BookControlHists(int var=0);
- //void SaveHists(const char* name="RF/TRD1/Digitizations/SDigiVar?",
- //Bool_t kSingleKey=kTRUE, const char* opt="RECREATE"); // *MENU*
private:
void Init() ;
Int_t fFirstEvent; // first event to process
Int_t fLastEvent; // last event to process
Float_t fSampling; // See AliEMCALGeometry
- //JLK
- // Control hists
- //Int_t fControlHists; //!
- //TList *fHists; //!
- //JLK
- //ClassDef(AliEMCALSDigitizer,5) // description
ClassDef(AliEMCALSDigitizer,6) // description
};
// --- ROOT system ---
#include <TBrowser.h>
#include <TClonesArray.h>
-#include <TH1.h>
#include <TH2.h>
#include <TParticle.h>
#include <TROOT.h>
//______________________________________________________________________
AliEMCALv2::AliEMCALv2()
- : AliEMCALv1(),
- fHDe(0),
- fHNhits(0)
+ : AliEMCALv1()
{
// ctor
}
//______________________________________________________________________
AliEMCALv2::AliEMCALv2(const char *name, const char *title)
- : AliEMCALv1(name,title),
- fHDe(0),
- fHNhits(0)
+ : AliEMCALv1(name,title)
{
// Standard Creator.
fTimeCut = 30e-09;
fGeometry = GetGeometry();
- fHDe = fHNhits = 0;
- // if (gDebug>0){
- if (1){
- TH1::AddDirectory(0);
- fHDe = new TH1F("fHDe","De in EMCAL", 1000, 0., 10.);
- fHNhits = new TH1F("fHNhits","#hits in EMCAL", 2001, -0.5, 2000.5);
- fHistograms->Add(fHDe);
- fHistograms->Add(fHNhits);
- TH1::AddDirectory(1);
- }
}
//______________________________________________________________________
}
}
-//_________________________________________________________________
-void AliEMCALv2::FinishEvent()
-{
- // Calculate deposit energy and fill control histogram; 26-may-05
- static double de=0.;
- fHNhits->Fill(double(fHits->GetEntries()));
- de = GetDepositEnergy(0);
- if(fHDe) fHDe->Fill(de);
-}
-
-//_________________________________________________________________
-Double_t AliEMCALv2::GetDepositEnergy(int print)
-{
- // 23-mar-05 - for testing
- if(fHits == 0) return 0.;
- AliEMCALHit *hit=0;
- Double_t de=0.;
- for(int ih=0; ih<fHits->GetEntries(); ih++) {
- hit = (AliEMCALHit*)fHits->UncheckedAt(ih);
- de += hit->GetEnergy();
- }
- if(print>0) {
- cout<<"AliEMCALv2::GetDepositEnergy() : fHits "<<fHits<<endl;
- printf(" #hits %i de %f \n", fHits->GetEntries(), de);
- if(print>1) {
- printf(" #primary particles %i\n", gAlice->GetHeader()->GetNprimary());
- }
- }
- return de;
-}
-
//___________________________________________________________
void AliEMCALv2::Browse(TBrowser* b)
{
//*-- Author: Aleksei Pavlinov
// --- ROOT system ---
-class TClonesArray;
-class TLorentzVector;
-class TFile;
-class TH1F;
-
class TBrowser;
// --- AliRoot header files ---
Int_t id, Float_t *hits, Float_t *p);
virtual void StepManager(void) ;
- virtual void FinishEvent();
// Gives the version number
virtual Int_t IsVersion(void) const {return 2;}
virtual const TString Version(void)const {return TString("v2");}
- // 23-mar-05
- virtual Double_t GetDepositEnergy(int print=1); // *MENU*
+
// 30-aug-04
virtual void Browse(TBrowser* b);
// drawing
void TestIndexTransition(int pri=0, int idmax=0); // *MENU*
protected:
- TH1F* fHDe; //!
- TH1F* fHNhits; //!
private:
AliEMCALv2(const AliEMCALv2 & emcal);
AliEMCALv2 & operator = (const AliEMCALv2 & /*rvalue*/);
- ClassDef(AliEMCALv2,1) //Implementation of EMCAL manager class to produce hits in a Shish-Kebab
+ ClassDef(AliEMCALv2,2) //Implementation of EMCAL manager class to produce hits in a Shish-Kebab
};
#include "TParticle.h"
#include "TVirtualMC.h"
#include "TBrowser.h"
-#include "TH1.h"
#include "TH2.h"
#include <cassert>
fTimeCut = 30e-09;
fGeometry = GetGeometry();
- fHDe = fHNhits = 0;
- // if (gDebug>0){
- if (1){
- TH1::AddDirectory(0);
- fHDe = new TH1F("fHDe","De in EMCAL", 1000, 0., 1.);
- fHNhits = new TH1F("fHNhits","#hits in EMCAL", 1001, -0.5, 1000.5);
- int nbin = 77*2;
- fHDeDz = new TH1F("fHDeDz","Longitudinal shower profile", 6*nbin, 0.0, 0.16*nbin);
-
- fHistograms->Add(fHDe);
- fHistograms->Add(fHNhits);
- fHistograms->Add(fHDeDz);
- TH1::AddDirectory(1);
- }
}
//______________________________________________________________________
AliEMCALv3::AliEMCALv3(const AliEMCALv3 & emcal):AliEMCALv1(emcal)
{
fGeometry = emcal.fGeometry;
- fHDe = emcal.fHDe;
- fHNhits = emcal.fHNhits;
- fHDeDz = emcal.fHDeDz;
-
}
//______________________________________________________________________
}
}
-void AliEMCALv3::FinishEvent()
-{ // 26-may-05
- static double de=0.;
- fHNhits->Fill(double(fHits->GetEntries()));
- de = GetDepositEnergy(0);
- if(fHDe) fHDe->Fill(de);
-}
-
-Double_t AliEMCALv3::GetDepositEnergy(int print)
-{ // 23-mar-05 - for testing
- // cout<<"AliEMCALv3::GetDepositEnergy() : fHits "<<fHits<<endl;
- if(fHits == 0) return 0.;
- static AliEMCALHitv1 *hit=0;
- static Double_t de=0., deHit=0., zhl=0.0, zShift = 0.16*77;
- for(int ih=0; ih<fHits->GetEntries(); ih++) {
- hit = (AliEMCALHitv1*)fHits->UncheckedAt(ih);
- deHit = (double)hit->GetEnergy();
- zhl = (double)hit->Z() + zShift;
- if(zhl>2.*zShift) zhl = 2*zShift - 0.002;
- de += deHit;
- if(fHDeDz) fHDeDz->Fill(zhl, deHit);
- }
- if(print>0) {
- printf(" #hits %i de %f \n", fHits->GetEntries(), de);
- if(print>1) {
- printf(" #primary particles %i\n", gAlice->GetHeader()->GetNprimary());
- }
- }
- return de;
-}
-
+//_____________________________________________
void AliEMCALv3::Browse(TBrowser* b)
{
TObject::Browse(b);
//*--
//*-- Author: Aleksei Pavlinov
-class TClonesArray;
-class TLorentzVector;
-class TFile;
-class TH1F;
-
class AliEMCALGeometry;
// --- AliRoot header files ---
Int_t id, Float_t *hits, Float_t *p);
virtual void StepManager(void) ;
- virtual void FinishEvent();
// Gives the version number
virtual Int_t IsVersion(void) const {return 3;}
Fatal("operator =", "not implemented") ;
return *this;}
- virtual Double_t GetDepositEnergy(int print=1); // *MENU*
virtual void Browse(TBrowser* b);
AliEMCALGeometry* fGeometry; //!
- TH1F* fHDe; //!
- TH1F* fHNhits; //!
- TH1F* fHDeDz; //!
- ClassDef(AliEMCALv3,0) //Implementation of EMCAL manager class to produce hits in a Shish-Kebab
+ ClassDef(AliEMCALv3,1) //Implementation of EMCAL manager class to produce hits in a Shish-Kebab
};