X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWGGA%2FEMCALJetTasks%2FAliAnalysisTaskRho.cxx;h=815bf4704691796c2c27567b20da512c70751f62;hb=9148fa8984cfa9dd128d8aeba1d9cce576c7950b;hp=54e564ed9cbb9eab81f3afffd31905f294ec4f76;hpb=c60e0a211965a410f83f6e1be4c0a59e1a19a088;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWGGA/EMCALJetTasks/AliAnalysisTaskRho.cxx b/PWGGA/EMCALJetTasks/AliAnalysisTaskRho.cxx index 54e564ed9cb..815bf470469 100644 --- a/PWGGA/EMCALJetTasks/AliAnalysisTaskRho.cxx +++ b/PWGGA/EMCALJetTasks/AliAnalysisTaskRho.cxx @@ -1,42 +1,47 @@ -// $Id: $ +// $Id$ // -// Calculation of rho +// Calculation of rho from a collection of jets. +// If scale function is given the scaled rho will be exported +// with the name as "fRhoName".apppend("_Scaled"). // // Authors: R.Reed, S.Aiola -#include -#include -#include +#include "AliAnalysisTaskRho.h" + #include #include +#include +#include +#include #include -#include "AliLog.h" #include "AliAnalysisManager.h" -#include "AliVEventHandler.h" #include "AliCentrality.h" #include "AliEmcalJet.h" +#include "AliLog.h" +#include "AliRhoParameter.h" #include "AliVCluster.h" - -#include "AliAnalysisTaskRho.h" +#include "AliVEventHandler.h" ClassImp(AliAnalysisTaskRho) //________________________________________________________________________ AliAnalysisTaskRho::AliAnalysisTaskRho() : AliAnalysisTaskRhoBase(), - fTracksName("tracks"), - fJetsName("KtJets"), - fClustersName("caloClusters"), + fTracksName(), + fJetsName(), fRhoScaledName(""), fPhiMin(0), fPhiMax(0), fEtaMin(0), fEtaMax(0), fAreaCut(0), + fAreaEmcCut(0), fNExclLeadJets(0), fScaleFunction(0), fCreateHisto(kFALSE), + fTracks(0), + fJets(0), fOutputList(0), fHistCentrality(0), fHistJetPt(0), @@ -55,23 +60,26 @@ AliAnalysisTaskRho::AliAnalysisTaskRho() : fHistNjetvsNtrack(0), fRhoScaled(0) { - // Constructor + // Constructor. } + //________________________________________________________________________ -AliAnalysisTaskRho::AliAnalysisTaskRho(const char *name) : +AliAnalysisTaskRho::AliAnalysisTaskRho(const char *name, Bool_t histo) : AliAnalysisTaskRhoBase(name), fTracksName("tracks"), fJetsName("KtJets"), - fClustersName("caloClusters"), fRhoScaledName(""), fPhiMin(0), - fPhiMax(0), - fEtaMin(0), - fEtaMax(0), - fAreaCut(0), - fNExclLeadJets(0), + fPhiMax(TMath::TwoPi()), + fEtaMin(-0.5), + fEtaMax(+0.5), + fAreaCut(0.01), + fAreaEmcCut(0), + fNExclLeadJets(1), fScaleFunction(0), - fCreateHisto(kFALSE), + fCreateHisto(histo), + fTracks(0), + fJets(0), fOutputList(0), fHistCentrality(0), fHistJetPt(0), @@ -90,29 +98,32 @@ AliAnalysisTaskRho::AliAnalysisTaskRho(const char *name) : fHistNjetvsNtrack(0), fRhoScaled(0) { - // Constructor + // Constructor. - DefineOutput(1, TList::Class()); + if (fCreateHisto) + DefineOutput(1, TList::Class()); } //________________________________________________________________________ void AliAnalysisTaskRho::UserCreateOutputObjects() { + // User create output objects, called at the beginning of the analysis. + AliAnalysisTaskRhoBase::UserCreateOutputObjects(); - fRhoScaledName = fRhoName; - fRhoScaledName += "_Scaled"; - fRhoScaled = new TParameter(fRhoScaledName, 0); + if (fScaleFunction) { + fRhoScaledName = fRhoName; + fRhoScaledName += "_Scaled"; + fRhoScaled = new AliRhoParameter(fRhoScaledName, 0); + } + + if (!fCreateHisto) + return; OpenFile(1); fOutputList = new TList(); fOutputList->SetOwner(); - if (!fCreateHisto) { - PostData(1, fOutputList); - return; - } - fHistCentrality = new TH1F("Centrality", "Centrality", 101, -1, 100); fHistRhovsCent = new TH2F("RhovsCent", "RhovsCent", 101, -1, 100, 500, 0, 500); fHistDeltaRhovsCent = new TH2F("DeltaRhovsCent", "DetlaRhovsCent", 101, -1, 100, 500, -250, 250); @@ -149,83 +160,55 @@ void AliAnalysisTaskRho::UserCreateOutputObjects() fOutputList->Add(fHistNjetvsNtrack); PostData(1, fOutputList); - } -//________________________________________________________________________ -Double_t AliAnalysisTaskRho::GetScaleFactor(Double_t cent) -{ - Double_t scale = 1; - if (fScaleFunction) - scale = fScaleFunction->Eval(cent); - return scale; -} - - //________________________________________________________________________ void AliAnalysisTaskRho::UserExec(Option_t *) { // Main loop, called for each event. - - AliAnalysisTaskRhoBase::UserExec(""); - - fRho->SetVal(-1); - // add rho to event if not yet there - if (!(InputEvent()->FindListObject(fRhoScaledName))) { - new(fRhoScaled) TParameter(fRhoScaledName, -1); - InputEvent()->AddObject(fRhoScaled); - } - else { - fRhoScaled->SetVal(-1); + if (!fIsInit) { + ExecOnce(); + fIsInit = 1; } - // optimization in case autobranch loading is off - AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); - if (fTracksName == "Tracks") - am->LoadBranch("Tracks"); + if (!fJets) + return; + + fRho->SetVal(-1); + if (fRhoScaled) + fRhoScaled->SetVal(-1); - TClonesArray *jets = 0; - TClonesArray *tracks = 0; + DetermineCent(); - tracks = dynamic_cast(InputEvent()->FindListObject(fTracksName)); - if (!tracks) { - AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() )); - return; - } - - jets = dynamic_cast(InputEvent()->FindListObject(fJetsName)); - if (!jets) { - AliError(Form("Pointer to jets %s == 0", fJetsName.Data() )); - return; - } - - if (fCreateHisto) - fHistCentrality->Fill(fCent); + const Int_t Njets = fJets->GetEntries(); - const Int_t Ntracks = tracks->GetEntries(); - const Int_t Njets = jets->GetEntries(); + Int_t maxJetIds[] = {-1, -1}; + Float_t maxJetPts[] = { 0, 0}; - Int_t maxJetIds[] = {-1, -1}; - Int_t maxJetPts[] = {0, 0}; if (fNExclLeadJets > 0) { - for (Int_t ij = 0; ij < Njets; ij++) { - - AliEmcalJet *jet = static_cast(jets->At(ij)); - + for (Int_t ij = 0; ij < Njets; ++ij) { + AliEmcalJet *jet = static_cast(fJets->At(ij)); if (!jet) { - AliError(Form("Could not receive jet %d", ij)); + AliError(Form("%s: Could not receive jet %d", GetName(), ij)); continue; } - + + if (jet->Area() < fAreaCut) + continue; + if (jet->AreaEmc()Phi() < fPhiMin) || (jet->Phi() > fPhiMax)) + continue; + if ((jet->Eta() < fEtaMin) || (jet->Eta() > fEtaMax)) + continue; + if (jet->Pt() > maxJetPts[0]) { maxJetPts[1] = maxJetPts[0]; - maxJetIds[1] = maxJetPts[0]; + maxJetIds[1] = maxJetIds[0]; maxJetPts[0] = jet->Pt(); maxJetIds[0] = ij; - } - - if (jet->Pt() > maxJetPts[1]) { + } else if (jet->Pt() > maxJetPts[1]) { maxJetPts[1] = jet->Pt(); maxJetIds[1] = ij; } @@ -241,32 +224,32 @@ void AliAnalysisTaskRho::UserExec(Option_t *) // push all jets within selected acceptance into stack for (Int_t iJets = 0; iJets < Njets; ++iJets) { - + // exlcuding lead jets if (iJets == maxJetIds[0] || iJets == maxJetIds[1]) continue; - AliEmcalJet *jet = static_cast(jets->At(iJets)); - - if (!jet) - continue; + AliEmcalJet *jet = static_cast(fJets->At(iJets)); + if (!jet) { + AliError(Form("%s: Could not receive jet %d", GetName(), iJets)); + continue; + } - // applying some other cuts if (jet->Area() < fAreaCut) continue; + if (jet->AreaEmc()Phi() < fPhiMin) || (jet->Phi() > fPhiMax)) continue; if ((jet->Eta() < fEtaMin) || (jet->Eta() > fEtaMax)) continue; - if (jet->Area() == 0) - continue; rhovec[NjetAcc] = jet->Pt() / jet->Area(); - - NjetAcc++; + ++NjetAcc; if (fCreateHisto) { // filling histograms + const Int_t Ntracks = fTracks->GetEntriesFast(); fHistJetPt->Fill(jet->Pt()); fHistJetArea->Fill(jet->Area()); fHistJetPtvsCent->Fill(fCent, jet->Pt()); @@ -275,44 +258,78 @@ void AliAnalysisTaskRho::UserExec(Option_t *) fHistJetAreavsNtrack->Fill(Ntracks, jet->Area()); } } - - if (fCreateHisto) { - fHistNjetvsCent->Fill(fCent, NjetAcc); - fHistNjetvsNtrack->Fill(Ntracks, NjetAcc); - } - - Double_t scale = GetScaleFactor(fCent); - Double_t rhochem = GetRhoFactor(fCent); - Double_t rho0 = -1; - - if (NjetAcc > 0){ + if (NjetAcc > 0) { //find median value - rho0 = TMath::Median(NjetAcc, rhovec); - - Double_t rhoScaled = rho0 * scale; - + Double_t rho0 = TMath::Median(NjetAcc, rhovec); fRho->SetVal(rho0); - fRhoScaled->SetVal(rhoScaled); + Double_t rhoScaled = rho0; + if (fRhoScaled) { + rhoScaled *= GetScaleFactor(fCent); + fRhoScaled->SetVal(rhoScaled); + } if (fCreateHisto) { // filling other histograms + Double_t rho = GetRhoFactor(fCent); + const Int_t Ntracks = fTracks->GetEntriesFast(); fHistRhovsCent->Fill(fCent, rho0); - fHistDeltaRhovsCent->Fill(fCent, rho0 - rhochem); - fHistDeltaRhoScalevsCent->Fill(fCent, rhoScaled - rhochem); + fHistDeltaRhovsCent->Fill(fCent, rho0 - rho); + fHistDeltaRhoScalevsCent->Fill(fCent, rhoScaled - rho); fHistRhovsNtrack->Fill(Ntracks, rho0); - fHistDeltaRhovsNtrack->Fill(Ntracks, rho0 - rhochem); - fHistDeltaRhoScalevsNtrack->Fill(Ntracks, rhoScaled - rhochem); + fHistDeltaRhovsNtrack->Fill(Ntracks, rho0 - rho); + fHistDeltaRhoScalevsNtrack->Fill(Ntracks, rhoScaled - rho); } } - if (fCreateHisto) + if (fCreateHisto) { + const Int_t Ntracks = fTracks->GetEntriesFast(); + fHistCentrality->Fill(fCent); + fHistNjetvsCent->Fill(fCent, NjetAcc); + fHistNjetvsNtrack->Fill(Ntracks, NjetAcc); PostData(1, fOutputList); + } } +//________________________________________________________________________ +void AliAnalysisTaskRho::ExecOnce() +{ + // Initialize some settings that need to be determined in UserExec. + + AliAnalysisTaskRhoBase::ExecOnce(); + + if (fRhoScaled) { + // add rho to event if not yet there + if (!(InputEvent()->FindListObject(fRhoScaledName))) { + InputEvent()->AddObject(fRhoScaled); + } else { + AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fRhoScaledName.Data())); + return; + } + } + + if (fCreateHisto) { + fTracks = dynamic_cast(InputEvent()->FindListObject(fTracksName)); + if (!fTracks) { + AliError(Form("%s: Pointer to tracks %s == 0", GetName(), fTracksName.Data() )); + return; + } + } + + fJets = dynamic_cast(InputEvent()->FindListObject(fJetsName)); + if (!fJets) { + AliError(Form("%s: Pointer to jets %s == 0", GetName(), fJetsName.Data() )); + return; + } +} //________________________________________________________________________ -void AliAnalysisTaskRho::Terminate(Option_t *) +Double_t AliAnalysisTaskRho::GetScaleFactor(Double_t cent) { + // Get scale factor. + Double_t scale = 1; + if (fScaleFunction) + scale = fScaleFunction->Eval(cent); + return scale; }