--- /dev/null
+
+/**************************************************************************
+ * 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. *
+ **************************************************************************/
+
+//======================================================================
+// *** November 2006
+// Author: magali.estienne@ires.in2p3.fr
+// 1) Define 2 grids and take (eta,phi) from the grid or use the grid for the TPC and
+// EtaPhiFromIndex and TowerIndexFromEtaPhi for the particles in EMCAL acceptance
+// 2 options are possible : fGrid==0, work only with full TPC acceptance (for the moment)
+// fGrid==1, work with a part of the TPC acceptance
+// 2) Try to implement 2 full grids for TPC and EMCal separately and to merge them
+// 3) Need to include Dead-zone -> Wait for exact positions in the new detector geometry
+// Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
+//======================================================================
+// ***September 2006
+// TTask : Fill Unit Array for the Tracks information
+// Called by ESD reader for jet analysis
+// Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
+//======================================================================
+// *** July 2006
+// 1) When the tracks are in the EMCal acceptance, the functions EtaPhiFromIndex
+// and TowerIndexFromEtaPhi in the AliEMCALGeometry class are used to extract the
+// index or the eta, phi position of a grid.
+// 2) Define a grid for TPC
+// Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
+//======================================================================
+// ***July 2006
+// Fill Unit Array class
+// Class used by AliJetESDReader to fill a UnitArray from the information extracted
+// from the particle tracks
+// Author: magali.estienne@ires.in2p3.fr
+//======================================================================
+
+
+// --- Standard library ---
+#include <Riostream.h>
+
+// --- ROOT system ---
+#include <TSystem.h>
+#include <TLorentzVector.h>
+#include <TVector3.h>
+//#include "Math/Vector3D.h"
+//#include "Math/Vector3Dfwd.h"
+#include "TTask.h"
+#include <TGeoManager.h>
+#include <TMatrixD.h>
+#include <TArrayD.h>
+
+// --- AliRoot header files ---
+#include "AliJetFinder.h"
+#include "AliJetReaderHeader.h"
+#include "AliJetReader.h"
+#include "AliJetESDReader.h"
+#include "AliJetESDReaderHeader.h"
+#include "AliESD.h"
+#include "AliEMCALGeometry.h"
+#include "AliJetUnitArray.h"
+#include "AliJetFillUnitArrayTracks.h"
+#include "AliJetGrid.h"
+
+ClassImp(AliJetFillUnitArrayTracks)
+
+AliEMCALGeometry *AliJetFillUnitArrayTracks::fGeom=0;
+
+//_____________________________________________________________________________
+AliJetFillUnitArrayTracks::AliJetFillUnitArrayTracks()
+ : TTask("AliJetFillUnitArrayTracks","Fill Unit Array with tpc/its and emcal information")
+{
+ // constructor
+ fNIn = 0;
+ fOpt = 0;
+ fDebug = 0;
+ fNphi = 0;
+ fNeta = 0;
+ fPhiMin = 0;
+ fPhiMax = 0;
+ fEtaMin = 0;
+ fEtaMax = 0;
+ fEtaBinInTPCAcc = 0;
+ fPhiBinInTPCAcc = 0;
+ fEtaBinInEMCalAcc = 0;
+ fPhiBinInEMCalAcc = 0;
+ fNbinPhi = 0;
+}
+
+//____________________________________________________________________________
+void AliJetFillUnitArrayTracks::SetEMCALGeometry()
+{
+ // Set EMCAL geometry information
+ fGeom = AliEMCALGeometry::GetInstance();
+ if (fGeom == 0)
+ fGeom = AliEMCALGeometry::GetInstance("SHISH_77_TRD1_2X2_FINAL_110DEG","EMCAL");
+ if(fDebug>1) printf("\n EMCAL Geometry setted ! \n");
+}
+
+//____________________________________________________________________________
+void AliJetFillUnitArrayTracks::InitParameters()
+{
+ fHCorrection = 0; // For hadron correction
+ fHadCorr = 0; // For hadron correction
+ fNumUnits = fGeom->GetNCells(); // Number of towers in EMCAL
+ fDebug = fReaderHeader->GetDebug();
+
+ fEtaMinCal = fGeom->GetArm1EtaMin();
+ fEtaMaxCal = fGeom->GetArm1EtaMax();
+ fPhiMinCal = fGeom->GetArm1PhiMin();
+ fPhiMaxCal = fGeom->GetArm1PhiMax()-10.; // A verifier quelle doit etre la derniere valeur !
+
+ if(fDebug>30){
+ cout << "fEtaMinCal : " << fEtaMinCal << endl;
+ cout << "fEtaMaxCal : " << fEtaMaxCal << endl;
+ cout << "fPhiMinCal : " << fPhiMinCal << endl;
+ cout << "fPhiMaxCal : " << fPhiMaxCal << endl;
+ }
+
+ fTPCGrid->GetAccParam(fNphi,fNeta,fPhiMin,
+ fPhiMax,fEtaMin,fEtaMax);
+ fTPCGrid->GetBinParam(fPhiBinInTPCAcc,fEtaBinInTPCAcc,
+ fPhiBinInEMCalAcc,fEtaBinInEMCalAcc,fNbinPhi);
+
+ fEta = fTPCGrid->GetArrayEta();
+ fPhi = fTPCGrid->GetArrayPhi();
+ fIndex = fTPCGrid->GetIndexObject();
+
+ if(fDebug>20){
+ for(Int_t i=0; i<fNphi+1; i++) cout << "phi[" << i << "] : " << (*fPhi)[i] << endl;
+ for(Int_t i=0; i<fNeta+1; i++) cout << "eta[" << i << "] : " << (*fEta)[i] << endl;
+
+ for(Int_t i=0; i<fNphi+1; i++)
+ for(Int_t j=0; j<fNeta+1; j++) {cout << "fIndex[" << i << "," << j << "] : " <<
+ (*fIndex)(i,j) << endl; }
+ }
+ if(fDebug>1) printf("\n Parameters initiated ! \n");
+}
+
+//_____________________________________________________________________________
+AliJetFillUnitArrayTracks::~AliJetFillUnitArrayTracks()
+{
+ // destructor
+}
+
+//_____________________________________________________________________________
+void AliJetFillUnitArrayTracks::Exec(Option_t* option)
+{
+ //
+ // Main method.
+ // Explain
+
+ fDebug = fReaderHeader->GetDebug();
+ if(fDebug>1) printf("In AliJetFillUnitArrayTracks !");
+ if(fDebug>3) printf("\nfGeom->GetEntries() = %d\n", fGeom->GetNCells());
+ // Set EMCal Geometry
+ SetEMCALGeometry();
+ // Set parameters
+ InitParameters();
+
+ TClonesArray *lvArray = fMomentumArray; // Correct checked !
+ Int_t nInT = lvArray->GetEntries(); // Correct checked !
+ Float_t ptMin = fReaderHeader->GetPtCut(); // Correct checked !
+
+ // sflag -> not yet implemented !!!!
+
+ if(fDebug>3) cout << "nInT : " << nInT << endl;
+
+ if (nInT == 0) return;
+
+ // local arrays for input
+ Float_t* enT = new Float_t[nInT];
+ Float_t* ptT = new Float_t[nInT];
+ Float_t* etaT = new Float_t[nInT];
+ Float_t* phiT = new Float_t[nInT];
+ Float_t* thetaT = new Float_t[nInT];
+
+ Int_t trackInEmcalAcc = 0;
+ Int_t trackInTpcAcc = 0;
+ Int_t trackInTpcAccOnly = 0;
+
+ Int_t nElements = fTPCGrid->GetNEntries();
+ Int_t nElements2 = fEMCalGrid->GetNEntries();
+ fGrid = fTPCGrid->GetGridType();
+
+ if(fDebug>3){
+ cout << "nElements : " << nElements << endl;
+ cout << "nElements2 : " << nElements2 << endl;
+ cout << "fNumUnits : " << fNumUnits << endl;
+ cout << "sum : " << fNumUnits+nElements << endl;
+ }
+
+ // Set energy exactly to zero
+ if(fGrid==0)
+ for(Int_t k=0; k<nElements; k++)
+ fUnitArray[k].SetUnitEnergy(0.);
+
+ if(fGrid==1)
+ for(Int_t k=0; k<fNumUnits+nElements; k++)
+ fUnitArray[k].SetUnitEnergy(0.);
+
+ // load input vectors
+ for (Int_t i = 0; i < nInT; i++)
+ {
+ TLorentzVector *lv = (TLorentzVector*) lvArray->At(i);
+ enT[i] = lv->Energy();
+ ptT[i] = lv->Pt();
+ etaT[i] = lv->Eta();
+ phiT[i] = ((lv->Phi() < 0) ? (lv->Phi()) + 2. * TMath::Pi() : lv->Phi());
+ thetaT[i] = 2.0*TMath::ATan(TMath::Exp(-etaT[i]));
+
+ if(fDebug>20){
+ cout << "enT[" << i << "] : " << enT[i] << endl;
+ cout << "ptT[" << i << "] : " << ptT[i] << endl;
+ cout << "etaT[" << i << "] : " << etaT[i] << endl;
+ cout << "phiT[" << i << "] : " << phiT[i] << endl;
+ cout << "thetaT[" << i << "] : " << thetaT[i] << endl;
+ cout << "fEtaMinCal : " << fEtaMinCal << ", fEtaMaxCal : " << fEtaMaxCal << endl;
+ cout << "fPhiMinCal : " << fPhiMinCal << ", fPhiMaxCal : " << fPhiMaxCal << endl;
+ cout << "fEtaMin : " << fEtaMin << ", fEtaMax : " << fEtaMax << endl;
+ cout << "fPhiMin : " << fPhiMin << ", fPhiMax : " << fPhiMax << endl;
+ }
+
+ if(fGrid==0)
+ {
+ // For the moment, only TPC filled from its grid in its total acceptance
+ if(fDebug>2)
+ cout << "In total TPC acceptance +++++++++++++++++++++++++++++++++++++++++++" << endl;
+
+ trackInTpcAccOnly += 1;
+
+ Int_t idTPC = fTPCGrid->GetIndex(phiT[i],etaT[i]);
+
+ Float_t unitEnergy = 0.;
+ unitEnergy = fUnitArray[idTPC].GetUnitEnergy();
+
+ if(unitEnergy > 0. && fDebug >10){
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "TPC id : " << idTPC << endl;
+ cout << "etaT[" << i << "] : " << etaT[i] << endl;
+ cout << "phiT[" << i << "] : " << phiT[i] << endl;
+ cout << "unitEnergy in TPC acceptance : " << unitEnergy << endl;
+ cout << "fUnitArray[idTPC].GetUnitEnergy(): " <<
+ fUnitArray[idTPC].GetUnitEnergy() << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ }
+
+ // Fill energy in TPC acceptance
+ fUnitArray[idTPC].SetUnitEnergy(unitEnergy + ptT[i]);
+
+ if(fDebug > 10){
+ cout << "ptT[" << i << "] : " << ptT[i] << endl;
+ cout << "unitEnergy in TPC acceptance after : " <<
+ fUnitArray[idTPC].GetUnitEnergy() << endl;
+ }
+
+ // Pt cut flag
+ if(fUnitArray[idTPC].GetUnitEnergy()<ptMin)
+ fUnitArray[idTPC].SetUnitCutFlag(kPtSmaller);
+ else fUnitArray[idTPC].SetUnitCutFlag(kPtHigher);
+ // Detector flag
+ fUnitArray[idTPC].SetUnitDetectorFlag(kTpc);
+ }
+
+ if(fGrid==1)
+ {
+ // Fill track information in EMCAL acceptance
+ if((etaT[i] >= fEtaMin && etaT[i] <= fEtaMax) &&
+ (phiT[i] >= fPhiMin && phiT[i] <= fPhiMax))// &&
+ // GetCutFlag(i) == 1)
+ // ptT[i] > ptMin)
+ {
+ trackInEmcalAcc += 1;
+
+ if(fDebug>20){
+ cout << "before : " << endl;
+ cout << "etaT[i] : " << etaT[i] << endl;
+ cout << "phiT[i] : " << phiT[i]*180/TMath::Pi() << endl;
+ }
+
+ // This function should be modified soon
+// Int_t towerID = fGeom->TowerIndexFromEtaPhi(etaT[i],180.0/TMath::Pi()*phiT[i]); // Obsolete
+ Int_t towerID = fGeom->TowerIndexFromEtaPhi2(etaT[i],180.0/TMath::Pi()*phiT[i]); // Mine modified
+// Int_t towerID = fEMCalGrid->GetIndexFromEtaPhi(phiT[i],etaT[i]); // Using an EMCal grid -> calculated
+// Int_t towerID = fEMCalGrid->GetIndex(phiT[i],etaT[i]); // Using an EMCal grid -> tabulated (faster)
+
+ Float_t unitEnergy = fUnitArray[towerID].GetUnitEnergy();
+
+ if(fDebug>20) {
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "after : " << endl;
+ cout << "towerID : " << towerID << endl;
+ cout << "etaT[i] : " << etaT[i] << endl;
+ cout << "phiT[i](rad) : " << phiT[i] << endl;
+ cout << "phiT[i] : " << phiT[i]*180/TMath::Pi() << endl;
+ cout << "unitEnergy in emcal acceptance : " << unitEnergy << endl;
+ cout << "fHadCorr : " << fHadCorr << endl;
+ cout << "fHCorrection : " << fHCorrection << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ }
+
+ //OLD WAY: //Do Hadron Correction
+ if (fHCorrection != 0)
+ {
+ Float_t hCEnergy = fHadCorr->GetEnergy(enT[i], (Double_t)etaT[i]);
+ unitEnergy -= hCEnergy*TMath::Sin(thetaT[i]);
+ if(fDebug>20){
+ cout << "Inside loop for hadron correction ==========================" << endl;
+ cout << "enT[" << i << "] : " << enT[i] << endl;
+ cout << "ptT[" << i << "] : " << ptT[i] << endl;
+ cout << "etaT[" << i << "] : " << etaT[i] << endl;
+ cout << "Energy correction : " << hCEnergy << endl;
+ cout << "unitEnergy : " << unitEnergy << endl;
+ cout << "fUnitArray[towerID].GetUnitEnergy() : " <<
+ fUnitArray[towerID].GetUnitEnergy() << endl;
+ }
+ } //end Hadron Correction loop
+
+ fUnitArray[towerID].SetUnitEnergy(unitEnergy + ptT[i]);
+
+ // Put a pt cut flag
+ if(fUnitArray[towerID].GetUnitEnergy()<ptMin)
+ fUnitArray[towerID].SetUnitCutFlag(kPtSmaller);
+ else fUnitArray[towerID].SetUnitCutFlag(kPtHigher);
+ // Detector flag
+ fUnitArray[towerID].SetUnitDetectorFlag(kTpc);
+
+ if(fDebug>10){
+ cout << "After pT filled ===============================================" << endl;
+ cout << "PtCut : " << ptMin << endl;
+ cout << "ptT[" << i << "] : " << ptT[i] << endl;
+ cout << "unitEnergy : " << unitEnergy << endl;
+ cout << "fUnitArray[towerID].GetUnitEnergy() : " << fUnitArray[towerID].GetUnitEnergy() << endl;
+ cout << "fUnitArray[towerID].GetUnitCutFlag() : " << fUnitArray[towerID].GetUnitCutFlag() << endl;
+ cout << "fUnitArray[towerID].GetUnitEta() : " << fUnitArray[towerID].GetUnitEta() << endl;
+ cout << "fUnitArray[towerID].GetUnitPhi() : " << fUnitArray[towerID].GetUnitPhi() << endl;
+ }
+ } // end loop on EMCal acceptance cut + tracks quality
+ else{
+ // Outside EMCal acceptance
+ if(fDebug>2)
+ cout << "Outside EMCal acceptance +++++++++++++++++++++++++++++++++++++++++++" << endl;
+
+ trackInTpcAcc += 1;
+
+ // Int_t idTPC = GetIndexFromEtaPhi(etaT[i],phiT[i]);
+ Int_t idTPC = fTPCGrid->GetIndex(phiT[i],etaT[i]);
+
+ Float_t unitEnergy2 = 0.;
+ unitEnergy2 = fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy(); // check if fNumUnits or fNumUnits-1
+
+ if(unitEnergy2 > 0. && fDebug >10){
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "TPC id : " << idTPC << endl;
+ cout << "Total id : " << fNumUnits-1+idTPC << endl;
+ cout << "etaT[" << i << "] : " << etaT[i] << endl;
+ cout << "phiT[" << i << "] : " << phiT[i] << endl;
+ cout << "unitEnergy outside emcal acceptance : " << unitEnergy2 << endl;
+ cout << "fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy(): " <<
+ fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy() << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
+ }
+
+ // Fill energy outside emcal acceptance
+ fUnitArray[fNumUnits-1+idTPC].SetUnitEnergy(unitEnergy2 + ptT[i]);
+
+ // Pt cut flag
+ if(fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy()<ptMin)
+ fUnitArray[fNumUnits-1+idTPC].SetUnitCutFlag(kPtSmaller);
+ else fUnitArray[fNumUnits-1+idTPC].SetUnitCutFlag(kPtHigher);
+ // Detector flag
+ fUnitArray[fNumUnits-1+idTPC].SetUnitDetectorFlag(kTpc);
+ }
+ } // end fGrid==1
+ } // end loop on entries (tpc tracks)
+
+ if(fDebug>3){
+ printf("Number of tracks in EMCAL acceptance: %d\n", trackInEmcalAcc);
+ printf("Number of tracks in TPC acceptance: %d\n", trackInTpcAcc);
+ }
+
+ if(fGrid==0) fNIn = trackInTpcAccOnly;
+ if(fGrid==1) fNIn = trackInEmcalAcc+trackInTpcAcc;
+
+ fOpt = fReaderHeader->GetDetector();
+
+ if(fDebug>1) printf("fOpt in FillUnitArrayFromTPCTracks : %d\n", fOpt);
+
+ if(fOpt==1 || option=="tpc") // if only TPC
+ { //Set all unit flags, Eta, Phi
+
+ if(fGrid==0)
+ {
+ if(fDebug>1) printf("In if(fOpt==1 || option==tpc)\n");
+ for(Int_t i=0; i<nElements; i++)
+ {
+ Float_t eta = 0.;
+ Float_t phi = 0.;
+ if (fDebug>10) Info("FillUnitArray","Setting all units outside jets \n");
+ fUnitArray[i].SetUnitFlag(kOutJet); // returns 0, 1, 2...
+ fUnitArray[i].SetUnitEntries(nElements);
+ fUnitArray[i].SetUnitID(i);
+ fTPCGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID(),phi,eta);
+ fUnitArray[i].SetUnitEta(eta);
+ fUnitArray[i].SetUnitPhi(phi);
+ fUnitArray[i].SetUnitDeta(fTPCGrid->GetDeta());
+ fUnitArray[i].SetUnitDphi(fTPCGrid->GetDphi());
+ }
+ }
+
+ if(fGrid==1)
+ {
+ if(fDebug>1) printf("In if(fOpt==1 || option==tpc)\n");
+ for(Int_t i=0; i<fNumUnits+nElements; i++)
+ {
+ Float_t eta = 0.;
+ Float_t phi = 0.;
+ if (fDebug>10) Info("FillUnitArray","Setting all units outside jets \n");
+ //Set all units to be outside a jet initially
+ fUnitArray[i].SetUnitFlag(kOutJet); // returns 0, 1, 2...
+ fUnitArray[i].SetUnitEntries(fNumUnits+nElements);
+ fUnitArray[i].SetUnitID(i);
+
+ if(fUnitArray[i].GetUnitID()<fNumUnits)
+ {
+ // fGeom->EtaPhiFromIndex2(fUnitArray[i].GetUnitID(), eta, phi); // My function in HEADPythia63
+ fGeom->EtaPhiFromIndex(fUnitArray[i].GetUnitID(), eta, phi); // From EMCal geometry
+ phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
+ // fEMCalGrid->GetEtaPhiFromIndex2(i,phi,eta); // My function from Grid
+ // fEMCalGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID(),phi,eta); // My function from Grid
+ fUnitArray[i].SetUnitEta(eta);
+ //fUnitArray[i].SetUnitPhi(phi*TMath::Pi()/180.0);
+ fUnitArray[i].SetUnitPhi(phi);
+ fUnitArray[i].SetUnitDeta(fEMCalGrid->GetDeta());
+ fUnitArray[i].SetUnitDphi(fEMCalGrid->GetDphi());
+ }
+ else {
+ fTPCGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID()+1-fNumUnits,phi,eta);
+ fUnitArray[i].SetUnitEta(eta);
+ fUnitArray[i].SetUnitPhi(phi);
+ fUnitArray[i].SetUnitDeta(fTPCGrid->GetDeta());
+ fUnitArray[i].SetUnitDphi(fTPCGrid->GetDphi());
+
+ if(fDebug>10)
+ {
+ if(fUnitArray[i].GetUnitEnergy()!=0.){
+ cout << "(fUnitArray[" << i << "].GetUnitID()+1-fNumUnits : " <<
+ fUnitArray[i].GetUnitID()+1-fNumUnits << endl;
+ cout << "(fUnitArray[" << i << "].GetUnitEnergy() : " <<
+ fUnitArray[i].GetUnitEnergy() << endl;
+ cout << "(fUnitArray[" << i << "].GetUnitEta() : " <<
+ fUnitArray[i].GetUnitEta() << endl;
+ cout << "(fUnitArray[" << i << "].GetUnitPhi() : " <<
+ fUnitArray[i].GetUnitPhi() << endl;
+ }
+ }
+ }
+ fUnitArray[i].SetUnitClusterID(0);
+ }//end loop over all units in array (same as all towers in EMCAL)
+ }
+
+ if(fDebug>20)
+ {
+ for(Int_t i=0; i<fNumUnits+nElements; i++)
+ {
+ if(fUnitArray[i].GetUnitEnergy()!=0.){
+ cout << "######################################################### " << endl;
+ cout << "Final UnitArray filled with energy != 0" << i << endl;
+ cout << "Pointeur UnitArray : " << fUnitArray << " ID : " <<
+ fUnitArray[i].GetUnitID() << " Energy : " << fUnitArray[i].GetUnitEnergy() <<
+ " Eta : " << fUnitArray[i].GetUnitEta() << " Phi : " << fUnitArray[i].GetUnitPhi() << endl;
+ }
+ }
+ }
+
+ } // end if(fOpt==1 || option=="tpc")
+
+ delete[] enT;
+ delete[] ptT;
+ delete[] etaT;
+ delete[] phiT;
+ delete[] thetaT;
+
+}
+
+//__________________________________________________________
+void AliJetFillUnitArrayTracks::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
+{
+ for(Int_t j=0; j<fNphi+1; j++) {
+ for(Int_t i=0; i<fNeta+1; i++) {
+
+ // TPC grid only
+ //-------------------------------------
+ if(fGrid==0) {
+ if(j*(fNeta+1)+i == index) {
+ eta = fEta2->At(i);
+ phi = fPhi2->At(j);
+ }
+ }
+
+ // TPC-EMCAL grid
+ //-------------------------------------
+ Int_t ii = 0;
+ if(i==0) ii = 0;
+ if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
+ if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
+
+ if(fGrid==1) {
+ if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
+ eta = fEta2->At(i);
+ phi = fPhi2->At(j);
+ }
+
+ if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
+ ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
+ if(ii==0) {Int_t ind = 0; eta = fEta2->At(ind);}
+ else eta = fEta2->At(i);
+ phi = fPhi2->At(j);
+ }
+
+ if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) && ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))+(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
+ eta = fEta2->At(i);
+ phi = fPhi2->At(j);
+ }
+ }
+ }
+ }
+}
+
+
+
+
+
+
+
+
+
+
+
+
--- /dev/null
+#ifndef ALIJETFILLUNITARRAYTRACKS_H
+#define ALIJETFILLUNITARRAYTRACKS_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//---------------------------------------------------------------------
+// Jet Fill Unit Array
+// Called by ESD Reader for jet analysis
+// Author: Magali Estienne (magali.estienne@subatech.in2p3.fr)
+//---------------------------------------------------------------------
+
+#ifndef ROOT_TTask
+#include "TTask.h"
+#endif
+
+#include <TMatrixD.h>
+#include <TArrayD.h>
+
+class AliEMCALGeometry;
+class AliJetHadronCorrection;
+class AliJetReader;
+class AliJetESDReader;
+class TClonesArray;
+class AliJetUnitArray;
+//class AliJetReaderHeader;
+class AliJetReader;
+class AliJetGrid;
+
+class AliJetFillUnitArrayTracks : public TTask
+{
+ public:
+ AliJetFillUnitArrayTracks();
+ virtual ~AliJetFillUnitArrayTracks();
+
+ // Setter
+ void SetReaderHeader(AliJetReaderHeader *readerHeader) {fReaderHeader = readerHeader;}
+ void SetMomentumArray(TClonesArray *momentumArray) {fMomentumArray = momentumArray;}
+ void SetUnitArray(AliJetUnitArray *unitArray) {fUnitArray = unitArray;}
+ void SetHadCorrection(Int_t flag = 1) {fHCorrection = flag;}
+ void SetHadCorrector(AliJetHadronCorrectionv1* corr) {fHadCorr = corr;}
+ void SetTPCGrid(AliJetGrid *grid) {fTPCGrid = grid;}
+ void SetEMCalGrid(AliJetGrid *grid) {fEMCalGrid = grid;}
+ void SetGrid(Double_t phiMin,Double_t phiMax,Double_t etaMin,Double_t etaMax);
+
+ // Getter
+ AliJetUnitArray* GetUnitArray() {return fUnitArray;}
+ // Int_t GetIndexFromEtaPhi(Double_t eta,Double_t phi) const;
+ void GetEtaPhiFromIndex(Int_t index,Float_t &eta,Float_t &phi);
+ Int_t GetNeta() {return fNeta;}
+ Int_t GetNphi() {return fNphi;}
+
+ void Exec(Option_t*);
+
+ protected:
+ Int_t fNumUnits; // Number of units in the unit object array (same as num towers in EMCAL)
+ Float_t fEtaMinCal; // Define EMCal acceptance in Eta
+ Float_t fEtaMaxCal; // Define EMCal acceptance in Eta
+ Float_t fPhiMinCal; // Define EMCal acceptance in Phi
+ Float_t fPhiMaxCal; // Define EMCal acceptance in Phi
+ AliJetHadronCorrectionv1 *fHadCorr; // Pointer to Hadron Correction Object
+ Int_t fHCorrection; // Hadron correction flag
+ Int_t fNIn; // Number of Array filled in UnitArray
+ Int_t fOpt; // Detector to be used for jet reconstruction
+ Int_t fDebug; // Debug option
+
+ AliJetReaderHeader *fReaderHeader; // ReaderHeader
+ TClonesArray *fMomentumArray; // MomentumArray
+ AliJetUnitArray *fUnitArray; // UnitArray
+ AliJetGrid *fTPCGrid; // Define filled grid
+ AliJetGrid *fEMCalGrid; // Define filled grid
+
+ // geometry info
+ static AliEMCALGeometry *fGeom; //!
+
+ Int_t fNphi; // number of points in the grid: phi
+ Int_t fNeta; // " eta
+ TArrayD* fPhi2; // grid points in phi
+ TArrayD* fEta2; // grid points in eta
+ TArrayD* fPhi; // grid points in phi
+ TArrayD* fEta; // grid points in eta
+ TMatrixD* fIndex; // grid points in (phi,eta)
+ TMatrixD* fParams; // matrix of parameters in the grid points
+ Int_t fGrid; // Select the grid acceptance you want to fill
+ // 0 = TPC acceptance, 1 = TPC-EMCal acceptance
+ Float_t fPhiMin;
+ Float_t fPhiMax;
+ Float_t fEtaMin;
+ Float_t fEtaMax;
+ Int_t fEtaBinInTPCAcc;
+ Int_t fPhiBinInTPCAcc;
+ Int_t fEtaBinInEMCalAcc;
+ Int_t fPhiBinInEMCalAcc;
+ Int_t fNbinPhi;
+
+ private:
+ void SetEMCALGeometry();
+ void InitParameters();
+
+ // Int_t fEvent;
+
+ ClassDef(AliJetFillUnitArrayTracks,1) // Fill Unit Array with tpc and/or emcal information
+};
+
+#endif
--- /dev/null
+/**************************************************************************
+ * 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. *
+ **************************************************************************/
+
+
+/* $Id$ */
+
+//=========================================================================
+// Modified class for JETAN
+// Some flags have been changed
+// Some are not used for the moment
+// Author: magali.estienne@ires.in2p3.fr
+//=========================================================================
+// Enumerated types for use in JetFinder classes
+//
+//*-- Author: Mark Horner (LBL/UCT)
+//
+
+#ifndef ALIJETFINDERBGCALCTYPE_T
+#define ALIJETFINDERBGCALCTYPE_T
+
+ typedef enum { kRatio, kCone, kConstant
+ } AliJetFinderBGCalcType_t;
+#endif
+
+#ifndef ALIJETFINDERRESETTYPE_T
+#define ALIJETFINDERRESETTYPE_T
+
+ typedef enum { kResetData, kResetTracks, kResetDigits, kResetParameters,
+ kResetAll, kResetPartons, kResetParticles, kResetJets
+ } AliJetFinderResetType_t;
+#endif
+
+#ifndef ALIJETFINDERTRACKTYPE_T
+#define ALIJETFINDERTRACKTYPE_T
+ typedef enum { kAllP, kEM, kCharged, kNeutral, kHadron, kChargedHadron, kNoTracks, kEMChargedPi0, kNoNeutronNeutrinoKlong
+ } AliJetFinderTrackType_t;
+#endif
+
+#ifndef ALIJETFINDERSMEARINGTYPE_T
+#define ALIJETFINDERSMEARINGTYPE_T
+ typedef enum { kSmear, kEfficiency , kSmearEffic, kPerfectTracks
+ } AliJetFinderSmearingType_t;
+#endif
+
+#ifndef ALIJETFINDEREMCALTYPE_T
+#define ALIJETFINDEREMCALTYPE_T
+typedef enum { kHits, kTimeCut,kNoHits, kDigits, kClusters
+ } AliJetFinderEmcalType_t;
+#endif
+
+#ifndef ALIJETFINDERFILETYPE_T
+#define ALIJETFINDERFILETYPE_T
+ typedef enum { kHijing,kPythia,kData
+ } AliJetFinderFileType_t;
+#endif
+
+#ifndef ALIJETFINDERUA1UNITFLAGTYPE_T
+#define ALIJETFINDERUA1UNITFLAGTYPE_T
+ typedef enum { kInCurrentJet, kInJet, kOutJet, kBelowMinEt
+ } AliJetFinderUnitFlagType_t;
+#endif
+
+#ifndef ALIJETFINDERUA1UNITCUTFLAGTYPE_T
+#define ALIJETFINDERUA1UNITCUTFLAGTYPE_T
+ typedef enum { kPtSmaller, kPtHigher
+ } AliJetFinderUnitCutFlagType_t;
+#endif
+
+#ifndef ALIJETFINDERUA1UNITSIGNALFLAGTYPE_T
+#define ALIJETFINDERUA1UNITSIGNALFLAGTYPE_T
+ typedef enum { kGood, kBad
+ } AliJetFinderUnitSignalFlagType_t;
+#endif
+
+#ifndef ALIJETFINDERUNITDETECTORFLAGTYPE_T
+#define ALIJETFINDERUNITDETECTORFLAGTYPE_T
+ typedef enum { kTpc, kEmcal, kAll
+ } AliJetFinderUnitDetectorFlagType_t;
+#endif
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 2001-2002, 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. *
+ **************************************************************************/
+
+/* $Id: AliJetGrid.cxx,v 1.0 07/09/2006 */
+
+//=========================================================================
+// *** 07 September 2006
+// This class allows to fill a really general grid in (eta,phi)
+// Two types of grid can be setted :
+// if fGrid==0 -> the complete acceptance of a rectangular detector acceptance is filled
+// if fGrid==1 -> the previous grid - an other rectangular detector acceptance is filled
+// A (eta,phi(rad)) position can be extracted from an index value without looping over all entries
+// An index position can be extracted from a (eta,phi(rad)) position without looping over all entries
+//
+// How to run it :
+// > AliJetGrid *grid = new AliJetGrid((NbinPhi-1),(NbinEta-1),PhiMin,PhiMax,EtaMin,EtaMax);
+// > grid->SetGridType(...); // 0 or 1 for the moment
+// > grid->SetMatrixIndexes();
+// > grid->SetIndexIJ();
+//
+// Author : magali.estienne@ires.in2p3.fr
+//=========================================================================
+
+// Standard headers
+#include <Riostream.h>
+// Root headers
+#include <TMatrixD.h>
+#include <TArrayD.h>
+#include <TArrayI.h>
+// AliRoot headers
+#include "AliJetGrid.h"
+
+ClassImp(AliJetGrid)
+
+//__________________________________________________________
+AliJetGrid::AliJetGrid() {
+
+ // Default constructor
+
+ fNphi = 0;
+ fNeta = 0;
+ fPhi = 0;
+ fEta = 0;
+ fIndex = 0x0;
+ fPhiMin = 0.; // acceptance removed inside
+ fPhiMax = 0.; // acceptance removed inside
+ fEtaMin = 0.; // acceptance removed inside
+ fEtaMax = 0.; // acceptance removed inside
+ fMinPhi = 0.; // total acceptance
+ fMaxPhi = 0.; // total acceptance
+ fMinEta = 0.; // total acceptance
+ fMaxEta = 0.; // total acceptance
+ fEtaBinInTPCAcc = 0;
+ fPhiBinInTPCAcc = 0;
+ fEtaBinInEMCalAcc = 0;
+ fPhiBinInEMCalAcc = 0;
+ fNbinPhi = 0;
+
+ fDebug = 1;
+}
+
+//__________________________________________________________
+AliJetGrid::AliJetGrid(Int_t nphi,Int_t neta,Double_t phiMin,Double_t phiMax,Double_t etaMin,Double_t etaMax) {
+
+ // Standard constructor
+
+ fNphi = nphi;
+ fNeta = neta;
+ fPhiMin = 0.; // rad - acceptance removed inside
+ fPhiMax = 0.; // rad - acceptance removed inside
+ fEtaMin = 0.; // acceptance removed inside
+ fEtaMax = 0.; // acceptance removed inside
+ fNbinPhi = 0;
+
+ fMaxPhi = phiMax; // total acceptance
+ fMinPhi = phiMin; // total acceptance
+ fMaxEta = etaMax; // total acceptance
+ fMinEta = etaMin; // total acceptance
+
+ fPhi = new TArrayD(fNphi+1);
+ fEta = new TArrayD(fNeta+1);
+ fIndexI = new TArrayI((fNeta+1)*(fNphi+1)+1);
+ fIndexJ = new TArrayI((fNeta+1)*(fNphi+1)+1);
+
+ for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = (phiMax-phiMin)/fNphi*i+phiMin;
+ for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = (etaMax-etaMin)/fNeta*i+etaMin;
+
+ if(fDebug > 3){
+ for(Int_t i=0; i<(fNphi+1); i++) cout << (*fPhi)[i] << endl;
+ for(Int_t i=0; i<(fNeta+1); i++) cout << (*fEta)[i] << endl;
+ }
+
+ fIndex = new TMatrixD(fNphi+1,fNeta+1);
+
+}
+
+//__________________________________________________________
+AliJetGrid::AliJetGrid(const AliJetGrid& grid):TNamed(grid) {
+
+ // Copy constructor
+
+ fNphi = grid.fNphi;
+ fNeta = grid.fNeta;
+ fPhiMin = grid.fPhiMin;
+ fPhiMax = grid.fPhiMax;
+ fEtaMin = grid.fEtaMin;
+ fEtaMax = grid.fEtaMax;
+ fEtaBinInTPCAcc = grid.fEtaBinInTPCAcc;
+ fPhiBinInTPCAcc = grid.fPhiBinInTPCAcc;
+ fEtaBinInEMCalAcc = grid.fEtaBinInEMCalAcc;
+ fPhiBinInEMCalAcc = grid.fPhiBinInEMCalAcc;
+ fNbinPhi = grid.fNbinPhi;
+ fMaxPhi = grid.fMaxPhi;
+ fMinPhi = grid.fMinPhi;
+ fMaxEta = grid.fMaxEta;
+ fMinEta = grid.fMinEta;
+
+ fPhi = new TArrayD(fNphi+1);
+ for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = grid.fPhi->At(i);
+ fEta = new TArrayD(fNeta+1);
+ for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = grid.fEta->At(i);
+
+ fIndex = new TMatrixD(fNphi+1,fNeta+1);
+ for(Int_t i=0; i<fNphi+1; i++) {
+ for(Int_t j=0; j<fNeta+1; j++) (*fIndex)(i,j)=(*grid.fIndex)(i,j);
+ }
+}
+
+//__________________________________________________________
+AliJetGrid::~AliJetGrid() {
+
+ // Destructor
+ delete fPhi;
+ delete fEta;
+ delete fIndexI;
+ delete fIndexJ;
+ delete fIndex;
+}
+
+//__________________________________________________________
+void AliJetGrid::InitParams(Double_t phiMinCal,Double_t phiMaxCal,Double_t etaMinCal,Double_t etaMaxCal)
+{ // To set initial parameters
+
+ fPhiMin = phiMinCal; // rad
+ fPhiMax = phiMaxCal; // rad
+ fEtaMin = etaMinCal;
+ fEtaMax = etaMaxCal;
+ fNbinPhi = static_cast<int>(fPhiMin/(2*TMath::Pi()/fNphi));
+
+ // Define some binning numbers
+ if(fGrid==0){
+ for(Int_t i=0; i<fNphi+1; i++) fPhiBinInTPCAcc++;
+ fPhiBinInEMCalAcc = 0;
+
+ for(Int_t i=0; i<fNeta+1; i++) fEtaBinInTPCAcc++;
+ fEtaBinInEMCalAcc = 0;
+ }
+
+ if(fGrid==1){
+ for(Int_t i=0; i<fNphi+1; i++) {
+ fPhiBinInTPCAcc++;
+ if(fPhi->At(i) >= fPhiMin &&
+ fPhi->At(i) <= fPhiMax)
+ fPhiBinInEMCalAcc++;
+ }
+ for(Int_t i=0; i<fNeta+1; i++) {
+ fEtaBinInTPCAcc++;
+ if((fEta->At(i) >= fEtaMin) &&
+ (fEta->At(i) <= fEtaMax))
+ fEtaBinInEMCalAcc++;
+ }
+ }
+
+}
+
+//__________________________________________________________
+TArrayD* AliJetGrid::GetArrayEta()
+{ // Returns an array with the eta points
+
+ return fEta;
+
+}
+
+//__________________________________________________________
+TArrayD* AliJetGrid::GetArrayPhi()
+{ // Returns an array with the phi points
+
+ return fPhi;
+
+}
+
+//__________________________________________________________
+TMatrixD* AliJetGrid::GetIndexObject()
+{ // Returns a pointer to the matrix
+
+ return fIndex;
+
+}
+
+//__________________________________________________________
+void AliJetGrid::GetAccParam(Int_t &nphi, Int_t &neta, Float_t &minphi, Float_t &maxphi,
+ Float_t &mineta, Float_t &maxeta)
+{ // Returns EMCAL acceptance initially setted
+
+ nphi = fNphi;
+ neta = fNeta;
+ minphi = fPhiMin;
+ maxphi = fPhiMax;
+ mineta = fEtaMin;
+ maxeta = fEtaMax;
+}
+
+//__________________________________________________________
+void AliJetGrid::GetBinParam(Int_t &phibintpc, Int_t &etabintpc,
+ Int_t &phibinemc, Int_t &etabinemc, Int_t &nbinphi)
+{ // Returns number of bins in TPC and EMCAL geometry
+
+ etabintpc = fEtaBinInTPCAcc;
+ phibintpc = fPhiBinInTPCAcc;
+ etabinemc = fEtaBinInEMCalAcc;
+ phibinemc = fPhiBinInEMCalAcc;
+ nbinphi = fNbinPhi;
+}
+
+//__________________________________________________________
+Int_t AliJetGrid::GetIndexFromEtaPhi(Double_t phi,Double_t eta) const
+{ // Tells the index value of a corresponding (eta,phi) real position
+ // Loop over all entries -> takes time.
+ // Used one time at the begining to fill the grids
+
+ /* this is how bins are numbered
+
+ in all TPC or in TPC - EMCAL
+ ... ... ... ..
+ ---------------
+ ... ... . 10 | | | 11
+ ---+---+--- ---------------
+ ^ 6 | 7 | 8 or 8 | | | 9
+ | ---+---+--- ---------------
+ 3 | 4 | 5 4 | 5 | 6 | 7
+ phi ---+---+--- ---------------
+ 0 | 1 | 2 0 | 1 | 2 | 3
+
+
+ eta ->
+ */
+
+ Int_t etaBin=0,phiBin=0,absID=0;
+ Int_t etaBin2=0,etaBin3=0;
+
+ // Fill all the grid in eta/phi (all TPC acceptance)
+ //-----------------------------------------------------
+ if(fGrid==0){
+ if(eta <= fEta->At(0)) {
+ etaBin = 0;
+ } else if(eta >= fEta->At(fNeta)) {
+ etaBin = fNeta;
+ } else {
+ for(Int_t i=0; i<fNeta+1; i++) {
+ if(eta < fEta->At(i)) {
+ etaBin = i-1;
+ break;
+ }
+ }
+ }
+ if(phi <= fPhi->At(0)) {
+ phiBin = 0;
+ } else if(phi >= fPhi->At(fNphi)) {
+ phiBin = fNphi;
+ } else {
+ for(Int_t i=0; i<fNphi+1; i++) {
+ if(phi < fPhi->At(i)) {
+ phiBin = i-1;
+ break;
+ }
+ }
+ }
+
+ // Calculate absolute id
+ absID = phiBin*(fNeta+1) + etaBin;
+
+ }
+
+ // Fill the grid but do not count id in EMCal acceptance
+ //------------------------------------------------------
+ if((eta >= fEtaMin && eta <= fEtaMax) &&
+ (phi >= fPhiMin && phi <= fPhiMax)){
+ etaBin = etaBin2 = etaBin3 = -1;
+ phiBin = -1;
+ }
+
+ if(fGrid == 1){
+ if(phi<fPhiMin) {
+ if(eta <= fEta->At(0)) {
+ etaBin = 0;
+ } else if(eta >= fEta->At(fNeta)) {
+ etaBin = fNeta;
+ } else {
+ for(Int_t i=0; i<fNeta+1; i++) {
+ if(eta < fEta->At(i)) {
+ etaBin = i-1;
+ break;
+ }
+ }
+ }
+ }
+ else {
+ if((phi>=fPhiMin) && (phi<=fPhiMax)) {
+ if(eta <= fEta->At(0)) {
+ etaBin2 = 0;
+ } else if(eta >= fEta->At(fNeta)) {
+ etaBin2 = fNeta-fEtaBinInEMCalAcc;
+ } else {
+ for(Int_t i=0; i<fNeta+1; i++) {
+ if(eta < fEta->At(i) && ((fEta->At(0)<eta && eta<fEtaMin) ||
+ (fEtaMax<eta && eta<fEta->At(fNeta)))) {
+ // cout << "i : " << i << endl;
+ if(eta<fEtaMin)
+ etaBin2 = i-1;
+ else etaBin2 = i-1-fEtaBinInEMCalAcc;
+ break;
+ }
+ }
+ }
+ }
+ else {
+ if(eta <= fEta->At(0)) {
+ etaBin3 = 0;
+ } else if(eta >= fEta->At(fNeta)) {
+ etaBin3 = fNeta;
+ } else {
+ for(Int_t i=0; i<fNeta+1; i++) {
+ if(eta < fEta->At(i)) {
+ etaBin3 = i-1;
+ break;
+ }
+ }
+ }
+ }
+ }
+
+ if(phi <= fPhi->At(0)) {
+ phiBin = 0;
+ } else if(phi >= fPhi->At(fNphi)) {
+ phiBin = fNphi;
+ } else {
+ for(Int_t i=0; i<fNphi+1; i++) {
+ if(phi < fPhi->At(i)) {
+ phiBin = i-1;
+ break;
+ }
+ }
+ }
+
+ // Calculate absID
+ //-------------------------------------------------------
+ if(phi<fPhiMin)
+ absID = phiBin*(fNeta+1) + etaBin;
+ if(phi>=fPhiMin && phi<=fPhiMax) {
+ if(eta >= fEtaMin && eta <= fEtaMax) absID = -1;
+ else{
+ absID = (fNbinPhi+1)*(fNeta+1)
+ + (phiBin-fNbinPhi-1)*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
+ + etaBin2;
+ }
+ }
+ if(phi>fPhiMax)
+ absID = (fNbinPhi+1)*(fNeta+1)
+ + fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
+ + (phiBin-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)
+ + etaBin3;
+
+ } // END OPTION==1
+
+ return absID;
+
+}
+
+//__________________________________________________________
+void AliJetGrid::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
+{ // Get (eta,phi) position for a given index BUT loop over all entries (takes time)
+
+ for(Int_t j=0; j<fNphi+1; j++) {
+ for(Int_t i=0; i<fNeta+1; i++) {
+
+ // TPC grid only
+ //-------------------------------------
+ if(fGrid==0) {
+ if(j*(fNeta+1)+i == index) {
+ eta = fEta->At(i);
+ phi = fPhi->At(j);
+ }
+ }
+
+ // TPC-EMCAL grid
+ //-------------------------------------
+ Int_t ii = 0;
+ if(i==0) ii = 0;
+ if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
+ if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
+
+ if(fGrid==1) {
+ if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
+ eta = fEta->At(i);
+ phi = fPhi->At(j);
+ }
+
+ if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
+ ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
+ if(ii==0) {Int_t ind = 0; eta = fEta->At(ind);}
+ else eta = fEta->At(i);
+ phi = fPhi->At(j);
+ }
+
+ if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) &&
+ ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
+ +(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
+ eta = fEta->At(i);
+ phi = fPhi->At(j);
+ }
+ }
+ }
+ }
+}
+
+//__________________________________________________________
+Int_t AliJetGrid::GetIndex(Double_t phi, Double_t eta)
+{ // Get index value for a (eta,phi) position - Direct value
+
+ Int_t ieta = GetIndexJFromEta(eta);
+ Int_t iphi = GetIndexIFromPhi(phi);
+
+ Int_t index = GetMatrixIndex(iphi,ieta);
+
+ if(fDebug>10){
+ cout << "(phi,eta) : " << phi << ", " << eta << endl;
+ cout << "index : " << index << endl;
+ }
+ return index;
+
+}
+
+//__________________________________________________________
+Int_t AliJetGrid::GetIndexJFromEta(Double_t eta)
+{ // Get eta id
+ // Eta discretized
+
+ Int_t idEta =0;
+ Double_t temp = (eta+fMaxEta)/(fMaxEta-fMinEta)*fNeta;
+ if(fDebug>20)
+ {
+ cout << "eta : " << eta << endl;
+ cout << "fMaxEta : " << fMaxEta << endl;
+ cout << "fMinEta : " << fMinEta << endl;
+ cout << "fNeta : " << fNeta << endl;
+ cout << "temp eta before cast : " << temp << endl;
+ }
+ idEta = static_cast<Int_t>(temp+0.5);
+ if(fDebug>20) cout << "temp eta after cast : " << idEta << endl;
+ return idEta;
+}
+//__________________________________________________________
+Int_t AliJetGrid::GetIndexIFromPhi(Double_t phi)
+{ // Get phi id
+ // Phi discretized
+
+ Int_t idPhi = 0;
+ Double_t temp = 0.;
+ if(fMinPhi==0) temp = phi/(fMaxPhi-fMinPhi)*fNphi;
+ else temp = (phi-fMinPhi)/(fMaxPhi-fMinPhi)*fNphi;
+
+ if(fDebug>20)
+ {
+ cout << "phi : " << phi << endl;
+ cout << "fMaxPhi : " << fMaxPhi << endl;
+ cout << "fMinPhi : " << fMinPhi << endl;
+ cout << "fNphi : " << fNphi << endl;
+ cout << "temp phi before cast : " << temp << endl;
+ }
+ idPhi = static_cast<Int_t>(temp+0.5);
+ if(fDebug>20) cout << "temp phi after cast : " << idPhi << endl;
+ return idPhi;
+
+
+}
+
+//__________________________________________________________
+void AliJetGrid::SetMatrixIndex(Int_t i,Double_t par)
+{ // Allows to set parameters using only one index (if fGrid==0) !!
+ // Not used !
+
+ Int_t iphi = (Int_t)i/fNeta;
+ Int_t ieta = i-iphi*fNeta;
+ SetMatrixIndex(iphi,ieta,par);
+
+ return;
+}
+
+//__________________________________________________________
+void AliJetGrid::SetMatrixIndexes()
+{ // Fill the final matrix object with the corresponding index in eta/phi
+
+ for(Int_t i=0; i<fNphi+1; i++){
+ for(Int_t j=0; j<fNeta+1; j++){
+ (*fIndex)(i,j) = GetIndexFromEtaPhi(fPhi->At(i),fEta->At(j))+1;
+ if(fDebug>2){
+ cout << "(*fIndex)(" << i << "," << j << ") : " << (*fIndex)(i,j) <<
+ ", phi : " << fPhi->At(i) << ", eta : " << fEta->At(j) << endl;
+ }
+ }
+ }
+ printf("##########################################################\n");
+ printf("TMatrix object filled !\n");
+ printf("Size of the object : phi x eta = (fNphi+1) x (fNeta+1) = %d\n",(fNphi+1)*(fNeta+1));
+ printf("##########################################################\n");
+}
+
+//__________________________________________________________
+void AliJetGrid::SetIndexIJ()
+{ //
+
+ for(Int_t i=0; i<fNphi+1; i++){
+ for(Int_t j=0; j<fNeta+1; j++){
+ Int_t id = static_cast<Int_t>((*fIndex)(i,j));
+
+ if(id!=-1)
+ {
+ (*fIndexI)[id] = i;
+ (*fIndexJ)[id] = j;
+ }
+ }
+ }
+
+ printf("##########################################################\n");
+ printf(" In SetIndexIJ - Grid indexes setted !\n");
+ printf("##########################################################\n");
+}
+
+//__________________________________________________________
+void AliJetGrid::GetIJFromIndex(Int_t index, Int_t i, Int_t j)
+{ // returns i position id of eta and j position id of phi for a given grid index
+ i = (*fIndexI)[index];
+ j = (*fIndexJ)[index];
+}
+
+//__________________________________________________________
+void AliJetGrid::GetEtaPhiFromIndex2(Int_t index, Float_t &phi, Float_t &eta)
+{ // returns eta, phi values for a given grid index
+
+ phi = fPhi->At((*fIndexI)[index]);
+ eta = fEta->At((*fIndexJ)[index]);
+}
+
+//__________________________________________________________
+Int_t AliJetGrid::GetNEntries()
+{ // Returns the number of entries of the grid
+
+ if(fDebug>20){
+ cout << "fMaxPhi : " << fMaxPhi << endl;
+ cout << "fMaxEta : " << fMaxEta << endl;
+ }
+
+ Int_t indexNum = GetIndex(fMaxPhi,fMaxEta);
+ if(fDebug>20) cout << "indexNum : " << indexNum << endl;
+ return indexNum;
+
+}
+
+//__________________________________________________________
+Int_t AliJetGrid::GetNEntries2()
+{ // Returns the number of entries of the grid
+
+ Int_t indexNum = GetIndex(fMaxPhi-1.,fMaxEta-0.5);
+ if(fDebug>20) cout << "indexNum : " << indexNum << endl;
+ return indexNum;
+
+}
+
+
+
+
+
+
--- /dev/null
+#ifndef ALIJETGRID_H
+#define ALIJETGRID_H
+/* Copyright(c) 2001-2002, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+// Class description :
+//
+// Author : Magali Estienne, IPHC Strasbourg - e-mail: magali.estienne@ires.in2p3.fr
+//
+// --- Standard library ---
+#include <Riostream.h>
+//-- Root headers ---
+#include <TNamed.h>
+#include <TMatrixD.h>
+#include <TArrayD.h>
+#include <TArrayI.h>
+//-------------------
+
+class AliJetGrid : public TNamed {
+
+ public:
+
+ AliJetGrid();
+ AliJetGrid(Int_t nphi,Int_t neta,Double_t phiMin,Double_t etaMin,Double_t phiMax,Double_t etaMax);
+ AliJetGrid(const AliJetGrid& grid);
+ virtual ~AliJetGrid();
+
+ // Getter
+ TArrayD* GetArrayEta();
+ TArrayD* GetArrayPhi();
+ TMatrixD* GetIndexObject();
+ void InitParams(Double_t phiMinCal,Double_t phiMaxCal,Double_t etaMinCal,Double_t etaMaxCal);
+ Int_t GetIndexFromEtaPhi(Double_t phi,Double_t eta) const;
+ void GetEtaPhiFromIndex(Int_t index,Float_t &eta,Float_t &phi);
+ Int_t GetBinsEta() const {return fNeta+1;}
+ Int_t GetBinsPhi() const {return fNphi+1;}
+ Int_t GetTotBins() const {return (fNeta+1)*(fNphi+1);}
+ Int_t GetPointsEta() const {return fNeta;}
+ Int_t GetPointsPhi() const {return fNphi;}
+ Int_t GetTotPoints() const {return fNeta*fNphi;}
+ Int_t GetMatrixIndex(Int_t iphi,Int_t ieta) {return static_cast<Int_t>((*fIndex)(iphi,ieta));}
+ Int_t GetIndexIFromPhi(Double_t phi);
+ Int_t GetIndexJFromEta(Double_t eta);
+ Int_t GetIndex(Double_t phi,Double_t eta);
+ void GetAccParam(Int_t &nphi, Int_t &neta, Float_t &minphi,
+ Float_t &maxphi, Float_t &mineta, Float_t &maxeta);
+ void GetBinParam(Int_t &phibintpc, Int_t &etabintpc,
+ Int_t &phibinemc, Int_t &etabinemc, Int_t &nbinphi);
+ void GetIJFromIndex(Int_t index, Int_t i, Int_t j);
+ void GetEtaPhiFromIndex2(Int_t index, Float_t &phi, Float_t &eta);
+ Int_t GetNEntries();
+ Int_t GetNEntries2();
+ Int_t GetDeta() {return static_cast<Int_t>((fEtaMax-fEtaMin)/fNeta);
+ if(fDebug>21) cout << "static_cast<Int_t>((fEtaMax-fEtaMin)/fNeta) : " <<
+ static_cast<Int_t>((fEtaMax-fEtaMin)/fNeta);}
+ Int_t GetDphi() {return static_cast<Int_t>((fPhiMax-fPhiMin)/fNphi);
+ if(fDebug>21) cout << "static_cast<Int_t>((fPhiMax-fPhiMin)/fNphi) : " <<
+ static_cast<Int_t>((fPhiMax-fPhiMin)/fNphi);}
+ Int_t GetGridType() {return fGrid;}
+
+ // Setter
+ void SetEtaRange(Double_t etaMin, Double_t etaMax) {fEtaMin = etaMin; fEtaMax = etaMax;}
+ void SetPhiRange(Double_t phiMin, Double_t phiMax) {fPhiMin = phiMin; fPhiMax = phiMax;}
+ void SetNeta(Int_t neta) {fNeta = neta;}
+ void SetNphi(Int_t nphi) {fNphi = nphi;}
+ void SetMatrixIndexes();
+ void SetMatrixIndex(Int_t i,Double_t par);
+ void SetMatrixIndex(Int_t iphi,Int_t ieta,Double_t par) {
+ (*fIndex)(iphi,ieta)=par; return; }
+ void SetGridType(Int_t type) {fGrid = type;}
+ void SetIndexIJ();
+
+ private:
+ Int_t fGrid; // Close the type of grid you want to fill
+ // 0 = grid in tpc acceptance, 1 = grid in (tpc-emcal) acceptance
+ Int_t fNphi; // number of points in the grid in phi
+ Int_t fNeta; // " eta
+ TArrayD* fPhi; // grid points in phi
+ TArrayD* fEta; // grid points in eta
+ TMatrixD* fIndex; // matrix of indexes in the grid points
+ TArrayI* fIndexI; // grid points in phi
+ TArrayI* fIndexJ; // grid points in eta
+ Double_t fPhiMin; // phi acceptance min
+ Double_t fPhiMax; // phi acceptance max
+ Double_t fEtaMin; // eta acceptance min
+ Double_t fEtaMax; // eta acceptance max
+ Int_t fEtaBinInTPCAcc; // number of points in TPC acceptance in eta
+ Int_t fPhiBinInTPCAcc; // number of points in TPC acceptance in phi
+ Int_t fEtaBinInEMCalAcc; // number of points in EMCal acceptance in eta
+ Int_t fPhiBinInEMCalAcc; // number of points in EMCal acceptance in phi
+ Int_t fNbinEta;
+ Int_t fNbinPhi;
+ Double_t fMaxPhi;
+ Double_t fMinPhi;
+ Double_t fMaxEta;
+ Double_t fMinEta;
+ Int_t fDebug;
+
+ ClassDef(AliJetGrid,1) // Parameters used by AliTPCtrackerParam
+};
+
+
+#endif
+
+
+
+
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2002, 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. *
+ **************************************************************************/
+
+/* $Id$ */
+
+//===============================================================
+// To be modified for hadron correction using particle ID and track merging
+// Author : magali.estienne@ires.in2p3.fr
+//===============================================================
+//--
+//--
+//-- Hadron correction base class
+//--
+//--
+//--
+
+// --- AliRoot header files ---
+#include "AliJetHadronCorrection.h"
+
+ClassImp(AliJetHadronCorrection)
+
+AliJetHadronCorrection::AliJetHadronCorrection(const char *name,const char *title)
+:TNamed(name,title) { }
--- /dev/null
+#ifndef ALIJETHADRONCORRECTION_H
+#define ALIJETHADRONCORRECTION_H
+/* Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+/* $Id$ */
+
+//_________________________________________________________________________
+//
+//*-- Author: Aleksei Pavlinov (WSU)
+// This pure abstract class which defines only interface
+//
+#include "TNamed.h"
+
+class AliJetHadronCorrection : public TNamed {
+
+ public:
+ AliJetHadronCorrection(const char *name="name", const char *title="title");
+ virtual ~AliJetHadronCorrection() {}
+
+ // Add for particle
+ virtual Double_t GetEnergy(Double_t pmom, Double_t eta, Int_t gid)=0;
+
+ ClassDef(AliJetHadronCorrection,1) // Hadron correction for EMC (abstract class)
+};
+
+#endif // ALIJETHADRONCORRECTION_H
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2002, 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. *
+ **************************************************************************/
+
+
+/* $Id$ */
+
+//===============================================================
+// To be modified for hadron correction using particle ID and track merging
+// Author : magali.estienne@ires.in2p3.fr
+//===============================================================
+
+
+// --- AliRoot header files ---
+#include "AliJetHadronCorrectionv0.h"
+
+const Int_t maxVariant = 8; // size eta grid
+const Int_t nVec = 10; // size momentum grid
+const Int_t nPol = 4; // number coefficients of polinom
+static Double_t etaGrid[maxVariant]={ 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.67};
+static Double_t pGrid[nVec]={ 0.2, 0.5, 1.0, 2.0, 3.0, 5.0,10.0,15.0,25.0,40.0};
+// c[][] - first index for eta, second for momentum
+static Double_t c[maxVariant][nPol] ={
+ {1.305705e-01, 3.725653e-01, -1.219962e-02, 1.806235e-04},
+ {1.296153e-01, 3.827408e-01, -1.238640e-02, 1.822804e-04},
+ {1.337690e-01, 3.797454e-01, -1.245227e-02, 1.848243e-04},
+ {1.395796e-01, 3.623994e-01, -9.196803e-03, 1.243278e-04},
+ {1.457184e-01, 3.753655e-01, -1.035324e-02, 1.473447e-04},
+ {1.329164e-01, 4.219044e-01, -1.310515e-02, 1.948883e-04},
+ {8.136581e-02, 4.646087e-01, -1.531917e-02, 2.274749e-04},
+ {1.119836e-01, 4.262497e-01, -1.160125e-02, 1.628738e-04} };
+
+ClassImp(AliJetHadronCorrectionv0)
+
+AliJetHadronCorrectionv0* AliJetHadronCorrectionv0::fHadrCorr = 0;
+
+AliJetHadronCorrectionv0::AliJetHadronCorrectionv0(const char *name,const char *title)
+ :AliJetHadronCorrection(name, title)
+{
+ fHadrCorr = this;
+}
+
+AliJetHadronCorrectionv0*
+AliJetHadronCorrectionv0::Instance()
+{
+ fHadrCorr = new AliJetHadronCorrectionv0();
+ return fHadrCorr;
+}
+
+Double_t
+AliJetHadronCorrectionv0::GetEnergy(Double_t pmom, Double_t eta, Int_t /*gid*/)
+{
+ Int_t iEta=0; // index
+ Double_t etaw = TMath::Abs(eta);
+ if(etaw > etaGrid[maxVariant-1]) etaw = etaGrid[maxVariant-1];
+ for(Int_t i=0; i<maxVariant; i++) if(eta>=etaGrid[i]) {iEta = i; break;}
+
+ Double_t e[2], y, pw = pmom;
+ if(pmom > pGrid[nVec-1]) pw = pGrid[nVec-1];
+ for(Int_t i=0; i<2; i++){ // e for two eta value
+ e[i] = c[iEta][0];
+ y = 1.;
+ for(Int_t j=1; j<nPol; j++){
+ y *= pw;
+ e[i] += c[iEta][j]*y;
+ }
+ if(i==0) iEta ++;
+ }
+
+ Double_t deta = etaGrid[iEta] - etaGrid[iEta-1];
+ Double_t a = (e[1] - e[0])/deta; // slope
+ Double_t energy = e[0] + a*(eta-etaGrid[iEta-1]);
+ return energy;
+}
--- /dev/null
+#ifndef ALIJETHADRONCORRECTIONV0_H
+#define ALIJETHADRONCORRECTIONV0_H
+/* Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+/* $Id$ */
+
+//_________________________________________________________________________
+//
+//*-- Author: Aleksei Pavlinov (WSU)
+#include "AliJetHadronCorrection.h"
+
+class AliJetHadronCorrectionv0: public AliJetHadronCorrection {
+
+ private:
+ static AliJetHadronCorrectionv0* fHadrCorr;
+
+ protected:
+ AliJetHadronCorrectionv0(const char *name="HadronCorrectionv0", const char *title="title");
+
+ public:
+ static AliJetHadronCorrectionv0* Instance();
+ virtual Double_t GetEnergy(Double_t pmom, Double_t eta, Int_t gid);
+ Double_t GetEnergy(Double_t pmom, Double_t eta)
+ {return GetEnergy(pmom,eta,7);}
+
+ virtual ~AliJetHadronCorrectionv0() {}
+
+ ClassDef(AliJetHadronCorrectionv0,1) // Hadron correction for EMC (version for MDC)
+};
+
+#endif // ALIJETHADRONCORRECTIONV0_H
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2002, 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. *
+ **************************************************************************/
+
+//===============================================================
+// To be modified for hadron correction using particle ID and track merging
+// Author : magali.estienne@ires.in2p3.fr
+//===============================================================
+// Author : Mark Horner (LBL/UCT)
+
+// --- Standard library ---
+#include "Riostream.h"
+#include "TMath.h"
+
+// --- AliRoot header files ---
+//#include "../EMCAL/AliJetGeometry.h"
+#include "AliEMCALGeometry.h"
+#include "AliJetHadronCorrectionv1.h"
+
+ClassImp(AliJetHadronCorrectionv1)
+
+Double_t AliJetHadronCorrectionv1::fgParLookup[HCPARAMETERS][HCPARAMETERSETS] =
+{
+ {-2.82271e-4 , -2.39954e-4},
+ {2.50796e-2 , 2.07172e-2},
+ {1.02861e-3 , 1.48576e-3},
+ {2.11539e-2 , -1.38473-2},
+ {2.27003e-2 , 2.78252e-2},
+ {1.65078e-6 , 1.51821e-6}
+};
+
+AliJetHadronCorrectionv1* AliJetHadronCorrectionv1::fgHadrCorr = 0;
+
+void AliJetHadronCorrectionv1::SetGeometry(AliEMCALGeometry *geometry)
+{
+ // Initialise EMCAL geometry
+ if (!geometry)
+ {
+ SetParameters();
+ fSamplingFraction = geometry->GetSampling();
+ cout<<"Setting the sampling fraction to :"<<fSamplingFraction<<endl;
+ }else
+ {
+ SetParameters(geometry->GetName());
+ fSamplingFraction = geometry->GetSampling();
+ cout<<"Setting the sampling fraction to :"<<fSamplingFraction<<endl;
+ }
+ return;
+}
+
+/*
+###########################################################################
+###########################################################################
+
+This will have to be modified with the inclusion of the new EMCAL geometry
+That means I guess a new study of the hadron correction...
+
+Geometry : SHISH...
+
+###########################################################################
+###########################################################################
+*/
+
+
+
+
+void AliJetHadronCorrectionv1::SetGeometry(TString name,Double_t fs)
+{
+ // Initialise EMCAL geometry
+ cout << "Setting sampling fraction to "<<fSamplingFraction<<endl;
+ fSamplingFraction = fs;
+ if ( name == "" ||
+ name == "EMCAL_5655_21" ||
+ name == "EMCALArch1a" ||
+ name == "EMCALArch1aN" ||
+ name == "G56_2_55_19" ||
+ name == "G56_2_55_19_104_14" )
+ { // set parameters to this hadron correction
+ cout<<"HC parameters!"<<endl;
+ for (Int_t i=0;i<6;i++)
+ {
+ fPar[i] = fgParLookup[i][0];
+ cout <<fPar[i]<<endl;
+ }
+ }else if( name == "EMCAL_6564_21" ||
+ name == "G65_2_64_19" )
+ {
+ cout<<"HC parameters!"<<endl;
+ for (Int_t i=0;i<6;i++)
+ {
+ fPar[i] = fgParLookup[i][1];
+ cout <<fPar[i]<<endl;
+ }
+ }else
+ {
+ printf("Geometry not defined in hadron correction\n");
+ }
+
+}
+
+
+AliJetHadronCorrectionv1::AliJetHadronCorrectionv1(const char *name,const char *title)
+ :AliJetHadronCorrection(name, title)
+{
+ fgHadrCorr = this;
+}
+
+/*
+AliJetHadronCorrectionv1::AliJetHadronCorrectionv1(const char *name,const char *title,AliJetGeometry *geometry)
+{
+
+ fHadrCorr = this;
+ SetGeometry(geometry);
+
+}
+*/
+
+AliJetHadronCorrectionv1*
+AliJetHadronCorrectionv1::Instance()
+{
+ // return pointer to global instance. Instantiate if needed
+ if (! fgHadrCorr) new AliJetHadronCorrectionv1();
+ return fgHadrCorr;
+}
+
+Double_t
+AliJetHadronCorrectionv1::GetEnergy(Double_t pmom, Double_t eta, Int_t /*gid*/)
+{
+ // Return parametrised energy response
+ Double_t etai = TMath::Abs(eta);
+ Double_t value = fPar[5]*pmom*pmom*pmom+ fPar[0]*pmom*pmom+fPar[1]*pmom +fPar[2]*pmom*etai +fPar[3]*etai + fPar[4];
+ return fSamplingFraction*value;
+
+}
--- /dev/null
+#ifndef ALIJETHADRONCORRECTIONV1_H
+#define ALIJETHADRONCORRECTIONV1_H
+/* Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//
+//*-- Author: Mark Horner (LBL/UCT)
+//
+#include "AliJetHadronCorrection.h"
+
+#define HCPARAMETERS 6
+#define HCPARAMETERSETS 2
+
+class AliEMCALGeometry;
+
+class AliJetHadronCorrectionv1 : public AliJetHadronCorrection
+{
+ public:
+ static AliJetHadronCorrectionv1* Instance();
+ virtual Double_t GetEnergy(Double_t pmom, Double_t eta, Int_t gid);
+ Double_t GetEnergy(Double_t pmom, Double_t eta) {return GetEnergy(pmom,eta,7);}
+
+ void SetGeometry(TString name, Double_t fs = 1.);
+ virtual ~AliJetHadronCorrectionv1() {}
+
+ protected:
+ AliJetHadronCorrectionv1(const char *name="HadronCorrectionv1", const char *title="Hadron Correction");
+ // AliJetHadronCorrectionv1(const char *name="HadronCorrectionv1", const char *title="Hadron Correction",AliJetGeometry *geometry = NULL);
+ void SetGeometry(AliEMCALGeometry *geometry);
+
+ private:
+ void SetParameters(TString name = "") {Warning("SetParameter","Dummy method with argument %s",name.Data());}
+
+ static AliJetHadronCorrectionv1* fgHadrCorr; // Pointer to global instance (singleton)
+ static Double_t fgParLookup[HCPARAMETERS][HCPARAMETERSETS]; // Global array with parameters for hadronic response
+ Double_t fPar[6];
+ Float_t fSamplingFraction; // Sampling fraction
+
+ ClassDef(AliJetHadronCorrectionv1,2) // Hadron correction for EMC (version for MDC)
+};
+
+#endif // ALIJETHADRONCORRECTIONV1_H
--- /dev/null
+/**************************************************************************
+ * 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. *
+ **************************************************************************/
+
+
+/* $Id$ */
+
+//_________________________________________________________________________
+// Unit used by UA1 algorithm
+// --
+//*-- Author: Sarah Blyth (LBL/UCT)
+// --
+// Revised Version for JETAN - 30/03/2006
+// -- Magali Estienne (IReS)
+
+#include "AliJetUnitArray.h"
+
+
+ClassImp(AliJetUnitArray)
+
+AliJetUnitArray::AliJetUnitArray()
+{
+ // Default constructor
+ fUnitEnergy = 0.0;
+ fUnitEta = 0.0;
+ fUnitPhi = 0.0;
+ fUnitID = 0;
+ fUnitFlag = kOutJet;
+}
+
+AliJetUnitArray::~AliJetUnitArray()
+{
+ // Destructor
+}
+
+Bool_t AliJetUnitArray::operator>(AliJetUnitArray unit) const
+{
+ // Greater than operator used by sort
+ if( fUnitEnergy > unit.GetUnitEnergy())
+ return kTRUE;
+ else
+ return kFALSE;
+}
+
+Bool_t AliJetUnitArray::operator<( AliJetUnitArray unit) const
+{
+ // Less than operator used by sort
+ if( fUnitEnergy < unit.GetUnitEnergy())
+ return kTRUE;
+ else
+ return kFALSE;
+}
+
+Bool_t AliJetUnitArray::operator==( AliJetUnitArray unit) const
+{
+ // equality operator used by sort
+ if( fUnitEnergy == unit.GetUnitEnergy())
+ return kTRUE;
+ else
+ return kFALSE;
+}
--- /dev/null
+#ifndef ALIUA1JETFINDERUNIT_H
+#define ALIUA1JETFINDERUNIT_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * * * See cxx source for full Copyright notice */
+
+/* $Id$ */
+
+//======================================================================
+// Revised Version for JETAN - 08/11/2006
+// Functions added
+// Author: magali.estienne@ires.in2p3.fr
+//======================================================================
+// Unit used by UA1 algorithm
+//
+//*-- Author: Sarah Blyth (LBL/UCT)
+
+#include <TObject.h>
+#include "AliJetFinderTypes.h"
+
+class AliJetUnitArray : public TObject
+{
+ public:
+ AliJetUnitArray();
+ ~AliJetUnitArray();
+
+ // Setter
+ void SetUnitEnergy(Float_t energy) {fUnitEnergy = energy;}
+ void SetUnitEta(Float_t eta) {fUnitEta = eta;}
+ void SetUnitPhi(Float_t phi) {fUnitPhi = phi;}
+ void SetUnitDeta(Float_t deta) {fUnitDeta = deta;}
+ void SetUnitDphi(Float_t dphi) {fUnitDphi = dphi;}
+ void SetUnitID(Int_t id) {fUnitID = id;}
+ void SetUnitEntries(Int_t num) {fUnitNum = num;}
+ void SetUnitClusterID(Int_t id) {fUnitClusterID = id;}
+ void SetUnitFlag(AliJetFinderUnitFlagType_t flag)
+ {
+ fUnitFlag = flag;
+ }
+ void SetUnitCutFlag(AliJetFinderUnitCutFlagType_t cutFlag)
+ {
+ fUnitCutFlag = cutFlag;
+ }
+ void SetUnitSignalFlag(AliJetFinderUnitSignalFlagType_t signalFlag)
+ {
+ fUnitSignalFlag = signalFlag;
+ }
+ void SetUnitDetectorFlag(AliJetFinderUnitDetectorFlagType_t detectorflag)
+ {
+ fUnitDetectorFlag = detectorflag;
+ }
+
+ // Getter
+ Float_t GetUnitEnergy() const {return fUnitEnergy;}
+ Float_t GetUnitEta() const {return fUnitEta;}
+ Float_t GetUnitPhi() const {return fUnitPhi;}
+ Float_t GetUnitDeta() const {return fUnitDeta;}
+ Float_t GetUnitDphi() const {return fUnitDphi;}
+ Int_t GetUnitID() const {return fUnitID;}
+ Int_t GetUnitEntries() const {return fUnitNum;}
+ Int_t GetUnitClusterID() const {return fUnitClusterID;}
+ AliJetFinderUnitFlagType_t GetUnitFlag() const
+ {
+ return fUnitFlag;
+ }
+ AliJetFinderUnitCutFlagType_t GetUnitCutFlag() const
+ {
+ return fUnitCutFlag;
+ }
+ AliJetFinderUnitSignalFlagType_t GetUnitSignalFlag() const
+ {
+ return fUnitSignalFlag;
+ }
+ AliJetFinderUnitDetectorFlagType_t GetUnitDetectorFlag() const
+ {
+ return fUnitDetectorFlag;
+ }
+
+ Bool_t operator> ( AliJetUnitArray unit1) const;
+ Bool_t operator< ( AliJetUnitArray unit1) const;
+ Bool_t operator== ( AliJetUnitArray unit1) const;
+
+ protected:
+ Float_t fUnitEnergy; // Energy of the unit
+ Float_t fUnitEta; // Eta of the unit
+ Float_t fUnitPhi; // Phi of the unit
+ Float_t fUnitDeta; // Delta Eta of the unit
+ Float_t fUnitDphi; // Delta Phi of the unit
+ Int_t fUnitID; // ID of the unit
+ Int_t fUnitNum; // number of units
+ Int_t fUnitClusterID; // ID of the unit
+ AliJetFinderUnitFlagType_t fUnitFlag; //Flag of the unit
+ AliJetFinderUnitCutFlagType_t fUnitCutFlag; //Flag of the unit
+ AliJetFinderUnitSignalFlagType_t fUnitSignalFlag; //Flag of the unit
+ AliJetFinderUnitDetectorFlagType_t fUnitDetectorFlag; // Detector flag of the unit
+
+ ClassDef(AliJetUnitArray,3)
+};
+
+#endif
+