[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / AliForwardMCFlowTaskQC.cxx

//
// Calculate flow on MC data in the forward and central regions using the Q cumulants method.
//
// Inputs:
//  - AliAODEvent
//
// Outputs:
//  - AnalysisResults.root
//
#include "AliForwardMCFlowTaskQC.h"
#include "AliAODMCParticle.h"
#include "AliAODMCHeader.h"
#include "TGraph.h"
#include "TF1.h"
#include "TProfile2D.h"
#include "AliAODEvent.h"
#include "AliAODForwardMult.h"
#include "AliAODCentralMult.h"

ClassImp(AliForwardMCFlowTaskQC)
#if 0
;
#endif
//_____________________________________________________________________
AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC() 
  : AliForwardFlowTaskQC(), 
    fBinsFMDTR(),          // List of FMDTR analysis objects
    fBinsSPDTR(),          // List of SPDTR analysis objects
    fBinsMC(),             // List of MC truth analysis objects
    fdNdedpMC(),           // MC truth d^2N/detadphi histogram
    fWeights(),            // Flow weights
    fImpactParToCent(),    // Impact parameter to centrality graph
    fUseImpactPar(0),      // Use impact par for centrality
    fAddFlow(0),           // Add flow to MC truth
    fAddType(0),           // Add type of flow to MC truth
    fAddOrder(0)           // Add order of flow to MC truth        
{} 
  //
  // Default Constructor
  //
//_____________________________________________________________________
AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC(const char* name) 
  : AliForwardFlowTaskQC(name),
    fBinsFMDTR(),          // List of FMDTR analysis objects
    fBinsSPDTR(),          // List of SPDTR analysis objects
    fBinsMC(),             // List of MC truth analysis objects
    fdNdedpMC(),           // MC truth d^2N/detadphi histogram
    fWeights(),            // Flow weights
    fImpactParToCent(),    // Impact parameter to centrality graph
    fUseImpactPar(0),      // Use impact par for centrality
    fAddFlow(0),           // Add flow to MC truth
    fAddType(0),           // Add type of flow to MC truth
    fAddOrder(0)           // Add order of flow to MC truth        
{ 
  // 
  // Constructor
  // Parameters:
  //  name: Name of task
  //
  fdNdedpMC = TH2D("fdNdedpMC", "fdNdedpMC", 
		   48, -6., 6., 200, 0., 2.*TMath::Pi());
  fdNdedpMC.Sumw2();
  
  //  Double_t impactParam[] = {0.,1.75,4.225,5.965,7.765,9.215,10.46,
  //                            11.565,12.575,13.515,16.679};
  //  Double_t centrality[] = {0.,2.5,7.5,15,25,35,45,55,65,75,90};
  Double_t impactParam[] = { 0.00,  3.72,  5.23,  7.31,  8.88, 10.20, 
			    11.38, 12.47, 13.50, 14.51, 16.679};
  Double_t centrality[]  = { 0.,    5.,   10.,   20.,   30.,   40., 
			    50.,   60.,   70.,   80.,  100.};

  Int_t nPoints = sizeof(impactParam)/sizeof(Double_t);
  fImpactParToCent = new TGraph(nPoints, impactParam, centrality);

}
//_____________________________________________________________________
AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC(const AliForwardMCFlowTaskQC& o) 
  : AliForwardFlowTaskQC(o), 
    fBinsFMDTR(),                          // List of FMDTR analysis objects
    fBinsSPDTR(),                          // List of SPDTR analysis objects
    fBinsMC(),                             // List of MC truth analysis objects
    fdNdedpMC(o.fdNdedpMC),                // MC truth d^2N/detadphi histogram
    fWeights(o.fWeights),                  // Flow weights
    fImpactParToCent(o.fImpactParToCent),  // Impact parameter to centrality
    fUseImpactPar(o.fUseImpactPar),        // Use impact par for centrality
    fAddFlow(o.fAddFlow),                  // Add flow to MC truth
    fAddType(o.fAddType),                  // Add type of flow to MC truth
    fAddOrder(o.fAddOrder)                 // Add order of flow to MC truth
{
  //
  // Copy Constructor
  //
} 
//_____________________________________________________________________
AliForwardMCFlowTaskQC&
AliForwardMCFlowTaskQC::operator=(const AliForwardMCFlowTaskQC& o)
{
  // 
  // Assignment operator
  // Parameters:
  //    o Object to copy from 
  //
  if (&o == this) return *this;
  fdNdedpMC        = o.fdNdedpMC;
  fWeights         = o.fWeights;
  fImpactParToCent = o.fImpactParToCent;
  fUseImpactPar    = o.fUseImpactPar;
  fAddFlow         = o.fAddFlow;
  fAddType         = o.fAddType;
  fAddOrder        = o.fAddOrder;
  return *this;
}
//_____________________________________________________________________
void AliForwardMCFlowTaskQC::InitVertexBins()
{
  //
  // Initiate VertexBin lists
  //
  AliForwardFlowTaskQC::InitVertexBins();

  for(UShort_t n = 2; n <= 6; n++) {
    if (!fv[n]) continue;
      for (Int_t v = 1; v <= fVtxAxis->GetNbins(); v++) {
      Int_t vL = Int_t(fVtxAxis->GetBinLowEdge(v));
      Int_t vH = Int_t(fVtxAxis->GetBinUpEdge(v));
      fBinsFMDTR.Add(new VertexBin(vL, vH, n, "FMDTR", (fgDispVtx ? kFALSE : kTRUE), fFMDCut));
      fBinsSPDTR.Add(new VertexBin(vL, vH, n, "SPDTR", kTRUE, fSPDCut));
      fBinsMC.Add(new VertexBin(vL, vH, n, "MC", kTRUE));
    }
  }
}
//_____________________________________________________________________
void AliForwardMCFlowTaskQC::InitHists()
{
  //
  // Initiate diagnostics hists and add to outputlist
  //
  AliForwardFlowTaskQC::InitHists();

  TIter nextFMDTR(&fBinsFMDTR);
  VertexBin* bin = 0;
  while ((bin = static_cast<VertexBin*>(nextFMDTR()))) {
    bin->AddOutput(fSumList);
  }
  TIter nextSPDTR(&fBinsSPDTR);
  while ((bin = static_cast<VertexBin*>(nextSPDTR()))) {
    bin->AddOutput(fSumList);
  }
  TIter nextMC(&fBinsMC);
  while ((bin = static_cast<VertexBin*>(nextMC()))) {
    bin->AddOutput(fSumList);
  }

  TList* wList = new TList();
  wList->SetName("FlowWeights");
  fWeights.Init(wList);
  fSumList->Add(wList);

}
//_____________________________________________________________________
Bool_t AliForwardMCFlowTaskQC::Analyze() 
{
  // 
  // Load FMD and SPD MC objects from aod tree and call the cumulants 
  // calculation for the correct vertexbin
  //
  if (!AliForwardFlowTaskQC::Analyze()) return kFALSE;

  // Run analysis on trackrefs from FMD and SPD
  const AliAODForwardMult* aodfmult = 
    static_cast<AliAODForwardMult*>(fAOD->FindListObject("ForwardMC"));
  const AliAODCentralMult* aodcmult = 
    static_cast<AliAODCentralMult*>(fAOD->FindListObject("CentralClustersMC"));
  Int_t vtx = fVtxAxis->FindBin(fVtx)-1;
  
  // if objects are present, get histograms
  if (aodfmult) {
    const TH2D& fmdTRdNdetadphi = aodfmult->GetHistogram();
    FillVtxBinList(fBinsFMDTR, fmdTRdNdetadphi, vtx);
  }
  if (aodcmult) {
    const TH2D& spdTRdNdetadphi = aodcmult->GetHistogram();
    FillVtxBinList(fBinsSPDTR, spdTRdNdetadphi, vtx);
  }

  // Run analysis on MC branch
  if (!LoopAODMC()) return kFALSE;
  FillVtxBinList(fBinsMC, fdNdedpMC, vtx);

  return kTRUE;
}
//_____________________________________________________________________
void AliForwardMCFlowTaskQC::Finalize() 
{
  //
  // Finalize command, called by Terminate()
  //
  AliForwardFlowTaskQC::Finalize();

  EndVtxBinList(fBinsFMDTR);
  EndVtxBinList(fBinsSPDTR);
  EndVtxBinList(fBinsMC);

  return;
}
//_____________________________________________________________________
Bool_t AliForwardMCFlowTaskQC::CheckTrigger(const AliAODForwardMult* aodfm) const 
{
  //
  // Function to look for a trigger string in the event.
  //
  // Parameters: 
  //  AliAODForwardMult: forward mult object with trigger and vertex info
  //
  // Returns true if B trigger is present - for some reason this is the one we use in MC
  //
  return aodfm->IsTriggerBits(AliAODForwardMult::kB);
}
// _____________________________________________________________________
Bool_t AliForwardMCFlowTaskQC::GetCentrality(const AliAODForwardMult* aodfm)
{
  // 
  // Function to use centrality parametrization from impact parameter
  // if flag is not set call AliForwardFlowTaskQC::GetCentrality
  //
  // Parameters:
  //  AliAODForwardMult: forward mult object with trigger and vertex info
  //
  // Returns true when centrality is set.
  //
  if (fUseImpactPar) {
    fCent = GetCentFromB();
//    if (fCent >  30) fCent = 42;
    if (fCent != -1) return kTRUE;
  }
  return AliForwardFlowTaskQC::GetCentrality(aodfm);
}
//_____________________________________________________________________
Bool_t AliForwardMCFlowTaskQC::LoopAODMC()  
{
  // 
  // Loop over AliAODParticle branch and fill d^2N/detadphi-histograms.
  // Add flow if set to do so in AddTask function
  fdNdedpMC.Reset();

  //retreive MC particles from event
  TClonesArray* mcArray = 
    static_cast<TClonesArray*>(fAOD->FindListObject(
                               AliAODMCParticle::StdBranchName()));
  if(!mcArray){
    AliWarning("No MC array found in AOD. Try making it again.");
    return kFALSE;
  }

  Double_t rp = 0;
  AliAODMCHeader* header = 
    dynamic_cast<AliAODMCHeader*>(fAOD->FindListObject(
                                  AliAODMCHeader::StdBranchName()));
  if (!header) 
    AliWarning("No header file found.");
  else 
    rp = header->GetReactionPlaneAngle();

  Int_t ntracks = mcArray->GetEntriesFast();

  UShort_t flowFlags = 0;
  if (fAddFlow.Length() > 1) {
    if (fAddFlow.Contains("pt"))   flowFlags |= AliForwardFlowWeights::kPt;
    if (fAddFlow.Contains("pid"))  flowFlags |= AliForwardFlowWeights::kPt;
    if (fAddFlow.Contains("eta"))  flowFlags |= AliForwardFlowWeights::kEta;
    if (fAddFlow.Contains("cent")) flowFlags |= AliForwardFlowWeights::kCent;
  }
  Double_t b = -1.;
  b = header->GetImpactParameter();

  // Track loop: chek how many particles will be accepted
  Double_t weight = 0;
  for (Int_t it = 0; it < ntracks; it++) {
    weight = 1;
    AliAODMCParticle* particle = (AliAODMCParticle*)mcArray->At(it);
    if (!particle) {
      AliError(Form("Could not receive track %d", it));
      continue;
    }
    if (!particle->IsPhysicalPrimary()) continue;
    if (particle->Charge() == 0) continue;
    Double_t pT = particle->Pt();
    Double_t eta = particle->Eta();
    Double_t phi = particle->Phi();
    if (eta > -4. && eta < 5.) {
      // Add flow if it is in the argument
      if (flowFlags != 0) 
/*	weight = fWeights.CalcWeight(eta, pT, phi, particle->PdgCode(), 
				     rp, fCent, fAddType, fAddOrder, 
				     flowFlags) + 1;*/
        weight = fWeights.CalcWeight(eta, pT, phi, particle->PdgCode(),
                                    rp, b); 

//        Printf("%f", weight);
        fdNdedpMC.Fill(eta, phi, weight);
    }
  }
  Int_t sBin = fdNdedpMC.GetXaxis()->FindBin(-4.);
  Int_t eBin = fdNdedpMC.GetXaxis()->FindBin(+5.);
  for ( ; sBin < eBin; sBin++) fdNdedpMC.SetBinContent(sBin, 0, 1.);

  return kTRUE;
}
//_____________________________________________________________________
Double_t AliForwardMCFlowTaskQC::GetCentFromB() const
{
  //
  // Get centrality from MC impact parameter.
  //
  Double_t cent = -1.;
  Double_t b = -1.;
  AliAODMCHeader* header = 
    static_cast<AliAODMCHeader*>(fAOD->FindListObject(AliAODMCHeader::
						      StdBranchName()));
  if (!header) return cent;
  b = header->GetImpactParameter();

  cent = fImpactParToCent->Eval(b);

  return cent;
}
//_____________________________________________________________________
void AliForwardMCFlowTaskQC::PrintFlowSetup() const  
{
  //
  // Print the setup of the flow task
  //
  Printf("AliForwardMCFlowTaskQC::Print");
  Printf("Number of bins in vertex axis:\t%d", fVtxAxis->GetNbins());
  Printf("Range of vertex axis         :\t[%3.1f,%3.1f]", 
			  fVtxAxis->GetXmin(), fVtxAxis->GetXmax());
  printf("Doing flow analysis for      :\t");
  for (Int_t n  = 2; n <= 6; n++) if (fv[n]) printf("v%d ", n);
  printf("\n");
  Printf("Displaced vertex flag:       :\t%s", (fgDispVtx ? "true" : "false"));
  Printf("FMD sigma cut:               :\t%f", fFMDCut);
  Printf("SPD sigma cut:               :\t%f", fSPDCut);
}
//_____________________________________________________________________
//
//
// EOF
Commit	Line	Data
2b556440	1	//
	2	// Calculate flow on MC data in the forward and central regions using the Q cumulants method.
	3	//
	4	// Inputs:
	5	// - AliAODEvent
	6	//
	7	// Outputs:
	8	// - AnalysisResults.root
	9	//
	10	#include "AliForwardMCFlowTaskQC.h"
	11	#include "AliAODMCParticle.h"
	12	#include "AliAODMCHeader.h"
	13	#include "TGraph.h"
	14	#include "TF1.h"
	15	#include "TProfile2D.h"
	16	#include "AliAODEvent.h"
	17	#include "AliAODForwardMult.h"
	18	#include "AliAODCentralMult.h"
	19
	20	ClassImp(AliForwardMCFlowTaskQC)
936b0a6c	21	#if 0
	22	;
	23	#endif
2b556440	24	//_____________________________________________________________________
936b0a6c	25	AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC()
	26	: AliForwardFlowTaskQC(),
	27	fBinsFMDTR(), // List of FMDTR analysis objects
	28	fBinsSPDTR(), // List of SPDTR analysis objects
	29	fBinsMC(), // List of MC truth analysis objects
	30	fdNdedpMC(), // MC truth d^2N/detadphi histogram
	31	fWeights(), // Flow weights
	32	fImpactParToCent(), // Impact parameter to centrality graph
	33	fUseImpactPar(0), // Use impact par for centrality
	34	fAddFlow(0), // Add flow to MC truth
	35	fAddType(0), // Add type of flow to MC truth
	36	fAddOrder(0) // Add order of flow to MC truth
	37	{}
2b556440	38	//
	39	// Default Constructor
	40	//
	41	//_____________________________________________________________________
936b0a6c	42	AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC(const char* name)
	43	: AliForwardFlowTaskQC(name),
	44	fBinsFMDTR(), // List of FMDTR analysis objects
	45	fBinsSPDTR(), // List of SPDTR analysis objects
	46	fBinsMC(), // List of MC truth analysis objects
	47	fdNdedpMC(), // MC truth d^2N/detadphi histogram
	48	fWeights(), // Flow weights
	49	fImpactParToCent(), // Impact parameter to centrality graph
	50	fUseImpactPar(0), // Use impact par for centrality
	51	fAddFlow(0), // Add flow to MC truth
	52	fAddType(0), // Add type of flow to MC truth
	53	fAddOrder(0) // Add order of flow to MC truth
2b556440	54	{
	55	//
	56	// Constructor
	57	// Parameters:
	58	// name: Name of task
	59	//
936b0a6c	60	fdNdedpMC = TH2D("fdNdedpMC", "fdNdedpMC",
936b0a6c	61	48, -6., 6., 200, 0., 2.*TMath::Pi());
2b556440	62	fdNdedpMC.Sumw2();
2b556440	63
936b0a6c	64	// Double_t impactParam[] = {0.,1.75,4.225,5.965,7.765,9.215,10.46,
	65	// 11.565,12.575,13.515,16.679};
	66	// Double_t centrality[] = {0.,2.5,7.5,15,25,35,45,55,65,75,90};
	67	Double_t impactParam[] = { 0.00, 3.72, 5.23, 7.31, 8.88, 10.20,
	68	11.38, 12.47, 13.50, 14.51, 16.679};
	69	Double_t centrality[] = { 0., 5., 10., 20., 30., 40.,
	70	50., 60., 70., 80., 100.};
2b556440	71
	72	Int_t nPoints = sizeof(impactParam)/sizeof(Double_t);
	73	fImpactParToCent = new TGraph(nPoints, impactParam, centrality);
	74
	75	}
	76	//_____________________________________________________________________
936b0a6c	77	AliForwardMCFlowTaskQC::AliForwardMCFlowTaskQC(const AliForwardMCFlowTaskQC& o)
	78	: AliForwardFlowTaskQC(o),
	79	fBinsFMDTR(), // List of FMDTR analysis objects
	80	fBinsSPDTR(), // List of SPDTR analysis objects
	81	fBinsMC(), // List of MC truth analysis objects
	82	fdNdedpMC(o.fdNdedpMC), // MC truth d^2N/detadphi histogram
	83	fWeights(o.fWeights), // Flow weights
	84	fImpactParToCent(o.fImpactParToCent), // Impact parameter to centrality
	85	fUseImpactPar(o.fUseImpactPar), // Use impact par for centrality
	86	fAddFlow(o.fAddFlow), // Add flow to MC truth
	87	fAddType(o.fAddType), // Add type of flow to MC truth
	88	fAddOrder(o.fAddOrder) // Add order of flow to MC truth
	89	{
2b556440	90	//
	91	// Copy Constructor
	92	//
936b0a6c	93	}
2b556440	94	//_____________________________________________________________________
	95	AliForwardMCFlowTaskQC&
	96	AliForwardMCFlowTaskQC::operator=(const AliForwardMCFlowTaskQC& o)
	97	{
	98	//
	99	// Assignment operator
	100	// Parameters:
	101	// o Object to copy from
	102	//
eeb6c967	103	if (&o == this) return *this;
2b556440	104	fdNdedpMC = o.fdNdedpMC;
936b0a6c	105	fWeights = o.fWeights;
2b556440	106	fImpactParToCent = o.fImpactParToCent;
4b5b52b7	107	fUseImpactPar = o.fUseImpactPar;
2b556440	108	fAddFlow = o.fAddFlow;
	109	fAddType = o.fAddType;
	110	fAddOrder = o.fAddOrder;
	111	return *this;
	112	}
	113	//_____________________________________________________________________
	114	void AliForwardMCFlowTaskQC::InitVertexBins()
	115	{
	116	//
	117	// Initiate VertexBin lists
	118	//
	119	AliForwardFlowTaskQC::InitVertexBins();
	120
d420e249	121	for(UShort_t n = 2; n <= 6; n++) {
2b556440	122	if (!fv[n]) continue;
d420e249	123	for (Int_t v = 1; v <= fVtxAxis->GetNbins(); v++) {
a7989168	124	Int_t vL = Int_t(fVtxAxis->GetBinLowEdge(v));
	125	Int_t vH = Int_t(fVtxAxis->GetBinUpEdge(v));
	126	fBinsFMDTR.Add(new VertexBin(vL, vH, n, "FMDTR", (fgDispVtx ? kFALSE : kTRUE), fFMDCut));
	127	fBinsSPDTR.Add(new VertexBin(vL, vH, n, "SPDTR", kTRUE, fSPDCut));
03f52e9d	128	fBinsMC.Add(new VertexBin(vL, vH, n, "MC", kTRUE));
2b556440	129	}
2b556440	130	}
2b556440	131	}
	132	//_____________________________________________________________________
	133	void AliForwardMCFlowTaskQC::InitHists()
	134	{
	135	//
	136	// Initiate diagnostics hists and add to outputlist
	137	//
	138	AliForwardFlowTaskQC::InitHists();
	139
	140	TIter nextFMDTR(&fBinsFMDTR);
	141	VertexBin* bin = 0;
	142	while ((bin = static_cast<VertexBin*>(nextFMDTR()))) {
	143	bin->AddOutput(fSumList);
	144	}
	145	TIter nextSPDTR(&fBinsSPDTR);
	146	while ((bin = static_cast<VertexBin*>(nextSPDTR()))) {
	147	bin->AddOutput(fSumList);
	148	}
	149	TIter nextMC(&fBinsMC);
	150	while ((bin = static_cast<VertexBin*>(nextMC()))) {
	151	bin->AddOutput(fSumList);
	152	}
	153
3112bd03	154	TList* wList = new TList();
	155	wList->SetName("FlowWeights");
	156	fWeights.Init(wList);
	157	fSumList->Add(wList);
	158
2b556440	159	}
	160	//_____________________________________________________________________
	161	Bool_t AliForwardMCFlowTaskQC::Analyze()
	162	{
	163	//
	164	// Load FMD and SPD MC objects from aod tree and call the cumulants
	165	// calculation for the correct vertexbin
	166	//
	167	if (!AliForwardFlowTaskQC::Analyze()) return kFALSE;
	168
	169	// Run analysis on trackrefs from FMD and SPD
936b0a6c	170	const AliAODForwardMult* aodfmult =
	171	static_cast<AliAODForwardMult*>(fAOD->FindListObject("ForwardMC"));
	172	const AliAODCentralMult* aodcmult =
	173	static_cast<AliAODCentralMult*>(fAOD->FindListObject("CentralClustersMC"));
d420e249	174	Int_t vtx = fVtxAxis->FindBin(fVtx)-1;
2b556440	175
4b5b52b7	176	// if objects are present, get histograms
d420e249	177	if (aodfmult) {
	178	const TH2D& fmdTRdNdetadphi = aodfmult->GetHistogram();
	179	FillVtxBinList(fBinsFMDTR, fmdTRdNdetadphi, vtx);
	180	}
	181	if (aodcmult) {
	182	const TH2D& spdTRdNdetadphi = aodcmult->GetHistogram();
	183	FillVtxBinList(fBinsSPDTR, spdTRdNdetadphi, vtx);
	184	}
2b556440	185
	186	// Run analysis on MC branch
	187	if (!LoopAODMC()) return kFALSE;
d420e249	188	FillVtxBinList(fBinsMC, fdNdedpMC, vtx);
2b556440	189
	190	return kTRUE;
	191	}
	192	//_____________________________________________________________________
	193	void AliForwardMCFlowTaskQC::Finalize()
	194	{
	195	//
	196	// Finalize command, called by Terminate()
	197	//
	198	AliForwardFlowTaskQC::Finalize();
	199
4b5b52b7	200	EndVtxBinList(fBinsFMDTR);
	201	EndVtxBinList(fBinsSPDTR);
	202	EndVtxBinList(fBinsMC);
2b556440	203
4b5b52b7	204	return;
	205	}
	206	//_____________________________________________________________________
d420e249	207	Bool_t AliForwardMCFlowTaskQC::CheckTrigger(const AliAODForwardMult* aodfm) const
	208	{
	209	//
	210	// Function to look for a trigger string in the event.
	211	//
	212	// Parameters:
	213	// AliAODForwardMult: forward mult object with trigger and vertex info
	214	//
	215	// Returns true if B trigger is present - for some reason this is the one we use in MC
	216	//
	217	return aodfm->IsTriggerBits(AliAODForwardMult::kB);
	218	}
	219	// _____________________________________________________________________
4b5b52b7	220	Bool_t AliForwardMCFlowTaskQC::GetCentrality(const AliAODForwardMult* aodfm)
	221	{
	222	//
	223	// Function to use centrality parametrization from impact parameter
	224	// if flag is not set call AliForwardFlowTaskQC::GetCentrality
	225	//
	226	// Parameters:
	227	// AliAODForwardMult: forward mult object with trigger and vertex info
	228	//
	229	// Returns true when centrality is set.
	230	//
	231	if (fUseImpactPar) {
	232	fCent = GetCentFromB();
d420e249	233	// if (fCent > 30) fCent = 42;
d420e249	234	if (fCent != -1) return kTRUE;
2b556440	235	}
d420e249	236	return AliForwardFlowTaskQC::GetCentrality(aodfm);
2b556440	237	}
	238	//_____________________________________________________________________
	239	Bool_t AliForwardMCFlowTaskQC::LoopAODMC()
	240	{
	241	//
	242	// Loop over AliAODParticle branch and fill d^2N/detadphi-histograms.
	243	// Add flow if set to do so in AddTask function
	244	fdNdedpMC.Reset();
	245
	246	//retreive MC particles from event
936b0a6c	247	TClonesArray* mcArray =
d420e249	248	static_cast<TClonesArray*>(fAOD->FindListObject(
d420e249	249	AliAODMCParticle::StdBranchName()));
2b556440	250	if(!mcArray){
d420e249	251	AliWarning("No MC array found in AOD. Try making it again.");
2b556440	252	return kFALSE;
	253	}
	254
	255	Double_t rp = 0;
936b0a6c	256	AliAODMCHeader* header =
d420e249	257	dynamic_cast<AliAODMCHeader*>(fAOD->FindListObject(
d420e249	258	AliAODMCHeader::StdBranchName()));
936b0a6c	259	if (!header)
2b556440	260	AliWarning("No header file found.");
936b0a6c	261	else
2b556440	262	rp = header->GetReactionPlaneAngle();
2b556440	263
	264	Int_t ntracks = mcArray->GetEntriesFast();
	265
936b0a6c	266	UShort_t flowFlags = 0;
	267	if (fAddFlow.Length() > 1) {
	268	if (fAddFlow.Contains("pt")) flowFlags \|= AliForwardFlowWeights::kPt;
	269	if (fAddFlow.Contains("pid")) flowFlags \|= AliForwardFlowWeights::kPt;
	270	if (fAddFlow.Contains("eta")) flowFlags \|= AliForwardFlowWeights::kEta;
	271	if (fAddFlow.Contains("cent")) flowFlags \|= AliForwardFlowWeights::kCent;
	272	}
03f52e9d	273	Double_t b = -1.;
03f52e9d	274	b = header->GetImpactParameter();
936b0a6c	275
2b556440	276	// Track loop: chek how many particles will be accepted
	277	Double_t weight = 0;
	278	for (Int_t it = 0; it < ntracks; it++) {
	279	weight = 1;
	280	AliAODMCParticle* particle = (AliAODMCParticle*)mcArray->At(it);
	281	if (!particle) {
	282	AliError(Form("Could not receive track %d", it));
	283	continue;
	284	}
	285	if (!particle->IsPhysicalPrimary()) continue;
	286	if (particle->Charge() == 0) continue;
	287	Double_t pT = particle->Pt();
	288	Double_t eta = particle->Eta();
	289	Double_t phi = particle->Phi();
4b5b52b7	290	if (eta > -4. && eta < 5.) {
2b556440	291	// Add flow if it is in the argument
936b0a6c	292	if (flowFlags != 0)
03f52e9d	293	/* weight = fWeights.CalcWeight(eta, pT, phi, particle->PdgCode(),
936b0a6c	294	rp, fCent, fAddType, fAddOrder,
03f52e9d	295	flowFlags) + 1;*/
	296	weight = fWeights.CalcWeight(eta, pT, phi, particle->PdgCode(),
	297	rp, b);
	298
	299	// Printf("%f", weight);
	300	fdNdedpMC.Fill(eta, phi, weight);
2b556440	301	}
2b556440	302	}
d420e249	303	Int_t sBin = fdNdedpMC.GetXaxis()->FindBin(-4.);
	304	Int_t eBin = fdNdedpMC.GetXaxis()->FindBin(+5.);
	305	for ( ; sBin < eBin; sBin++) fdNdedpMC.SetBinContent(sBin, 0, 1.);
2b556440	306
	307	return kTRUE;
	308	}
	309	//_____________________________________________________________________
2b556440	310	Double_t AliForwardMCFlowTaskQC::GetCentFromB() const
	311	{
	312	//
	313	// Get centrality from MC impact parameter.
2b556440	314	//
	315	Double_t cent = -1.;
	316	Double_t b = -1.;
936b0a6c	317	AliAODMCHeader* header =
	318	static_cast<AliAODMCHeader*>(fAOD->FindListObject(AliAODMCHeader::
	319	StdBranchName()));
2b556440	320	if (!header) return cent;
	321	b = header->GetImpactParameter();
	322
	323	cent = fImpactParToCent->Eval(b);
	324
	325	return cent;
	326	}
	327	//_____________________________________________________________________
d420e249	328	void AliForwardMCFlowTaskQC::PrintFlowSetup() const
	329	{
	330	//
	331	// Print the setup of the flow task
	332	//
	333	Printf("AliForwardMCFlowTaskQC::Print");
	334	Printf("Number of bins in vertex axis:\t%d", fVtxAxis->GetNbins());
	335	Printf("Range of vertex axis :\t[%3.1f,%3.1f]",
	336	fVtxAxis->GetXmin(), fVtxAxis->GetXmax());
	337	printf("Doing flow analysis for :\t");
	338	for (Int_t n = 2; n <= 6; n++) if (fv[n]) printf("v%d ", n);
	339	printf("\n");
	340	Printf("Displaced vertex flag: :\t%s", (fgDispVtx ? "true" : "false"));
	341	Printf("FMD sigma cut: :\t%f", fFMDCut);
	342	Printf("SPD sigma cut: :\t%f", fSPDCut);
	343	}
	344	//_____________________________________________________________________
2b556440	345	//
	346	//
	347	// EOF
	348