+
+}
+
+//-----------------------------------------------------------------------
+void AliFlowEventSimple::Get2Qsub(AliFlowVector* Qarray, Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
+{
+
+ // calculate Q-vector in harmonic n without weights (default harmonic n=2)
+ Double_t dQX = 0.;
+ Double_t dQY = 0.;
+
+ Int_t iOrder = n;
+ Double_t sumOfWeights = 0.;
+ Double_t dPhi = 0.;
+ Double_t dPt = 0.;
+ Double_t dEta = 0.;
+ Double_t dWeight = 1.;
+
+ AliFlowTrackSimple* pTrack = NULL;
+
+ Int_t iNbinsPhiSub0 = 0;
+ Int_t iNbinsPhiSub1 = 0;
+ Double_t dBinWidthPt = 0.;
+ Double_t dPtMin = 0.;
+ Double_t dBinWidthEta= 0.;
+ Double_t dEtaMin = 0.;
+
+ Double_t dWphi = 1.; // weight Phi
+ Double_t dWpt = 1.; // weight Pt
+ Double_t dWeta = 1.; // weight Eta
+
+ TH1F* phiWeightsSub0 = NULL;
+ TH1F* phiWeightsSub1 = NULL;
+ TH1D* ptWeights = NULL;
+ TH1D* etaWeights = NULL;
+
+ if(weightsList)
+ {
+ if(usePhiWeights)
+ {
+ phiWeightsSub0 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub0"));
+ if(phiWeightsSub0) {
+ iNbinsPhiSub0 = phiWeightsSub0->GetNbinsX();
+ }
+ phiWeightsSub1 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub1"));
+ if(phiWeightsSub1) {
+ iNbinsPhiSub1 = phiWeightsSub1->GetNbinsX();
+ }
+ }
+ if(usePtWeights)
+ {
+ ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
+ if(ptWeights)
+ {
+ dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
+ dPtMin = (ptWeights->GetXaxis())->GetXmin();
+ }
+ }
+ if(useEtaWeights)
+ {
+ etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
+ if(etaWeights)
+ {
+ dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
+ dEtaMin = (etaWeights->GetXaxis())->GetXmin();
+ }
+ }
+ } // end of if(weightsList)
+
+ //loop over the two subevents
+ for (Int_t s=0; s<2; s++)
+ {
+ // loop over tracks
+ for(Int_t i=0; i<fNumberOfTracks; i++)
+ {
+ pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
+ if(pTrack)
+ {
+ if(pTrack->InRPSelection())
+ {
+ if (pTrack->InSubevent(s))
+ {
+ dPhi = pTrack->Phi();
+ dPt = pTrack->Pt();
+ dEta = pTrack->Eta();
+ dWeight = pTrack->Weight();
+
+ // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
+ //subevent 0
+ if(s == 0) {
+ if(phiWeightsSub0 && iNbinsPhiSub0) {
+ Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub0/TMath::TwoPi()));
+ //use the phi value at the center of the bin
+ dPhi = phiWeightsSub0->GetBinCenter(phiBin);
+ dWphi = phiWeightsSub0->GetBinContent(phiBin);
+ }
+ }
+ //subevent 1
+ else if (s == 1) {
+ if(phiWeightsSub1 && iNbinsPhiSub1) {
+ Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub1/TMath::TwoPi()));
+ //use the phi value at the center of the bin
+ dPhi = phiWeightsSub1->GetBinCenter(phiBin);
+ dWphi = phiWeightsSub1->GetBinContent(phiBin);
+ }
+ }
+
+ // determine v'(pt) weight:
+ if(ptWeights && dBinWidthPt)
+ {
+ dWpt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
+ }
+
+ // determine v'(eta) weight:
+ if(etaWeights && dBinWidthEta)
+ {
+ dWeta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
+ }
+
+ // building up the weighted Q-vector:
+ dQX += dWeight*dWphi*dWpt*dWeta*TMath::Cos(iOrder*dPhi);
+ dQY += dWeight*dWphi*dWpt*dWeta*TMath::Sin(iOrder*dPhi);
+
+ // weighted multiplicity:
+ sumOfWeights+=dWeight*dWphi*dWpt*dWeta;
+
+ } // end of subevent
+ } // end of if (pTrack->InRPSelection())
+ } // end of if (pTrack)
+ else
+ {
+ cerr << "no particle!!!"<<endl;
+ }
+ } // loop over particles
+ Qarray[s].Set(dQX,dQY);
+ Qarray[s].SetMult(sumOfWeights);
+ //reset
+ sumOfWeights = 0.;
+ dQX = 0.;
+ dQY = 0.;
+ }
+
+}
+
+
+//-----------------------------------------------------------------------
+void AliFlowEventSimple::Print(Option_t *option) const
+{
+ // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
+ // ===============================================
+ // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
+ printf( "Class.Print Name = %s, #tracks= %d, Number of RPs= %d, MC EventPlaneAngle= %f\n",
+ GetName(),fNumberOfTracks, fNumberOfRPs, fMCReactionPlaneAngle );
+
+ TString optionstr(option);
+ if (!optionstr.Contains("all")) return;
+ if (fTrackCollection)
+ {
+ fTrackCollection->Print(option);
+ }
+ else
+ {
+ printf( "Empty track collection \n");
+ }
+}
+
+//-----------------------------------------------------------------------
+void AliFlowEventSimple::Browse(TBrowser *b)
+{
+ //browse in TBrowser
+ if (!b) return;
+ if (!fNumberOfTracksWrap)
+ {
+ fNumberOfTracksWrap = new TParameter<int>("fNumberOfTracks", fNumberOfTracks);
+ b->Add(fNumberOfTracksWrap);
+ }
+ if (!fNumberOfRPsWrap)
+ {
+ fNumberOfRPsWrap = new TParameter<int>("fNumberOfRPs", fNumberOfRPs);
+ b->Add(fNumberOfRPsWrap);
+ }
+ if (!fMCReactionPlaneAngleWrap)
+ {
+ fMCReactionPlaneAngleWrap = new TParameter<double>(" fMCReactionPlaneAngle", fMCReactionPlaneAngle);
+ b->Add( fMCReactionPlaneAngleWrap);
+ }
+ if (fTrackCollection) b->Add(fTrackCollection,"AliFlowTracksSimple");
+}
+
+//-----------------------------------------------------------------------
+AliFlowEventSimple::AliFlowEventSimple( TTree* inputTree,
+ const AliFlowTrackSimpleCuts* rpCuts,
+ const AliFlowTrackSimpleCuts* poiCuts):
+ fTrackCollection(NULL),
+ fReferenceMultiplicity(0),
+ fNumberOfTracks(0),
+ fNumberOfRPs(0),
+ fMCReactionPlaneAngle(0.),
+ fMCReactionPlaneAngleIsSet(kFALSE),
+ fAfterBurnerPrecision(0.001),
+ fUserModified(kFALSE),
+ fNumberOfTracksWrap(NULL),
+ fNumberOfRPsWrap(NULL),
+ fMCReactionPlaneAngleWrap(NULL)
+{
+ //constructor, fills the event from a TTree of kinematic.root files
+ //applies RP and POI cuts, tags the tracks
+
+ Int_t numberOfInputTracks = inputTree->GetEntries() ;
+ fTrackCollection = new TObjArray(numberOfInputTracks/2);
+
+ TParticle* pParticle = new TParticle();
+ inputTree->SetBranchAddress("Particles",&pParticle);
+
+ Int_t iSelParticlesPOI = 0;
+
+ for (Int_t i=0; i<numberOfInputTracks; i++)
+ {
+ inputTree->GetEntry(i); //get input particle
+
+ if (!pParticle) continue; //no particle
+ if (!pParticle->IsPrimary()) continue;
+
+ Bool_t rpOK = rpCuts->PassesCuts(pParticle);
+ Bool_t poiOK = poiCuts->PassesCuts(pParticle);
+
+ if (rpOK || poiOK)
+ {
+ AliFlowTrackSimple* pTrack = new AliFlowTrackSimple(pParticle);
+
+ //marking the particles used for int. flow:
+ if(rpOK)
+ {
+ pTrack->SetForRPSelection(kTRUE);
+ fNumberOfRPs++;
+ }
+ //marking the particles used for diff. flow:
+ if(poiOK)
+ {
+ pTrack->SetForPOISelection(kTRUE);
+ iSelParticlesPOI++;
+ }
+ //adding a particles which were used either for int. or diff. flow to the list
+ AddTrack(pTrack);
+ }
+ }//for i
+ delete pParticle;