[u/mrichter/AliRoot.git] / JETAN / AliJetESDFillUnitArrayEMCalDigits.cxx


/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
 
//---------------------------------------------------------------------
// Fill Unit Array 
// Called by ESD reader for jet analysis
// Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
//---------------------------------------------------------------------

// --- ROOT system ---
#include <TLorentzVector.h>
#include <TMath.h>

// --- AliRoot header files ---
#include "AliESDCaloCluster.h"
#include "AliESDCaloCells.h"
#include "AliJetUnitArray.h"
#include "AliJetESDFillUnitArrayEMCalDigits.h"
// #include "AliEMCALCalibData.h"
// #include "AliCDBManager.h"
// // class AliCDBStorage;
// #include "AliCDBEntry.h"

// --- ROOT system ---
class TSystem;
class TGeoManager;
class TArrayS;

// --- AliRoot header files ---
class AliJetFinder;
class AliJetReader;
class AliJetESDReader;
// class AliEMCALCalibData;
// class AliCDBManager;
// class AliCDBStorage;
// class AliCDBEntry;

ClassImp(AliJetESDFillUnitArrayEMCalDigits)

//_____________________________________________________________________________
AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits()
  : AliJetFillUnitArray(),
    fESD(0),
    fNIn(0),
    fOpt(0),
    fCluster(0),
    fDebug(0),
    fNCEMCAL(0),
    fNCPHOS(0),
    fNCCalo(0),
    fApplyElectronCorrection(kFALSE),
    fApplyFractionHadronicCorrection(kFALSE),
    fFractionHadronicCorrection(0.3),
    fClus(0x0),
    fNDigitEmcal(0),
    fNDigitEmcalCut(0)
{
  // constructor
}

//_____________________________________________________________________________
AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits(AliESDEvent *esd)
  : AliJetFillUnitArray(),
    fESD(esd),
    fNIn(0),
    fOpt(0),
    fCluster(0),
    fDebug(0),
    fNCEMCAL(0),
    fNCPHOS(0),
    fNCCalo(0),
    fApplyElectronCorrection(kFALSE),
    fApplyFractionHadronicCorrection(kFALSE),
    fFractionHadronicCorrection(0.3),
    fClus(0x0),
    fNDigitEmcal(0),
    fNDigitEmcalCut(0)
{
  // constructor
}

//_____________________________________________________________________________
AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits(const AliJetESDFillUnitArrayEMCalDigits &det)
  : AliJetFillUnitArray(det),
    fESD(det.fESD),
    fNIn(det.fNIn),
    fOpt(det.fOpt),
    fCluster(det.fCluster),
    fDebug(det.fDebug),
    fNCEMCAL(det.fNCEMCAL),
    fNCPHOS(det.fNCPHOS),
    fNCCalo(det.fNCCalo),
    fApplyElectronCorrection(det.fApplyElectronCorrection),
    fApplyFractionHadronicCorrection(det.fApplyFractionHadronicCorrection),
    fFractionHadronicCorrection(det.fFractionHadronicCorrection),
    fClus(det.fClus),
    fNDigitEmcal(det.fNDigitEmcal),
    fNDigitEmcalCut(det.fNDigitEmcalCut)
{
  // Copy constructor
}

//_____________________________________________________________________________
AliJetESDFillUnitArrayEMCalDigits& AliJetESDFillUnitArrayEMCalDigits::operator=(const AliJetESDFillUnitArrayEMCalDigits& other)
{
  // Assignment

  fESD = other.fESD;
  fNIn = other.fNIn;
  fOpt = other.fOpt;
  fCluster = other.fCluster;
  fDebug = other.fDebug;
  fNCEMCAL = other.fNCEMCAL;
  fNCPHOS = other.fNCPHOS;
  fNCCalo = other.fNCCalo;
  fApplyElectronCorrection = other.fApplyElectronCorrection;
  fApplyFractionHadronicCorrection = other.fApplyFractionHadronicCorrection;
  fFractionHadronicCorrection = other.fFractionHadronicCorrection;
  fClus = other.fClus;
  fNDigitEmcal = other.fNDigitEmcal;
  fNDigitEmcalCut = other.fNDigitEmcalCut;

  return (*this);

}

//_____________________________________________________________________________
AliJetESDFillUnitArrayEMCalDigits::~AliJetESDFillUnitArrayEMCalDigits()
{
  // destructor
}

//_____________________________________________________________________________
void AliJetESDFillUnitArrayEMCalDigits::Exec(Option_t* const /*option*/)
{
  //
  // Main method.
  // Fill the unit array with the neutral particle information from the EMCal cells in ESD
  //

  fDebug = fReaderHeader->GetDebug();
  fOpt = fReaderHeader->GetDetector();
  fCluster = fReaderHeader->GetCluster();

  // Init parameters
  InitParameters();

//(not used ?)  Int_t   nDigitTot      = 0;
  Int_t   goodDigit      = 0;
//(not used ?)  Int_t   beg            = 0;
//(not used ?)  Int_t   end            = 0; 
  Int_t   index          = 0;
//(not used ?)  Int_t   count          = 0;

//(not used ?)  Int_t nRefEnt = fRefArray->GetEntries();

  if(!fCluster) { // Keep all digit information
    // Loop over all cell information
    //------------------------------------------------------------------
    AliESDCaloCells &cells= *(fESD->GetEMCALCells());
    Int_t ncell = cells.GetNumberOfCells() ;
//(not used ?)    Int_t type = cells.GetType();

    for (Int_t icell=  0; icell <  ncell; icell++) {
      Int_t   digitID   = cells.GetCellNumber(icell);
      Float_t digitAmp  = cells.GetAmplitude(icell);
//(not used ?)      Float_t digitTime = cells.GetTime(icell);
      Float_t digitEn   = digitAmp*0.0153; // Last correct
      //  Float_t digitEn = Calibrate(digitAmp,digitID);

      Float_t etaD=-10., phiD=-10.;
      fGeom->EtaPhiFromIndex(digitID,etaD,phiD); 
      //  fEMCalGrid->GetEtaPhiFromIndex2(digitID,phiD,etaD);

      phiD = ((phiD < 0) ? phiD + 2.* TMath::Pi() : phiD);
      
      Float_t digitEt = digitEn*TMath::Abs(TMath::Sin(EtaToTheta(etaD)));

      AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(digitID);

      if(uArray->GetUnitEnergy() == 0.) goodDigit++;
      uArray->SetUnitTrackID(digitID);

      Float_t unitEnergy = 0.;
      Bool_t ok = kFALSE;
      unitEnergy = uArray->GetUnitEnergy();

      if(unitEnergy==0){
	if(!fProcId){
	  new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
	  fProcId = kTRUE;
	}
	fRefArray->Add(uArray);
	fNDigitEmcal++;
	ok = kTRUE;
      }

      // Detector flag
      if(unitEnergy>0.)
	uArray->SetUnitDetectorFlag(kAll);
      else uArray->SetUnitDetectorFlag(kEmcal);
      
      uArray->SetUnitEnergy(unitEnergy+digitEt);

      uArray->SetUnitCutFlag(kPtHigher);

      // To be modified !!!
      uArray->SetUnitSignalFlag(kGood);

      // This is for jet multiplicity
      uArray->SetUnitClusterID(index);
	
      if(fDebug > 1) printf("goodDigit : %d\n", goodDigit);

    } // End loop over cells
  } // end if !fCluster
  else { // Keep digit information from clusterization

    // Loop over calo clusters
    //------------------------------------------------------------------

    //select EMCAL clusters only
    TRefArray * caloClusters  = new TRefArray();
    fESD->GetEMCALClusters(caloClusters);

    // Total number of EMCAL cluster
    Int_t nclus = caloClusters->GetEntries() ;
    Int_t beg   = 0;
    Float_t pos[3] ;

    // Get reconstructed vertex position
    Double_t vertexPosition[3] ;
    fESD->GetVertex()->GetXYZ(vertexPosition) ;

    // Get CaloCells
    AliESDCaloCells &cells= *(fESD->GetEMCALCells());
//(not used ?)    Int_t ncell = cells.GetNumberOfCells() ;
//(not used ?)    Int_t type = cells.GetType();

    for(Int_t j = beg; j < nclus; j++) { // loop over clusters
      // Retrieve cluster from esd
      fClus = (AliESDCaloCluster *) caloClusters->At(j) ;

      // Get the cluster info
//(not used ?)      Float_t energy         = fClus->E() ;
//(not used ?)      Int_t   iprim          = fClus->GetLabel();

      fClus->GetPosition(pos) ;
      TVector3 vpos(pos[0],pos[1],pos[2]) ;
      TLorentzVector p ;
      fClus->GetMomentum(p,vertexPosition);

      Int_t     digMult = fClus->GetNCells() ;
      UShort_t *digID   = fClus->GetCellsAbsId() ;
//(not used ?)      Double_t *digAmpFrac = fClus->GetCellsAmplitudeFraction() ;
      Int_t     trackIndex = fClus->GetTrackMatched();

      // Do double-counted electron correction 
      if (fApplyElectronCorrection != 0 && trackIndex !=-1 )
	{
          // The electron correction go there
          // Under construction !!!!
	}  // End of Electron correction

      // Get CaloCells of cluster and fill the unitArray
      for(Int_t i = 0; i < digMult ; i++){
        Int_t    digitID   = digID[i]; // or clus->GetCellNumber(i) ;
//(not used ?)        Double_t digitAmpFrac = digAmpFrac[i];
        Float_t  digitAmp     = cells.GetCellAmplitude(digitID) ;
//(not used ?)        Float_t  digitTime    = cells.GetCellTime(digitID);

        // Calibration for an energy in GeV
        Float_t digitEn = digitAmp*0.0153;

        Float_t etaD=-10., phiD=-10.;
        fGeom->EtaPhiFromIndex(digitID,etaD,phiD);
        //  fEMCalGrid->GetEtaPhiFromIndex2(digitID,phiD,etaD);

        phiD = ((phiD < 0) ? phiD + 2.* TMath::Pi() : phiD);

        Float_t digitEt = digitEn*TMath::Abs(TMath::Sin(EtaToTheta(etaD)));

	cout << "Digit " << i << ", eta: " << etaD << ", phi: " << phiD << endl;

	AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(digitID);
	if(uArray->GetUnitEnergy() == 0.) goodDigit++;
	uArray->SetUnitTrackID(digitID);

	Float_t unitEnergy = 0.;
	Bool_t ok = kFALSE;
	unitEnergy = uArray->GetUnitEnergy();

	if(unitEnergy==0){
	  if(!fProcId){
	    new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
	    fProcId = kTRUE;
	  }
	  fRefArray->Add(uArray);
	  fNDigitEmcal++;
	  ok = kTRUE;
	}

	// Detector flag
	if(unitEnergy>0.)
	  uArray->SetUnitDetectorFlag(kAll);
	else uArray->SetUnitDetectorFlag(kEmcal);
      
	uArray->SetUnitEnergy(unitEnergy+digitEt);

	uArray->SetUnitCutFlag(kPtHigher);

	// To be modified !!!
	uArray->SetUnitSignalFlag(kGood);

	// This is for jet multiplicity
	uArray->SetUnitClusterID(index);
	
	if(fDebug > 12) printf("goodDigit : %d\n", goodDigit);
      
      } // End loop over cells
    } // End loop over clusters
  } // end else

  fNIn += goodDigit;

  if(fDebug>1)
    {
      printf("End of digits %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%");
      printf("goodDigit : %d\n", goodDigit);
    }
}

/*
//____________________________________________________________________________
void AliJetESDFillUnitArrayEMCalDigits::GetCalibrationParameters()
{
  // Set calibration parameters:
  // if calibration database exists, they are read from database,
  // otherwise, they are taken from digitizer.
  //
  // It is a user responsilibity to open CDB before reconstruction,
  // for example:
  // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");

  //Check if calibration is stored in data base

  if(!fCalibData && (AliCDBManager::Instance()->IsDefaultStorageSet()))
    {
      AliCDBEntry *entry = (AliCDBEntry*)
        AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
      if (entry) fCalibData =  (AliEMCALCalibData*) entry->GetObject();
    }

  if(!fCalibData)
    printf("************* Calibration parameters not found in CDB! ****************");
//    AliFatal("Calibration parameters not found in CDB!");


}

//____________________________________________________________________________
Float_t  AliJetESDFillUnitArrayEMCalDigits::Calibrate(Int_t amp, Int_t AbsId)
{

  // Convert digitized amplitude into energy.
  // Calibration parameters are taken from calibration data base for raw data,
  // or from digitizer parameters for simulated data.

  if(fCalibData){

    if (fGeom==0)
      printf("************* Did not get geometry from EMCALLoader ***************");
//      AliFatal("Did not get geometry from EMCALLoader") ;

    Int_t iSupMod = -1;
    Int_t nModule  = -1;
    Int_t nIphi   = -1;
    Int_t nIeta   = -1;
    Int_t iphi    = -1;
    Int_t ieta    = -1;

    Bool_t bCell = fGeom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
    if(!bCell) {
      // fGeom->PrintGeometry();
      Error("Calibrate()"," Wrong cell id number : %i", AbsId);
      assert(0);
    }

    fGeom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);

    fADCchannelECA  = fCalibData->GetADCchannel (iSupMod,ieta,iphi);
    fADCpedestalECA = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);

   return -fADCpedestalECA + amp * fADCchannelECA ;

  }
  else //Return energy with default parameters if calibration is not available
    return -fADCpedestalECA + amp * fADCchannelECA ;

}
*/
Commit	Line	Data
be6e5811	1
	2	/**************************************************************************
	3	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	4	* *
	5	* Author: The ALICE Off-line Project. *
	6	* Contributors are mentioned in the code where appropriate. *
	7	* *
	8	* Permission to use, copy, modify and distribute this software and its *
	9	* documentation strictly for non-commercial purposes is hereby granted *
	10	* without fee, provided that the above copyright notice appears in all *
	11	* copies and that both the copyright notice and this permission notice *
	12	* appear in the supporting documentation. The authors make no claims *
	13	* about the suitability of this software for any purpose. It is *
	14	* provided "as is" without express or implied warranty. *
	15	**************************************************************************/
	16
	17	//---------------------------------------------------------------------
	18	// Fill Unit Array
	19	// Called by ESD reader for jet analysis
	20	// Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
	21	//---------------------------------------------------------------------
	22
	23	// --- ROOT system ---
	24	#include <TLorentzVector.h>
	25	#include <TMath.h>
	26
	27	// --- AliRoot header files ---
	28	#include "AliESDCaloCluster.h"
	29	#include "AliESDCaloCells.h"
	30	#include "AliJetUnitArray.h"
	31	#include "AliJetESDFillUnitArrayEMCalDigits.h"
	32	// #include "AliEMCALCalibData.h"
	33	// #include "AliCDBManager.h"
	34	// // class AliCDBStorage;
	35	// #include "AliCDBEntry.h"
	36
	37	// --- ROOT system ---
	38	class TSystem;
	39	class TGeoManager;
	40	class TArrayS;
	41
	42	// --- AliRoot header files ---
	43	class AliJetFinder;
	44	class AliJetReader;
	45	class AliJetESDReader;
	46	// class AliEMCALCalibData;
	47	// class AliCDBManager;
	48	// class AliCDBStorage;
	49	// class AliCDBEntry;
	50
	51	ClassImp(AliJetESDFillUnitArrayEMCalDigits)
	52
	53	//_____________________________________________________________________________
	54	AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits()
	55	: AliJetFillUnitArray(),
	56	fESD(0),
	57	fNIn(0),
	58	fOpt(0),
	59	fCluster(0),
	60	fDebug(0),
	61	fNCEMCAL(0),
	62	fNCPHOS(0),
	63	fNCCalo(0),
	64	fApplyElectronCorrection(kFALSE),
65	fApplyFractionHadronicCorrection(kFALSE),
66	fFractionHadronicCorrection(0.3),
67	fClus(0x0),
68	fNDigitEmcal(0),
69	fNDigitEmcalCut(0)
70	{
71	// constructor
72	}
73
74	//_____________________________________________________________________________
75	AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits(AliESDEvent *esd)
76	: AliJetFillUnitArray(),
77	fESD(esd),
78	fNIn(0),
79	fOpt(0),
80	fCluster(0),
81	fDebug(0),
82	fNCEMCAL(0),
83	fNCPHOS(0),
84	fNCCalo(0),
85	fApplyElectronCorrection(kFALSE),
86	fApplyFractionHadronicCorrection(kFALSE),
87	fFractionHadronicCorrection(0.3),
88	fClus(0x0),
89	fNDigitEmcal(0),
90	fNDigitEmcalCut(0)
91	{
92	// constructor
93	}
94
95	//_____________________________________________________________________________
96	AliJetESDFillUnitArrayEMCalDigits::AliJetESDFillUnitArrayEMCalDigits(const AliJetESDFillUnitArrayEMCalDigits &det)
97	: AliJetFillUnitArray(det),
98	fESD(det.fESD),
99	fNIn(det.fNIn),
100	fOpt(det.fOpt),
101	fCluster(det.fCluster),
102	fDebug(det.fDebug),
103	fNCEMCAL(det.fNCEMCAL),
104	fNCPHOS(det.fNCPHOS),
105	fNCCalo(det.fNCCalo),
106	fApplyElectronCorrection(det.fApplyElectronCorrection),
107	fApplyFractionHadronicCorrection(det.fApplyFractionHadronicCorrection),
108	fFractionHadronicCorrection(det.fFractionHadronicCorrection),
109	fClus(det.fClus),
110	fNDigitEmcal(det.fNDigitEmcal),
111	fNDigitEmcalCut(det.fNDigitEmcalCut)
112	{
113	// Copy constructor
114	}
115
116	//_____________________________________________________________________________
117	AliJetESDFillUnitArrayEMCalDigits& AliJetESDFillUnitArrayEMCalDigits::operator=(const AliJetESDFillUnitArrayEMCalDigits& other)
118	{
119	// Assignment
120
121	fESD = other.fESD;
122	fNIn = other.fNIn;
123	fOpt = other.fOpt;
124	fCluster = other.fCluster;
125	fDebug = other.fDebug;
126	fNCEMCAL = other.fNCEMCAL;
127	fNCPHOS = other.fNCPHOS;
128	fNCCalo = other.fNCCalo;
129	fApplyElectronCorrection = other.fApplyElectronCorrection;
130	fApplyFractionHadronicCorrection = other.fApplyFractionHadronicCorrection;
131	fFractionHadronicCorrection = other.fFractionHadronicCorrection;
132	fClus = other.fClus;
133	fNDigitEmcal = other.fNDigitEmcal;
134	fNDigitEmcalCut = other.fNDigitEmcalCut;
135
136	return (*this);
137
138	}
139
140	//_____________________________________________________________________________
141	AliJetESDFillUnitArrayEMCalDigits::~AliJetESDFillUnitArrayEMCalDigits()
142	{
143	// destructor
144	}
145
146	//_____________________________________________________________________________
147	void AliJetESDFillUnitArrayEMCalDigits::Exec(Option_t* const /option/)
148	{
149	//
150	// Main method.
151	// Fill the unit array with the neutral particle information from the EMCal cells in ESD
152	//
153
154	fDebug = fReaderHeader->GetDebug();
155	fOpt = fReaderHeader->GetDetector();
156	fCluster = fReaderHeader->GetCluster();
157
158	// Init parameters
159	InitParameters();
160
161	//(not used ?) Int_t nDigitTot = 0;
162	Int_t goodDigit = 0;
163	//(not used ?) Int_t beg = 0;
164	//(not used ?) Int_t end = 0;
165	Int_t index = 0;
166	//(not used ?) Int_t count = 0;
167
168	//(not used ?) Int_t nRefEnt = fRefArray->GetEntries();
169
170	if(!fCluster) { // Keep all digit information
171	// Loop over all cell information
172	//------------------------------------------------------------------
173	AliESDCaloCells &cells= *(fESD->GetEMCALCells());
174	Int_t ncell = cells.GetNumberOfCells() ;
175	//(not used ?) Int_t type = cells.GetType();
176
177	for (Int_t icell= 0; icell < ncell; icell++) {
178	Int_t digitID = cells.GetCellNumber(icell);
179	Float_t digitAmp = cells.GetAmplitude(icell);
180	//(not used ?) Float_t digitTime = cells.GetTime(icell);
181	Float_t digitEn = digitAmp*0.0153; // Last correct
182	// Float_t digitEn = Calibrate(digitAmp,digitID);
183
184	Float_t etaD=-10., phiD=-10.;
185	fGeom->EtaPhiFromIndex(digitID,etaD,phiD);
186	// fEMCalGrid->GetEtaPhiFromIndex2(digitID,phiD,etaD);
187
188	phiD = ((phiD < 0) ? phiD + 2.* TMath::Pi() : phiD);
189
190	Float_t digitEt = digitEn*TMath::Abs(TMath::Sin(EtaToTheta(etaD)));
191
192	AliJetUnitArray uArray = (AliJetUnitArray)fUnitArray->At(digitID);
193
194	if(uArray->GetUnitEnergy() == 0.) goodDigit++;
195	uArray->SetUnitTrackID(digitID);
196
197	Float_t unitEnergy = 0.;
198	Bool_t ok = kFALSE;
199	unitEnergy = uArray->GetUnitEnergy();
200
201	if(unitEnergy==0){
202	if(!fProcId){
203	new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
204	fProcId = kTRUE;
205	}
206	fRefArray->Add(uArray);
207	fNDigitEmcal++;
208	ok = kTRUE;
209	}
210
211	// Detector flag
212	if(unitEnergy>0.)
213	uArray->SetUnitDetectorFlag(kAll);
214	else uArray->SetUnitDetectorFlag(kEmcal);
215
216	uArray->SetUnitEnergy(unitEnergy+digitEt);
217
218	uArray->SetUnitCutFlag(kPtHigher);
219
220	// To be modified !!!
221	uArray->SetUnitSignalFlag(kGood);
222
223	// This is for jet multiplicity
224	uArray->SetUnitClusterID(index);
225
226	if(fDebug > 1) printf("goodDigit : %d\n", goodDigit);
227
228	} // End loop over cells
229	} // end if !fCluster
230	else { // Keep digit information from clusterization
231
232	// Loop over calo clusters
233	//------------------------------------------------------------------
234
235	//select EMCAL clusters only
236	TRefArray * caloClusters = new TRefArray();
237	fESD->GetEMCALClusters(caloClusters);
238
239	// Total number of EMCAL cluster
240	Int_t nclus = caloClusters->GetEntries() ;
241	Int_t beg = 0;
242	Float_t pos[3] ;
243
244	// Get reconstructed vertex position
245	Double_t vertexPosition[3] ;
246	fESD->GetVertex()->GetXYZ(vertexPosition) ;
247
248	// Get CaloCells
249	AliESDCaloCells &cells= *(fESD->GetEMCALCells());
250	//(not used ?) Int_t ncell = cells.GetNumberOfCells() ;
251	//(not used ?) Int_t type = cells.GetType();
252
253	for(Int_t j = beg; j < nclus; j++) { // loop over clusters
254	// Retrieve cluster from esd
08c28025	255	fClus = (AliESDCaloCluster *) caloClusters->At(j) ;
be6e5811	256
	257	// Get the cluster info
	258	//(not used ?) Float_t energy = fClus->E() ;
	259	//(not used ?) Int_t iprim = fClus->GetLabel();
	260
	261	fClus->GetPosition(pos) ;
	262	TVector3 vpos(pos[0],pos[1],pos[2]) ;
	263	TLorentzVector p ;
	264	fClus->GetMomentum(p,vertexPosition);
	265
	266	Int_t digMult = fClus->GetNCells() ;
	267	UShort_t *digID = fClus->GetCellsAbsId() ;
	268	//(not used ?) Double_t *digAmpFrac = fClus->GetCellsAmplitudeFraction() ;
	269	Int_t trackIndex = fClus->GetTrackMatched();
	270
	271	// Do double-counted electron correction
	272	if (fApplyElectronCorrection != 0 && trackIndex !=-1 )
	273	{
	274	// The electron correction go there
	275	// Under construction !!!!
	276	} // End of Electron correction
	277
	278	// Get CaloCells of cluster and fill the unitArray
	279	for(Int_t i = 0; i < digMult ; i++){
	280	Int_t digitID = digID[i]; // or clus->GetCellNumber(i) ;
	281	//(not used ?) Double_t digitAmpFrac = digAmpFrac[i];
	282	Float_t digitAmp = cells.GetCellAmplitude(digitID) ;
	283	//(not used ?) Float_t digitTime = cells.GetCellTime(digitID);
	284
	285	// Calibration for an energy in GeV
	286	Float_t digitEn = digitAmp*0.0153;
	287
	288	Float_t etaD=-10., phiD=-10.;
	289	fGeom->EtaPhiFromIndex(digitID,etaD,phiD);
	290	// fEMCalGrid->GetEtaPhiFromIndex2(digitID,phiD,etaD);
	291
	292	phiD = ((phiD < 0) ? phiD + 2.* TMath::Pi() : phiD);
	293
	294	Float_t digitEt = digitEn*TMath::Abs(TMath::Sin(EtaToTheta(etaD)));
	295
	296	cout << "Digit " << i << ", eta: " << etaD << ", phi: " << phiD << endl;
	297
	298	AliJetUnitArray uArray = (AliJetUnitArray)fUnitArray->At(digitID);
	299	if(uArray->GetUnitEnergy() == 0.) goodDigit++;
	300	uArray->SetUnitTrackID(digitID);
	301
	302	Float_t unitEnergy = 0.;
	303	Bool_t ok = kFALSE;
	304	unitEnergy = uArray->GetUnitEnergy();
	305
	306	if(unitEnergy==0){
	307	if(!fProcId){
	308	new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
	309	fProcId = kTRUE;
	310	}
	311	fRefArray->Add(uArray);
	312	fNDigitEmcal++;
	313	ok = kTRUE;
	314	}
	315
	316	// Detector flag
	317	if(unitEnergy>0.)
	318	uArray->SetUnitDetectorFlag(kAll);
	319	else uArray->SetUnitDetectorFlag(kEmcal);
320
321	uArray->SetUnitEnergy(unitEnergy+digitEt);
322
323	uArray->SetUnitCutFlag(kPtHigher);
324
325	// To be modified !!!
326	uArray->SetUnitSignalFlag(kGood);
327
328	// This is for jet multiplicity
329	uArray->SetUnitClusterID(index);
330
331	if(fDebug > 12) printf("goodDigit : %d\n", goodDigit);
332
333	} // End loop over cells
334	} // End loop over clusters
335	} // end else
336
337	fNIn += goodDigit;
338
339	if(fDebug>1)
340	{
341	printf("End of digits %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%");
342	printf("goodDigit : %d\n", goodDigit);
343	}
344	}
345
346	/*
347	//____________________________________________________________________________
348	void AliJetESDFillUnitArrayEMCalDigits::GetCalibrationParameters()
349	{
350	// Set calibration parameters:
351	// if calibration database exists, they are read from database,
352	// otherwise, they are taken from digitizer.
353	//
354	// It is a user responsilibity to open CDB before reconstruction,
355	// for example:
356	// AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
357
358	//Check if calibration is stored in data base
359
360	if(!fCalibData && (AliCDBManager::Instance()->IsDefaultStorageSet()))
361	{
362	AliCDBEntry entry = (AliCDBEntry)
363	AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
364	if (entry) fCalibData = (AliEMCALCalibData*) entry->GetObject();
365	}
366
367	if(!fCalibData)
368	printf("*********** Calibration parameters not found in CDB! **************");
369	// AliFatal("Calibration parameters not found in CDB!");
370
371
372	}
373
374	//____________________________________________________________________________
375	Float_t AliJetESDFillUnitArrayEMCalDigits::Calibrate(Int_t amp, Int_t AbsId)
376	{
377
378	// Convert digitized amplitude into energy.
379	// Calibration parameters are taken from calibration data base for raw data,
380	// or from digitizer parameters for simulated data.
381
382	if(fCalibData){
383
384	if (fGeom==0)
385	printf("*********** Did not get geometry from EMCALLoader *************");
386	// AliFatal("Did not get geometry from EMCALLoader") ;
387
388	Int_t iSupMod = -1;
389	Int_t nModule = -1;
390	Int_t nIphi = -1;
391	Int_t nIeta = -1;
392	Int_t iphi = -1;
393	Int_t ieta = -1;
394
395	Bool_t bCell = fGeom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
396	if(!bCell) {
397	// fGeom->PrintGeometry();
398	Error("Calibrate()"," Wrong cell id number : %i", AbsId);
399	assert(0);
400	}
401
402	fGeom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
403
404	fADCchannelECA = fCalibData->GetADCchannel (iSupMod,ieta,iphi);
405	fADCpedestalECA = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
406
407	return -fADCpedestalECA + amp * fADCchannelECA ;
408
409	}
410	else //Return energy with default parameters if calibration is not available
411	return -fADCpedestalECA + amp * fADCchannelECA ;
412
413	}
414	*/
415