ROOT logo
/*************************************************************************
* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
*                                                                        *
* Author: The ALICE Off-line Project.                                    *
* Contributors are mentioned in the code where appropriate.              *
*                                                                        *
* Permission to use, copy, modify and distribute this software and its   *
* documentation strictly for non-commercial purposes is hereby granted   *
* without fee, provided that the above copyright notice appears in all   *
* copies and that both the copyright notice and this permission notice   *
* appear in the supporting documentation. The authors make no claims     *
* about the suitability of this software for any purpose. It is          *
* provided "as is" without express or implied warranty.                  * 
**************************************************************************/

#define ALIFLOWANALYSISWITHSCALARPRODUCT_CXX
 
#include "TFile.h"      
#include "TList.h"
#include "TMath.h"
#include "TString.h"
#include "TProfile.h"
#include "TVector2.h"
#include "TH1D.h"
#include "TH1F.h"
#include "TH2D.h"

#include "AliFlowCommonConstants.h"
#include "AliFlowEventSimple.h"
#include "AliFlowVector.h"
#include "AliFlowTrackSimple.h"
#include "AliFlowCommonHist.h"
#include "AliFlowCommonHistResults.h"
#include "AliFlowAnalysisWithScalarProduct.h"

//////////////////////////////////////////////////////////////////////////////
// AliFlowAnalysisWithScalarProduct:
// Description: Maker to analyze Flow from the Scalar Product method.
// authors: Naomi van del Kolk (kolk@nikhef.nl)
//          Ante Bilandzic (anteb@nikhef.nl)
// mods:    Carlos Perez (cperez@nikhef.nl)
//////////////////////////////////////////////////////////////////////////////

ClassImp(AliFlowAnalysisWithScalarProduct)

//-----------------------------------------------------------------------
AliFlowAnalysisWithScalarProduct::AliFlowAnalysisWithScalarProduct():
fDebug(0),
fMinimalBook(kFALSE),
fUsePhiWeights(0),
fApplyCorrectionForNUA(0),
fHarmonic(2),
fNormalizationType(1),
fTotalQvector(3),
fWeightsList(NULL),
fHistList(NULL),
fHistProConfig(NULL),
fHistProQaQbNorm(NULL),
fHistSumOfWeights(NULL),
fHistProNUAq(NULL),
fHistProQNorm(NULL),
fHistProQaQb(NULL),
fHistProQaQbM(NULL),
fHistMaMb(NULL),
fHistQNormQaQbNorm(NULL),
fHistQaNormMa(NULL),
fHistQbNormMb(NULL),
fResolution(NULL),
fHistQaQb(NULL),
fHistQaQbCos(NULL),
fHistNumberOfSubtractedDaughters(NULL),
fCommonHists(NULL),
fCommonHistsuQ(NULL),
fCommonHistsRes(NULL)
{
  //ctor
  for(int i=0; i!=2; ++i) {
    fPhiWeightsSub[i] = NULL;
    for(int j=0; j!=2; ++j) {
      fHistProUQ[i][j] = NULL;
      fHistProUQQaQb[i][j] = NULL;
      fHistProNUAu[i][j][0] = NULL;
      fHistProNUAu[i][j][1] = NULL;
      for(int k=0; k!=3; ++k)
        fHistSumOfWeightsu[i][j][k] = NULL;
    }
  }
}
//-----------------------------------------------------------------------
AliFlowAnalysisWithScalarProduct::~AliFlowAnalysisWithScalarProduct() 
{
  //destructor
  delete fWeightsList;
  delete fHistList;
}
//-----------------------------------------------------------------------
void AliFlowAnalysisWithScalarProduct::Init() {
  //Define all histograms
  //printf("---Analysis with the Scalar Product Method--- Init\n");
  //printf("--- fNormalizationType %d ---\n", fNormalizationType);
  //save old value and prevent histograms from being added to directory
  //to avoid name clashes in case multiple analaysis objects are used
  //in an analysis
  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);
 
  fHistList = new TList();
  fHistList->SetName("cobjSP");
  fHistList->SetOwner();

  TList *uQRelated = new TList();
  uQRelated->SetName("uQ");
  uQRelated->SetOwner();

  TList *nuaRelated = new TList();
  nuaRelated->SetName("NUA");
  nuaRelated->SetOwner();

  TList *errorRelated = new TList();
  errorRelated->SetName("error");
  errorRelated->SetOwner();

  TList *tQARelated = new TList();
  tQARelated->SetName("QA");
  tQARelated->SetOwner();

  fCommonHists = new AliFlowCommonHist("AliFlowCommonHist_SP","AliFlowCommonHist",fMinimalBook);
  (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic); // store harmonic 
  fHistList->Add(fCommonHists);

  //  fCommonHistsuQ = new AliFlowCommonHist("AliFlowCommonHist_uQ");
  //  (fCommonHistsuQ->GetHarmonic())->Fill(0.5,fHarmonic); // store harmonic 
  //  fHistList->Add(fCommonHistsuQ);

  fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResults_SP","",fHarmonic);
  fHistList->Add(fCommonHistsRes);

  fHistProConfig = new TProfile("FlowPro_Flags_SP","Flow_Flags_SP",4,0.5,4.5,"s");
  fHistProConfig->GetXaxis()->SetBinLabel(1,"fApplyCorrectionForNUA");
  fHistProConfig->GetXaxis()->SetBinLabel(2,"fNormalizationType");
  fHistProConfig->GetXaxis()->SetBinLabel(3,"fUsePhiWeights");
  fHistProConfig->GetXaxis()->SetBinLabel(4,"fHarmonic");
  fHistProConfig->Fill(1,fApplyCorrectionForNUA);
  fHistProConfig->Fill(2,fNormalizationType);
  fHistProConfig->Fill(3,fUsePhiWeights);
  fHistProConfig->Fill(4,fHarmonic);
  fHistList->Add(fHistProConfig);

  fHistProQaQbNorm = new TProfile("FlowPro_QaQbNorm_SP","FlowPro_QaQbNorm_SP", 1, 0.5, 1.5,"s");
  fHistProQaQbNorm->SetYTitle("<QaQb/Na/Nb>");
  errorRelated->Add(fHistProQaQbNorm);

  fHistProNUAq = new TProfile("FlowPro_NUAq_SP","FlowPro_NUAq_SP", 6, 0.5, 6.5,"s");
  fHistProNUAq->GetXaxis()->SetBinLabel( 1,"<<sin(#Phi_{a})>>");
  fHistProNUAq->GetXaxis()->SetBinLabel( 2,"<<cos(#Phi_{a})>>");
  fHistProNUAq->GetXaxis()->SetBinLabel( 3,"<<sin(#Phi_{b})>>");
  fHistProNUAq->GetXaxis()->SetBinLabel( 4,"<<cos(#Phi_{b})>>");
  fHistProNUAq->GetXaxis()->SetBinLabel( 5,"<<sin(#Phi_{t})>>");
  fHistProNUAq->GetXaxis()->SetBinLabel( 6,"<<cos(#Phi_{t})>>");
  nuaRelated->Add(fHistProNUAq);

  fHistSumOfWeights = new TH1D("Flow_SumOfWeights_SP","Flow_SumOfWeights_SP",2,0.5, 2.5);
  fHistSumOfWeights->GetXaxis()->SetBinLabel( 1,"#Sigma Na*Nb");
  fHistSumOfWeights->GetXaxis()->SetBinLabel( 2,"#Sigma (Na*Nb)^2");
  errorRelated->Add(fHistSumOfWeights);

  TString sPOI[2] = {"RP","POI"}; // backward compatibility
  TString sEta[2] = {"Pt","eta"}; // backward compatibility
  TString sTitle[2] = {"p_{T} [GeV]","#eta"};
  TString sWeights[3] = {"uQ","uQuQ","uQQaQb"};
  Int_t iNbins[2];
  Double_t dMin[2], dMax[2];
  iNbins[0] = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
  iNbins[1] = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
  dMin[0]   = AliFlowCommonConstants::GetMaster()->GetPtMin();
  dMin[1]   = AliFlowCommonConstants::GetMaster()->GetEtaMin();
  dMax[0]   = AliFlowCommonConstants::GetMaster()->GetPtMax();
  dMax[1]   = AliFlowCommonConstants::GetMaster()->GetEtaMax();
  for(Int_t iPOI=0; iPOI!=2; ++iPOI) 
    for(Int_t iSpace=0; iSpace!=2; ++iSpace) {
    // uQ
    fHistProUQ[iPOI][iSpace] = new TProfile( Form( "FlowPro_UQ_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                       Form( "FlowPro_UQ%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                       iNbins[iSpace], dMin[iSpace], dMax[iSpace], "s");
    fHistProUQ[iPOI][iSpace]->SetXTitle(sTitle[iSpace].Data());
    fHistProUQ[iPOI][iSpace]->SetYTitle("<uQ>");
    uQRelated->Add(fHistProUQ[iPOI][iSpace]);

    // NUAu
    fHistProNUAu[iPOI][iSpace][0] = new TProfile( Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           iNbins[iSpace], dMin[iSpace], dMax[iSpace]);
    fHistProNUAu[iPOI][iSpace][0]->SetXTitle(sTitle[iSpace].Data());
    nuaRelated->Add(fHistProNUAu[iPOI][iSpace][0]);
    fHistProNUAu[iPOI][iSpace][1] = new TProfile( Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           iNbins[iSpace], dMin[iSpace], dMax[iSpace]);
    fHistProNUAu[iPOI][iSpace][1]->SetXTitle(sTitle[iSpace].Data());
    nuaRelated->Add(fHistProNUAu[iPOI][iSpace][1]);

    // uQ QaQb
    fHistProUQQaQb[iPOI][iSpace] = new TProfile( Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ),
                                           iNbins[iSpace], dMin[iSpace], dMax[iSpace]);
    fHistProUQQaQb[iPOI][iSpace]->SetXTitle(sTitle[iSpace].Data());
    fHistProUQQaQb[iPOI][iSpace]-> SetYTitle("<Qu QaQb>");
    errorRelated->Add(fHistProUQQaQb[iPOI][iSpace]);

    // uWeights
    for(Int_t i=0; i!=3; ++i) {
      fHistSumOfWeightsu[iPOI][iSpace][i] = new TH1D( Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()),
                                                      Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()),
                                                      iNbins[iSpace], dMin[iSpace], dMax[iSpace]);
      fHistSumOfWeightsu[iPOI][iSpace][i]->SetYTitle(sWeights[i].Data());
      fHistSumOfWeightsu[iPOI][iSpace][i]->SetXTitle(sTitle[iSpace].Data());
      errorRelated->Add(fHistSumOfWeightsu[iPOI][iSpace][i]);
    }
  }
  //weights
  if(fUsePhiWeights) {
    if(!fWeightsList) {
      printf( "WARNING: fWeightsList is NULL in the Scalar Product method.\n" );
      exit(0);
    }
    fPhiWeightsSub[0] = dynamic_cast<TH1F*>
                        (fWeightsList->FindObject("phi_weights_sub0"));
    if(!fPhiWeightsSub[0]) {
      printf( "WARNING: phi_weights_sub0 not found in the Scalar Product method.\n" );
      exit(0);  
    }
    nuaRelated->Add( fPhiWeightsSub[0] );
    fPhiWeightsSub[1] = dynamic_cast<TH1F*>
                      (fWeightsList->FindObject("phi_weights_sub1"));
    if(!fPhiWeightsSub[1]) {
      printf( "WARNING: phi_weights_sub1 not found in the Scalar Product method.\n" );
      exit(0);  
    }
    nuaRelated->Add( fPhiWeightsSub[1] );
  }

  if(!fMinimalBook) {
    fHistProQNorm = new TProfile("FlowPro_QNorm_SP","FlowPro_QNorm_SP",       1,0.5,1.5,"s");
    fHistProQNorm->SetYTitle("<|Qa+Qb|>");
    tQARelated->Add(fHistProQNorm);

    fHistProQaQb  = new TProfile("FlowPro_QaQb_SP","FlowPro_QaQb_SP",         1,0.5,1.5,"s");
    fHistProQaQb->SetYTitle("<QaQb>");
    tQARelated->Add(fHistProQaQb);

    fHistProQaQbM = new TProfile("FlowPro_QaQbvsM_SP","FlowPro_QaQbvsM_SP",1000,0.0,10000);
    fHistProQaQbM->SetYTitle("<QaQb>");
    fHistProQaQbM->SetXTitle("M");
    fHistProQaQbM->Sumw2();
    tQARelated->Add(fHistProQaQbM);

    fHistMaMb = new TH2D("Flow_MavsMb_SP","Flow_MavsMb_SP",100,0.,100.,100,0.,100.);
    fHistMaMb->SetYTitle("Ma");
    fHistMaMb->SetXTitle("Mb");
    tQARelated->Add(fHistMaMb);

    fHistQNormQaQbNorm = new TH2D("Flow_QNormvsQaQbNorm_SP","Flow_QNormvsQaQbNorm_SP",88,-1.1,1.1,22,0.,1.1);
    fHistQNormQaQbNorm->SetYTitle("|Q/Mq|");
    fHistQNormQaQbNorm->SetXTitle("QaQb/MaMb");
    tQARelated->Add(fHistQNormQaQbNorm);

    fHistQaNormMa = new TH2D("Flow_QaNormvsMa_SP","Flow_QaNormvsMa_SP",100,0.,100.,22,0.,1.1);
    fHistQaNormMa->SetYTitle("|Qa/Ma|");
    fHistQaNormMa->SetXTitle("Ma");
    tQARelated->Add(fHistQaNormMa);

    fHistQbNormMb = new TH2D("Flow_QbNormvsMb_SP","Flow_QbNormvsMb_SP",100,0.,100.,22,0.,1.1);
    fHistQbNormMb->SetYTitle("|Qb/Mb|");
    fHistQbNormMb->SetXTitle("Mb");
    tQARelated->Add(fHistQbNormMb);

    fResolution = new TH1D("Flow_resolution_SP","Flow_resolution_SP",100,-1.0,1.0);
    fResolution->SetYTitle("dN/d(Cos2(#phi_a - #phi_b))");
    fResolution->SetXTitle("Cos2(#phi_a - #phi_b)");
    tQARelated->Add(fResolution);

    fHistQaQb = new TH1D("Flow_QaQb_SP","Flow_QaQb_SP",200,-100.,100.);
    fHistQaQb->SetYTitle("dN/dQaQb");
    fHistQaQb->SetXTitle("dQaQb");
    tQARelated->Add(fHistQaQb);

    fHistQaQbCos = new TH1D("Flow_QaQbCos_SP","Flow_QaQbCos_SP",63,0.,TMath::Pi());
    fHistQaQbCos->SetYTitle("dN/d#phi");
    fHistQaQbCos->SetXTitle("#phi");
    tQARelated->Add(fHistQaQbCos);
    
    fHistNumberOfSubtractedDaughters = new TH1I("NumberOfSubtractedDaughtersInQ",";daughters;counts;Number of daughters subtracted from Q",5,0.,5.);
    tQARelated->Add(fHistNumberOfSubtractedDaughters);
  }

  fHistList->Add(uQRelated);
  fHistList->Add(nuaRelated);
  fHistList->Add(errorRelated);
  fHistList->Add(tQARelated);

  TH1::AddDirectory(oldHistAddStatus);
}

//-----------------------------------------------------------------------
void AliFlowAnalysisWithScalarProduct::Make(AliFlowEventSimple* anEvent) {
  // Scalar Product method
  if (!anEvent) return; // for coverity

  // Get Q vectors for the subevents
  AliFlowVector* vQarray = new AliFlowVector[2];
  if (fUsePhiWeights)
    anEvent->Get2Qsub(vQarray,fHarmonic,fWeightsList,kTRUE);
  else
    anEvent->Get2Qsub(vQarray,fHarmonic);
  // Subevent a
  AliFlowVector vQa = vQarray[0];
  // Subevent b
  AliFlowVector vQb = vQarray[1];
  delete [] vQarray;

  Double_t dMa = vQa.GetMult();
  if( dMa < 2 ) return;
  Double_t dMb = vQb.GetMult();
  if( dMb < 2 ) return;
  //fill control histograms
  if (fUsePhiWeights) {
    fCommonHists->FillControlHistograms(anEvent,fWeightsList,kTRUE);
  } else {
    fCommonHists->FillControlHistograms(anEvent);
  }

  //Normalizing: weight the Q vectors for the subevents
  Double_t dNa = fNormalizationType ? dMa: vQa.Mod(); // SP corresponds to true
  Double_t dNb = fNormalizationType ? dMb: vQb.Mod(); // SP corresponds to true
  Double_t dWa = fNormalizationType ? dMa: 1; // SP corresponds to true
  Double_t dWb = fNormalizationType ? dMb: 1; // SP corresponds to true

  //Scalar product of the two subevents vectors
  Double_t dQaQb = (vQa*vQb);

  //printf("==============\n");
  //printf("vQa: { %f, %f }\n",vQa.X(),vQa.Y());
  //printf("QaQb/|Qa||Qb|: %f\n",dQaQb/vQa.Mod()/vQb.Mod());
  //printf("QaQb/|Ma||Mb|: %f\n",dQaQb/dMa/dMb);
  
  //      01    10     11   <===   fTotalQVector
  // Q ?= Qa or Qb or QaQb
  AliFlowVector vQm;
  vQm.Set(0.0,0.0);
  Double_t dNq=0;
  if( (fTotalQvector%2)>0 ) { // 01 or 11
    vQm += vQa;
    dNq += dMa;
  }
  if( fTotalQvector>1 ) { // 10 or 11
    vQm += vQb;
    dNq += dMb;
  }
  Double_t dWq = fNormalizationType ? dNq: 1; // SP corresponds to true
  dNq = fNormalizationType ? dNq: vQm.Mod(); // SP corresponds to true

  //fill some EP control histograms
  fHistProQaQbNorm->Fill(1.,dQaQb/dNa/dNb,dWa*dWb);  //Fill (QaQb/NaNb) with weight (WaWb).
  //needed for the error calculation:
  fHistSumOfWeights -> Fill(1.,dWa*dWb);
  fHistSumOfWeights -> Fill(2.,pow(dWa*dWb,2.));
  //needed for correcting non-uniform acceptance: 
  fHistProNUAq->Fill(1.,vQa.Y()/dNa,dWa); // to get <<sin(phi_a)>>
  fHistProNUAq->Fill(2.,vQa.X()/dNa,dWa); // to get <<cos(phi_a)>>
  fHistProNUAq->Fill(3.,vQb.Y()/dNb,dWb); // to get <<sin(phi_b)>>
  fHistProNUAq->Fill(4.,vQb.X()/dNb,dWb); // to get <<cos(phi_b)>>
  fHistProNUAq->Fill(5.,vQm.Y()/dNq,dWq);
  fHistProNUAq->Fill(6.,vQm.X()/dNq,dWq);

  //loop over the tracks of the event
  AliFlowTrackSimple*   pTrack = NULL; 
  Int_t iNumberOfTracks = anEvent->NumberOfTracks(); 
  for (Int_t i=0;i<iNumberOfTracks;i++) {
    pTrack = anEvent->GetTrack(i) ; 
    if (!pTrack) continue;
    Double_t dPhi = pTrack->Phi();
    Double_t dPt  = pTrack->Pt();
    Double_t dEta = pTrack->Eta();

    //calculate vU
    TVector2 vU;
    Double_t dUX = TMath::Cos(fHarmonic*dPhi);
    Double_t dUY = TMath::Sin(fHarmonic*dPhi);
    vU.Set(dUX,dUY);

    //      01    10     11   <===   fTotalQVector
    // Q ?= Qa or Qb or QaQb
    vQm.Set(0.0,0.0); //start the loop fresh
    Double_t dMq=0;
    if( (fTotalQvector%2)>0 ) { // 01 or 11
      vQm += vQa;
      dMq += dMa;
    }
    if( fTotalQvector>1 ) { // 10 or 11
      vQm += vQb;
      dMq += dMb;
    }

    //remove track if in subevent
    for(Int_t inSubEvent=0; inSubEvent<2; ++inSubEvent) {
      if( !pTrack->InSubevent( inSubEvent ) )
        continue;
      if(inSubEvent==0)
        if( (fTotalQvector%2)!=1 )
          continue;
      if(inSubEvent==1)
        if( fTotalQvector<2 )
          continue;
      Double_t dW=1;
      //Double_t dPhiCenter = dPhi;
      //if phi weights are used
      if(fUsePhiWeights && fPhiWeightsSub[inSubEvent]) {
        Int_t iNBinsPhiSub = fPhiWeightsSub[inSubEvent]->GetNbinsX();
        Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNBinsPhiSub/TMath::TwoPi()));
        dW = fPhiWeightsSub[inSubEvent]->GetBinContent(phiBin);
        //dPhiCenter = fPhiWeightsSub[inSubEvent]->GetBinCenter(phiBin);
      }
      //Double_t dQmX = vQm.X() - dW*(pTrack->Weight())* TMath::Cos(fHarmonic*dPhiCenter);
      //Double_t dQmY = vQm.Y() - dW*(pTrack->Weight())* TMath::Sin(fHarmonic*dPhiCenter);
      //vQm.Set(dQmX,dQmY);
      
      //subtrack the track from the Q vector, but only if it was used to construct this
      //Q vector: i.e. check wether it has the same tags and is in the same subevent
      //this is especially important for the daughters (as for the mother it is already checked)
      Int_t numberOfsubtractedDaughters=vQm.SubtractTrackWithDaughters(pTrack,dW);
      
      if(!fMinimalBook) {
        fHistNumberOfSubtractedDaughters->Fill(numberOfsubtractedDaughters);
      }

      dMq = dMq-dW*pTrack->Weight();
    }
    dNq = fNormalizationType ? dMq : vQm.Mod();
    dWq = fNormalizationType ? dMq : 1;

    //Filling QA (for compatibility with previous version)
    if(!fMinimalBook) {
      fHistProQaQb->Fill(1,vQa*vQb,1);
      fHistProQaQbM->Fill(anEvent->GetNumberOfRPs()+0.5,(vQa*vQb)/dMa/dMb,dMa*dMb);
      fHistQaQbCos->Fill(TMath::ACos((vQa*vQb)/vQa.Mod()/vQb.Mod()));
      fResolution->Fill( TMath::Cos( vQa.Phi()-vQb.Phi() ) );
      fHistQaQb->Fill(vQa*vQb);
      fHistMaMb->Fill(dMb,dMa);
      fHistProQNorm->Fill(1,vQm.Mod()/dMq,dMq);
      fHistQNormQaQbNorm->Fill((vQa*vQb)/dMa/dMb,vQm.Mod()/dMq);
      fHistQaNormMa->Fill(dMa,vQa.Mod()/dMa);
      fHistQbNormMb->Fill(dMb,vQb.Mod()/dMb);
    }

    Double_t dUQ = vU*vQm;

    //fill the profile histograms
    for(Int_t iPOI=0; iPOI!=2; ++iPOI) {
      if( (iPOI==0)&&(!pTrack->InRPSelection()) )
        continue;
      if( (iPOI==1)&&(!pTrack->InPOISelection(fPOItype)) )
        continue;
      fHistProUQ[iPOI][0]->Fill(dPt ,dUQ/dNq,dWq); //Fill (uQ/Nq') with weight (Nq')
      fHistProUQ[iPOI][1]->Fill(dEta,dUQ/dNq,dWq); //Fill (uQ/Nq') with weight (Nq')
      //needed for the error calculation:
      fHistProUQQaQb[iPOI][0]-> Fill(dPt ,dUQ/dNq*dQaQb/dNa/dNb,dWq*dWa*dWb); //Fill [Qu/Nq']*[QaQb/NaNb] with weight (Nq')NaNb
      fHistProUQQaQb[iPOI][1]-> Fill(dEta,dUQ/dNq*dQaQb/dNa/dNb,dWq*dWa*dWb); //Fill [Qu/Nq']*[QaQb/NaNb] with weight (Nq')NaNb
      fHistSumOfWeightsu[iPOI][0][0]->Fill(dPt ,dWq);        // sum of Nq'     
      fHistSumOfWeightsu[iPOI][0][1]->Fill(dPt ,pow(dWq,2.));// sum of Nq'^2     
      fHistSumOfWeightsu[iPOI][0][2]->Fill(dPt ,dWq*dWa*dWb);// sum of Nq'*Na*Nb     
      fHistSumOfWeightsu[iPOI][1][0]->Fill(dEta,dWq);        // sum of Nq'     
      fHistSumOfWeightsu[iPOI][1][1]->Fill(dEta,pow(dWq,2.));// sum of Nq'^2     
      fHistSumOfWeightsu[iPOI][1][2]->Fill(dEta,dNq*dWa*dWb);// sum of Nq'*Na*Nb
      //NUA:
      fHistProNUAu[iPOI][0][0]->Fill(dPt,dUY,1.); //sin u
      fHistProNUAu[iPOI][0][1]->Fill(dPt,dUX,1.); //cos u
      fHistProNUAu[iPOI][1][0]->Fill(dEta,dUY,1.); //sin u
      fHistProNUAu[iPOI][1][1]->Fill(dEta,dUX,1.); //cos u
    }
  }//loop over tracks

}

//--------------------------------------------------------------------  
void AliFlowAnalysisWithScalarProduct::GetOutputHistograms(TList *outputListHistos){
  //get pointers to all output histograms (called before Finish())
  fHistList = outputListHistos;

  fCommonHists = (AliFlowCommonHist*) fHistList->FindObject("AliFlowCommonHist_SP");
  //  fCommonHistsuQ = (AliFlowCommonHist*) fHistList->FindObject("AliFlowCommonHist_uQ");
  fCommonHistsRes = (AliFlowCommonHistResults*) fHistList->FindObject("AliFlowCommonHistResults_SP");
  fHistProConfig = (TProfile*) fHistList->FindObject("FlowPro_Flags_SP");
  if(!fHistProConfig) printf("Error loading fHistProConfig\n");
  TList *uQ = (TList*) fHistList->FindObject("uQ");
  TList *nua = (TList*) fHistList->FindObject("NUA");
  TList *error = (TList*) fHistList->FindObject("error");

  fHistProQaQbNorm = (TProfile*) error->FindObject("FlowPro_QaQbNorm_SP");
  if(!fHistProQaQbNorm) printf("Error loading fHistProQaQbNorm\n");
  fHistProNUAq = (TProfile*) nua->FindObject("FlowPro_NUAq_SP");
  if(!fHistProNUAq) printf("Error loading fHistProNUAq\n");
  fHistSumOfWeights = (TH1D*) error->FindObject("Flow_SumOfWeights_SP");
  if(!fHistSumOfWeights) printf("Error loading fHistSumOfWeights\n");

  TString sPOI[2] = {"RP","POI"};
  TString sEta[2] = {"Pt","eta"};
  TString sWeights[3] = {"uQ","uQuQ","uQQaQb"};
  for(Int_t iPOI=0; iPOI!=2; ++iPOI) for(Int_t iSpace=0; iSpace!=2; ++iSpace) {
    fHistProUQ[iPOI][iSpace] = (TProfile*) uQ->FindObject( Form( "FlowPro_UQ_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) );
    if(!fHistProUQ[iPOI][iSpace]) printf("Error loading fHistProUQ[%d][%d]\n",iPOI,iSpace);
    fHistProNUAu[iPOI][iSpace][0] = (TProfile*) nua->FindObject( Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) );
    if(!fHistProNUAu[iPOI][iSpace][0]) printf("Error loading fHistProNUAu[%d][%d][0]\n",iPOI,iSpace);
    fHistProNUAu[iPOI][iSpace][1] = (TProfile*) nua->FindObject( Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) );
    if(!fHistProNUAu[iPOI][iSpace][1]) printf("Error loading fHistProNUAu[%d][%d][1]\n",iPOI,iSpace);
    fHistProUQQaQb[iPOI][iSpace] = (TProfile*) error->FindObject( Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) );
    for(Int_t i=0; i!=3; ++i){
      fHistSumOfWeightsu[iPOI][iSpace][i] = (TH1D*) error->FindObject( Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()) );
      if(!fHistSumOfWeightsu[iPOI][iSpace][i]) printf("Error loading fHistSumOfWeightsu[%d][%d][%d]\n",iPOI,iSpace,i);
    }
  }
  if(fHistProConfig) {
    fApplyCorrectionForNUA = (Int_t) fHistProConfig->GetBinContent(1);
    fNormalizationType  = (Int_t) fHistProConfig->GetBinContent(2);
    fUsePhiWeights = (Int_t) fHistProConfig->GetBinContent(3);
    fHarmonic = (Int_t) fHistProConfig->GetBinContent(4);
  }
}            

//--------------------------------------------------------------------            
void AliFlowAnalysisWithScalarProduct::Finish() {
  //calculate flow and fill the AliFlowCommonHistResults
  printf("AliFlowAnalysisWithScalarProduct::Finish()\n");
  
  // access harmonic:
  fApplyCorrectionForNUA = (Int_t)(fHistProConfig->GetBinContent(1));
  fNormalizationType = (Int_t)(fHistProConfig->GetBinContent(2));
  fHarmonic = (Int_t)(fHistProConfig->GetBinContent(4));
  
  printf("*************************************\n");
  printf("*************************************\n");
  printf("      Integrated flow from           \n");
  printf("         Scalar Product              \n\n");
  if(!fNormalizationType)
    printf("          (BehaveAsEP)               \n\n");
  
  //Calculate reference flow
  //----------------------------------
  //weighted average over (QaQb/NaNb) with weight (NaNb)
  Double_t dEntriesQaQb = fHistProQaQbNorm->GetEntries();
  if( dEntriesQaQb < 1 )
    return;
  Double_t dQaQb  = fHistProQaQbNorm->GetBinContent(1);
  if( dQaQb < 0 )
    return;
  Double_t dSpreadQaQb = fHistProQaQbNorm->GetBinError(1);

  //NUA qq:
  Double_t dImQa = fHistProNUAq->GetBinContent(1);
  Double_t dReQa = fHistProNUAq->GetBinContent(2);
  Double_t dImQb = fHistProNUAq->GetBinContent(3);
  Double_t dReQb = fHistProNUAq->GetBinContent(4);
  if(fApplyCorrectionForNUA) 
    dQaQb = dQaQb - dImQa*dImQb - dReQa*dReQb; 
  printf("QaQb = %f +- %f\n", dQaQb, (dSpreadQaQb/TMath::Sqrt(dEntriesQaQb)) );
  Double_t dV = TMath::Sqrt(dQaQb);

  printf("ResSub = %f\n", dV );
  printf("fTotalQvector %d \n",fTotalQvector);
  if(!fNormalizationType) {
    if(fTotalQvector>2) {
      dV = ComputeResolution( TMath::Sqrt2()*FindXi(dV,1e-6) );
      printf("An estimate of the event plane resolution is: %f\n", dV );
    }
  }
  printf("ResTot = %f\n", dV );
  //statistical error of dQaQb: 
  //  statistical error = term1 * spread * term2:
  //  term1 = sqrt{sum_{i=1}^{N} w^2}/(sum_{i=1}^{N} w)
  //  term2 = 1/sqrt(1-term1^2) 
  Double_t dSumOfLinearWeights = fHistSumOfWeights->GetBinContent(1);
  Double_t dSumOfQuadraticWeights = fHistSumOfWeights->GetBinContent(2);
  Double_t dTerm1 = 0.;
  Double_t dTerm2 = 0.;
  if(dSumOfLinearWeights)
    dTerm1 = pow(dSumOfQuadraticWeights,0.5)/dSumOfLinearWeights;
  if(1.-pow(dTerm1,2.) > 0.)
    dTerm2 = 1./pow(1-pow(dTerm1,2.),0.5);
  Double_t dStatErrorQaQb = dTerm1 * dSpreadQaQb * dTerm2;
  Double_t dVerr = 0.;
  if(dQaQb > 0.)
    dVerr = (1./(2.*pow(dQaQb,0.5)))*dStatErrorQaQb;
  fCommonHistsRes->FillIntegratedFlow(dV,dVerr);
  printf("v%d(subevents) = %f +- %f\n",fHarmonic,dV,dVerr);

  Int_t iNbins[2];
  iNbins[0] = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
  iNbins[1] = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
   
  //Calculate differential flow and integrated flow (RP, POI)
  //---------------------------------------------------------
  //v as a function of eta for RP selection
  for(Int_t iRFPorPOI=0; iRFPorPOI != 2; ++iRFPorPOI)
  for(Int_t iPTorETA=0; iPTorETA != 2; ++iPTorETA)
  for(Int_t b=1; b != iNbins[iPTorETA]+1; ++b) {
    Double_t duQpro = fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b);
    if(fApplyCorrectionForNUA)
      duQpro = duQpro 
	- fHistProNUAu[iRFPorPOI][iPTorETA][1]->GetBinContent(b)*fHistProNUAq->GetBinContent(6)
	- fHistProNUAu[iRFPorPOI][iPTorETA][0]->GetBinContent(b)*fHistProNUAq->GetBinContent(5);
    Double_t dv2pro = -999.;
    if( TMath::Abs(dV!=0.) ) { dv2pro = duQpro/dV; }
    //calculate the statistical error
    Double_t dv2ProErr = CalculateStatisticalError(iRFPorPOI, iPTorETA, b, dStatErrorQaQb);
    if( (iRFPorPOI==0)&&(iPTorETA==0) )
      fCommonHistsRes->FillDifferentialFlowPtRP(  b, dv2pro, dv2ProErr);   
    if( (iRFPorPOI==0)&&(iPTorETA==1) )
      fCommonHistsRes->FillDifferentialFlowEtaRP( b, dv2pro, dv2ProErr);   
    if( (iRFPorPOI==1)&&(iPTorETA==0) )
      fCommonHistsRes->FillDifferentialFlowPtPOI( b, dv2pro, dv2ProErr);   
    if( (iRFPorPOI==1)&&(iPTorETA==1) )
      fCommonHistsRes->FillDifferentialFlowEtaPOI(b, dv2pro, dv2ProErr);
    //printf("POI %d | PT %d >>> %f +- %f\n",iRFPorPOI,iPTorETA,dv2pro,dv2ProErr);
  }
  
  printf("\n");
  printf("*************************************\n");
  printf("*************************************\n");
}

//-----------------------------------------------------------------------
void AliFlowAnalysisWithScalarProduct::WriteHistograms(TDirectoryFile *outputFileName) const {
 //store the final results in output .root file
 outputFileName->Add(fHistList);
 outputFileName->Write(outputFileName->GetName(), TObject::kSingleKey);
}

//--------------------------------------------------------------------            
Double_t AliFlowAnalysisWithScalarProduct::CalculateStatisticalError(Int_t iRFPorPOI, Int_t iPTorETA, Int_t b, Double_t aStatErrorQaQb) const {
  //calculate the statistical error for differential flow for bin b
  Double_t duQproSpread = fHistProUQ[iRFPorPOI][iPTorETA]->GetBinError(b);
  Double_t sumOfMq = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][0]->GetBinContent(b);
  Double_t sumOfMqSquared = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][1]->GetBinContent(b);
  Double_t dQaQb = fHistProQaQbNorm->GetBinContent(1);

  //non-isotropic terms:  
  if(fApplyCorrectionForNUA) {
    Double_t dImQa = fHistProNUAq->GetBinContent(1);  // <<sin(phi_a)>>
    Double_t dReQa = fHistProNUAq->GetBinContent(2);  // <<cos(phi_a)>>
    Double_t dImQb = fHistProNUAq->GetBinContent(3);  // <<sin(phi_b)>>
    Double_t dReQb = fHistProNUAq->GetBinContent(4);  // <<cos(phi_b)>>
    dQaQb = dQaQb - dImQa*dImQb - dReQa*dReQb; 
  }

  Double_t dTerm1 = 0.;
  Double_t dTerm2 = 0.;
  if(sumOfMq) {
    dTerm1 = (pow(sumOfMqSquared,0.5)/sumOfMq);
  } 
  if(1.-pow(dTerm1,2.)>0.) {
    dTerm2 = 1./pow(1.-pow(dTerm1,2.),0.5); 
  }
  Double_t duQproErr = dTerm1*duQproSpread*dTerm2;
  // covariances:
  Double_t dTerm1Cov = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][2]->GetBinContent(b);
  Double_t dTerm2Cov = fHistSumOfWeights->GetBinContent(1);
  Double_t dTerm3Cov = sumOfMq;
  Double_t dWeightedCovariance = 0.;
  if(dTerm2Cov*dTerm3Cov>0.) {
    Double_t dDenominator = 1.-dTerm1Cov/(dTerm2Cov*dTerm3Cov);
    Double_t dPrefactor = dTerm1Cov/(dTerm2Cov*dTerm3Cov);
    if(dDenominator!=0) {
      Double_t dCovariance = ( fHistProUQQaQb[iRFPorPOI][iPTorETA]->GetBinContent(b)-dQaQb*fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b))/dDenominator;
      dWeightedCovariance = dCovariance*dPrefactor; 
    }
  }
  Double_t dv2ProErr = 0.; // final statitical error: 
  if(dQaQb>0.) {
    Double_t dv2ProErrorSquared = (1./4.)*pow(dQaQb,-3.)*
      (pow(fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b),2.)*pow(aStatErrorQaQb,2.)
       + 4.*pow(dQaQb,2.)*pow(duQproErr,2.)
       - 4.*dQaQb*fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b)*dWeightedCovariance);
    if(dv2ProErrorSquared>0.) dv2ProErr = pow(dv2ProErrorSquared,0.5);
  }
  return dv2ProErr;
}

Double_t AliFlowAnalysisWithScalarProduct::ComputeResolution( Double_t x ) const {
  // Computes resolution for Event Plane method
  if(x > 51.3) {
    printf("Warning: Estimation of total resolution might be WRONG. Please check!");
    return 0.99981;
  }
  Double_t a = x*x/4;
  Double_t b = TMath::Exp(-a)*TMath::BesselI0(a)+TMath::Exp(-a)*TMath::BesselI1(a);
  return TMath::Sqrt(TMath::PiOver2())/2*x*b;
}

Double_t AliFlowAnalysisWithScalarProduct::FindXi( Double_t res, Double_t prec ) const {
  // Computes x(res) for Event Plane method
  if(res > 0.99981) {
    printf("Warning: Resolution for subEvent is high. You reached the precision limit.");
    return 51.3;
  }
  int nSteps =0;
  Double_t xtmp=0, xmin=0, xmax=51.3, rtmp=0, delta=2*prec;
  while( delta > prec ) {
    xtmp = 0.5*(xmin+xmax);
    rtmp = ComputeResolution(xtmp);
    delta = TMath::Abs( res-rtmp );
    if(rtmp>res) xmax = xtmp;
    if(rtmp<res) xmin = xtmp;
    nSteps++;
  }
  return xtmp;
}

 AliFlowAnalysisWithScalarProduct.cxx:1
 AliFlowAnalysisWithScalarProduct.cxx:2
 AliFlowAnalysisWithScalarProduct.cxx:3
 AliFlowAnalysisWithScalarProduct.cxx:4
 AliFlowAnalysisWithScalarProduct.cxx:5
 AliFlowAnalysisWithScalarProduct.cxx:6
 AliFlowAnalysisWithScalarProduct.cxx:7
 AliFlowAnalysisWithScalarProduct.cxx:8
 AliFlowAnalysisWithScalarProduct.cxx:9
 AliFlowAnalysisWithScalarProduct.cxx:10
 AliFlowAnalysisWithScalarProduct.cxx:11
 AliFlowAnalysisWithScalarProduct.cxx:12
 AliFlowAnalysisWithScalarProduct.cxx:13
 AliFlowAnalysisWithScalarProduct.cxx:14
 AliFlowAnalysisWithScalarProduct.cxx:15
 AliFlowAnalysisWithScalarProduct.cxx:16
 AliFlowAnalysisWithScalarProduct.cxx:17
 AliFlowAnalysisWithScalarProduct.cxx:18
 AliFlowAnalysisWithScalarProduct.cxx:19
 AliFlowAnalysisWithScalarProduct.cxx:20
 AliFlowAnalysisWithScalarProduct.cxx:21
 AliFlowAnalysisWithScalarProduct.cxx:22
 AliFlowAnalysisWithScalarProduct.cxx:23
 AliFlowAnalysisWithScalarProduct.cxx:24
 AliFlowAnalysisWithScalarProduct.cxx:25
 AliFlowAnalysisWithScalarProduct.cxx:26
 AliFlowAnalysisWithScalarProduct.cxx:27
 AliFlowAnalysisWithScalarProduct.cxx:28
 AliFlowAnalysisWithScalarProduct.cxx:29
 AliFlowAnalysisWithScalarProduct.cxx:30
 AliFlowAnalysisWithScalarProduct.cxx:31
 AliFlowAnalysisWithScalarProduct.cxx:32
 AliFlowAnalysisWithScalarProduct.cxx:33
 AliFlowAnalysisWithScalarProduct.cxx:34
 AliFlowAnalysisWithScalarProduct.cxx:35
 AliFlowAnalysisWithScalarProduct.cxx:36
 AliFlowAnalysisWithScalarProduct.cxx:37
 AliFlowAnalysisWithScalarProduct.cxx:38
 AliFlowAnalysisWithScalarProduct.cxx:39
 AliFlowAnalysisWithScalarProduct.cxx:40
 AliFlowAnalysisWithScalarProduct.cxx:41
 AliFlowAnalysisWithScalarProduct.cxx:42
 AliFlowAnalysisWithScalarProduct.cxx:43
 AliFlowAnalysisWithScalarProduct.cxx:44
 AliFlowAnalysisWithScalarProduct.cxx:45
 AliFlowAnalysisWithScalarProduct.cxx:46
 AliFlowAnalysisWithScalarProduct.cxx:47
 AliFlowAnalysisWithScalarProduct.cxx:48
 AliFlowAnalysisWithScalarProduct.cxx:49
 AliFlowAnalysisWithScalarProduct.cxx:50
 AliFlowAnalysisWithScalarProduct.cxx:51
 AliFlowAnalysisWithScalarProduct.cxx:52
 AliFlowAnalysisWithScalarProduct.cxx:53
 AliFlowAnalysisWithScalarProduct.cxx:54
 AliFlowAnalysisWithScalarProduct.cxx:55
 AliFlowAnalysisWithScalarProduct.cxx:56
 AliFlowAnalysisWithScalarProduct.cxx:57
 AliFlowAnalysisWithScalarProduct.cxx:58
 AliFlowAnalysisWithScalarProduct.cxx:59
 AliFlowAnalysisWithScalarProduct.cxx:60
 AliFlowAnalysisWithScalarProduct.cxx:61
 AliFlowAnalysisWithScalarProduct.cxx:62
 AliFlowAnalysisWithScalarProduct.cxx:63
 AliFlowAnalysisWithScalarProduct.cxx:64
 AliFlowAnalysisWithScalarProduct.cxx:65
 AliFlowAnalysisWithScalarProduct.cxx:66
 AliFlowAnalysisWithScalarProduct.cxx:67
 AliFlowAnalysisWithScalarProduct.cxx:68
 AliFlowAnalysisWithScalarProduct.cxx:69
 AliFlowAnalysisWithScalarProduct.cxx:70
 AliFlowAnalysisWithScalarProduct.cxx:71
 AliFlowAnalysisWithScalarProduct.cxx:72
 AliFlowAnalysisWithScalarProduct.cxx:73
 AliFlowAnalysisWithScalarProduct.cxx:74
 AliFlowAnalysisWithScalarProduct.cxx:75
 AliFlowAnalysisWithScalarProduct.cxx:76
 AliFlowAnalysisWithScalarProduct.cxx:77
 AliFlowAnalysisWithScalarProduct.cxx:78
 AliFlowAnalysisWithScalarProduct.cxx:79
 AliFlowAnalysisWithScalarProduct.cxx:80
 AliFlowAnalysisWithScalarProduct.cxx:81
 AliFlowAnalysisWithScalarProduct.cxx:82
 AliFlowAnalysisWithScalarProduct.cxx:83
 AliFlowAnalysisWithScalarProduct.cxx:84
 AliFlowAnalysisWithScalarProduct.cxx:85
 AliFlowAnalysisWithScalarProduct.cxx:86
 AliFlowAnalysisWithScalarProduct.cxx:87
 AliFlowAnalysisWithScalarProduct.cxx:88
 AliFlowAnalysisWithScalarProduct.cxx:89
 AliFlowAnalysisWithScalarProduct.cxx:90
 AliFlowAnalysisWithScalarProduct.cxx:91
 AliFlowAnalysisWithScalarProduct.cxx:92
 AliFlowAnalysisWithScalarProduct.cxx:93
 AliFlowAnalysisWithScalarProduct.cxx:94
 AliFlowAnalysisWithScalarProduct.cxx:95
 AliFlowAnalysisWithScalarProduct.cxx:96
 AliFlowAnalysisWithScalarProduct.cxx:97
 AliFlowAnalysisWithScalarProduct.cxx:98
 AliFlowAnalysisWithScalarProduct.cxx:99
 AliFlowAnalysisWithScalarProduct.cxx:100
 AliFlowAnalysisWithScalarProduct.cxx:101
 AliFlowAnalysisWithScalarProduct.cxx:102
 AliFlowAnalysisWithScalarProduct.cxx:103
 AliFlowAnalysisWithScalarProduct.cxx:104
 AliFlowAnalysisWithScalarProduct.cxx:105
 AliFlowAnalysisWithScalarProduct.cxx:106
 AliFlowAnalysisWithScalarProduct.cxx:107
 AliFlowAnalysisWithScalarProduct.cxx:108
 AliFlowAnalysisWithScalarProduct.cxx:109
 AliFlowAnalysisWithScalarProduct.cxx:110
 AliFlowAnalysisWithScalarProduct.cxx:111
 AliFlowAnalysisWithScalarProduct.cxx:112
 AliFlowAnalysisWithScalarProduct.cxx:113
 AliFlowAnalysisWithScalarProduct.cxx:114
 AliFlowAnalysisWithScalarProduct.cxx:115
 AliFlowAnalysisWithScalarProduct.cxx:116
 AliFlowAnalysisWithScalarProduct.cxx:117
 AliFlowAnalysisWithScalarProduct.cxx:118
 AliFlowAnalysisWithScalarProduct.cxx:119
 AliFlowAnalysisWithScalarProduct.cxx:120
 AliFlowAnalysisWithScalarProduct.cxx:121
 AliFlowAnalysisWithScalarProduct.cxx:122
 AliFlowAnalysisWithScalarProduct.cxx:123
 AliFlowAnalysisWithScalarProduct.cxx:124
 AliFlowAnalysisWithScalarProduct.cxx:125
 AliFlowAnalysisWithScalarProduct.cxx:126
 AliFlowAnalysisWithScalarProduct.cxx:127
 AliFlowAnalysisWithScalarProduct.cxx:128
 AliFlowAnalysisWithScalarProduct.cxx:129
 AliFlowAnalysisWithScalarProduct.cxx:130
 AliFlowAnalysisWithScalarProduct.cxx:131
 AliFlowAnalysisWithScalarProduct.cxx:132
 AliFlowAnalysisWithScalarProduct.cxx:133
 AliFlowAnalysisWithScalarProduct.cxx:134
 AliFlowAnalysisWithScalarProduct.cxx:135
 AliFlowAnalysisWithScalarProduct.cxx:136
 AliFlowAnalysisWithScalarProduct.cxx:137
 AliFlowAnalysisWithScalarProduct.cxx:138
 AliFlowAnalysisWithScalarProduct.cxx:139
 AliFlowAnalysisWithScalarProduct.cxx:140
 AliFlowAnalysisWithScalarProduct.cxx:141
 AliFlowAnalysisWithScalarProduct.cxx:142
 AliFlowAnalysisWithScalarProduct.cxx:143
 AliFlowAnalysisWithScalarProduct.cxx:144
 AliFlowAnalysisWithScalarProduct.cxx:145
 AliFlowAnalysisWithScalarProduct.cxx:146
 AliFlowAnalysisWithScalarProduct.cxx:147
 AliFlowAnalysisWithScalarProduct.cxx:148
 AliFlowAnalysisWithScalarProduct.cxx:149
 AliFlowAnalysisWithScalarProduct.cxx:150
 AliFlowAnalysisWithScalarProduct.cxx:151
 AliFlowAnalysisWithScalarProduct.cxx:152
 AliFlowAnalysisWithScalarProduct.cxx:153
 AliFlowAnalysisWithScalarProduct.cxx:154
 AliFlowAnalysisWithScalarProduct.cxx:155
 AliFlowAnalysisWithScalarProduct.cxx:156
 AliFlowAnalysisWithScalarProduct.cxx:157
 AliFlowAnalysisWithScalarProduct.cxx:158
 AliFlowAnalysisWithScalarProduct.cxx:159
 AliFlowAnalysisWithScalarProduct.cxx:160
 AliFlowAnalysisWithScalarProduct.cxx:161
 AliFlowAnalysisWithScalarProduct.cxx:162
 AliFlowAnalysisWithScalarProduct.cxx:163
 AliFlowAnalysisWithScalarProduct.cxx:164
 AliFlowAnalysisWithScalarProduct.cxx:165
 AliFlowAnalysisWithScalarProduct.cxx:166
 AliFlowAnalysisWithScalarProduct.cxx:167
 AliFlowAnalysisWithScalarProduct.cxx:168
 AliFlowAnalysisWithScalarProduct.cxx:169
 AliFlowAnalysisWithScalarProduct.cxx:170
 AliFlowAnalysisWithScalarProduct.cxx:171
 AliFlowAnalysisWithScalarProduct.cxx:172
 AliFlowAnalysisWithScalarProduct.cxx:173
 AliFlowAnalysisWithScalarProduct.cxx:174
 AliFlowAnalysisWithScalarProduct.cxx:175
 AliFlowAnalysisWithScalarProduct.cxx:176
 AliFlowAnalysisWithScalarProduct.cxx:177
 AliFlowAnalysisWithScalarProduct.cxx:178
 AliFlowAnalysisWithScalarProduct.cxx:179
 AliFlowAnalysisWithScalarProduct.cxx:180
 AliFlowAnalysisWithScalarProduct.cxx:181
 AliFlowAnalysisWithScalarProduct.cxx:182
 AliFlowAnalysisWithScalarProduct.cxx:183
 AliFlowAnalysisWithScalarProduct.cxx:184
 AliFlowAnalysisWithScalarProduct.cxx:185
 AliFlowAnalysisWithScalarProduct.cxx:186
 AliFlowAnalysisWithScalarProduct.cxx:187
 AliFlowAnalysisWithScalarProduct.cxx:188
 AliFlowAnalysisWithScalarProduct.cxx:189
 AliFlowAnalysisWithScalarProduct.cxx:190
 AliFlowAnalysisWithScalarProduct.cxx:191
 AliFlowAnalysisWithScalarProduct.cxx:192
 AliFlowAnalysisWithScalarProduct.cxx:193
 AliFlowAnalysisWithScalarProduct.cxx:194
 AliFlowAnalysisWithScalarProduct.cxx:195
 AliFlowAnalysisWithScalarProduct.cxx:196
 AliFlowAnalysisWithScalarProduct.cxx:197
 AliFlowAnalysisWithScalarProduct.cxx:198
 AliFlowAnalysisWithScalarProduct.cxx:199
 AliFlowAnalysisWithScalarProduct.cxx:200
 AliFlowAnalysisWithScalarProduct.cxx:201
 AliFlowAnalysisWithScalarProduct.cxx:202
 AliFlowAnalysisWithScalarProduct.cxx:203
 AliFlowAnalysisWithScalarProduct.cxx:204
 AliFlowAnalysisWithScalarProduct.cxx:205
 AliFlowAnalysisWithScalarProduct.cxx:206
 AliFlowAnalysisWithScalarProduct.cxx:207
 AliFlowAnalysisWithScalarProduct.cxx:208
 AliFlowAnalysisWithScalarProduct.cxx:209
 AliFlowAnalysisWithScalarProduct.cxx:210
 AliFlowAnalysisWithScalarProduct.cxx:211
 AliFlowAnalysisWithScalarProduct.cxx:212
 AliFlowAnalysisWithScalarProduct.cxx:213
 AliFlowAnalysisWithScalarProduct.cxx:214
 AliFlowAnalysisWithScalarProduct.cxx:215
 AliFlowAnalysisWithScalarProduct.cxx:216
 AliFlowAnalysisWithScalarProduct.cxx:217
 AliFlowAnalysisWithScalarProduct.cxx:218
 AliFlowAnalysisWithScalarProduct.cxx:219
 AliFlowAnalysisWithScalarProduct.cxx:220
 AliFlowAnalysisWithScalarProduct.cxx:221
 AliFlowAnalysisWithScalarProduct.cxx:222
 AliFlowAnalysisWithScalarProduct.cxx:223
 AliFlowAnalysisWithScalarProduct.cxx:224
 AliFlowAnalysisWithScalarProduct.cxx:225
 AliFlowAnalysisWithScalarProduct.cxx:226
 AliFlowAnalysisWithScalarProduct.cxx:227
 AliFlowAnalysisWithScalarProduct.cxx:228
 AliFlowAnalysisWithScalarProduct.cxx:229
 AliFlowAnalysisWithScalarProduct.cxx:230
 AliFlowAnalysisWithScalarProduct.cxx:231
 AliFlowAnalysisWithScalarProduct.cxx:232
 AliFlowAnalysisWithScalarProduct.cxx:233
 AliFlowAnalysisWithScalarProduct.cxx:234
 AliFlowAnalysisWithScalarProduct.cxx:235
 AliFlowAnalysisWithScalarProduct.cxx:236
 AliFlowAnalysisWithScalarProduct.cxx:237
 AliFlowAnalysisWithScalarProduct.cxx:238
 AliFlowAnalysisWithScalarProduct.cxx:239
 AliFlowAnalysisWithScalarProduct.cxx:240
 AliFlowAnalysisWithScalarProduct.cxx:241
 AliFlowAnalysisWithScalarProduct.cxx:242
 AliFlowAnalysisWithScalarProduct.cxx:243
 AliFlowAnalysisWithScalarProduct.cxx:244
 AliFlowAnalysisWithScalarProduct.cxx:245
 AliFlowAnalysisWithScalarProduct.cxx:246
 AliFlowAnalysisWithScalarProduct.cxx:247
 AliFlowAnalysisWithScalarProduct.cxx:248
 AliFlowAnalysisWithScalarProduct.cxx:249
 AliFlowAnalysisWithScalarProduct.cxx:250
 AliFlowAnalysisWithScalarProduct.cxx:251
 AliFlowAnalysisWithScalarProduct.cxx:252
 AliFlowAnalysisWithScalarProduct.cxx:253
 AliFlowAnalysisWithScalarProduct.cxx:254
 AliFlowAnalysisWithScalarProduct.cxx:255
 AliFlowAnalysisWithScalarProduct.cxx:256
 AliFlowAnalysisWithScalarProduct.cxx:257
 AliFlowAnalysisWithScalarProduct.cxx:258
 AliFlowAnalysisWithScalarProduct.cxx:259
 AliFlowAnalysisWithScalarProduct.cxx:260
 AliFlowAnalysisWithScalarProduct.cxx:261
 AliFlowAnalysisWithScalarProduct.cxx:262
 AliFlowAnalysisWithScalarProduct.cxx:263
 AliFlowAnalysisWithScalarProduct.cxx:264
 AliFlowAnalysisWithScalarProduct.cxx:265
 AliFlowAnalysisWithScalarProduct.cxx:266
 AliFlowAnalysisWithScalarProduct.cxx:267
 AliFlowAnalysisWithScalarProduct.cxx:268
 AliFlowAnalysisWithScalarProduct.cxx:269
 AliFlowAnalysisWithScalarProduct.cxx:270
 AliFlowAnalysisWithScalarProduct.cxx:271
 AliFlowAnalysisWithScalarProduct.cxx:272
 AliFlowAnalysisWithScalarProduct.cxx:273
 AliFlowAnalysisWithScalarProduct.cxx:274
 AliFlowAnalysisWithScalarProduct.cxx:275
 AliFlowAnalysisWithScalarProduct.cxx:276
 AliFlowAnalysisWithScalarProduct.cxx:277
 AliFlowAnalysisWithScalarProduct.cxx:278
 AliFlowAnalysisWithScalarProduct.cxx:279
 AliFlowAnalysisWithScalarProduct.cxx:280
 AliFlowAnalysisWithScalarProduct.cxx:281
 AliFlowAnalysisWithScalarProduct.cxx:282
 AliFlowAnalysisWithScalarProduct.cxx:283
 AliFlowAnalysisWithScalarProduct.cxx:284
 AliFlowAnalysisWithScalarProduct.cxx:285
 AliFlowAnalysisWithScalarProduct.cxx:286
 AliFlowAnalysisWithScalarProduct.cxx:287
 AliFlowAnalysisWithScalarProduct.cxx:288
 AliFlowAnalysisWithScalarProduct.cxx:289
 AliFlowAnalysisWithScalarProduct.cxx:290
 AliFlowAnalysisWithScalarProduct.cxx:291
 AliFlowAnalysisWithScalarProduct.cxx:292
 AliFlowAnalysisWithScalarProduct.cxx:293
 AliFlowAnalysisWithScalarProduct.cxx:294
 AliFlowAnalysisWithScalarProduct.cxx:295
 AliFlowAnalysisWithScalarProduct.cxx:296
 AliFlowAnalysisWithScalarProduct.cxx:297
 AliFlowAnalysisWithScalarProduct.cxx:298
 AliFlowAnalysisWithScalarProduct.cxx:299
 AliFlowAnalysisWithScalarProduct.cxx:300
 AliFlowAnalysisWithScalarProduct.cxx:301
 AliFlowAnalysisWithScalarProduct.cxx:302
 AliFlowAnalysisWithScalarProduct.cxx:303
 AliFlowAnalysisWithScalarProduct.cxx:304
 AliFlowAnalysisWithScalarProduct.cxx:305
 AliFlowAnalysisWithScalarProduct.cxx:306
 AliFlowAnalysisWithScalarProduct.cxx:307
 AliFlowAnalysisWithScalarProduct.cxx:308
 AliFlowAnalysisWithScalarProduct.cxx:309
 AliFlowAnalysisWithScalarProduct.cxx:310
 AliFlowAnalysisWithScalarProduct.cxx:311
 AliFlowAnalysisWithScalarProduct.cxx:312
 AliFlowAnalysisWithScalarProduct.cxx:313
 AliFlowAnalysisWithScalarProduct.cxx:314
 AliFlowAnalysisWithScalarProduct.cxx:315
 AliFlowAnalysisWithScalarProduct.cxx:316
 AliFlowAnalysisWithScalarProduct.cxx:317
 AliFlowAnalysisWithScalarProduct.cxx:318
 AliFlowAnalysisWithScalarProduct.cxx:319
 AliFlowAnalysisWithScalarProduct.cxx:320
 AliFlowAnalysisWithScalarProduct.cxx:321
 AliFlowAnalysisWithScalarProduct.cxx:322
 AliFlowAnalysisWithScalarProduct.cxx:323
 AliFlowAnalysisWithScalarProduct.cxx:324
 AliFlowAnalysisWithScalarProduct.cxx:325
 AliFlowAnalysisWithScalarProduct.cxx:326
 AliFlowAnalysisWithScalarProduct.cxx:327
 AliFlowAnalysisWithScalarProduct.cxx:328
 AliFlowAnalysisWithScalarProduct.cxx:329
 AliFlowAnalysisWithScalarProduct.cxx:330
 AliFlowAnalysisWithScalarProduct.cxx:331
 AliFlowAnalysisWithScalarProduct.cxx:332
 AliFlowAnalysisWithScalarProduct.cxx:333
 AliFlowAnalysisWithScalarProduct.cxx:334
 AliFlowAnalysisWithScalarProduct.cxx:335
 AliFlowAnalysisWithScalarProduct.cxx:336
 AliFlowAnalysisWithScalarProduct.cxx:337
 AliFlowAnalysisWithScalarProduct.cxx:338
 AliFlowAnalysisWithScalarProduct.cxx:339
 AliFlowAnalysisWithScalarProduct.cxx:340
 AliFlowAnalysisWithScalarProduct.cxx:341
 AliFlowAnalysisWithScalarProduct.cxx:342
 AliFlowAnalysisWithScalarProduct.cxx:343
 AliFlowAnalysisWithScalarProduct.cxx:344
 AliFlowAnalysisWithScalarProduct.cxx:345
 AliFlowAnalysisWithScalarProduct.cxx:346
 AliFlowAnalysisWithScalarProduct.cxx:347
 AliFlowAnalysisWithScalarProduct.cxx:348
 AliFlowAnalysisWithScalarProduct.cxx:349
 AliFlowAnalysisWithScalarProduct.cxx:350
 AliFlowAnalysisWithScalarProduct.cxx:351
 AliFlowAnalysisWithScalarProduct.cxx:352
 AliFlowAnalysisWithScalarProduct.cxx:353
 AliFlowAnalysisWithScalarProduct.cxx:354
 AliFlowAnalysisWithScalarProduct.cxx:355
 AliFlowAnalysisWithScalarProduct.cxx:356
 AliFlowAnalysisWithScalarProduct.cxx:357
 AliFlowAnalysisWithScalarProduct.cxx:358
 AliFlowAnalysisWithScalarProduct.cxx:359
 AliFlowAnalysisWithScalarProduct.cxx:360
 AliFlowAnalysisWithScalarProduct.cxx:361
 AliFlowAnalysisWithScalarProduct.cxx:362
 AliFlowAnalysisWithScalarProduct.cxx:363
 AliFlowAnalysisWithScalarProduct.cxx:364
 AliFlowAnalysisWithScalarProduct.cxx:365
 AliFlowAnalysisWithScalarProduct.cxx:366
 AliFlowAnalysisWithScalarProduct.cxx:367
 AliFlowAnalysisWithScalarProduct.cxx:368
 AliFlowAnalysisWithScalarProduct.cxx:369
 AliFlowAnalysisWithScalarProduct.cxx:370
 AliFlowAnalysisWithScalarProduct.cxx:371
 AliFlowAnalysisWithScalarProduct.cxx:372
 AliFlowAnalysisWithScalarProduct.cxx:373
 AliFlowAnalysisWithScalarProduct.cxx:374
 AliFlowAnalysisWithScalarProduct.cxx:375
 AliFlowAnalysisWithScalarProduct.cxx:376
 AliFlowAnalysisWithScalarProduct.cxx:377
 AliFlowAnalysisWithScalarProduct.cxx:378
 AliFlowAnalysisWithScalarProduct.cxx:379
 AliFlowAnalysisWithScalarProduct.cxx:380
 AliFlowAnalysisWithScalarProduct.cxx:381
 AliFlowAnalysisWithScalarProduct.cxx:382
 AliFlowAnalysisWithScalarProduct.cxx:383
 AliFlowAnalysisWithScalarProduct.cxx:384
 AliFlowAnalysisWithScalarProduct.cxx:385
 AliFlowAnalysisWithScalarProduct.cxx:386
 AliFlowAnalysisWithScalarProduct.cxx:387
 AliFlowAnalysisWithScalarProduct.cxx:388
 AliFlowAnalysisWithScalarProduct.cxx:389
 AliFlowAnalysisWithScalarProduct.cxx:390
 AliFlowAnalysisWithScalarProduct.cxx:391
 AliFlowAnalysisWithScalarProduct.cxx:392
 AliFlowAnalysisWithScalarProduct.cxx:393
 AliFlowAnalysisWithScalarProduct.cxx:394
 AliFlowAnalysisWithScalarProduct.cxx:395
 AliFlowAnalysisWithScalarProduct.cxx:396
 AliFlowAnalysisWithScalarProduct.cxx:397
 AliFlowAnalysisWithScalarProduct.cxx:398
 AliFlowAnalysisWithScalarProduct.cxx:399
 AliFlowAnalysisWithScalarProduct.cxx:400
 AliFlowAnalysisWithScalarProduct.cxx:401
 AliFlowAnalysisWithScalarProduct.cxx:402
 AliFlowAnalysisWithScalarProduct.cxx:403
 AliFlowAnalysisWithScalarProduct.cxx:404
 AliFlowAnalysisWithScalarProduct.cxx:405
 AliFlowAnalysisWithScalarProduct.cxx:406
 AliFlowAnalysisWithScalarProduct.cxx:407
 AliFlowAnalysisWithScalarProduct.cxx:408
 AliFlowAnalysisWithScalarProduct.cxx:409
 AliFlowAnalysisWithScalarProduct.cxx:410
 AliFlowAnalysisWithScalarProduct.cxx:411
 AliFlowAnalysisWithScalarProduct.cxx:412
 AliFlowAnalysisWithScalarProduct.cxx:413
 AliFlowAnalysisWithScalarProduct.cxx:414
 AliFlowAnalysisWithScalarProduct.cxx:415
 AliFlowAnalysisWithScalarProduct.cxx:416
 AliFlowAnalysisWithScalarProduct.cxx:417
 AliFlowAnalysisWithScalarProduct.cxx:418
 AliFlowAnalysisWithScalarProduct.cxx:419
 AliFlowAnalysisWithScalarProduct.cxx:420
 AliFlowAnalysisWithScalarProduct.cxx:421
 AliFlowAnalysisWithScalarProduct.cxx:422
 AliFlowAnalysisWithScalarProduct.cxx:423
 AliFlowAnalysisWithScalarProduct.cxx:424
 AliFlowAnalysisWithScalarProduct.cxx:425
 AliFlowAnalysisWithScalarProduct.cxx:426
 AliFlowAnalysisWithScalarProduct.cxx:427
 AliFlowAnalysisWithScalarProduct.cxx:428
 AliFlowAnalysisWithScalarProduct.cxx:429
 AliFlowAnalysisWithScalarProduct.cxx:430
 AliFlowAnalysisWithScalarProduct.cxx:431
 AliFlowAnalysisWithScalarProduct.cxx:432
 AliFlowAnalysisWithScalarProduct.cxx:433
 AliFlowAnalysisWithScalarProduct.cxx:434
 AliFlowAnalysisWithScalarProduct.cxx:435
 AliFlowAnalysisWithScalarProduct.cxx:436
 AliFlowAnalysisWithScalarProduct.cxx:437
 AliFlowAnalysisWithScalarProduct.cxx:438
 AliFlowAnalysisWithScalarProduct.cxx:439
 AliFlowAnalysisWithScalarProduct.cxx:440
 AliFlowAnalysisWithScalarProduct.cxx:441
 AliFlowAnalysisWithScalarProduct.cxx:442
 AliFlowAnalysisWithScalarProduct.cxx:443
 AliFlowAnalysisWithScalarProduct.cxx:444
 AliFlowAnalysisWithScalarProduct.cxx:445
 AliFlowAnalysisWithScalarProduct.cxx:446
 AliFlowAnalysisWithScalarProduct.cxx:447
 AliFlowAnalysisWithScalarProduct.cxx:448
 AliFlowAnalysisWithScalarProduct.cxx:449
 AliFlowAnalysisWithScalarProduct.cxx:450
 AliFlowAnalysisWithScalarProduct.cxx:451
 AliFlowAnalysisWithScalarProduct.cxx:452
 AliFlowAnalysisWithScalarProduct.cxx:453
 AliFlowAnalysisWithScalarProduct.cxx:454
 AliFlowAnalysisWithScalarProduct.cxx:455
 AliFlowAnalysisWithScalarProduct.cxx:456
 AliFlowAnalysisWithScalarProduct.cxx:457
 AliFlowAnalysisWithScalarProduct.cxx:458
 AliFlowAnalysisWithScalarProduct.cxx:459
 AliFlowAnalysisWithScalarProduct.cxx:460
 AliFlowAnalysisWithScalarProduct.cxx:461
 AliFlowAnalysisWithScalarProduct.cxx:462
 AliFlowAnalysisWithScalarProduct.cxx:463
 AliFlowAnalysisWithScalarProduct.cxx:464
 AliFlowAnalysisWithScalarProduct.cxx:465
 AliFlowAnalysisWithScalarProduct.cxx:466
 AliFlowAnalysisWithScalarProduct.cxx:467
 AliFlowAnalysisWithScalarProduct.cxx:468
 AliFlowAnalysisWithScalarProduct.cxx:469
 AliFlowAnalysisWithScalarProduct.cxx:470
 AliFlowAnalysisWithScalarProduct.cxx:471
 AliFlowAnalysisWithScalarProduct.cxx:472
 AliFlowAnalysisWithScalarProduct.cxx:473
 AliFlowAnalysisWithScalarProduct.cxx:474
 AliFlowAnalysisWithScalarProduct.cxx:475
 AliFlowAnalysisWithScalarProduct.cxx:476
 AliFlowAnalysisWithScalarProduct.cxx:477
 AliFlowAnalysisWithScalarProduct.cxx:478
 AliFlowAnalysisWithScalarProduct.cxx:479
 AliFlowAnalysisWithScalarProduct.cxx:480
 AliFlowAnalysisWithScalarProduct.cxx:481
 AliFlowAnalysisWithScalarProduct.cxx:482
 AliFlowAnalysisWithScalarProduct.cxx:483
 AliFlowAnalysisWithScalarProduct.cxx:484
 AliFlowAnalysisWithScalarProduct.cxx:485
 AliFlowAnalysisWithScalarProduct.cxx:486
 AliFlowAnalysisWithScalarProduct.cxx:487
 AliFlowAnalysisWithScalarProduct.cxx:488
 AliFlowAnalysisWithScalarProduct.cxx:489
 AliFlowAnalysisWithScalarProduct.cxx:490
 AliFlowAnalysisWithScalarProduct.cxx:491
 AliFlowAnalysisWithScalarProduct.cxx:492
 AliFlowAnalysisWithScalarProduct.cxx:493
 AliFlowAnalysisWithScalarProduct.cxx:494
 AliFlowAnalysisWithScalarProduct.cxx:495
 AliFlowAnalysisWithScalarProduct.cxx:496
 AliFlowAnalysisWithScalarProduct.cxx:497
 AliFlowAnalysisWithScalarProduct.cxx:498
 AliFlowAnalysisWithScalarProduct.cxx:499
 AliFlowAnalysisWithScalarProduct.cxx:500
 AliFlowAnalysisWithScalarProduct.cxx:501
 AliFlowAnalysisWithScalarProduct.cxx:502
 AliFlowAnalysisWithScalarProduct.cxx:503
 AliFlowAnalysisWithScalarProduct.cxx:504
 AliFlowAnalysisWithScalarProduct.cxx:505
 AliFlowAnalysisWithScalarProduct.cxx:506
 AliFlowAnalysisWithScalarProduct.cxx:507
 AliFlowAnalysisWithScalarProduct.cxx:508
 AliFlowAnalysisWithScalarProduct.cxx:509
 AliFlowAnalysisWithScalarProduct.cxx:510
 AliFlowAnalysisWithScalarProduct.cxx:511
 AliFlowAnalysisWithScalarProduct.cxx:512
 AliFlowAnalysisWithScalarProduct.cxx:513
 AliFlowAnalysisWithScalarProduct.cxx:514
 AliFlowAnalysisWithScalarProduct.cxx:515
 AliFlowAnalysisWithScalarProduct.cxx:516
 AliFlowAnalysisWithScalarProduct.cxx:517
 AliFlowAnalysisWithScalarProduct.cxx:518
 AliFlowAnalysisWithScalarProduct.cxx:519
 AliFlowAnalysisWithScalarProduct.cxx:520
 AliFlowAnalysisWithScalarProduct.cxx:521
 AliFlowAnalysisWithScalarProduct.cxx:522
 AliFlowAnalysisWithScalarProduct.cxx:523
 AliFlowAnalysisWithScalarProduct.cxx:524
 AliFlowAnalysisWithScalarProduct.cxx:525
 AliFlowAnalysisWithScalarProduct.cxx:526
 AliFlowAnalysisWithScalarProduct.cxx:527
 AliFlowAnalysisWithScalarProduct.cxx:528
 AliFlowAnalysisWithScalarProduct.cxx:529
 AliFlowAnalysisWithScalarProduct.cxx:530
 AliFlowAnalysisWithScalarProduct.cxx:531
 AliFlowAnalysisWithScalarProduct.cxx:532
 AliFlowAnalysisWithScalarProduct.cxx:533
 AliFlowAnalysisWithScalarProduct.cxx:534
 AliFlowAnalysisWithScalarProduct.cxx:535
 AliFlowAnalysisWithScalarProduct.cxx:536
 AliFlowAnalysisWithScalarProduct.cxx:537
 AliFlowAnalysisWithScalarProduct.cxx:538
 AliFlowAnalysisWithScalarProduct.cxx:539
 AliFlowAnalysisWithScalarProduct.cxx:540
 AliFlowAnalysisWithScalarProduct.cxx:541
 AliFlowAnalysisWithScalarProduct.cxx:542
 AliFlowAnalysisWithScalarProduct.cxx:543
 AliFlowAnalysisWithScalarProduct.cxx:544
 AliFlowAnalysisWithScalarProduct.cxx:545
 AliFlowAnalysisWithScalarProduct.cxx:546
 AliFlowAnalysisWithScalarProduct.cxx:547
 AliFlowAnalysisWithScalarProduct.cxx:548
 AliFlowAnalysisWithScalarProduct.cxx:549
 AliFlowAnalysisWithScalarProduct.cxx:550
 AliFlowAnalysisWithScalarProduct.cxx:551
 AliFlowAnalysisWithScalarProduct.cxx:552
 AliFlowAnalysisWithScalarProduct.cxx:553
 AliFlowAnalysisWithScalarProduct.cxx:554
 AliFlowAnalysisWithScalarProduct.cxx:555
 AliFlowAnalysisWithScalarProduct.cxx:556
 AliFlowAnalysisWithScalarProduct.cxx:557
 AliFlowAnalysisWithScalarProduct.cxx:558
 AliFlowAnalysisWithScalarProduct.cxx:559
 AliFlowAnalysisWithScalarProduct.cxx:560
 AliFlowAnalysisWithScalarProduct.cxx:561
 AliFlowAnalysisWithScalarProduct.cxx:562
 AliFlowAnalysisWithScalarProduct.cxx:563
 AliFlowAnalysisWithScalarProduct.cxx:564
 AliFlowAnalysisWithScalarProduct.cxx:565
 AliFlowAnalysisWithScalarProduct.cxx:566
 AliFlowAnalysisWithScalarProduct.cxx:567
 AliFlowAnalysisWithScalarProduct.cxx:568
 AliFlowAnalysisWithScalarProduct.cxx:569
 AliFlowAnalysisWithScalarProduct.cxx:570
 AliFlowAnalysisWithScalarProduct.cxx:571
 AliFlowAnalysisWithScalarProduct.cxx:572
 AliFlowAnalysisWithScalarProduct.cxx:573
 AliFlowAnalysisWithScalarProduct.cxx:574
 AliFlowAnalysisWithScalarProduct.cxx:575
 AliFlowAnalysisWithScalarProduct.cxx:576
 AliFlowAnalysisWithScalarProduct.cxx:577
 AliFlowAnalysisWithScalarProduct.cxx:578
 AliFlowAnalysisWithScalarProduct.cxx:579
 AliFlowAnalysisWithScalarProduct.cxx:580
 AliFlowAnalysisWithScalarProduct.cxx:581
 AliFlowAnalysisWithScalarProduct.cxx:582
 AliFlowAnalysisWithScalarProduct.cxx:583
 AliFlowAnalysisWithScalarProduct.cxx:584
 AliFlowAnalysisWithScalarProduct.cxx:585
 AliFlowAnalysisWithScalarProduct.cxx:586
 AliFlowAnalysisWithScalarProduct.cxx:587
 AliFlowAnalysisWithScalarProduct.cxx:588
 AliFlowAnalysisWithScalarProduct.cxx:589
 AliFlowAnalysisWithScalarProduct.cxx:590
 AliFlowAnalysisWithScalarProduct.cxx:591
 AliFlowAnalysisWithScalarProduct.cxx:592
 AliFlowAnalysisWithScalarProduct.cxx:593
 AliFlowAnalysisWithScalarProduct.cxx:594
 AliFlowAnalysisWithScalarProduct.cxx:595
 AliFlowAnalysisWithScalarProduct.cxx:596
 AliFlowAnalysisWithScalarProduct.cxx:597
 AliFlowAnalysisWithScalarProduct.cxx:598
 AliFlowAnalysisWithScalarProduct.cxx:599
 AliFlowAnalysisWithScalarProduct.cxx:600
 AliFlowAnalysisWithScalarProduct.cxx:601
 AliFlowAnalysisWithScalarProduct.cxx:602
 AliFlowAnalysisWithScalarProduct.cxx:603
 AliFlowAnalysisWithScalarProduct.cxx:604
 AliFlowAnalysisWithScalarProduct.cxx:605
 AliFlowAnalysisWithScalarProduct.cxx:606
 AliFlowAnalysisWithScalarProduct.cxx:607
 AliFlowAnalysisWithScalarProduct.cxx:608
 AliFlowAnalysisWithScalarProduct.cxx:609
 AliFlowAnalysisWithScalarProduct.cxx:610
 AliFlowAnalysisWithScalarProduct.cxx:611
 AliFlowAnalysisWithScalarProduct.cxx:612
 AliFlowAnalysisWithScalarProduct.cxx:613
 AliFlowAnalysisWithScalarProduct.cxx:614
 AliFlowAnalysisWithScalarProduct.cxx:615
 AliFlowAnalysisWithScalarProduct.cxx:616
 AliFlowAnalysisWithScalarProduct.cxx:617
 AliFlowAnalysisWithScalarProduct.cxx:618
 AliFlowAnalysisWithScalarProduct.cxx:619
 AliFlowAnalysisWithScalarProduct.cxx:620
 AliFlowAnalysisWithScalarProduct.cxx:621
 AliFlowAnalysisWithScalarProduct.cxx:622
 AliFlowAnalysisWithScalarProduct.cxx:623
 AliFlowAnalysisWithScalarProduct.cxx:624
 AliFlowAnalysisWithScalarProduct.cxx:625
 AliFlowAnalysisWithScalarProduct.cxx:626
 AliFlowAnalysisWithScalarProduct.cxx:627
 AliFlowAnalysisWithScalarProduct.cxx:628
 AliFlowAnalysisWithScalarProduct.cxx:629
 AliFlowAnalysisWithScalarProduct.cxx:630
 AliFlowAnalysisWithScalarProduct.cxx:631
 AliFlowAnalysisWithScalarProduct.cxx:632
 AliFlowAnalysisWithScalarProduct.cxx:633
 AliFlowAnalysisWithScalarProduct.cxx:634
 AliFlowAnalysisWithScalarProduct.cxx:635
 AliFlowAnalysisWithScalarProduct.cxx:636
 AliFlowAnalysisWithScalarProduct.cxx:637
 AliFlowAnalysisWithScalarProduct.cxx:638
 AliFlowAnalysisWithScalarProduct.cxx:639
 AliFlowAnalysisWithScalarProduct.cxx:640
 AliFlowAnalysisWithScalarProduct.cxx:641
 AliFlowAnalysisWithScalarProduct.cxx:642
 AliFlowAnalysisWithScalarProduct.cxx:643
 AliFlowAnalysisWithScalarProduct.cxx:644
 AliFlowAnalysisWithScalarProduct.cxx:645
 AliFlowAnalysisWithScalarProduct.cxx:646
 AliFlowAnalysisWithScalarProduct.cxx:647
 AliFlowAnalysisWithScalarProduct.cxx:648
 AliFlowAnalysisWithScalarProduct.cxx:649
 AliFlowAnalysisWithScalarProduct.cxx:650
 AliFlowAnalysisWithScalarProduct.cxx:651
 AliFlowAnalysisWithScalarProduct.cxx:652
 AliFlowAnalysisWithScalarProduct.cxx:653
 AliFlowAnalysisWithScalarProduct.cxx:654
 AliFlowAnalysisWithScalarProduct.cxx:655
 AliFlowAnalysisWithScalarProduct.cxx:656
 AliFlowAnalysisWithScalarProduct.cxx:657
 AliFlowAnalysisWithScalarProduct.cxx:658
 AliFlowAnalysisWithScalarProduct.cxx:659
 AliFlowAnalysisWithScalarProduct.cxx:660
 AliFlowAnalysisWithScalarProduct.cxx:661
 AliFlowAnalysisWithScalarProduct.cxx:662
 AliFlowAnalysisWithScalarProduct.cxx:663
 AliFlowAnalysisWithScalarProduct.cxx:664
 AliFlowAnalysisWithScalarProduct.cxx:665
 AliFlowAnalysisWithScalarProduct.cxx:666
 AliFlowAnalysisWithScalarProduct.cxx:667
 AliFlowAnalysisWithScalarProduct.cxx:668
 AliFlowAnalysisWithScalarProduct.cxx:669
 AliFlowAnalysisWithScalarProduct.cxx:670
 AliFlowAnalysisWithScalarProduct.cxx:671
 AliFlowAnalysisWithScalarProduct.cxx:672
 AliFlowAnalysisWithScalarProduct.cxx:673
 AliFlowAnalysisWithScalarProduct.cxx:674
 AliFlowAnalysisWithScalarProduct.cxx:675
 AliFlowAnalysisWithScalarProduct.cxx:676
 AliFlowAnalysisWithScalarProduct.cxx:677
 AliFlowAnalysisWithScalarProduct.cxx:678
 AliFlowAnalysisWithScalarProduct.cxx:679
 AliFlowAnalysisWithScalarProduct.cxx:680
 AliFlowAnalysisWithScalarProduct.cxx:681
 AliFlowAnalysisWithScalarProduct.cxx:682
 AliFlowAnalysisWithScalarProduct.cxx:683
 AliFlowAnalysisWithScalarProduct.cxx:684
 AliFlowAnalysisWithScalarProduct.cxx:685
 AliFlowAnalysisWithScalarProduct.cxx:686
 AliFlowAnalysisWithScalarProduct.cxx:687
 AliFlowAnalysisWithScalarProduct.cxx:688
 AliFlowAnalysisWithScalarProduct.cxx:689
 AliFlowAnalysisWithScalarProduct.cxx:690
 AliFlowAnalysisWithScalarProduct.cxx:691
 AliFlowAnalysisWithScalarProduct.cxx:692
 AliFlowAnalysisWithScalarProduct.cxx:693
 AliFlowAnalysisWithScalarProduct.cxx:694
 AliFlowAnalysisWithScalarProduct.cxx:695
 AliFlowAnalysisWithScalarProduct.cxx:696
 AliFlowAnalysisWithScalarProduct.cxx:697
 AliFlowAnalysisWithScalarProduct.cxx:698
 AliFlowAnalysisWithScalarProduct.cxx:699
 AliFlowAnalysisWithScalarProduct.cxx:700
 AliFlowAnalysisWithScalarProduct.cxx:701
 AliFlowAnalysisWithScalarProduct.cxx:702
 AliFlowAnalysisWithScalarProduct.cxx:703
 AliFlowAnalysisWithScalarProduct.cxx:704
 AliFlowAnalysisWithScalarProduct.cxx:705
 AliFlowAnalysisWithScalarProduct.cxx:706
 AliFlowAnalysisWithScalarProduct.cxx:707
 AliFlowAnalysisWithScalarProduct.cxx:708
 AliFlowAnalysisWithScalarProduct.cxx:709
 AliFlowAnalysisWithScalarProduct.cxx:710
 AliFlowAnalysisWithScalarProduct.cxx:711
 AliFlowAnalysisWithScalarProduct.cxx:712
 AliFlowAnalysisWithScalarProduct.cxx:713
 AliFlowAnalysisWithScalarProduct.cxx:714
 AliFlowAnalysisWithScalarProduct.cxx:715