#include "TList.h"
#include "AliESDpid.h"
#include <iostream>
+#include "TH3F.h"
#include "TH2F.h"
+#include "TH2I.h"
+#include "TH1I.h"
+#include "TFile.h"
#include "AliAnalysisHadEtCorrections.h"
+#include "AliAnalysisEtSelector.h"
#include "AliLog.h"
#include "AliCentrality.h"
-
-
+#include "AliPHOSGeoUtils.h"
+#include "AliPHOSGeometry.h"
+#include "AliAnalysisEtRecEffCorrection.h"
+#include "AliESDpid.h"
using namespace std;
AliAnalysisEt()
,fCorrections(0)
,fPidCut(0)
- ,fClusterType(0)
+ ,nChargedHadronsMeasured(0)
+ ,nChargedHadronsTotal(0)
,fHistChargedPionEnergyDeposit(0)
,fHistProtonEnergyDeposit(0)
,fHistAntiProtonEnergyDeposit(0)
,fHistChargedKaonEnergyDeposit(0)
,fHistMuonEnergyDeposit(0)
- ,fHistRemovedEnergy(0)
+ ,fHistRemovedEnergy(0)
,fGeomCorrection(1.0)
- ,fEMinCorrection(1.0)
+ ,fEMinCorrection(1.0/0.687)
+ ,fRecEffCorrection(1.0)
+ ,fClusterPositionAccepted(0)
+ ,fClusterPositionAll(0)
+ ,fClusterPositionAcceptedEnergy(0)
+ ,fClusterPositionAllEnergy(0)
+ ,fClusterEnergy(0)
+ ,fClusterEnergyCent(0)
+ ,fClusterEnergyCentMatched(0)
+ ,fClusterEnergyCentNotMatched(0)
+ ,fClusterEt(0)
+ ,fHistChargedEnergyRemoved(0)
+ ,fHistNeutralEnergyRemoved(0)
+ ,fHistGammaEnergyAdded(0)
+ ,fHistMatchedTracksEvspTvsCent(0)
+ ,fHistMatchedTracksEvspTvsCentEffCorr(0)
+ ,fHistMatchedTracksEvspTvsCentEffTMCorr(0)
+ ,fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr(0)
+ ,fHistMatchedTracksEvspTvsCentEffTMCorr500MeV(0)
+ ,fHistFoundHadronsvsCent(0)
+ ,fHistNotFoundHadronsvsCent(0)
+ ,fHistFoundHadronsEtvsCent(0)
+ ,fHistNotFoundHadronsEtvsCent(0)
+ ,fHistFoundHadronsvsCent500MeV(0)
+ ,fHistNotFoundHadronsvsCent500MeV(0)
+ ,fHistFoundHadronsEtvsCent500MeV(0)
+ ,fHistNotFoundHadronsEtvsCent500MeV(0)
+ ,fHistNominalRawEt(0)
+ ,fHistNominalNonLinHighEt(0)
+ ,fHistNominalNonLinLowEt(0)
+ ,fHistNominalEffHighEt(0)
+ ,fHistNominalEffLowEt(0)
+ ,fHistTotRawEtEffCorr(0)
+ ,fHistTotRawEt(0)
+ ,fHistTotRawEtEffCorr500MeV(0)
+ ,fHistTotAllRawEt(0)
+ ,fHistTotAllRawEtEffCorr(0)
+ ,fHistNClustersPhosVsEmcal(0)
+ ,fHistClusterSizeVsCent(0)
+ ,fHistMatchedClusterSizeVsCent(0)
+ ,fHistTotAllRawEtVsTotalPt(0)
+ ,fHistTotAllRawEtVsTotalPtVsCent(0)
+ ,fHistTotMatchedRawEtVsTotalPtVsCent(0)
+ ,fHistPIDProtonsTrackMatchedDepositedVsNch(0)
+ ,fHistPIDAntiProtonsTrackMatchedDepositedVsNch(0)
+ ,fHistPIDProtonsTrackMatchedDepositedVsNcl(0)
+ ,fHistPIDAntiProtonsTrackMatchedDepositedVsNcl(0)
+ ,fHistPiKPTrackMatchedDepositedVsNch(0)
+ //,
+ ,fHistPIDProtonsTrackMatchedDepositedVsNchNoEff(0)
+ ,fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff(0)
+ ,fHistPIDProtonsTrackMatchedDepositedVsNclNoEff(0)
+ ,fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff(0)
+ ,fHistPiKPTrackMatchedDepositedVsNchNoEff(0)
+ ,fHistCentVsNchVsNclReco(0)
+ ,fHistRawSignalReco(0)
+ ,fHistEffCorrSignalReco(0)
+ ,fHistRecoRCorrVsPtVsCent(0)
{
}
AliAnalysisEtReconstructed::~AliAnalysisEtReconstructed()
{//destructor
delete fCorrections;
- delete fHistChargedPionEnergyDeposit; /** Energy deposited in calorimeter by charged pions */
- delete fHistProtonEnergyDeposit; /** Energy deposited in calorimeter by protons */
- delete fHistAntiProtonEnergyDeposit; /** Energy deposited in calorimeter by anti-protons */
- delete fHistChargedKaonEnergyDeposit; /** Energy deposited in calorimeter by charged kaons */
+ delete fHistChargedPionEnergyDeposit; /** Energy deposited in calorimeter by charged pions */
+ delete fHistProtonEnergyDeposit; /** Energy deposited in calorimeter by protons */
+ delete fHistAntiProtonEnergyDeposit; /** Energy deposited in calorimeter by anti-protons */
+ delete fHistChargedKaonEnergyDeposit; /** Energy deposited in calorimeter by charged kaons */
delete fHistMuonEnergyDeposit; /** Energy deposited in calorimeter by muons */
delete fHistRemovedEnergy; // removed energy
-
+ delete fClusterPositionAccepted;
+ delete fClusterPositionAll;
+ delete fClusterPositionAcceptedEnergy;
+ delete fClusterPositionAllEnergy;
+ delete fClusterEnergy;
+ delete fClusterEnergyCent;
+ delete fClusterEnergyCentMatched;
+ delete fClusterEnergyCentNotMatched;
+ delete fClusterEt;
+ delete fHistChargedEnergyRemoved;
+ delete fHistNeutralEnergyRemoved;
+ delete fHistGammaEnergyAdded;
+ delete fHistMatchedTracksEvspTvsCent;
+ delete fHistMatchedTracksEvspTvsCentEffCorr;
+ delete fHistMatchedTracksEvspTvsCentEffTMCorr;
+ delete fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr;
+ delete fHistMatchedTracksEvspTvsCentEffTMCorr500MeV;
+ delete fHistFoundHadronsvsCent;
+ delete fHistNotFoundHadronsvsCent;
+ delete fHistFoundHadronsEtvsCent;
+ delete fHistNotFoundHadronsEtvsCent;
+ delete fHistFoundHadronsvsCent500MeV;
+ delete fHistNotFoundHadronsvsCent500MeV;
+ delete fHistFoundHadronsEtvsCent500MeV;
+ delete fHistNotFoundHadronsEtvsCent500MeV;
+ delete fHistNominalRawEt;
+ delete fHistNominalNonLinHighEt;
+ delete fHistNominalNonLinLowEt;
+ delete fHistNominalEffHighEt;
+ delete fHistNominalEffLowEt;
+ delete fHistTotRawEtEffCorr;
+ delete fHistTotRawEt;
+ delete fHistTotAllRawEt;
+ delete fHistTotAllRawEtEffCorr;
+ delete fHistTotRawEtEffCorr500MeV;
+ delete fHistNClustersPhosVsEmcal;
+ delete fHistClusterSizeVsCent;
+ delete fHistMatchedClusterSizeVsCent;
+ delete fHistTotAllRawEtVsTotalPt;
+ delete fHistTotAllRawEtVsTotalPtVsCent;
+ delete fHistTotMatchedRawEtVsTotalPtVsCent;
+ delete fHistPIDProtonsTrackMatchedDepositedVsNch;
+ delete fHistPIDAntiProtonsTrackMatchedDepositedVsNch;
+ delete fHistPIDProtonsTrackMatchedDepositedVsNcl;
+ delete fHistPIDAntiProtonsTrackMatchedDepositedVsNcl;
+ delete fHistPiKPTrackMatchedDepositedVsNch;
+ delete fHistPIDProtonsTrackMatchedDepositedVsNchNoEff;
+ delete fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff;
+ delete fHistPIDProtonsTrackMatchedDepositedVsNclNoEff;
+ delete fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff;
+ delete fHistPiKPTrackMatchedDepositedVsNchNoEff;
+ delete fHistCentVsNchVsNclReco;
+ delete fHistRawSignalReco;
+ delete fHistEffCorrSignalReco;
+ delete fHistRecoRCorrVsPtVsCent;
}
Int_t AliAnalysisEtReconstructed::AnalyseEvent(AliVEvent* ev)
{
+
+ //AliAnalysisEt::AnalyseEvent(ev);
// analyse ESD event
ResetEventValues();
-
if (!ev) {
AliFatal("ERROR: Event does not exist");
return 0;
AliFatal("ERROR: ESD Event does not exist");
return 0;
}
-
+ if(!fSelector){
+ AliFatal("ERROR: fSelector does not exist");
+ return 0;
+ }
+ fSelector->SetEvent(event);
+
Int_t cent = -1;
- if (fCentrality)
+ fCentrality = event->GetCentrality();
+ if (fCentrality && cent)
{
- cent = fCentrality->GetCentralityClass10("V0M");
- fCentClass = fCentrality->GetCentralityClass10("V0M");
+ cent = fCentrality->GetCentralityClass5("V0M");
+ fCentClass = fCentrality->GetCentralityClass5("V0M");
}
- Double_t protonMass = fgProtonMass;
-
- //for PID
- AliESDpid pID;
- pID.MakePID(event);
- TObjArray* list = fEsdtrackCutsTPC->GetAcceptedTracks(event);
-
- Int_t nGoodTracks = list->GetEntries();
- // printf("nGoodTracks %d nCaloClusters %d\n", nGoodTracks, event->GetNumberOfCaloClusters());
-
- for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
- {
- AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
- if (!track)
- {
- AliError(Form("ERROR: Could not get track %d", iTrack));
- continue;
- }
-
- fMultiplicity++;
-
-
- Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
- nSigmaPion = TMath::Abs(pID.NumberOfSigmasTPC(track,AliPID::kPion));
- nSigmaProton = TMath::Abs(pID.NumberOfSigmasTPC(track,AliPID::kProton));
- nSigmaKaon = TMath::Abs(pID.NumberOfSigmasTPC(track,AliPID::kKaon));
- nSigmaElectron = TMath::Abs(pID.NumberOfSigmasTPC(track,AliPID::kElectron));
- /*
- bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
- bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
- bool isKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
- bool isProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
- */
- Int_t nItsClusters = dynamic_cast<AliESDtrack*>(track)->GetNcls(0);
- Int_t nTPCClusters = dynamic_cast<AliESDtrack*>(track)->GetNcls(1);
+ //for PID
+ //AliESDpid *pID = new AliESDpid();
+ //pID->MakePID(event);
+ Float_t etPIDProtons = 0.0;
+ Float_t etPIDAntiProtons = 0.0;
+ Float_t etPiKPMatched = 0.0;
+ Float_t etPIDProtonsNoEff = 0.0;
+ Float_t etPIDAntiProtonsNoEff = 0.0;
+ Float_t etPiKPMatchedNoEff = 0.0;
+ Float_t multiplicity = fEsdtrackCutsTPC->GetReferenceMultiplicity(event,kTRUE);
- Float_t massPart = 0;
- const Double_t *pidWeights = track->PID();
- Int_t maxpid = -1;
- Double_t maxpidweight = 0;
-
- if (pidWeights)
- {
- for (Int_t p =0; p < AliPID::kSPECIES; p++)
- {
- if (pidWeights[p] > maxpidweight)
- {
- maxpidweight = pidWeights[p];
- maxpid = p;
- }
- }
- if (maxpid == AliPID::kProton)
- {
- //by definition of ET
- massPart = -protonMass*track->Charge();
- }
-
- }
-
- Double_t et = track->E() * TMath::Sin(track->Theta()) + massPart;
-
- fSparseTracks[0] = maxpid;
- fSparseTracks[1] = track->Charge();
- fSparseTracks[2] = track->M();
- fSparseTracks[3] = et;
- fSparseTracks[4] = track->Pt();
- fSparseTracks[5] = track->Eta();
- fSparseTracks[6] = cent;
- fSparseHistTracks->Fill(fSparseTracks);
- //printf("Rec track: iTrack %03d eta %4.3f phi %4.3f nITSCl %d nTPCCl %d\n", iTrack, track->Eta(), track->Phi(), nItsClusters, nTPCClusters); // tmp/debug printout
-
- if (TMath::Abs(track->Eta()) < fCuts->GetCommonEtaCut() && CheckGoodVertex(track) && nItsClusters > fCuts->GetReconstructedNItsClustersCut() && nTPCClusters > fCuts->GetReconstructedNTpcClustersCut() )
- {
- fTotChargedEt += et;
- fChargedMultiplicity++;
- if (maxpid != -1)
- {
- if (maxpid == AliPID::kProton)
- {
- fProtonEt += et;
- }
- if (maxpid == AliPID::kPion)
- {
- fPionEt += et;
- }
- if (maxpid == AliPID::kKaon)
- {
- fChargedKaonEt += et;
- }
- if (maxpid == AliPID::kMuon)
- {
- fMuonEt += et;
- }
- if (maxpid == AliPID::kElectron)
- {
- fElectronEt += et;
- }
- }
-
- if (TMath::Abs(track->Eta()) < fEtaCutAcc && track->Phi() < fPhiCutAccMax && track->Phi() > fPhiCutAccMin)
- {
- fTotChargedEtAcc += track->E()*TMath::Sin(track->Theta()) + massPart;
- if (maxpid != -1)
- {
- if (maxpid == AliPID::kProton)
- {
- fProtonEtAcc += et;
- }
- if (maxpid == AliPID::kPion)
- {
- fPionEtAcc += et;
- }
- if (maxpid == AliPID::kKaon)
- {
- fChargedKaonEtAcc += et;
- }
- if (maxpid == AliPID::kMuon)
- {
- fMuonEtAcc += et;
- }
- if (maxpid == AliPID::kElectron)
- {
- fElectronEtAcc += et;
- }
- }
-
- }
- }
-
- if (TrackHitsCalorimeter(track, event->GetMagneticField()))
- {
- Double_t phi = track->Phi();
- Double_t pt = track->Pt();
- // printf("Rec track hit: iTrack %03d phi %4.3f pt %4.3f\n", iTrack, phi, pt); // tmp/debug printout
- if (track->Charge() > 0) fHistPhivsPtPos->Fill(phi, pt);
- else fHistPhivsPtNeg->Fill(phi, pt);
- }
+ Float_t totalMatchedPt = 0.0;
+ Float_t totalPt = 0.0;
+ TObjArray* list = fEsdtrackCutsTPC->GetAcceptedTracks(event);
+ Int_t nGoodTracks = list->GetEntries();
+ for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++){
+ AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
+ if (!track)
+ {
+ Printf("ERROR: Could not get track %d", iTrack);
+ continue;
+ }
+ else{
+ totalPt +=track->Pt();
+ //pID->MakeITSPID(track);
+
+
+ }
}
- for (Int_t iCluster = 0; iCluster < event->GetNumberOfCaloClusters(); iCluster++)
+ //TRefArray *caloClusters = fSelector->GetClusters();//just gets the correct set of clusters - does not apply any cuts
+ //Float_t fClusterMult = caloClusters->GetEntries();
+
+ Float_t nominalRawEt = 0;
+ Float_t totEt500MeV = 0;
+ Float_t nonlinHighRawEt = 0;
+ Float_t nonlinLowRawEt = 0;
+ Float_t effHighRawEt = 0;
+ Float_t effLowRawEt = 0;
+ Float_t uncorrEt = 0;
+ Float_t rawSignal = 0;
+ Float_t effCorrSignal = 0;
+
+ nChargedHadronsMeasured = 0.0;
+ nChargedHadronsTotal = 0.0;
+ Float_t nChargedHadronsEtMeasured = 0.0;
+ Float_t nChargedHadronsEtTotal = 0.0;
+ Float_t nChargedHadronsMeasured500MeV = 0.0;
+ Float_t nChargedHadronsTotal500MeV = 0.0;
+ Float_t nChargedHadronsEtMeasured500MeV = 0.0;
+ Float_t nChargedHadronsEtTotal500MeV = 0.0;
+ Float_t fTotAllRawEt = 0.0;
+ Float_t fTotRawEt = 0.0;
+ Float_t fTotRawEtEffCorr = 0.0;
+ Float_t fTotAllRawEtEffCorr = 0.0;
+ Int_t nPhosClusters = 0;
+ Int_t nEmcalClusters = 0;
+
+
+ TRefArray *caloClusters = fSelector->GetClusters();
+ Int_t nCluster = caloClusters->GetEntries();
+
+ for (int iCluster = 0; iCluster < nCluster; iCluster++ )
{
- AliESDCaloCluster* cluster = event->GetCaloCluster(iCluster);
+ AliESDCaloCluster* cluster = ( AliESDCaloCluster* )caloClusters->At( iCluster );
if (!cluster)
{
AliError(Form("ERROR: Could not get cluster %d", iCluster));
continue;
}
- if (cluster->GetType() != fClusterType) continue;
-
- //if(cluster->GetTracksMatched() > 0)
-// printf("Rec Cluster: iCluster %03d E %4.3f type %.qd NCells %d, nmatched: %d, distance to closest: %f\n", iCluster, cluster->E(), (int)(cluster->GetType()), cluster->GetNCells(), cluster->GetNTracksMatched(), cluster->GetEmcCpvDistance()); // tmp/debug printout
-
-
- if (cluster->E() < fClusterEnergyCut) continue;
-
+ int x = 0;
+ fCutFlow->Fill(x++);
+ if(cluster->IsEMCAL()) nEmcalClusters++;
+ else nPhosClusters++;
+ if(!fSelector->IsDetectorCluster(*cluster)) continue;
+ fCutFlow->Fill(x++);
+ if(!fSelector->PassMinEnergyCut(*cluster)) continue;
+ fCutFlow->Fill(x++);
+ if (!fSelector->PassDistanceToBadChannelCut(*cluster)) continue;
+ fCutFlow->Fill(x++);
+ if (!fSelector->CutGeometricalAcceptance(*cluster)) continue;
+ //fCutFlow->Fill(x++);
Float_t pos[3];
cluster->GetPosition(pos);
TVector3 cp(pos);
+ fClusterPositionAll->Fill(cp.Phi(), cp.PseudoRapidity());
+ Float_t fReconstructedE = cluster->E();
+ fClusterPositionAllEnergy->Fill(cp.Phi(), cp.PseudoRapidity(),GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE);
+
+ //if(TMath::Abs(cp.Eta())> fCuts->fCuts->GetGeometryEmcalEtaAccCut() || cp.Phi() > fCuts->GetGeometryEmcalPhiAccMaxCut()*TMath::Pi()/180. || cp.Phi() > fCuts->GetGeometryEmcalPhiAccMinCut()*TMath::Pi()/180.) continue;//Do not accept if cluster is not in the acceptance
+ fTotAllRawEt += TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ fTotAllRawEtEffCorr +=GetCorrectionModification(*cluster,0,0,cent)* CorrectForReconstructionEfficiency(*cluster,cent);
+
+ fClusterEnergyCent->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE,cent);
+ Bool_t matched = kTRUE;//default to no track matched
+ Int_t trackMatchedIndex = cluster->GetTrackMatchedIndex();//find the index of the matched track
+ matched = !(fSelector->PassTrackMatchingCut(*cluster));//PassTrackMatchingCut is false if there is a matched track
+ if(matched){//if the track match is good (, is the track good?
+ if(trackMatchedIndex < 0) matched=kFALSE;//If the index is bad, don't count it
+ if(matched){
+ AliESDtrack *track = event->GetTrack(trackMatchedIndex);
+ //if this is a good track, accept track will return true. The track matched is good, so not track matched is false
+ matched = fEsdtrackCutsTPC->AcceptTrack(track);//If the track is bad, don't count it
+ if(matched){//if it is still matched see if the track p was less than the energy
+ Float_t rcorr = TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ fHistRecoRCorrVsPtVsCent->Fill(rcorr,track->Pt(), fCentClass);
+ if(fSelector->PassMinEnergyCut( (fReconstructedE - fsub* track->P())*TMath::Sin(cp.Theta()) )){//if more energy was deposited than the momentum of the track and more than one particle led to the cluster
+ // if(fReconstructedE - fsub* track->P() > 0.0){
+ fReconstructedE = fReconstructedE - fsub* track->P();
+ matched = kFALSE;
+ }
+ }
+ }
+ }
- Double_t distance = cluster->GetEmcCpvDistance();
- Int_t trackMatchedIndex = cluster->GetTrackMatchedIndex();
- if ( cluster->IsEMCAL() ) {
- distance = CalcTrackClusterDistance(pos, &trackMatchedIndex, event);
- }
-
- fSparseClusters[0] = 0;
- fSparseClusters[1] = 0;
- fSparseClusters[2] = 0;
- fSparseClusters[6] = 0;
- fSparseClusters[7] = 0;
- fSparseClusters[8] = 0;
- fSparseClusters[9] = cent;
- fSparseClusters[10] = 0;
-
-
- if (cluster->GetNTracksMatched() > 0 && trackMatchedIndex > -1)
- {
- AliVTrack *tmptrack = event->GetTrack(trackMatchedIndex);
- if (!tmptrack)
- {
- AliError("Error: track does not exist");
- return -1;
- }
- const Double_t *pidWeights = tmptrack->PID();
-
- Double_t maxpidweight = 0;
- Int_t maxpid = 0;
- Double_t massPart = 0;
- if (pidWeights)
- {
- for (Int_t p =0; p < AliPID::kSPECIES; p++)
- {
- if (pidWeights[p] > maxpidweight)
- {
- maxpidweight = pidWeights[p];
- maxpid = p;
- }
- }
- if (maxpid == AliPID::kProton)
- {
- //by definition of ET
- massPart = -protonMass*tmptrack->Charge();
- }
- }
- fSparseClusters[0] = maxpid;
- fSparseClusters[1] = tmptrack->Charge();
- fSparseClusters[2] = tmptrack->M();
- fSparseClusters[6] = tmptrack->E() * TMath::Sin(tmptrack->Theta()) + massPart;;
- fSparseClusters[7] = tmptrack->Pt();
- fSparseClusters[8] = tmptrack->Eta();
- }
-
- fSparseClusters[10] = distance;
-
- fHistTMDeltaR->Fill(distance);
- fHistTMDxDz->Fill(cluster->GetTrackDx(), cluster->GetTrackDz());
-
-// Float_t clusteret = cluster->E() * TMath::Sin(cp.Theta());
-
- Bool_t matched = false;
-
- if (cluster->IsEMCAL()) matched = distance < fTrackDistanceCut;
- else matched = (TMath::Abs(cluster->GetTrackDx()) < fTrackDxCut && TMath::Abs(cluster->GetTrackDz()) < fTrackDzCut);
if (matched)
{
+
if (cluster->GetNTracksMatched() > 0 && trackMatchedIndex>=0)
{
AliVTrack *track = event->GetTrack(trackMatchedIndex);
AliError("Error: track does not exist");
}
else {
+ totalMatchedPt +=track->Pt();
+ fClusterEnergyCentMatched->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE,cent);
+ fHistMatchedClusterSizeVsCent->Fill(cluster->GetNCells(),cent);
+
+ float eff = fTmCorrections->TrackMatchingEfficiency(track->Pt(),cent);
+ if(TMath::Abs(eff)<1e-5) eff = 1.0;
+ //cout<<"pt "<<track->Pt()<<" eff "<<eff<<" total "<<nChargedHadronsTotal<<endl;
+ nChargedHadronsMeasured++;
+ nChargedHadronsTotal += 1/eff;
+ Double_t effCorrEt = GetCorrectionModification(*cluster,0,0,cent) * CorrectForReconstructionEfficiency(*cluster,fReconstructedE,cent);
+ nChargedHadronsEtMeasured+= TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ //One efficiency is the gamma efficiency and the other is the track matching efficiency.
+ nChargedHadronsEtTotal+= 1/eff *TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ //cout<<"nFound "<<1<<" nFoundTotal "<<1/eff<<" etMeas "<<TMath::Sin(cp.Theta())*fReconstructedE<<" ET total "<< 1/eff *TMath::Sin(cp.Theta())*fReconstructedE<<endl;
+
+ Float_t nSigmaPion = fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion);
+ Float_t nSigmaProton = fPIDResponse->NumberOfSigmasTPC(track, AliPID::kProton);
+ bool isProton = (nSigmaPion>3.0 && nSigmaProton<3.0 && track->Pt()<0.9);
+ //cout<<"NSigmaProton "<<nSigmaProton<<endl;
+ etPiKPMatched += effCorrEt;
+ etPiKPMatchedNoEff +=TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ if(isProton){
+ if(track->Charge()>0){
+ etPIDProtons += effCorrEt;
+ etPIDProtonsNoEff +=TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ }
+ else{
+ etPIDAntiProtonsNoEff +=TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ etPIDAntiProtons += effCorrEt;
+ }
+ }
+ if(TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE>0.5){
+ nChargedHadronsMeasured500MeV++;
+ nChargedHadronsTotal500MeV += 1/eff;
+ nChargedHadronsEtMeasured500MeV+= TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ nChargedHadronsEtTotal500MeV+= 1/eff *TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ }
+ fHistMatchedTracksEvspTvsCent->Fill(track->P(),TMath::Sin(cp.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE,cent);
+ fHistMatchedTracksEvspTvsCentEffCorr->Fill(track->P(),effCorrEt,cent);
+ //Weighed by the number of tracks we didn't find
+ fHistMatchedTracksEvspTvsCentEffTMCorr->Fill(track->P(), effCorrEt,cent, (1/eff-1) );
+ if(cent<16 && cent>11){//centralities 60-80% where false track matches are low
+ for(int cbtest = 0; cbtest<20; cbtest++){//then we calculate the deposit matched to hadrons with different centrality bins' efficiencies
+ float efftest = fTmCorrections->TrackMatchingEfficiency(track->Pt(),cbtest);
+ if(TMath::Abs(efftest)<1e-5) efftest = 1.0;
+ Double_t effCorrEttest = GetCorrectionModification(*cluster,0,0,cent)*CorrectForReconstructionEfficiency(*cluster,fReconstructedE,cbtest);
+ fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr->Fill(track->P(), effCorrEttest,cbtest, (1/efftest-1) );
+ }
+ }
+ cluster->GetPosition(pos);
+ TVector3 p2(pos);
+ uncorrEt += TMath::Sin(p2.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ if(uncorrEt>=0.5) fHistMatchedTracksEvspTvsCentEffTMCorr500MeV->Fill(track->P(), effCorrEt,cent, (1/eff-1) );
const Double_t *pidWeights = track->PID();
Double_t maxpidweight = 0;
maxpid = p;
}
}
- if (fCuts->GetHistMakeTreeDeposit() && fTreeDeposit)
+ if (fCuts->GetHistMakeTreeDeposit() && fDepositTree)
{
- fEnergyDeposited = cluster->E();
- fEnergyTPC = track->E();
+ fEnergyDeposited =GetCorrectionModification(*cluster,0,0,cent)* fReconstructedE;
+ fMomentumTPC = track->P();
fCharge = track->Charge();
fParticlePid = maxpid;
fPidProb = maxpidweight;
}
else {
if (esdTrack) fTrackPassedCut = fEsdtrackCutsTPC->AcceptTrack(esdTrack);
- fTreeDeposit->Fill();
+ fDepositTree->Fill();
}
}
if (maxpidweight > fPidCut)
{
- Float_t dist = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
+ //Float_t dist = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
- Float_t theta = TMath::ATan(pos[2]/dist)+TMath::Pi()/2;
+ //Float_t theta = TMath::ATan(pos[2]/dist)+TMath::Pi()/2;
- Float_t et = cluster->E() * TMath::Sin(theta);
+ //Float_t et = fReconstructedE * TMath::Sin(theta);
if (maxpid == AliPID::kProton)
{
if (track->Charge() == 1)
{
- fBaryonEt += et;
- fHistProtonEnergyDeposit->Fill(cluster->E(), track->E());
+ fHistProtonEnergyDeposit->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE, track->E());
}
else if (track->Charge() == -1)
{
- fAntiBaryonEt += et;
- fHistAntiProtonEnergyDeposit->Fill(cluster->E(), track->E());
+ fHistAntiProtonEnergyDeposit->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE, track->E());
}
}
else if (maxpid == AliPID::kPion)
{
- fMesonEt += et;
- fHistChargedPionEnergyDeposit->Fill(cluster->E(), track->E());
+ fHistChargedPionEnergyDeposit->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE, track->E());
}
else if (maxpid == AliPID::kKaon)
{
- fMesonEt += et;
- fHistChargedKaonEnergyDeposit->Fill(cluster->E(), track->E());
+ fHistChargedKaonEnergyDeposit->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE, track->E());
}
else if (maxpid == AliPID::kMuon)
{
- fHistMuonEnergyDeposit->Fill(cluster->E(), track->E());
+ fHistMuonEnergyDeposit->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE, track->E());
}
}
}
}
//continue;
} // distance
- else
- {
- fSparseClusters[0] = AliPID::kPhoton;
- fSparseClusters[1] = 0;
-
- if (cluster->E() > fSingleCellEnergyCut && cluster->GetNCells() == fCuts->GetCommonSingleCell()) continue;
- if (cluster->E() < fClusterEnergyCut) continue;
- cluster->GetPosition(pos);
-
- // TODO: replace with TVector3, too lazy now...
-
- float dist = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
-
- float theta = TMath::ATan(pos[2]/dist)+TMath::Pi()/2;
- // float eta = TMath::Log(TMath::Abs( TMath::Tan( 0.5 * theta ) ) );
- fTotNeutralEt += cluster->E() * TMath::Sin(theta);
+ else{//these are clusters which were not track matched
+ fCutFlow->Fill(x++);
+ //std::cout << x++ << std::endl;
+
+ //if (fReconstructedE > fSingleCellEnergyCut && cluster->GetNCells() == fCuts->GetCommonSingleCell()) continue;
+ //if (fReconstructedE < fClusterEnergyCut) continue;
+ cluster->GetPosition(pos);
+
+ TVector3 p2(pos);
+
+ fClusterPositionAccepted->Fill(p2.Phi(), p2.PseudoRapidity());
+ fClusterPositionAcceptedEnergy->Fill(p2.Phi(), p2.PseudoRapidity(),GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE);
+ fClusterEnergy->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE);
+ fClusterEnergyCentNotMatched->Fill(GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE,cent);
+ fHistClusterSizeVsCent->Fill(cluster->GetNCells(),cent);
+ fClusterEt->Fill(TMath::Sin(p2.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE);
+ uncorrEt += TMath::Sin(p2.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ float myuncorrEt = TMath::Sin(p2.Theta())*GetCorrectionModification(*cluster,0,0,cent)*fReconstructedE;
+ fTotRawEt += myuncorrEt;
+
+ Double_t effCorrEt = CorrectForReconstructionEfficiency(*cluster,fReconstructedE,cent)*GetCorrectionModification(*cluster,0,0,cent);
+ rawSignal += myuncorrEt;
+ effCorrSignal +=effCorrEt;
+ //cout<<"cluster energy "<<fReconstructedE<<" eff corr Et "<<effCorrEt<<endl;
+ fTotRawEtEffCorr += effCorrEt;
+ fTotNeutralEt += effCorrEt;
+ nominalRawEt += effCorrEt;
+ if(myuncorrEt>=0.5){
+ totEt500MeV += effCorrEt;
+ //cout<<"test "<<myuncorrEt<<"> 0.5"<<endl;
+ }
+ else{
+ //cout<<"test "<<myuncorrEt<<"< 0.5"<<endl;
+ }
+ nonlinHighRawEt += effCorrEt*GetCorrectionModification(*cluster,1,0,cent);
+ nonlinLowRawEt += effCorrEt*GetCorrectionModification(*cluster,-1,0,cent);
+ effHighRawEt += effCorrEt*GetCorrectionModification(*cluster,0,1,cent);
+ effLowRawEt += effCorrEt*GetCorrectionModification(*cluster,0,-1,cent);
fNeutralMultiplicity++;
}
+ fMultiplicity++;
+ }
+
- cluster->GetPosition(pos);
-
- // TODO: replace with TVector3
-
- float dist = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
- float theta = TMath::ATan(pos[2]/dist)+TMath::Pi()/2;
- float eta = TMath::Log(TMath::Abs( TMath::Tan( 0.5 * theta ) ) );
- fSparseClusters[3] = cluster->E() * TMath::Sin(theta);
- fSparseClusters[4] = cluster->E();
- fSparseClusters[5] = eta;
+ fHistRawSignalReco->Fill(rawSignal);
+ fHistEffCorrSignalReco->Fill(effCorrSignal);
- fSparseHistClusters->Fill(fSparseClusters);
+ fHistNClustersPhosVsEmcal->Fill(nPhosClusters,nEmcalClusters,cent);
+ fChargedEnergyRemoved = GetChargedContribution(fNeutralMultiplicity);
+ fNeutralEnergyRemoved = GetNeutralContribution(fNeutralMultiplicity);
+ fHistChargedEnergyRemoved->Fill(fChargedEnergyRemoved, fNeutralMultiplicity);
+ fHistNeutralEnergyRemoved->Fill(fNeutralEnergyRemoved, fNeutralMultiplicity);
+
+ fGammaEnergyAdded = GetGammaContribution(fNeutralMultiplicity);
+ fHistGammaEnergyAdded->Fill(fGammaEnergyAdded, fNeutralMultiplicity);
- fMultiplicity++;
- }
- Double_t removedEnergy = GetChargedContribution(fNeutralMultiplicity) + GetNeutralContribution(fNeutralMultiplicity) - GetGammaContribution(fNeutralMultiplicity);
+ Double_t removedEnergy = GetChargedContribution(fNeutralMultiplicity) + GetNeutralContribution(fNeutralMultiplicity) + GetGammaContribution(fNeutralMultiplicity) + GetSecondaryContribution(fNeutralMultiplicity);
fHistRemovedEnergy->Fill(removedEnergy);
-// std::cout << "fTotNeutralEt: " << fTotNeutralEt << ", Contribution from non-removed charged: " << GetChargedContribution(fNeutralMultiplicity) << ", neutral: " << GetNeutralContribution(fNeutralMultiplicity) << ", gammas: " << GetGammaContribution(fNeutralMultiplicity) << ", multiplicity: " << fNeutralMultiplicity<< std::endl;
- fTotNeutralEt = fGeomCorrection * fEMinCorrection * fTotNeutralEt - removedEnergy;
- fTotNeutralEtAcc = fTotNeutralEt/fGeomCorrection;
+
+ fTotNeutralEtAcc = fTotNeutralEt;
+ //fHistTotRawEtEffCorr->Fill(fTotNeutralEt,cent);
+ fHistTotRawEtEffCorr->Fill(fTotRawEtEffCorr,cent);
+ fHistTotRawEt->Fill(fTotRawEt,cent);
+ fHistTotAllRawEt->Fill(fTotAllRawEt,cent);
+ fHistTotAllRawEtVsTotalPt->Fill(fTotAllRawEt,totalPt);
+ fHistTotAllRawEtVsTotalPtVsCent->Fill(fTotAllRawEt,totalPt,cent);
+ fHistTotMatchedRawEtVsTotalPtVsCent->Fill(fTotAllRawEt,totalMatchedPt,cent);
+ fHistTotAllRawEtEffCorr->Fill(fTotAllRawEtEffCorr,cent);
+ //cout<<"fTotAllRawEtEffCorr "<<fTotAllRawEtEffCorr<<" fTotAllRawEt "<<fTotAllRawEt<<" fTotRawEtEffCorr "<<fTotRawEtEffCorr<<"("<<fTotNeutralEt<<")"<<" fTotRawEt "<<fTotRawEt<<endl;
+ //cout<<"uncorr "<<uncorrEt<<" raw "<<nominalRawEt<<" tot raw "<<fTotNeutralEt;
+ fTotNeutralEt = fGeomCorrection * fEMinCorrection * (fTotNeutralEt - removedEnergy);
+ //cout<<" tot corr "<<fTotNeutralEt<<endl;
fTotEt = fTotChargedEt + fTotNeutralEt;
- fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;
- fSparseEt[0] = fTotEt;
- fSparseEt[1] = fTotNeutralEt;
- fSparseEt[2] = fTotChargedEtAcc;
- fSparseEt[3] = fMultiplicity;
- fSparseEt[4] = fNeutralMultiplicity;
- fSparseEt[5] = fChargedMultiplicity;
- fSparseEt[6] = cent;
- // Fill the histograms...
+// Fill the histograms...0
FillHistograms();
-
+ //std::cout << "fTotNeutralEt: " << fTotNeutralEt << ", Contribution from non-removed charged: " << GetChargedContribution(fNeutralMultiplicity) << ", neutral: " << GetNeutralContribution(fNeutralMultiplicity) << ", gammas: " << GetGammaContribution(fNeutralMultiplicity) << ", multiplicity: " << fNeutralMultiplicity<< std::endl;
+ //cout<<"cent "<<cent<<" cluster mult "<<fClusterMult<<" fTotNeutralEt "<<fTotNeutralEt<<" nominalRawEt "<<nominalRawEt<<endl;
+ fHistNominalRawEt->Fill(nominalRawEt,cent);
+ fHistTotRawEtEffCorr500MeV->Fill(totEt500MeV,cent);
+ fHistNominalNonLinHighEt->Fill(nonlinHighRawEt,cent);
+ fHistNominalNonLinLowEt->Fill(nonlinLowRawEt,cent);
+ fHistNominalEffHighEt->Fill(effHighRawEt,cent);
+ fHistNominalEffLowEt->Fill(effLowRawEt,cent);
+ fHistFoundHadronsvsCent->Fill(nChargedHadronsMeasured,cent);
+ fHistNotFoundHadronsvsCent->Fill(nChargedHadronsTotal-nChargedHadronsMeasured,cent);
+ fHistFoundHadronsEtvsCent->Fill(nChargedHadronsEtMeasured,cent);
+ fHistNotFoundHadronsEtvsCent->Fill(nChargedHadronsEtTotal-nChargedHadronsEtMeasured,cent);
+ //cout<<"found "<<nChargedHadronsMeasured<<" total "<<nChargedHadronsTotal<<" not found "<<nChargedHadronsTotal-nChargedHadronsMeasured<<" found "<< nChargedHadronsMeasured500MeV<<" not found "<< nChargedHadronsTotal500MeV-nChargedHadronsMeasured500MeV <<" total "<<nChargedHadronsTotal500MeV<<endl;
+ fHistFoundHadronsvsCent500MeV->Fill(nChargedHadronsMeasured500MeV,cent);
+ fHistNotFoundHadronsvsCent500MeV->Fill(nChargedHadronsTotal500MeV-nChargedHadronsMeasured500MeV,cent);
+ fHistFoundHadronsEtvsCent500MeV->Fill(nChargedHadronsEtMeasured500MeV,cent);
+ fHistNotFoundHadronsEtvsCent500MeV->Fill(nChargedHadronsEtTotal500MeV-nChargedHadronsEtMeasured500MeV,cent);
+// cout<<"Number of hadrons measured: "<<nChargedHadronsMeasured<<" Estimated total number of hadrons "<<nChargedHadronsTotal<<" ET in track matched hadrons "<<
+// nChargedHadronsEtMeasured;
+// if(nChargedHadronsMeasured>0)cout<<" ("<<nChargedHadronsEtMeasured/nChargedHadronsMeasured<<") ";
+// cout<<" ET in all hadrons ";
+// cout<<nChargedHadronsEtTotal;
+// if(nChargedHadronsTotal>0) cout<<" ("<<nChargedHadronsEtTotal/nChargedHadronsTotal<<") ";
+// cout<<endl;
+ fHistPIDProtonsTrackMatchedDepositedVsNch->Fill(etPIDProtons,multiplicity);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNch->Fill(etPIDAntiProtons,multiplicity);
+ fHistPIDProtonsTrackMatchedDepositedVsNcl->Fill(etPIDProtons,nCluster);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNcl->Fill(etPIDAntiProtons,nCluster);
+ fHistPIDProtonsTrackMatchedDepositedVsNchNoEff->Fill(etPIDProtonsNoEff,multiplicity);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff->Fill(etPIDAntiProtonsNoEff,multiplicity);
+ fHistPIDProtonsTrackMatchedDepositedVsNclNoEff->Fill(etPIDProtonsNoEff,nCluster);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff->Fill(etPIDAntiProtonsNoEff,nCluster);
+ fHistCentVsNchVsNclReco->Fill(cent,multiplicity,nCluster);
+ fHistPiKPTrackMatchedDepositedVsNch->Fill(etPiKPMatched,multiplicity);
+ fHistPiKPTrackMatchedDepositedVsNchNoEff->Fill(etPiKPMatchedNoEff,multiplicity);
+ //delete pID;
return 0;
}
bool AliAnalysisEtReconstructed::TrackHitsCalorimeter(AliVParticle* track, Double_t magField)
{ // propagate track to detector radius
-
if (!track) {
cout<<"Warning: track empty"<<endl;
return kFALSE;
Bool_t prop = esdTrack->PropagateTo(fDetectorRadius, magField);
// if (prop) Printf("Track propagated, eta: %f, phi: %f, pt: %f", esdTrack->Eta(), esdTrack->Phi(), esdTrack->Pt());
- return prop &&
- TMath::Abs(esdTrack->Eta()) < fEtaCutAcc &&
- esdTrack->Phi() > fPhiCutAccMin*TMath::Pi()/180. &&
- esdTrack->Phi() < fPhiCutAccMax*TMath::Pi()/180.;
+ return prop && fSelector->CutGeometricalAcceptance(*esdTrack);
}
void AliAnalysisEtReconstructed::FillOutputList(TList* list)
list->Add(fHistMuonEnergyDeposit);
list->Add(fHistRemovedEnergy);
+ list->Add(fClusterPositionAccepted);
+ list->Add(fClusterPositionAll);
+ list->Add(fClusterPositionAcceptedEnergy);
+ list->Add(fClusterPositionAllEnergy);
+ list->Add(fClusterEnergy);
+ list->Add(fClusterEnergyCent);
+ list->Add(fClusterEnergyCentMatched);
+ list->Add(fClusterEnergyCentNotMatched);
+ list->Add(fClusterEt);
+
+ list->Add(fHistChargedEnergyRemoved);
+ list->Add(fHistNeutralEnergyRemoved);
+ list->Add(fHistGammaEnergyAdded);
+ list->Add(fHistMatchedTracksEvspTvsCent);
+ list->Add(fHistMatchedTracksEvspTvsCentEffCorr);
+ list->Add(fHistMatchedTracksEvspTvsCentEffTMCorr);
+ list->Add(fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr);
+ list->Add(fHistMatchedTracksEvspTvsCentEffTMCorr500MeV);
+ list->Add(fHistFoundHadronsvsCent);
+ list->Add(fHistNotFoundHadronsvsCent);
+ list->Add(fHistFoundHadronsEtvsCent);
+ list->Add(fHistNotFoundHadronsEtvsCent);
+ list->Add(fHistFoundHadronsvsCent500MeV);
+ list->Add(fHistNotFoundHadronsvsCent500MeV);
+ list->Add(fHistFoundHadronsEtvsCent500MeV);
+ list->Add(fHistNotFoundHadronsEtvsCent500MeV);
+ list->Add(fHistNominalRawEt);
+ list->Add(fHistNominalNonLinHighEt);
+ list->Add(fHistNominalNonLinLowEt);
+ list->Add(fHistNominalEffHighEt);
+ list->Add(fHistNominalEffLowEt);
+ list->Add(fHistTotRawEtEffCorr);
+ list->Add(fHistTotRawEtEffCorr500MeV);
+ list->Add(fHistTotAllRawEtEffCorr);
+ list->Add(fHistTotRawEt);
+ list->Add(fHistTotAllRawEt);
+ list->Add(fHistNClustersPhosVsEmcal);
+ list->Add(fHistClusterSizeVsCent);
+ list->Add(fHistMatchedClusterSizeVsCent);
+ list->Add(fHistTotAllRawEtVsTotalPt);
+ list->Add(fHistTotAllRawEtVsTotalPtVsCent);
+ list->Add(fHistTotMatchedRawEtVsTotalPtVsCent);
+ list->Add(fHistPIDProtonsTrackMatchedDepositedVsNch);
+ list->Add(fHistPIDAntiProtonsTrackMatchedDepositedVsNch);
+ list->Add(fHistPIDProtonsTrackMatchedDepositedVsNcl);
+ list->Add(fHistPIDAntiProtonsTrackMatchedDepositedVsNcl);
+ list->Add(fHistPiKPTrackMatchedDepositedVsNch);
+ list->Add(fHistPIDProtonsTrackMatchedDepositedVsNchNoEff);
+ list->Add(fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff);
+ list->Add(fHistPIDProtonsTrackMatchedDepositedVsNclNoEff);
+ list->Add(fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff);
+ list->Add(fHistPiKPTrackMatchedDepositedVsNchNoEff);
+ list->Add(fHistCentVsNchVsNclReco);
+ list->Add(fHistRawSignalReco);
+ list->Add(fHistEffCorrSignalReco);
+ list->Add(fHistRecoRCorrVsPtVsCent);
}
void AliAnalysisEtReconstructed::CreateHistograms()
Double_t maxEt = 10;
// possibly change histogram limits
- if (fCuts) {
- nbinsEt = fCuts->GetHistNbinsParticleEt();
- minEt = fCuts->GetHistMinParticleEt();
- maxEt = fCuts->GetHistMaxParticleEt();
- }
+// if (fCuts) {
+// nbinsEt = fCuts->GetHistNbinsParticleEt();
+// minEt = fCuts->GetHistMinParticleEt();
+// maxEt = fCuts->GetHistMaxParticleEt();
+// }
TString histname;
histname = "fHistChargedPionEnergyDeposit" + fHistogramNameSuffix;
//fHistMuonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
//fHistMuonEnergyDeposit->SetYTitle("Energy of track");
-
+ histname = "fClusterPositionAccepted" + fHistogramNameSuffix;
+ fClusterPositionAccepted = new TH2D(histname.Data(), "Position of accepted neutral clusters",300, -TMath::Pi(),TMath::Pi(), 100, -0.7 , 0.7);
+ fClusterPositionAccepted->SetXTitle("#phi");
+ fClusterPositionAccepted->SetYTitle("#eta");
+
+ histname = "fClusterPositionAll" + fHistogramNameSuffix;
+ fClusterPositionAll = new TH2D(histname.Data(), "Position of accepted neutral clusters",300, -TMath::Pi(),TMath::Pi(), 100, -0.7 , 0.7);
+ fClusterPositionAll->SetXTitle("#phi");
+ fClusterPositionAll->SetYTitle("#eta");
+
+ histname = "fClusterPositionAcceptedEnergy" + fHistogramNameSuffix;
+ fClusterPositionAcceptedEnergy = new TH2D(histname.Data(), "Position of accepted neutral clusters",300, -TMath::Pi(),TMath::Pi(), 100, -0.7 , 0.7);
+ fClusterPositionAcceptedEnergy->SetXTitle("#phi");
+ fClusterPositionAcceptedEnergy->SetYTitle("#eta");
+
+ histname = "fClusterPositionAllEnergy" + fHistogramNameSuffix;
+ fClusterPositionAllEnergy = new TH2D(histname.Data(), "Position of accepted neutral clusters",300, -TMath::Pi(),TMath::Pi(), 100, -0.7 , 0.7);
+ fClusterPositionAllEnergy->SetXTitle("#phi");
+ fClusterPositionAllEnergy->SetYTitle("#eta");
+
+ histname = "fClusterEnergy" + fHistogramNameSuffix;
+ fClusterEnergy = new TH1F(histname.Data(), histname.Data(), 100, 0, 5);
+ fClusterEnergy->SetYTitle("Number of clusters");
+ fClusterEnergy->SetXTitle("Energy of cluster");
+
+ histname = "fClusterEnergyCent" + fHistogramNameSuffix;
+ fClusterEnergyCent = new TH2F(histname.Data(), histname.Data(), 100, 0, 5,20,-0.5,19.5);
+ fClusterEnergyCent->SetXTitle("Energy of cluster");
+ fClusterEnergyCent->SetYTitle("Centrality Bin");
+ fClusterEnergyCent->SetZTitle("Number of clusters");
+
+ histname = "fClusterEnergyCentMatched" + fHistogramNameSuffix;
+ fClusterEnergyCentMatched = new TH2F(histname.Data(), histname.Data(), 100, 0, 5,20,-0.5,19.5);
+ fClusterEnergyCentMatched->SetXTitle("Energy of cluster");
+ fClusterEnergyCentMatched->SetYTitle("Centrality Bin");
+ fClusterEnergyCentMatched->SetZTitle("Number of Clusters");
+
+ histname = "fClusterEnergyCentNotMatched" + fHistogramNameSuffix;
+ fClusterEnergyCentNotMatched = new TH2F(histname.Data(), histname.Data(), 100, 0, 5,20,-0.5,19.5);
+ fClusterEnergyCentNotMatched->SetXTitle("Energy of cluster");
+ fClusterEnergyCentNotMatched->SetYTitle("Centrality Bin");
+ fClusterEnergyCentNotMatched->SetZTitle("Number of clusters");
+
+ histname = "fClusterEt" + fHistogramNameSuffix;
+ fClusterEt = new TH1F(histname.Data(), histname.Data(), 100, 0, 5);
+ fClusterEt->SetXTitle("Number of clusters");
+ fClusterEt->SetYTitle("E_{T} of cluster");
+
+ histname = "fHistChargedEnergyRemoved" + fHistogramNameSuffix;
+ fHistChargedEnergyRemoved = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
+
+ histname = "fHistNeutralEnergyRemoved" + fHistogramNameSuffix;
+ fHistNeutralEnergyRemoved = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
+
+ histname = "fHistGammaEnergyAdded" + fHistogramNameSuffix;
+ fHistGammaEnergyAdded = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
+
+ fHistMatchedTracksEvspTvsCent = new TH3F("fHistMatchedTracksEvspTvsCent", "fHistMatchedTracksEvspTvsCent",100, 0, 3,100,0,3,20,-0.5,19.5);
+ fHistMatchedTracksEvspTvsCentEffCorr = new TH3F("fHistMatchedTracksEvspTvsCentEffCorr", "fHistMatchedTracksEvspTvsCentEffCorr",100, 0, 3,100,0,3,20,-0.5,19.5);
+ fHistMatchedTracksEvspTvsCentEffTMCorr = new TH3F("fHistMatchedTracksEvspTvsCentEffTMCorr", "fHistMatchedTracksEvspTvsCentEffTMCorr",100, 0, 3,100,0,3,20,-0.5,19.5);
+ fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr = new TH3F("fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr", "fHistPeripheralMatchedTracksEvspTvsCentEffTMCorr",100, 0, 3,100,0,3,20,-0.5,19.5);
+ fHistMatchedTracksEvspTvsCentEffTMCorr500MeV = new TH3F("fHistMatchedTracksEvspTvsCentEffTMCorr500MeV", "fHistMatchedTracksEvspTvsCentEffTMCorr500MeV",100, 0, 3,100,0,3,20,-0.5,19.5);
+
+ float max = 200;
+ if(fHistogramNameSuffix.Contains("P")){max = 100;}
+ fHistFoundHadronsvsCent = new TH2F("fHistFoundHadronsvsCent","fHistFoundHadronsvsCent",100,0,max,20,-0.5,19.5);
+ fHistNotFoundHadronsvsCent = new TH2F("fHistNotFoundHadronsvsCent","fHistNotFoundHadronsvsCent",100,0,max,20,-0.5,19.5);
+ fHistFoundHadronsEtvsCent = new TH2F("fHistFoundHadronsEtvsCent","fHistFoundHadronsEtvsCent",100,0,max,20,-0.5,19.5);
+ fHistNotFoundHadronsEtvsCent = new TH2F("fHistNotFoundHadronsEtvsCent","fHistNotFoundHadronsEtvsCent",100,0,max,20,-0.5,19.5);
+ fHistFoundHadronsvsCent500MeV = new TH2F("fHistFoundHadronsvsCent500MeV","fHistFoundHadronsvsCent500MeV",100,0,max,20,-0.5,19.5);
+ fHistNotFoundHadronsvsCent500MeV = new TH2F("fHistNotFoundHadronsvsCent500MeV","fHistNotFoundHadronsvsCent500MeV",100,0,max,20,-0.5,19.5);
+ fHistFoundHadronsEtvsCent500MeV = new TH2F("fHistFoundHadronsEtvsCent500MeV","fHistFoundHadronsEtvsCent500MeV",100,0,max,20,-0.5,19.5);
+ fHistNotFoundHadronsEtvsCent500MeV = new TH2F("fHistNotFoundHadronsEtvsCent500MeV","fHistNotFoundHadronsEtvsCent500MeV",100,0,max,20,-0.5,19.5);
+
+ fHistTotRawEtEffCorr = new TH2F("fHistTotRawEtEffCorr","fHistTotRawEtEffCorr",250,0,250,20,-0.5,19.5);
+ fHistTotRawEt = new TH2F("fHistTotRawEt","fHistTotRawEt",250,0,250,20,-0.5,19.5);
+ fHistTotRawEtEffCorr500MeV = new TH2F("fHistTotRawEtEffCorr500MeV","fHistTotRawEtEffCorr500MeV",250,0,250,20,-0.5,19.5);
+ fHistTotAllRawEt = new TH2F("fHistTotAllRawEt","fHistTotAllRawEt",250,0,250,20,-0.5,19.5);
+ fHistTotAllRawEtEffCorr = new TH2F("fHistTotAllRawEtEffCorr","fHistTotAllRawEtEffCorr",250,0,250,20,-0.5,19.5);
+ fHistNClustersPhosVsEmcal = new TH3F("fHistNClustersPhosVsEmcal","fHistNClustersPhosVsEmcal",50,0,50,250,0,250,20,-0.5,19);
+ fHistClusterSizeVsCent = new TH2F("fHistClusterSizeVsCent","fHistClusterSizeVsCent",10,0.5,10.5,20,-0.5,19.5);
+ fHistMatchedClusterSizeVsCent = new TH2F("fHistMatchedClusterSizeVsCent","fHistMatchedClusterSizeVsCent",10,0.5,10.5,20,-0.5,19.5);
+ fHistTotAllRawEtVsTotalPt = new TH2F("fHistTotAllRawEtVsTotalPt","fHistTotAllRawEtVsTotalPt",125,0,250,200,0,2000);
+ fHistTotAllRawEtVsTotalPtVsCent = new TH3F("fHistTotAllRawEtVsTotalPtVsCent","fHistTotAllRawEtVsTotalPtVsCent",125,0,250,200,0,2000,20,-0.5,19.5);
+ fHistTotMatchedRawEtVsTotalPtVsCent = new TH3F("fHistTotMatchedRawEtVsTotalPtVsCent","fHistTotMatchedRawEtVsTotalPtVsCent",250,0,250,100,0,200,20,-0.5,19.5);
+
+ maxEt = 500;
+ histname = "fHistNominalRawEt" + fHistogramNameSuffix;
+ fHistNominalRawEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
+ histname = "fHistNominalNonLinHighEt" + fHistogramNameSuffix;
+ fHistNominalNonLinHighEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
+ histname = "fHistNominalNonLinLowEt" + fHistogramNameSuffix;
+ fHistNominalNonLinLowEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
+ histname = "fHistNominalEffHighEt" + fHistogramNameSuffix;
+ fHistNominalEffHighEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
+ histname = "fHistNominalEffLowEt" + fHistogramNameSuffix;
+ fHistNominalEffLowEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
+
+ Float_t maxEtRange = 25;
+ Float_t maxEtRangeHigh = 125;
+ Float_t minEtRange = 0;
+ Int_t nbinsMult = 100;
+ Float_t maxMult = 3000;
+ Float_t minMult = 0;
+ Int_t nbinsCl = 250;
+ Float_t maxCl = 500;
+ Float_t minCl = 0;
+ fHistPIDProtonsTrackMatchedDepositedVsNch = new TH2F("fHistPIDProtonsTrackMatchedDepositedVsNch","PID'd protons deposited in calorimeter vs multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsMult,minMult,maxMult);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNch = new TH2F("fHistPIDAntiProtonsTrackMatchedDepositedVsNch","PID'd #bar{p} E_{T} deposited in calorimeter vs multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsMult,minMult,maxMult);
+ fHistPIDProtonsTrackMatchedDepositedVsNcl = new TH2F("fHistPIDProtonsTrackMatchedDepositedVsNcl","PID'd protons deposited in calorimeter vs cluster multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsCl,minCl,maxCl);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNcl = new TH2F("fHistPIDAntiProtonsTrackMatchedDepositedVsNcl","PID'd #bar{p} E_{T} deposited in calorimeter vs cluster multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsCl,minCl,maxCl);
+ fHistPiKPTrackMatchedDepositedVsNch = new TH2F("fHistPiKPTrackMatchedDepositedVsNch","PiKP track matched",nbinsEt,minEtRange,maxEtRangeHigh,nbinsMult,minMult,maxMult);
+
+ fHistPIDProtonsTrackMatchedDepositedVsNchNoEff = new TH2F("fHistPIDProtonsTrackMatchedDepositedVsNchNoEff","PID'd protons deposited in calorimeter vs multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsMult,minMult,maxMult);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff = new TH2F("fHistPIDAntiProtonsTrackMatchedDepositedVsNchNoEff","PID'd #bar{p} E_{T} deposited in calorimeter vs multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsMult,minMult,maxMult);
+ fHistPIDProtonsTrackMatchedDepositedVsNclNoEff = new TH2F("fHistPIDProtonsTrackMatchedDepositedVsNclNoEff","PID'd protons deposited in calorimeter vs cluster multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsCl,minCl,maxCl);
+ fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff = new TH2F("fHistPIDAntiProtonsTrackMatchedDepositedVsNclNoEff","PID'd #bar{p} E_{T} deposited in calorimeter vs cluster multiplicity",nbinsEt,minEtRange,maxEtRange,nbinsCl,minCl,maxCl);
+ fHistPiKPTrackMatchedDepositedVsNchNoEff = new TH2F("fHistPiKPTrackMatchedDepositedVsNchNoEff","PiKP track matched",nbinsEt,minEtRange,maxEtRangeHigh,nbinsMult,minMult,maxMult);
+
+
+ fHistCentVsNchVsNclReco = new TH3F("fHistCentVsNchVsNclReco","Cent bin vs Nch Vs NCl",20,-0.5,19.5,nbinsMult,minMult,maxMult,nbinsCl,minCl,maxCl);
+
+ fHistRawSignalReco = new TH1F("fHistRawSignalReco","fHistRawSignalReco",20,-0.5,19.5);
+ fHistEffCorrSignalReco = new TH1F("fHistEffCorrSignalReco","fHistEffCorrSignalReco",20,-0.5,19.5);
+ fHistRecoRCorrVsPtVsCent = new TH3F("fHistRecoRCorrVsPtVsCent","fHistRecoRCorrVsPtVsCent",72,0,2,50,0,10,20,-0.5,19.5);
}
-
-Double_t
-AliAnalysisEtReconstructed::CalcTrackClusterDistance(const Float_t clsPos[3],
- Int_t *trkMatchId,
- const AliESDEvent *event)
-{ // calculate distance between cluster and closest track
-
- Double_t trkPos[3] = {0,0,0};
-
- Int_t bestTrkMatchId = -1;
- Double_t distance = 9999; // init to a big number
-
- Double_t dist = 0;
- Double_t distX = 0, distY = 0, distZ = 0;
-
- for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++) {
- AliESDtrack *track = event->GetTrack(iTrack);
- if (!track) {
- AliError(Form("ERROR: Could not get track %d", iTrack));
- continue;
- }
-
- // check for approx. eta and phi range before we propagate..
- // TBD
-
- AliEMCALTrack *emctrack = new AliEMCALTrack(*track);
- if (!emctrack->PropagateToGlobal(clsPos[0],clsPos[1],clsPos[2],0.,0.) ) {
- continue;
- }
- emctrack->GetXYZ(trkPos);
- delete emctrack;
-
- distX = clsPos[0]-trkPos[0];
- distY = clsPos[1]-trkPos[1];
- distZ = clsPos[2]-trkPos[2];
- dist = TMath::Sqrt(distX*distX + distY*distY + distZ*distZ);
-
- if (dist < distance) {
- distance = dist;
- bestTrkMatchId = iTrack;
- }
- } // iTrack
-
- // printf("CalcTrackClusterDistance: bestTrkMatch %d origTrkMatch %d distance %f\n", bestTrkMatchId, *trkMatchId, distance);
- *trkMatchId = bestTrkMatchId;
- return distance;
+Double_t AliAnalysisEtReconstructed::ApplyModifiedCorrections(const AliESDCaloCluster& cluster,Int_t nonLinCorr, Int_t effCorr, Int_t cent)
+{
+ Float_t pos[3];
+ cluster.GetPosition(pos);
+ TVector3 cp(pos);
+ Double_t corrEnergy = fReCorrections->CorrectedEnergy(cluster.E(),cent);
+
+ Double_t factorNonLin = GetCorrectionModification(cluster, nonLinCorr,effCorr,cent);
+
+ cout<<"Warning: This function should not get called!"<<endl;
+ //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
+ return TMath::Sin(cp.Theta())*corrEnergy*factorNonLin;
}
+Double_t AliAnalysisEtReconstructed::GetCorrectionModification(const AliESDCaloCluster& cluster,Int_t nonLinCorr, Int_t effCorr, Int_t cent){//nonLinCorr 0 = nominal 1 = high -1 = low, effCorr 0 = nominal 1 = high -1 = low
+ if(nonLinCorr==0){
+ cout<<"Warning: This function should not get called!"<<endl;//this statement is basically here to avoid a compilation warning
+ }
+ if(effCorr==0){
+ cout<<"Warning: This function should not get called!"<<endl;//this statement is basically here to avoid a compilation warning
+ }
+ return cluster.E()*cent;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff