#include "AliCentrality.h"
#include "AliESDEvent.h"
#include "AliVTrack.h"
-#include "AliESDtrackCuts.h"
+#include "AliESDtrack.h"
#include "AliAODTrack.h"
+#include "AliAnalysisFilter.h"
+#include "AliAODMCHeader.h"
+#include "AliForwardFlowWeights.h"
ClassImp(AliForwardFlowTaskQC)
#if 0
fFMDCut(-1), // FMD sigma cut
fSPDCut(-1), // SPD sigma cut
fFlowFlags(0), // Flow flags
- fEtaGap(-1), // Eta gap value
+ fEtaGap(0.), // Eta gap value
fBinsForward(), // List with forward flow hists
fBinsCentral(), // List with central flow hists
fSumList(0), // Event sum list
fOutputList(0), // Result output list
fAOD(0), // AOD input event
- fESDTrackCuts(0), // ESD track cuts
+ fTrackCuts(0), // ESD track cuts
fMaxMoment(0), // Max flow moment
fVtx(1111), // Z vertex coordinate
fCent(-1), // Centrality
fCentAxis(), // Axis to control centrality/multiplicity binning
fFMDCut(-1), // FMD sigma cut
fSPDCut(-1), // SPD sigma cut
- fFlowFlags(kSymEta|kStdQC), // Flow flags
- fEtaGap(2.), // Eta gap value
+ fFlowFlags(0x0), // Flow flags
+ fEtaGap(0.), // Eta gap value
fBinsForward(), // List with forward flow hists
fBinsCentral(), // List with central flow hists
fSumList(0), // Event sum list
fOutputList(0), // Result output list
fAOD(0), // AOD input event
- fESDTrackCuts(0), // ESD track cuts
+ fTrackCuts(0), // ESD track cuts
fMaxMoment(4), // Max flow moment
fVtx(1111), // Z vertex coordinate
fCent(-1), // Centrality
fSumList(o.fSumList), // Event sum list
fOutputList(o.fOutputList), // Result output list
fAOD(o.fAOD), // AOD input event
- fESDTrackCuts(o.fESDTrackCuts), // ESD track cuts
+ fTrackCuts(o.fTrackCuts), // ESD track cuts
fMaxMoment(o.fMaxMoment), // Flow moments
fVtx(o.fVtx), // Z vertex coordinate
fCent(o.fCent), // Centrality
fSumList = o.fSumList;
fOutputList = o.fOutputList;
fAOD = o.fAOD;
- fESDTrackCuts = o.fESDTrackCuts;
+ fTrackCuts = o.fTrackCuts;
fMaxMoment = o.fMaxMoment;
fVtx = o.fVtx;
fCent = o.fCent;
//
InitVertexBins();
InitHists();
- if (fFlowFlags & kTPC) {
- fESDTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
- }
+ if ((fFlowFlags & kTracks) && !fTrackCuts) AliFatal("No track cuts set!");
PrintFlowSetup();
PostData(1, fSumList);
Int_t vH = Int_t(fVtxAxis->GetBinUpEdge(v));
if ((fFlowFlags & kFMD)) {
fBinsForward.Add(new VertexBin(vL, vH, fMaxMoment, "FMD", fFlowFlags, fFMDCut, fEtaGap));
- if (!(fFlowFlags & k3Cor)) fBinsCentral.Add(new VertexBin(vL, vH, fMaxMoment, "SPD-FMD", fFlowFlags|kNUAcorr, fSPDCut, fEtaGap));
+ if (!(fFlowFlags & k3Cor))
+ fBinsCentral.Add(new VertexBin(vL, vH, fMaxMoment, "SPD-FMD", fFlowFlags|kNUAcorr|kSPD, fSPDCut, fEtaGap));
}
else if ((fFlowFlags & kVZERO)) {
fBinsForward.Add(new VertexBin(vL, vH, fMaxMoment, "VZERO", fFlowFlags, 0, fEtaGap));
- if (!(fFlowFlags & k3Cor)) fBinsCentral.Add(new VertexBin(vL, vH, fMaxMoment, "SPD-VZERO", fFlowFlags|kNUAcorr, fSPDCut, fEtaGap));
+ if ((fFlowFlags & kEtaGap) && !(fFlowFlags & kTracks))
+ fBinsCentral.Add(new VertexBin(vL, vH, fMaxMoment, "SPD-VZERO", fFlowFlags|kNUAcorr|kSPD, fSPDCut, fEtaGap));
}
}
}
}
else if ((fFlowFlags & kVZERO) && !vzero) {
fHistEventSel->Fill(kNoForward);
- return kFALSE;
+// return kFALSE;
}
if (!aodcmult) fHistEventSel->Fill(kNoCentral);
Int_t nVtxBins = fVtxAxis->GetNbins();
while ((bin = static_cast<VertexBin*>(list.At(vtx+(nVtxBins*i))))) {
+ i++;
// If no tracks do things normally
- if (!(fFlowFlags & kTPC) && !bin->FillHists(h, fCent, kFillBoth|flags|kReset)) return;
+ if (!(fFlowFlags & kTracks) || (flags & kMC)) {
+ if (!bin->FillHists(h, fCent, kFillBoth|flags|kReset)) continue;
+ }
// if tracks things are more complicated
- else if ((fFlowFlags & kTPC)) {
- TObjArray* trList = GetTracks();
- if (!trList) return;
- Bool_t useEvent = bin->FillTracks(trList, kFillRef|kReset|flags);
- // If esd input trList is a new object owned by this task and should be cleaned up
- if (AliForwardUtil::CheckForAOD() == 2) delete trList;
- if (!useEvent) return;
- if (!bin->FillHists(h, fCent, kFillDiff|kReset|flags)) return;
+ else if ((fFlowFlags & kTracks)) {
+ if (!FillTracks(bin, kFillRef|kReset|flags)) continue;
+ if (!bin->FillHists(h, fCent, kFillDiff|kReset|flags)) continue;
}
bin->CumulantsAccumulate(fCent);
- i++;
}
return;
Int_t nVtxBins = fVtxAxis->GetNbins();
while ((bin = static_cast<VertexBin*>(list.At(vtx+(nVtxBins*i))))) {
- if (!bin->FillHists(href, fCent, kFillRef|flags|kReset)) return;
- bin->FillHists(hdiff, fCent, kFillDiff|kReset);
- bin->CumulantsAccumulate(fCent);
i++;
+ if (!(fFlowFlags & kTracks) || (flags & kMC)) {
+ if(!bin->FillHists(href, fCent, kFillRef|flags|kReset)) continue;
+ }
+ else if ((fFlowFlags & kTracks)) {
+ if (!FillTracks(bin, kFillRef|kReset|flags)) continue;
+ }
+ if (!bin->FillHists(hdiff, fCent, kFillDiff|kReset|flags)) continue;
+ bin->CumulantsAccumulate(fCent);
}
return;
TH2D& h = CombineHists(hcent, hfwd);
while ((bin = static_cast<VertexBin*>(list.At(vtx+(nVtxBins*i))))) {
- if (!bin->FillHists(h, fCent, kFillBoth|flags|kReset)) return;
- bin->CumulantsAccumulate3Cor(fCent);
i++;
+ if (!bin->FillHists(h, fCent, kFillBoth|flags|kReset)) continue;
+ bin->CumulantsAccumulate3Cor(fCent);
}
return;
return fHistdNdedp3Cor;
}
//_____________________________________________________________________
-TObjArray* AliForwardFlowTaskQC::GetTracks() const
+Bool_t AliForwardFlowTaskQC::FillTracks(VertexBin* bin, UShort_t mode) const
{
//
// Get TPC tracks to use for reference flow.
// Return: TObjArray with tracks
//
TObjArray* trList = 0;
- // Get input type
- UShort_t input = AliForwardUtil::CheckForAOD();
- switch (input) {
- // If AOD input, simply get the track array from the event
- case 1: trList = static_cast<TObjArray*>(fAOD->GetTracks());
- break;
- case 2: {
- // If ESD input get event, apply track cuts
- AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
- if (!esd) return 0;
- // Warning! trList is now a new array, we need to delete it after use
- // this is not a very good implementation!
- trList = fESDTrackCuts->GetAcceptedTracks(esd, kTRUE);
- break;
- }
- default: AliFatal("Neither ESD or AOD input. This should never happen");
- break;
- }
- return trList;
+ AliESDEvent* esdEv = 0;
+ if (AliForwardUtil::CheckForAOD() == 1) // AOD tracks
+ trList = static_cast<TObjArray*>(fAOD->GetTracks());
+ else
+ esdEv = dynamic_cast<AliESDEvent*>(InputEvent());
+
+ Bool_t useEvent = bin->FillTracks(trList, esdEv, fTrackCuts, mode);
+ return useEvent;
}
//_____________________________________________________________________
void AliForwardFlowTaskQC::Terminate(Option_t */*option*/)
fOutputList->SetName("Results");
fOutputList->SetOwner();
- if ((fFlowFlags & kEtaGap)) {
+ if ((fFlowFlags & kEtaGap) || (fFlowFlags & kTracks)) {
TParameter<Double_t>* etaGap = new TParameter<Double_t>("EtaGap", fEtaGap);
fOutputList->Add(etaGap);
}
// _____________________________________________________________________
AliVVZERO* AliForwardFlowTaskQC::GetVZERO() const
{
+ //
+ // Get VZERO object from ESD or AOD
+ //
+ // Return: VZERO data object
+ //
AliVVZERO* vzero = 0;
// Get input type
UShort_t input = AliForwardUtil::CheckForAOD();
if ((fFlags & kFMD)) nEtaBins = 24;
else if ((fFlags & kVZERO)) nEtaBins = 19;
}
- else if ((fFlags & kVZERO) && !fType.Contains("SPD")) nEtaBins = 11;
+ else if ((fFlags & kVZERO) && (fFlags & kEtaGap)) nEtaBins = 19;
+ else if ((fFlags & kVZERO)) nEtaBins = 11;
+
Double_t vzeroBins[12] = { -6, -3.7, -3.2, -2.7, -2.2, -1.7,
2.8, 3.4, 3.9, 4.5, 5.1, 6 };
Double_t vzeroBins2[20] = { -6, -3.7, -3.2, -2.7, -2.2, // VZERO
-2.0, -1.5, -1.0, -0.5 , 0., 0.5, 1.0, 1.5, 2.0, // SPD
2.8, 3.4, 3.9, 4.5, 5.1, 6 }; // VZERO
-
+
+ Int_t nRefBins = nEtaBins; // needs to be something as default
+ if ((fFlags & kStdQC)) {
+ if ((fFlags & kSymEta) && !((fFlags & kTracks) && (fFlags & kSPD))) nRefBins = 1;
+ else nRefBins = 2;
+ } else if ((fFlags & kEtaGap )) {
+ nRefBins = 2;
+ } else if ((fFlags & k3Cor)) {
+ nRefBins = 24;
+ }
+
// We initiate the reference histogram
fCumuRef = new TH2D(Form("%s_%d_%d%s_ref", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("%s_%d_%d%s_ref", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
- ((fFlags & k3Cor) ? 24 : ((fFlags & kEtaGap) || !(fFlags & kSymEta) ? 2 : 1)), -6., 6.,
+ nRefBins, -6., 6.,
nHBins, 0.5, nHBins+0.5);
if ((fFlags & kVZERO) && (fFlags & k3Cor)) fCumuRef->GetXaxis()->Set(nEtaBins, vzeroBins2);
SetupNUALabels(fCumuRef->GetYaxis());
list->Add(fCumuRef);
-
// We initiate the differential histogram
fCumuDiff = new TH2D(Form("%s_%d_%d%s_diff", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("%s_%d_%d%s_diff", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
nEtaBins, -6., 6., nHBins, 0.5, nHBins+0.5);
if ((fFlags & kVZERO)) {
- if ((fFlags & k3Cor)) fCumuDiff->GetXaxis()->Set(nEtaBins, vzeroBins2);
- else if (!fType.Contains("SPD")) fCumuDiff->GetXaxis()->Set(nEtaBins, vzeroBins);
+ if ((fFlags & kSPD) || (fFlags & k3Cor) || (fFlags & kEtaGap))
+ fCumuDiff->GetXaxis()->Set(nEtaBins, vzeroBins2);
+ else fCumuDiff->GetXaxis()->Set(nEtaBins, vzeroBins);
}
SetupNUALabels(fCumuDiff->GetYaxis());
list->Add(fCumuDiff);
// Initiate the ref cumulant sum histogram
cumuHist = new TH3D(Form("%sv%d_vertex_%d_%d%s_cumuRef", fType.Data(), n, fVzMin, fVzMax, GetQCType(fFlags)),
Form("%sv%d_vertex_%d_%d%s_cumuRef", fType.Data(), n, fVzMin, fVzMax, GetQCType(fFlags)),
- ((fFlags & k3Cor) ? 24 : ((fFlags & kEtaGap) || !(fFlags & kSymEta) ? 2 : 1)), -6., 6.,
+ nRefBins, -6., 6.,
cBins, 0., 100., nC2Bins, 0.5, nC2Bins+0.5);
if ((fFlags & kVZERO) && (fFlags & k3Cor)) cumuHist->GetXaxis()->Set(nEtaBins, vzeroBins2);
cumuHist->GetYaxis()->Set(cBins, centAxis->GetXbins()->GetArray());
Form("%sv%d_vertex_%d_%d%s_cumuDiff", fType.Data(), n, fVzMin, fVzMax, GetQCType(fFlags)),
nEtaBins, -6., 6., cBins, 0., 100., nC4Bins, 0.5, nC4Bins+0.5);
if ((fFlags & kVZERO)) {
- if ((fFlags & k3Cor)) cumuHist->GetXaxis()->Set(nEtaBins, vzeroBins2);
- else if (!fType.Contains("SPD")) cumuHist->GetXaxis()->Set(nEtaBins, vzeroBins);
+ if ((fFlags & kSPD) || (fFlags & k3Cor) || (fFlags & kEtaGap))
+ cumuHist->GetXaxis()->Set(nEtaBins, vzeroBins2);
+ else cumuHist->GetXaxis()->Set(nEtaBins, vzeroBins);
}
cumuHist->GetYaxis()->Set(cBins, centAxis->GetXbins()->GetArray());
fCumuHists.Add(cumuHist);
// Common NUA histograms
fCumuNUARef = new TH3D(Form("%s_vertex_%d_%d%s_cumuNUARef", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("%s_vertex_%d_%d%s_cumuNUARef", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
- ((fFlags & k3Cor) ? 24 : ((fFlags & kEtaGap) || !(fFlags & kSymEta) ? 2 : 1)), -6., 6.,
+ nRefBins, -6., 6.,
cBins, 0., 100., nHBins, 0.5, nHBins+0.5);
if ((fFlags & kVZERO) && (fFlags & k3Cor)) fCumuNUARef->GetXaxis()->Set(nEtaBins, vzeroBins2);
fCumuNUARef->GetYaxis()->Set(cBins, centAxis->GetXbins()->GetArray());
Form("%s_vertex_%d_%d%s_cumuNUADiff", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
nEtaBins, -6., 6., cBins, 0., 100., nHBins, 0.5, nHBins+0.5);
if ((fFlags & kVZERO)) {
- if ((fFlags & k3Cor)) fCumuNUADiff->GetXaxis()->Set(nEtaBins, vzeroBins2);
- else if (!fType.Contains("SPD")) fCumuNUADiff->GetXaxis()->Set(nEtaBins, vzeroBins);
+ if ((fFlags & kSPD) || (fFlags & k3Cor) || (fFlags & kEtaGap))
+ fCumuNUADiff->GetXaxis()->Set(nEtaBins, vzeroBins2);
+ else fCumuNUADiff->GetXaxis()->Set(nEtaBins, vzeroBins);
}
fCumuNUADiff->GetYaxis()->Set(cBins, centAxis->GetXbins()->GetArray());
SetupNUALabels(fCumuNUADiff->GetZaxis());
if (!dList) AliFatal("No diagnostics list found");
// Acceptance hist
- Double_t nPhiBins = ((fFlags & kFMD) ? 20 : 8);
+ Int_t nPhiBins = ((fFlags & kFMD) ? 20 : 8);
fdNdedpRefAcc = new TH2F(Form("h%sdNdedpRefAcc_%d_%d%s", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("%s reference flow acceptance map for %d cm < v_{z} < %d cm", fType.Data(), fVzMin, fVzMax),
- nEtaBins, -6, 6, nPhiBins, 0, TMath::TwoPi());
+ nEtaBins, -6., 6., nPhiBins, 0, TMath::TwoPi());
if ((fFlags & kVZERO)) {
- if ((fFlags & k3Cor)) fdNdedpRefAcc->GetXaxis()->Set(nEtaBins, vzeroBins2);
- else if (!fType.Contains("SPD")) fdNdedpRefAcc->GetXaxis()->Set(nEtaBins, vzeroBins);
+ if ((fFlags & kSPD) || (fFlags & k3Cor) || (fFlags & kEtaGap))
+ fdNdedpRefAcc->GetXaxis()->Set(nEtaBins, vzeroBins2);
+ else fdNdedpRefAcc->GetXaxis()->Set(nEtaBins, vzeroBins);
}
dList->Add(fdNdedpRefAcc);
fdNdedpDiffAcc = new TH2F(Form("h%sdNdedpDiffAcc_%d_%d%s", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("%s differential flow acceptance map for %d cm < v_{z} < %d cm", fType.Data(), fVzMin, fVzMax),
- nEtaBins, -6, 6, nPhiBins, 0, TMath::TwoPi());
+ nEtaBins, -6., 6., nPhiBins, 0, TMath::TwoPi());
if ((fFlags & kVZERO)) {
- if ((fFlags & k3Cor)) fdNdedpDiffAcc->GetXaxis()->Set(nEtaBins, vzeroBins2);
- else if (!fType.Contains("SPD")) fdNdedpDiffAcc->GetXaxis()->Set(nEtaBins, vzeroBins);
+ if ((fFlags & kSPD) || (fFlags & k3Cor) || (fFlags & kEtaGap))
+ fdNdedpDiffAcc->GetXaxis()->Set(nEtaBins, vzeroBins2);
+ else fdNdedpDiffAcc->GetXaxis()->Set(nEtaBins, vzeroBins);
}
dList->Add(fdNdedpDiffAcc);
fOutliers = new TH2F(Form("hOutliers_%s_%d_%d%s", fType.Data(), fVzMin, fVzMax, GetQCType(fFlags)),
Form("Maximum #sigma from mean N_{ch} pr. bin - %s, %d < v_{z} < %d",
fType.Data(), fVzMin, fVzMax),
- 20, 0., 100., 500, 0., (fType.Contains("MC") ? 15. : 5.));
+ 20, 0., 100., 500, 0., ((fFlags & kMC) ? 15. : 5.));
dList->Add(fOutliers);
return;
if ((mode & kFillRef)) fCumuRef->Reset();
if ((mode & kFillDiff)) fCumuDiff->Reset();
}
-
// Then we loop over the input and calculate sum cos(k*n*phi)
// and fill it in the reference and differential histograms
Int_t nBadBins = 0;
TMath::Abs(eta) < fEtaGap) break;
// Backward and forward eta gap break for reference flow
if ((fFlags & kEtaGap) && (mode & kFillRef) && TMath::Abs(eta) > TMath::Abs(limit)) break;
+ if ((fFlags & kStdQC) && (fFlags & kMC) && !(fFlags & kTracks)) {
+ if (!(fFlags & kSPD) && TMath::Abs(eta) < 1.75) break;
+ if ((fFlags & kSPD) && TMath::Abs(eta) > 2.00) break;
+ }
if (limit > 1e3) limit = dNdetadphi.GetXaxis()->GetBinLowEdge(etaBin);
continue;
} // End of phiBin == 0
nInAvg++;
if (weight == 0) continue;
if (weight > max) max = weight;
-
// Fill into Cos() and Sin() hists
- if ((mode & kFillRef)) {
+ if ((mode & kFillRef) && !((fFlags & kTracks) && (fFlags & kMC) && TMath::Abs(eta) > 0.75)) {
fCumuRef->Fill(eta, 0., weight);// mult goes in underflowbin - no visual, but not needed?
fdNdedpRefAcc->Fill(eta, phi, weight);
}
Double_t cosnPhi = weight*TMath::Cos(n*phi);
Double_t sinnPhi = weight*TMath::Sin(n*phi);
// fill ref
- if ((mode & kFillRef)) {
+ if ((mode & kFillRef) && !((fFlags & kTracks) && (fFlags & kMC) && TMath::Abs(eta) > 0.75)) {
fCumuRef->Fill(eta, cosBin, cosnPhi);
fCumuRef->Fill(eta, sinBin, sinnPhi);
}
return useEvent;
}
//_____________________________________________________________________
-Bool_t AliForwardFlowTaskQC::VertexBin::FillTracks(TObjArray* trList, UShort_t mode)
+Bool_t AliForwardFlowTaskQC::VertexBin::FillTracks(TObjArray* trList, AliESDEvent* esd,
+ AliAnalysisFilter* trFilter, UShort_t mode)
{
//
// Fill reference and differential eta-histograms
// mode: filling mode: kFillRef/kFillDiff/kFillBoth
//
if (!fCumuRef) AliFatal("You have not called AddOutput() - Terminating!");
+ if (!trList && !esd) {
+ AliError("FillTracks: No AOD track list or ESD event - something might be wrong!");
+ return kFALSE;
+ }
// Fist we reset histograms
if ((mode & kReset)) {
if ((mode & kFillDiff)) fCumuDiff->Reset();
}
- UShort_t input = AliForwardUtil::CheckForAOD();
// Then we loop over the input and calculate sum cos(k*n*phi)
// and fill it in the reference and differential histograms
- Int_t nTr = trList->GetEntries();
+ Int_t nTr = 0;
+ if (trList) nTr = trList->GetEntries();
+ if (esd) nTr = esd->GetNumberOfTracks();
if (nTr == 0) return kFALSE;
AliVTrack* tr = 0;
- AliAODTrack* aodTr = 0;
- // Cuts for AOD tracks (have already been applied to ESD tracks)
- const Double_t pTMin = 0.5, pTMax = 20., etaMin = -0.8, etaMax = 0.8, minNCl = 50;
+ // Cuts for AOD tracks (have already been applied to ESD tracks) - except dEdx
+// const tpcdEdx = 10;
for (Int_t i = 0; i < nTr; i++) { // track loop
- tr = (AliVTrack*)trList->At(i);
+ tr = (trList ? (AliVTrack*)trList->At(i) : (AliVTrack*)esd->GetTrack(i));
if (!tr) continue;
- if (input == 1) { // If AOD input
- // A dynamic cast would be more safe here, but this is faster...
- aodTr = (AliAODTrack*)tr;
+ if (esd) {
+ AliESDtrack* esdTr = (AliESDtrack*)tr;
+ if (!trFilter->IsSelected(esdTr)) continue;
+ }
+ else if (trList) { // If AOD input
+ Double_t pTMin = 0, pTMax = 0, etaMin = 0, etaMax = 0, minNCl = 0;
+ UInt_t bit = 0;
+ if ((fFlags & kTPC) == kTPC) pTMin = 0.2, pTMax = 5., etaMin = -0.8, etaMax = 0.8, minNCl = 70, bit = 128;
+ if ((fFlags & kHybrid) == kHybrid) pTMin = 0.2, pTMax = 5., etaMin = -0.8, etaMax = 0.8, minNCl = 70, bit = 272;
+
+ AliAODTrack* aodTr = (AliAODTrack*)tr;
if (aodTr->GetID() > -1) continue;
- if (!aodTr->TestFilterBit(128) || !aodTr->Pt() > pTMax || aodTr->Pt() < pTMin ||
+ if (!aodTr->TestFilterBit(bit) || aodTr->Pt() > pTMax || aodTr->Pt() < pTMin ||
aodTr->Eta() > etaMax || aodTr->Eta() < etaMin || aodTr->GetTPCNcls() < minNCl) continue;
}
+
+// if (tr->GetTPCsignal() < tpcdEdx) continue;
// Track accepted
Double_t eta = tr->Eta();
+ if (((fFlags & kSPD) || (fFlags & kEtaGap)) && TMath::Abs(eta) < fEtaGap) continue;
Double_t phi = tr->Phi();
+// Double_t pT = tr->Pt();
+// AliAODMCHeader* head = static_cast<AliAODMCHeader*>(aod->FindListObject(AliAODMCHeader::StdBranchName()));
+// Double_t rp = head->GetReactionPlaneAngle();
+// Double_t b = head->GetImpactParameter();
+ Double_t weight = 1.;//AliForwardFlowWeights::Instance().CalcWeight(eta, pT, phi, 0, rp, b);
+
if ((mode & kFillRef)) {
- fCumuRef->Fill(eta, 0.);// mult goes in underflowbin - no visual, but not needed?
- fdNdedpRefAcc->Fill(eta, phi);
+ fCumuRef->Fill(eta, 0., weight);// mult goes in underflowbin - no visual, but not needed?
+ fdNdedpRefAcc->Fill(eta, phi, weight);
}
if ((mode & kFillDiff)) {
- fCumuDiff->Fill(eta, 0);
- fdNdedpDiffAcc->Fill(eta, phi);
+ fCumuDiff->Fill(eta, 0., weight);
+ fdNdedpDiffAcc->Fill(eta, phi, weight);
}
for (Int_t n = 1; n <= 2*fMaxMoment; n++) {
Double_t cosBin = fCumuDiff->GetYaxis()->GetBinCenter(GetBinNumberCos(n));
Double_t sinBin = fCumuDiff->GetYaxis()->GetBinCenter(GetBinNumberSin(n));
- Double_t cosnPhi = TMath::Cos(n*phi);
- Double_t sinnPhi = TMath::Sin(n*phi);
+ Double_t cosnPhi = weight*TMath::Cos(n*phi);
+ Double_t sinnPhi = weight*TMath::Sin(n*phi);
// fill ref
if ((mode & kFillRef)) {
fCumuRef->Fill(eta, cosBin, cosnPhi);
for (Int_t etaBin = 1; etaBin <= fCumuRef->GetNbinsX(); etaBin++) {
Double_t eta = fCumuRef->GetXaxis()->GetBinCenter(etaBin);
if (fCumuRef->GetBinContent(etaBin, 0) == 0) continue;
+ if ((fFlags & kTracks) && (fFlags && kSPD) && !(fFlags & kEtaGap)) eta = -eta;
for (Int_t qBin = 0; qBin <= fCumuRef->GetNbinsY(); qBin++) {
fCumuNUARef->Fill(eta, cent, Double_t(qBin), fCumuRef->GetBinContent(etaBin, qBin));
}
Double_t sinPhi1Phi2 = 0, sinPhi1Phi2Phi3m = 0;
for (Int_t etaBin = 1; etaBin <= fCumuDiff->GetNbinsX(); etaBin++) {
Double_t eta = fCumuDiff->GetXaxis()->GetBinCenter(etaBin);
- Int_t refEtaBinA = fCumuRef->GetXaxis()->FindBin(eta);
- Int_t refEtaBinB = fCumuRef->GetXaxis()->FindBin(-eta);
+ Double_t refEta = eta;
+ if ((fFlags & kTracks) && (fFlags && kSPD) && !(fFlags & kEtaGap)) refEta = -eta;
+ Int_t refEtaBinA = fCumuRef->GetXaxis()->FindBin(refEta);
+ if ((fFlags & kEtaGap)) refEta = -eta;
+ Int_t refEtaBinB = fCumuRef->GetXaxis()->FindBin(refEta);
if (refEtaBinA != prevRefEtaBin) {
prevRefEtaBin = refEtaBinA;
// Reference flow
// Differential flow calculations for each eta bin is done:
// 2-particle differential flow
- if ((fFlags & kStdQC) && !(fFlags & kTPC)) {
+ if ((fFlags & kStdQC) && (!(fFlags & kTracks) || ((fFlags & kTracks) && (fFlags & kMC) && !(fFlags & kSPD) && TMath::Abs(eta) < 0.75))) {
mq = mp;
qnRe = pnRe;
qnIm = pnIm;
if (fDebug > 0) AliInfo(Form("%s - v_%d: centrality %3.1f:..", fType.Data(), n, cent));
for (Int_t etaBin = 1; etaBin <= cumuRef->GetNbinsX(); etaBin++) { // Eta loop begins
Double_t eta = cumuRef->GetXaxis()->GetBinCenter(etaBin);
- Int_t refEtaBinA = fCumuNUARef->GetXaxis()->FindBin(eta);
- Int_t refEtaBinB = fCumuNUARef->GetXaxis()->FindBin(-eta);
+ Double_t refEta = eta;
+ Int_t refEtaBinA = fCumuNUARef->GetXaxis()->FindBin(refEta);
+ if ((fFlags & kEtaGap)) refEta = -eta;
+ Int_t refEtaBinB = fCumuNUARef->GetXaxis()->FindBin(refEta);
// 2-particle reference flow
Double_t w2Two = cumuRef->GetBinContent(refEtaBinA, cBin, kW2Two);
Double_t w2 = cumuRef->GetBinContent(refEtaBinA, cBin, kW2);
if (fDebug > 0) AliInfo(Form("%s - v_%d: centrality %3.1f:..", fType.Data(), n, cent));
for (Int_t etaBin = 1; etaBin <= cumuDiff->GetNbinsX(); etaBin++) { // Eta loop begins
Double_t eta = cumuDiff->GetXaxis()->GetBinCenter(etaBin);
- Int_t refEtaBinA = fCumuNUARef->GetXaxis()->FindBin(eta);
- Int_t refEtaBinB = fCumuNUARef->GetXaxis()->FindBin(-eta);
+ Double_t refEta = eta;
+ Int_t refEtaBinA = fCumuNUARef->GetXaxis()->FindBin(refEta);
+ if ((fFlags & kEtaGap)) refEta = -eta;
+ Int_t refEtaBinB = fCumuNUARef->GetXaxis()->FindBin(refEta);
// Reference objects
Double_t w2 = cumuRef->GetBinContent(refEtaBinA, cBin, kW2);
sin2nPhiB += fCumuNUARef->GetBinContent(b, cBin, GetBinNumberSin(2*n));
multB += fCumuNUARef->GetBinContent(b, cBin, 0);
}
- if (multA == 0 || multB == 0) AliFatal("Empty NUA values for 3Cor!");
+ if (multA == 0 || multB == 0) {
+ AliWarning(Form("Empty NUA values for 3Cor! (%s)", cumuRef->GetName()));
+ continue;
+ }
cosP1nPhiA /= multA;
sinP1nPhiA /= multA;
cos2nPhiA /= multA;
qc2PrimeA -= cosP1nPsi*cosP1nPhiA + sinP1nPsi*sinP1nPhiA;
qc2PrimeB -= cosP1nPsi*cosP1nPhiB + sinP1nPsi*sinP1nPhiB; // Is this OK?
// Extra NUA term from 2n cosines and sines
- qc2PrimeA /= (1.-(cos2nPsi*cos2nPhiA + sin2nPsi*sin2nPhiA));
- qc2PrimeB /= (1.-(cos2nPsi*cos2nPhiB + sin2nPsi*sin2nPhiB));
+ if (cos2nPsi*cos2nPhiA + sin2nPsi*sin2nPhiA != 1.) qc2PrimeA /= (1.-(cos2nPsi*cos2nPhiA + sin2nPsi*sin2nPhiA));
+ if (cos2nPsi*cos2nPhiB + sin2nPsi*sin2nPhiB != 1.) qc2PrimeB /= (1.-(cos2nPsi*cos2nPhiB + sin2nPsi*sin2nPhiB));
}
if (!TMath::IsNaN(qc2PrimeA) && !TMath::IsNaN(qc2PrimeB) && qc2 != 0) {
if (qc2PrimeA*qc2PrimeB >= 0) {
for (Int_t n = 2; n <= fMaxMoment; n++) printf("v%d ", n);
printf("\n");
TString type = "Standard QC{2} and QC{4} calculations.";
- if ((fFlowFlags & kTPC)) type.ReplaceAll(".", " with TPC tracks for reference.");
if ((fFlowFlags & kEtaGap)) type = "QC{2} with a rapidity gap.";
if ((fFlowFlags & k3Cor)) type = "QC{2} with 3 correlators.";
+ if ((fFlowFlags & kTPC) == kTPC) type.ReplaceAll(".", " with TPC tracks for reference.");
+ if ((fFlowFlags & kHybrid) == kHybrid) type.ReplaceAll(".", " with hybrid tracks for reference.");
Printf("QC calculation type :\t%s", type.Data());
Printf("Symmetrize ref. flow wrt. eta = 0 :\t%s", ((fFlowFlags & kSymEta) ? "true" : "false"));
Printf("Apply NUA correction terms :\t%s", ((fFlowFlags & kNUAcorr) ? "true" : "false"));
Printf("Satellite vertex flag :\t%s", ((fFlowFlags & kSatVtx) ? "true" : "false"));
Printf("FMD sigma cut: :\t%f", fFMDCut);
Printf("SPD sigma cut: :\t%f", fSPDCut);
- if ((fFlowFlags & kEtaGap))
+ if ((fFlowFlags & kEtaGap) || (fFlowFlags & kTracks))
Printf("Eta gap: :\t%f", fEtaGap);
Printf("=======================================================");
}
else if ((flags & k3Cor)) type = "3Cor";
else type = "UNKNOWN";
if (prependUS) type.Prepend("_");
- if ((flags & kTPC)) type.Append("TPCTr");
+ if ((flags & kTPC) == kTPC) type.Append("TPCTr");
+ if ((flags & kHybrid) == kHybrid) type.Append("HybTr");
+
return type.Data();
}
//_____________________________________________________________________
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff