ROOT logo
/*************************************************************************
* Copyright(c) 1998-2009, 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.                  *
**************************************************************************/

/* $Id$ */

////////////////////////////////////////////////////////////////////////////
//  Class to perform pt-spectra (and ptMin-spectra) extraction of mothers //
//  particles starting from measured pt-spectra of daughter particles     //
//  that come from inclusive decays.                                      // 
//  E.g.: B->J/psi+X , B->e+X, B->D0+X, etc.                              //
//                                                                        //
//  In order to use this class, one first has to run a simulation         // 
//  (only kinematics) of the decay channel under study. The analysis      //
//  can be runned using the class AliAnalysisTaskPtMothFromPtDaugh        //  
//  which loops on events to create a TNtupla that stores just            // 
//  kinematics informations for mothers and daughters of the decay        //
//  under study (this is made in order to speed up).                      //
//                                                                        //
//  Therefore the standard inputs of this class are:                      //
//  (1) The TNtupla (created by the task using a TChain of galice.root)   //    
//  (2) pT-spectrum of the daughter particles                             //
//                                                                        //
//  Output would be the pT (and pTMin) spectrum of the mother, based      //
//  on the correction factors computed from the Kinematics.root files     //
//                                                                        //
//                                                                        //  
//  Authors: Giuseppe E. Bruno           &  Fiorella Fionda               //
//           (Giuseppe.Bruno@ba.infn.it)    (Fiorella.Fionda@ba.infn.it)  //
////////////////////////////////////////////////////////////////////////////

#include "TH1F.h"
#include "TNtuple.h"
#include "TFile.h"
#include "TParticle.h"
#include "TArrayI.h"

#include "AliPtMothFromPtDaugh.h"
#include "AliStack.h"
#include "AliLog.h"

ClassImp(AliPtMothFromPtDaugh)
//________________________________________________________________
AliPtMothFromPtDaugh::AliPtMothFromPtDaugh() :
  TNamed("AliPtMoth","AliPtMoth"),
  fDecayKine(0x0), 
  fWij(0x0),
  fFi(0x0),
  fWijMin(0x0),
  fFiMin(0x0),
  fHistoPtDaughter(0x0),
  fHistoPtMothers(0x0),
  fHistoPtMinMothers(0x0),
  fMothers(0x0),  
  fDaughter(0), 
  fyMothMax(0),
  fyMothMin(0),
  fyDaughMax(0),
  fyDaughMin(0),
  fUseEta(kFALSE),
  fAnalysisMode(kUserAnalysis)
  {
  //
  // Default constructor
  //
  }

//________________________________________________________________
AliPtMothFromPtDaugh::AliPtMothFromPtDaugh(const char* name, const char* title) :
  TNamed(name,title),
  fDecayKine(0x0),
  fWij(0x0),
  fFi(0x0),
  fWijMin(0x0),
  fFiMin(0x0),
  fHistoPtDaughter(0x0),
  fHistoPtMothers(0x0),
  fHistoPtMinMothers(0x0),
  fMothers(0x0),
  fDaughter(0),
  fyMothMax(0),
  fyMothMin(0),
  fyDaughMax(0),
  fyDaughMin(0),
  fUseEta(kFALSE),
  fAnalysisMode(kUserAnalysis)
  {
  //
  // Named constructor
  //
  }

//________________________________________________________________
AliPtMothFromPtDaugh::~AliPtMothFromPtDaugh()
  {
  //
  // Default destructor
  //
  if(fDecayKine) {delete fDecayKine;}
   fDecayKine=0;
   if(fMothers) {delete fMothers;}
   fMothers=0;
   if(fHistoPtMothers) { delete fHistoPtMothers; }
   fHistoPtMothers=0;
   if(fHistoPtMinMothers) { delete fHistoPtMinMothers;}
   fHistoPtMinMothers=0;
   if(fHistoPtDaughter) { delete fHistoPtDaughter; fHistoPtDaughter=0;}
  }

//______________________________________________________________________________________
AliPtMothFromPtDaugh::AliPtMothFromPtDaugh(const AliPtMothFromPtDaugh& extraction) :
  TNamed(extraction),
  fDecayKine(0),
  fWij(0),
  fFi(0),
  fWijMin(0),
  fFiMin(0),
  fHistoPtDaughter(0),
  fHistoPtMothers(0),
  fHistoPtMinMothers(0),
  fMothers(0),
  fDaughter(extraction.fDaughter),
  fyMothMax(extraction.fyMothMax),
  fyMothMin(extraction.fyMothMin),
  fyDaughMax(extraction.fyDaughMax),
  fyDaughMin(extraction.fyDaughMin),
  fUseEta(extraction.fUseEta),
  fAnalysisMode(extraction.fAnalysisMode)
    {
   // Copy constructor
   Int_t nbinsM=0;Int_t nbinsD=0;  Int_t nbinsMmin=0;
   if(extraction.fHistoPtDaughter) 
   {fHistoPtDaughter =(TH1F*)extraction.fHistoPtDaughter->Clone("fHistoPtDaughter_copy"); nbinsD = fHistoPtDaughter->GetNbinsX();}
   if(extraction.fHistoPtMothers) 
   {fHistoPtMothers = (TH1F*)extraction.fHistoPtMothers->Clone("fHistoPtMothers_copy"); nbinsM = fHistoPtMothers->GetNbinsX();}
   if(extraction.fHistoPtMinMothers) 
   {fHistoPtMinMothers = (TH1F*)extraction.fHistoPtMinMothers->Clone("fHistoPtMinMothers_copy"); nbinsMmin = fHistoPtMinMothers->GetNbinsX();}
 
   if(nbinsD>0 && nbinsM>0){
   if(extraction.fWij){
     fWij = new Double_t*[2*nbinsM]; 
     for(Int_t i=0;i<2*nbinsM;i++) *(fWij+i)=new Double_t[nbinsD];
        for(Int_t i=0;i<2*nbinsM;i++){
          for(Int_t j=0;j<nbinsD;j++){
             fWij[i][j]=extraction.fWij[i][j];
             } 
           } 
     } 
   
   if(extraction.fFi){
     fFi=new Double_t[2*nbinsM];
     for(Int_t i=0;i<2*nbinsM;i++) fFi[i]=extraction.fFi[i];
     }
   } // if nbinsD, nbinsM > 0 
  
  if(nbinsD>0 && nbinsMmin>0){
   if(extraction.fWijMin){ 
     fWijMin = new Double_t*[2*nbinsMmin];
     for(Int_t i=0;i<2*nbinsMmin;i++) *(fWijMin+i)=new Double_t[nbinsD];
        for(Int_t i=0;i<2*nbinsMmin;i++){
          for(Int_t j=0;j<nbinsD;j++){
             fWijMin[i][j]=extraction.fWijMin[i][j];
             }
           }
     }

   if(extraction.fFiMin){
     fFiMin=new Double_t[2*nbinsMmin];
     for(Int_t i=0;i<2*nbinsMmin;i++) fFiMin[i]=extraction.fFiMin[i];
     }

   } // if nbinsD, nbinsMmin > 0 
   
   if(extraction.fDecayKine) fDecayKine = (TNtuple*)extraction.fDecayKine->CloneTree(); 

   if(extraction.fMothers) fMothers = new TArrayI(*(extraction.fMothers));
  
   extraction.Copy(*this);
  }

//______________________________________________________________________________________
AliPtMothFromPtDaugh& AliPtMothFromPtDaugh::operator=(const AliPtMothFromPtDaugh &extraction)
  { 
    // operator assignment
    if (this!=&extraction) { 
    this->~AliPtMothFromPtDaugh();
    new(this) AliPtMothFromPtDaugh(extraction); 
    }
    return *this; 
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::CreateWeights()
  {
   // Set mothers and daughters pdg codes if not
   // Put a control if mothers histograms binning-size, rapidity (or pseudoapidity) range
   // are setting. Read daughter histogram (histName) from the file (pathFileName)
   // Initialize dimensions for correction factors

   DeleteWeights();
   if(!fMothers)  {AliError("Set pdg codes of mothers by SetPdgMothers!\n"); return kFALSE;}
   if(!fDaughter) {AliError("Set pdg code of daughter by SetPdgDaughter!\n"); return kFALSE;}  
   if(!fHistoPtDaughter) { AliError("Daughter histogram doesn't exist! \n"); return kFALSE;} 
   
   //Set Rapidity or Pseudorapidity range for mothers if not
   if(!fyMothMax || !fyMothMin ){ AliError("Set rapidity range or pseudoRapidity range for mothers: use SetYmothers(ymin,ymax) or SetEtaMothers(etamin,etamax)"); return kFALSE;}   
   if(!fyDaughMax || !fyDaughMin){ AliError("Set rapidity range or pseudoRapidity range for daughters:use SetYdaughters(ymin,ymax) or SetEtaDaughters(etamin,etamax)"); return kFALSE;}  
   if(!fHistoPtMothers) {AliError("Call method SetBinsPtMoth to set pT-histogram "); return kFALSE;}
   if(!fHistoPtMinMothers){AliError("Call method SetBinsPtMinMoth to set pTmin-histogram "); return kFALSE;}

   Int_t nbinsM=(fHistoPtMothers->GetNbinsX()+2);
   Int_t nbinsD=(fHistoPtDaughter->GetNbinsX()+2);
   Int_t nbinsMmin=(fHistoPtMinMothers->GetNbinsX()+2);

   //Create pointers for weights to reconstruct daughter and mothers pT-spectra 
   fWij=new Double_t*[2*nbinsM];
   for(Int_t i=0;i<2*nbinsM;i++) 
      {*(fWij+i)=new Double_t[nbinsD];}   
   fFi=new Double_t[2*nbinsM];   
    
   fWijMin=new Double_t*[2*nbinsMmin];
   for(Int_t i=0;i<2*nbinsMmin;i++) 
      {*(fWijMin+i)=new Double_t[nbinsD];}
   fFiMin=new Double_t[2*nbinsMmin];
   AliInfo(Form("Pt-mothers distribution: pt_min = %f  pt_max=%f  n_bins=%d \n",
           fHistoPtMothers->GetBinLowEdge(1),fHistoPtMothers->GetBinLowEdge(nbinsM-1),nbinsM-2));
   AliInfo(Form("PtMinimum-mothers distribution: pt_min = %f  pt_max=%f  n_bins=%d \n",
           fHistoPtMinMothers->GetBinLowEdge(1),fHistoPtMinMothers->GetBinLowEdge(nbinsMmin-1),nbinsMmin-2));
   AliInfo(Form("Pt-daughters distribution: pt_min = %f  pt_max=%f n_bins=%d \n",
           fHistoPtDaughter->GetBinLowEdge(1),fHistoPtDaughter->GetBinLowEdge(nbinsD-1),nbinsD-2));
   return kTRUE;
  }  

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::DeleteWeights()
  {
   //delete correction factors   
   //delete histogram of daughters
   if(!fHistoPtMothers || !fHistoPtMinMothers) {AliError("Mothers histograms don't exist! Cannot delete correction factors"); return;}
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;  
   if(fWij){
     for(Int_t i=0; i<(2*nbinsM); i++) delete fWij[i];
     delete [] fWij; fWij=0;
    }
   if(fFi) { delete fFi; fFi=0;} 
   if(fWijMin){
     for(Int_t i=0; i<(2*nbinsMmin); i++) delete fWijMin[i];
     delete [] fWijMin; fWijMin=0;
    }
   if(fFiMin) { delete fFiMin; fFiMin=0;}
 
   return;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::ReadHistoPtDaught(const TH1F *hist)
  {
   //Initialize daughter histograms with hist
   if(!hist) {AliError("Set correct histograms of daughter! It doesn't exist!\n"); return kFALSE;}
   if(fHistoPtDaughter) delete fHistoPtDaughter;
   fHistoPtDaughter = (TH1F*)hist->Clone();
   return kTRUE;
  }

//______________________________________________________________________________________
Double_t* AliPtMothFromPtDaugh::SetBinsSize(Double_t ptmin, Double_t ptmax,Int_t nbins,Double_t alpha)
  {
   // return a pointer of double which contains the binning size for (mothers) histograms:
   // alpha have to be > 0 and:  
   // alpha = 1 equal binning size
   // alpha < 1 increasing  " 
   // alpha > 1 decreasing  "
   if(ptmin<0 || ptmax<0 || nbins<=0 || alpha<=0) {AliError("Set correct bin-size: ptmin>=0, ptmax>=0, nbins>0, alpha>0! \n"); return 0;}
   Double_t *edgebin = new Double_t[nbins+1];
   Double_t ptmin1=TMath::Power(ptmin,alpha);
   Double_t ptmax1=TMath::Power(ptmax,alpha);
   Double_t size=(ptmax1-ptmin1)/nbins;
   for(Int_t i=0;i<nbins+1;i++) *(edgebin+i)=TMath::Power((ptmin1+i*size),(Double_t)1/alpha);
   return edgebin;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetBinsPtMoth(Double_t ptmin, Double_t ptmax,Int_t nbins,Double_t alpha)
  {
   // Set bin size for pt-spectrum of mothers using SetBinsSize:
   // alpha have to be > 0 and:
   // alpha = 1 equal binning size
   // alpha < 1 increasing  " 
   // alpha > 1 decreasing  "
   Double_t* edges = SetBinsSize(ptmin,ptmax,nbins,alpha);
   SetBinsPtMoth(nbins,edges);
   delete [] edges;
   return;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetBinsPtMoth(Int_t nbins,Double_t *edgeBins)
  {
   //set bin size given by the pointer edgeBins for pt-spectrum of mothers:
   //the dimension of the pointer edgeBins is nbins+1 and the points
   //has to be written in increasing order 
   if(nbins<0) {AliError("Numbers of bins should be > 0 !\n"); return;}
   if(fHistoPtMothers) {delete fHistoPtMothers; fHistoPtMothers=0;}
   fHistoPtMothers=new TH1F("fHistoPtMothers","Reconstructed p_{T}(Mothers)-spectrum",nbins,edgeBins);   
   return;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetBinsPtMinMoth(Int_t nbins,Double_t *edgeBins)
  {
   //set bin size given by the pointer edgeBins for ptMin-spectrum of mothers:
   //the dimension of the pointer edgeBins is nbins+1 and the points
   //has to be written in increasing order 
   if(nbins<0) {AliError("Numbers of bins should be > 0 !\n"); return;}
   if(fHistoPtMinMothers) {delete fHistoPtMinMothers; fHistoPtMinMothers=0;}
   fHistoPtMinMothers = new TH1F("fHistoPtMinMothers","Reconstructed p_{T}^{MIN}(Mothers)-spectrum",nbins,edgeBins);
   return;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetBinsPtMinMoth(Double_t ptmin, Double_t ptmax,Int_t nbins,Double_t alpha)
  {
   // Set bin size for ptMin-spectrum of mothers using SetBinsSize:
   // alpha have to be > 0 and:
   // alpha = 1 equal binning size
   // alpha < 1 increasing  " 
   // alpha > 1 decreasing  "
   Double_t* edges = SetBinsSize(ptmin,ptmax,nbins,alpha);
   SetBinsPtMinMoth(nbins,edges);
   delete [] edges;
   return;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetPdgMothersPrivate(Int_t n_mothers,Int_t *pdgM)
  {
   // Set pdg codes of mothers given by the pointer pdgM for the analysis. 
   // This is a private method. 
   if(fMothers) { delete fMothers; fMothers = 0; }
   fMothers = new TArrayI(n_mothers,pdgM); 
   return;
  }
//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetPdgMothers(Int_t n_mothers,Int_t *pdgM)
  {
   // Set user pdg codes of mothers: first check 
   // that the kUserAnalysis is the selected Analysis_mode. 
   // If not print out a message of error. 
   if(fAnalysisMode!=kUserAnalysis) {
     AliError("Nothing done: set the mode to  kUserAnalysis first");
     return;
   }
   //Set pdg codes of mothers given by the pointer pdgM for the analysis 
   SetPdgMothersPrivate(n_mothers,pdgM);
   return;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetBeautyMothers()
  {
   // Set pdg codes of beauty particles:
   // B-mesons (1-24) B-barions (25-59)
   Int_t pdgBeauty[]={511,521,513,523,10511,10521,10513,10523,20513,
      20523,515,525,531,10531,533,10533,205433,535,541,10541,543,10543,
      20543,545,5122,5112,5212,5222,5114,5214,5224,5132,5232,5312,5322,
      5314,5324,5332,5334,5142,5242,5412,5422,5414,5424,5342,5432,5434,
      5442,5444,5512,5522,5514,5524,5532,5534,5542,5544,5554};
   Int_t *pdgB=new Int_t[59];
   for(Int_t i=0;i<59;i++) pdgB[i]=pdgBeauty[i];
   SetPdgMothersPrivate(59,pdgB);
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::InitDefaultAnalysis()
 {
  // Set mothers and daughter pdg codes depending from the selected analysis. 
  // case kUserAnalysis: mothers and daughter are set by user (nothing to be done)
  // case kBtoJPSI: inclusive B-> J/Psi + X channels 
  // case kBtoEle: inclusive B-> e + X channels
  // case kBtoMuon: inclusive B-> mu + X channels
  // case kBtoD0: inclusive B-> D0 + X channels 
  
  switch(fAnalysisMode)
   {
   case kUserAnalysis:
    break; 
   case kBtoJPSI:
    SetBeautyMothers();
    fDaughter = 443;
    break;
   case kBtoEle: 
   SetBeautyMothers();
    fDaughter = 11;
    break;
   case kBtoMuon: 
    SetBeautyMothers();
    fDaughter = 13;
    break;
   case kBtoD0: 
    SetBeautyMothers();
    fDaughter = 421;
    break;
   }
  }

//______________________________________________________________________________________
Double_t* AliPtMothFromPtDaugh::GetBinsSize(const TH1F *hist,Int_t &n) const
  {
   // return the binning size of the histogram hist
   // n return the number of bins 
   Double_t* edges = (Double_t*)hist->GetXaxis()->GetXbins()->GetArray();
   n=hist->GetNbinsX();
   return edges;  
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::GetEtaMothers(Double_t &etaMin, Double_t &etaMax) const
  {
   // method to get the bounds of the pseudorapidity range 
   // for mothers. Return kTRUE if pseudorapidity is used and put 
   // pseudorapidity edges in etaMin and etaMax   
 
  if(fUseEta == kFALSE){
        AliError("You are using RAPIDITY range! \n"); 
        etaMin = 0.;
        etaMax = 0.;
        return kFALSE;
        }  
   etaMin = fyMothMin;
   etaMax = fyMothMax;  
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::GetEtaDaughter(Double_t &etaMin, Double_t &etaMax) const
  {
   // method to get the bounds of the pseudorapidity range 
   // for daughters. Return kTRUE if pseudorapidity is used and put 
   // pseudorapidity edges in etaMin and etaMax   

   if(fUseEta == kFALSE){
        AliError("You are using RAPIDITY range! \n"); 
        etaMin = 0.;
        etaMax = 0.; 
        return kFALSE;
        }
   etaMin = fyDaughMin;
   etaMax = fyDaughMax;
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::GetYMothers(Double_t &yMin, Double_t &yMax) const
  {
   // method to get the bounds of the rapidity range 
   // for mothers. Return kTRUE if rapidity is used and put 
   // rapidity edges in yMin and yMax   

   if(fUseEta == kTRUE){
        AliError("You are using PSEUDORAPIDITY range! \n"); 
        yMin = 0.;
        yMax = 0.;
        return kFALSE;
        }
   yMin = fyMothMin;
   yMax = fyMothMax;
   return kTRUE;
  }
 
//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::GetYDaughter(Double_t &yMin, Double_t &yMax) const
  {
   // method to get the bounds of the rapidity range 
   // for daughters. Return kTRUE if rapidity is used and put 
   // rapidity edges in yMin and yMax 

   if(fUseEta == kTRUE){
        AliError("You are using PSEUDORAPIDITY range! \n"); 
        yMin = 0.;
        yMax = 0.;
        return kFALSE;
        }
   yMin = fyDaughMin;
   yMax = fyDaughMax;
   return kTRUE;
  }


//______________________________________________________________________________________
void AliPtMothFromPtDaugh::SetPdgDaugh(Int_t pdgD)
  { 
   // Set pdg code for daughter particle. Check 
   // that the kUserAnalysis is the selected Analysis_mode. 
   // If not print out a message of error. 
   switch(fAnalysisMode)
    {
    case kUserAnalysis:
     fDaughter = pdgD;
     break;
    case kBtoJPSI:
     if(pdgD!= 443) {AliError("Nothing done: first change AnalysisMode to kUserAnalysis");}
     else {fDaughter = pdgD;} 
     break;
    case kBtoEle: 
     if(pdgD!= 11) {AliError("Nothing done: first change AnalysisMode to kUserAnalysis");}
     else {fDaughter = pdgD;} 
     break;
    case kBtoMuon: 
     if(pdgD!= 13) {AliError("Nothing done: first change AnalysisMode to kUserAnalysis");}
     else {fDaughter = pdgD;} 
     break; 
    case kBtoD0: 
     if(pdgD!= 421) {AliError("Nothing done: first change AnalysisMode to kUserAnalysis");}
     else {fDaughter = pdgD;}
     break;
    }
   return;
  } 

//______________________________________________________________________________________
Int_t* AliPtMothFromPtDaugh::GetPdgMothers(Int_t &n_mothers) const
  {
   // return the pointer to the array of pdg codes of mothers particles
   // if it exist. Put its dimension in n_mothers   
 
   if(!fMothers) {AliError("Mothers pdg are not defined! \n"); return 0x0;} 
   n_mothers = fMothers->GetSize(); 
   return fMothers->GetArray();
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetW(Int_t i,Int_t j) const
  {
   // Return value of correction factors Wij at the position i (pt-mothers bin index)-
   // j (pt-daughter bin index). Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Wij don't exist or the indices i or j are out of the variability range return 0 
    
   if(!fWij) {AliError("Correction factors Wij are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMothers) {AliError("mothers pt-histogram doesn't exist!\n"); return 0;}
   if(!fHistoPtDaughter) {AliError("daughters pt-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   if(i<0 || i>nbinsM-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow ",0,nbinsM-1,nbinsM-1)); return 0;} 
  
   if(j<0 || j>nbinsD-1) {AliError(Form("Index j out of range: %d =< j =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsD-1,nbinsD-1)); return 0;}  
   return fWij[i][j];
  } 

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetStatErrW(Int_t i,Int_t j) const
  {
   // Return value of statistical error on correction factors Wij at the position 
   // i (pt-mothers bin index)- j (pt-daughter bin index). Bin 0 is the underflow, 
   // bin nbins+1 the overflow. If Wij don't exist or the indices i or j are out of the 
   // variability range return 0 
  
   if(!fHistoPtMothers) {AliError("mothers pt-histogram doesn't exist!\n"); return 0;}
   if(!fHistoPtDaughter) {AliError("daughters pt-histogram doesn't exist!\n"); return 0;}
 
   if(!fWij) {AliError("Correction factors Wij are not been evaluated yet!\n"); return 0;}
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   if(i<0 || i>nbinsM-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsM-1,nbinsM-1)); return 0;}
   if(j<0 || j>nbinsD-1) {AliError(Form("Index j out of range: %d =< j =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsD-1,nbinsD-1)); return 0;}
   return fWij[i+nbinsM][j];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetF(Int_t i) const
  {
   // Return value of correction factors Fi at the position i (pt-mothers bin index).
   // Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Fi don't exist or the index i is out of the variability range return 0     
  
   if(!fFi) {AliError("Correction factors Fi are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMothers) {AliError("mothers pt-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   if(i<0 || i>nbinsM-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsM-1,nbinsM-1)); return 0;}
   return fFi[i];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetStatErrF(Int_t i) const
  {
   // Return statistical error on correction factors Fi at the position i (pt-mothers bin index).
   // Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Fi don't exist or the index i is out of the variability range return 0     

   if(!fFi) {AliError("Correction factors Fi are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMothers) {AliError("mothers pt-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   if(i<0 || i>nbinsM-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsM-1,nbinsM-1)); return 0;}

   return fFi[i+nbinsM];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetWmin(Int_t i,Int_t j) const
  {
   // Return value of correction factors Wij_min at the position i (ptMin-mothers bin index)-
   // j (pt-daughter bin index). Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Wij_min don't exist or the indices i or j are out of the variability range return 0  
    
   if(!fWijMin) {AliError("Correction factors Wij_min are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMinMothers) {AliError("mothers ptMin-histogram doesn't exist!\n"); return 0;}
   if(!fHistoPtDaughter) {AliError("daughters pt-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   if(i<0 || i>nbinsMmin-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsMmin-1,nbinsMmin-1)); return 0;}
   if(j<0 || j>nbinsD-1) {AliError(Form("Index j out of range: %d =< j =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsD-1,nbinsD-1)); return 0;}
   return fWijMin[i][j];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetStatErrWmin(Int_t i,Int_t j) const
  {
   // Return value of statistical error on correction factors Wij_min at the position 
   // i (ptMin-mothers bin index)- j (pt-daughter bin index). Bin 0 is the underflow, 
   // bin nbins+1 the overflow. If Wij_min don't exist or the indices i or j are out of the 
   // variability range return 0 

   if(!fWijMin) {AliError("Correction factors Wij_min are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMinMothers) {AliError("mothers ptMin-histogram doesn't exist!\n"); return 0;}
   if(!fHistoPtDaughter) {AliError("daughters pt-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   if(i<0 || i>nbinsMmin-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsMmin-1,nbinsMmin-1)); return 0;}
   if(j<0 || j>nbinsD-1) {AliError(Form("Index j out of range: %d =< j =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsD-1,nbinsD-1)); return 0;}
 
   return fWijMin[i+nbinsMmin][j];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetFmin(Int_t i) const
  {
   // Return value of correction factors Fi_min at the position i (ptMin-mothers bin index).
   // Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Fi_min don't exist or the index i is out of the variability range return 0     

   if(!fFiMin) {AliError("Correction factors Fi_min are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMinMothers) {AliError("mothers ptMin-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsMmin=fHistoPtMothers->GetNbinsX()+2;
   if(i<0 || i>nbinsMmin-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsMmin-1,nbinsMmin-1)); return 0;}
   return fFiMin[i];
  }

//______________________________________________________________________________________
Double_t AliPtMothFromPtDaugh::GetStatErrFmin(Int_t i) const
  {
   // Return statistical error on correction factors Fi_min at the position i (ptMin-mothers bin index).
   // Bin 0 is the underflow, bin nbins+1 the overflow. 
   // If Fi_min don't exist or the index i is out of the variability range return 0     

   if(!fFiMin) {AliError("Correction factors Fi_min are not been evaluated yet!\n"); return 0;}
   if(!fHistoPtMinMothers) {AliError("mothers ptMin-histogram doesn't exist!\n"); return 0;}
   Int_t nbinsMmin=fHistoPtMothers->GetNbinsX()+2;
   if(i<0 || i>nbinsMmin-1) {AliError(Form("Index i out of range: %d =< i =< %d. Bin 0 is the underflow - bin %d is the overflow.",0,nbinsMmin-1,nbinsMmin-1)); return 0;}

   return fFiMin[i+nbinsMmin];
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::IsMothers(Int_t pdgCode)
  {
   //return kTRUE if pdgCode is in the list of pdg codes of mothers
   Int_t dim = fMothers->GetSize();
   for(Int_t i=0;i<dim;i++) { 
    Int_t pdgMother = (Int_t)fMothers->GetAt(i); 
    if(pdgCode == pdgMother) return kTRUE; }
   return kFALSE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::IsSelectedDaugh(const TParticle *part, Int_t &labelDaugh, AliStack *const stack)
  {
   // return kTRUE if particle part has the selected daughter
   // if yes put the label of the track in labelDaugh  
   TParticle *daugh;
   Int_t nDg=part->GetNDaughters();
   if(nDg<=0) {AliError("I have no daugh!\n"); return kFALSE;}
   for(Int_t i=part->GetFirstDaughter();i<part->GetLastDaughter()+1;i++)
     {
      daugh=stack->Particle(i);
      if(TMath::Abs(daugh->GetPdgCode())==fDaughter) {
        labelDaugh=i;
        return kTRUE;
      }
     }
   return kFALSE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::Rapidity(const TParticle *particle, Double_t &y)
  {
   // Evaluated rapidity of particle  and put it in y. Return kFALSE if
   // cannot compute rapidity 
   y=-999;
   if(particle->Energy()-particle->Pz()<=0) return kFALSE;
   y=0.5*TMath::Log((particle->Energy()+particle->Pz())/(particle->Energy()-particle->Pz()));
   return kTRUE;
  }

//______________________________________________________________________________________
Int_t AliPtMothFromPtDaugh::GiveBinIndex(Double_t Ptpart,const TH1F *ptHist) const
  {
   // Return the bin index of pt respect to the binning size of ptHist
   // bin 0 is the underflow - nbins+1 is the overflow  
   Int_t nbins=ptHist->GetNbinsX();
   Int_t index=0;
   for(Int_t i=1;i<nbins+2;i++)
     {
      if(Ptpart<(ptHist->GetBinLowEdge(i)))
      {
      index=i-1;  
      break;
      }
     }
   return index;
  }

//______________________________________________________________________________________
Bool_t  AliPtMothFromPtDaugh::CutDaugh(Double_t yD,Double_t ptD)
  {
   // put a control for rapidity yD and transverse momentum ptD of daughter
   // return kTRUE if  fyDaughMin < yD < fyDaughMax and ptMinDaugh < ptD < ptMaxDaugh   
   Double_t ptMinDaugh = fHistoPtDaughter->GetBinLowEdge(1);
   Double_t ptMaxDaugh = fHistoPtDaughter->GetBinLowEdge(fHistoPtDaughter->GetNbinsX());
   if( yD < fyDaughMin || yD > fyDaughMax ) return kFALSE;
   if( ptD > ptMaxDaugh || ptD < ptMinDaugh ) return kFALSE;
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluateWij()
  {
   // Evaluate correction factors using to extract the ptRaw and
   // ptMinRaw distributions. Statistical errors on those are computed too
  
   if(!fHistoPtMothers || !fHistoPtDaughter || !fHistoPtMinMothers) 
   {AliError("Control mother and daughter histograms!\n"); return kFALSE;}

   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Float_t pdgM, pdgD, pxM,pyM,pzM,yM,pxD,pyD,pzD,yD,etaM,etaD,cutVarM,cutVarD;
   Double_t* entries=new Double_t[nbinsD];
   for(Int_t ii=0;ii<nbinsD;ii++) {entries[ii]=0.;}
   Int_t i,j,iMin;
   if(!fDecayKine) 
   {
    AliError("TNtupla is not defined!\n"); 
    delete [] entries;
    return kFALSE;
   }
   fDecayKine->SetBranchAddress("pdgM",&pdgM);
   fDecayKine->SetBranchAddress("pxM",&pxM);
   fDecayKine->SetBranchAddress("pyM",&pyM);
   fDecayKine->SetBranchAddress("pzM",&pzM);
   fDecayKine->SetBranchAddress("yM",&yM);
   fDecayKine->SetBranchAddress("etaM",&etaM);
   fDecayKine->SetBranchAddress("pdgD",&pdgD);
   fDecayKine->SetBranchAddress("pxD",&pxD);
   fDecayKine->SetBranchAddress("pyD",&pyD);
   fDecayKine->SetBranchAddress("pzD",&pzD);
   fDecayKine->SetBranchAddress("yD",&yD);
   fDecayKine->SetBranchAddress("etaD",&etaD);
   Double_t ptD,ptM;
   // Initialize correction factors for pT and pTmin if those exist
   if(!fWij)
    {
     AliError("Correction factors Wij have not been created!\n"); 
     delete [] entries;
     return kFALSE;
    }

   if(!fWijMin)
    {
     AliError("Correction factors Wij_min have not been created!\n"); 
     delete [] entries;
     return kFALSE;
    }
   for(Int_t ii=0;ii<(2*nbinsM);ii++){
     for(Int_t jj=0;jj<nbinsD;jj++){
        fWij[ii][jj]=0;
        }
      }
   for(Int_t ii=0;ii<(2*nbinsMmin);ii++){
     for(Int_t jj=0;jj<nbinsD;jj++){
       fWijMin[ii][jj]=0;
      }
    }
   Int_t nentries = (Int_t)fDecayKine->GetEntries();
   Int_t fNcurrent=0;
   Int_t nb = (Int_t)fDecayKine->GetEvent(fNcurrent);
   for (Int_t iev=0; iev<nentries; iev++){
   // check if rapidity or pseudorapidity range is set
   if(fUseEta == kFALSE) {cutVarD = yD; cutVarM = yM;}
   else {cutVarD = etaD; cutVarM = etaM;}
   ptD=TMath::Sqrt(pxD*pxD+pyD*pyD);
   ptM=TMath::Sqrt(pxM*pxM+pyM*pyM);
   pdgM = TMath::Abs(pdgM);
   if((TMath::Abs(pdgD))==fDaughter && IsMothers((Int_t)pdgM))
      {
       j=GiveBinIndex(ptD,fHistoPtDaughter);
       i=GiveBinIndex(ptM,fHistoPtMothers);
       iMin=GiveBinIndex(ptM,fHistoPtMinMothers);
       if(!CutDaugh(cutVarD,ptD)) { fNcurrent++; nb = (Int_t)fDecayKine->GetEvent(fNcurrent); continue;}
         if(cutVarM>fyMothMin && cutVarM<fyMothMax){
         fWij[i][j]+=1.;
         for(Int_t k=0;k<iMin+1;k++) {fWijMin[k][j]+=1.;}
          }
         entries[j]++;
         }
    fNcurrent++;
    nb = (Int_t)fDecayKine->GetEvent(fNcurrent);
   }
  for(Int_t jj=0;jj<nbinsD;jj++){
      for(Int_t ii=0;ii<nbinsM;ii++){
         if(entries[jj]>0){
           fWij[ii][jj]=fWij[ii][jj]/entries[jj];
           // evaluate statistical errors on fWij
           fWij[ii+nbinsM][jj]=(TMath::Sqrt(fWij[ii][jj]*(1-(fWij[ii][jj]/entries[jj]))))/entries[jj];
           }
         else{
           // if there are no entries in the bin-j of daughter distribution 
           // set factor = -1 and error = 999
           fWij[ii][jj]=-1;
           fWij[ii+nbinsM][jj]=999;
           }
         }
      for(Int_t ii=0;ii<nbinsMmin;ii++){
         if(entries[jj]>0){
            fWijMin[ii][jj]=fWijMin[ii][jj]/entries[jj];
            //evaluate statistical errors on fWijMin
            fWijMin[ii+nbinsMmin][jj] = (TMath::Sqrt(fWijMin[ii][jj]*(1-(fWijMin[ii][jj]/entries[jj]))))/entries[jj]; 
            }
             else{
            //if there are no entries set correction factor = -1 and error = -999
            fWijMin[ii][jj]=-1;
            fWijMin[ii+nbinsMmin][jj]=999;
           }
        }
    }
   delete [] entries;
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluateFi()
  {  
   // Evaluate acceptance correction factors that are applied on the 
   // raw distributions. Statistical errors on those are computed too  

   if(!fHistoPtMothers || !fHistoPtMinMothers)
   {AliError("Control mother histograms!\n"); return kFALSE;}

   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Double_t* entries=new Double_t[nbinsM];
   Double_t* entries1=new Double_t[nbinsMmin];
   if(!fFi)
     {
      AliError("Correction factors Fi have not been created!\n"); 
      delete [] entries;
      delete [] entries1;
      return kFALSE;
     }

   if(!fFiMin)
     {
      AliError("Correction factors Fi_min have not been created!\n"); 
      delete [] entries;
      delete [] entries1;
      return kFALSE;
     }

   //initialize the correction factor for pT and pTmin
   for(Int_t ii=0;ii<nbinsM;ii++){
       fFi[ii]=0.; //for correction factors
       fFi[ii+nbinsM]=0.; //for statistical error on correction factors
       entries[ii]=0.;
      }
   for(Int_t ii=0;ii<nbinsMmin;ii++){
       entries1[ii]=0.;
       fFiMin[ii]=0.; //for correction factors
       fFiMin[ii+nbinsMmin]=0.; //for statistical error on correction factors       
      }
   Float_t pdgM, pdgD, pxM,pyM,pzM,yM,pxD,pyD,pzD,yD,etaM,etaD,cutVarD,cutVarM;
   Int_t i,iMin;
   if(!fDecayKine) {
     AliError("TNtupla is not defined!\n"); 
     delete [] entries;
     delete [] entries1;
     return kFALSE;
   } 
   fDecayKine->SetBranchAddress("pdgM",&pdgM);
   fDecayKine->SetBranchAddress("pxM",&pxM);
   fDecayKine->SetBranchAddress("pyM",&pyM);
   fDecayKine->SetBranchAddress("pzM",&pzM);
   fDecayKine->SetBranchAddress("yM",&yM);
   fDecayKine->SetBranchAddress("etaM",&etaM);
   fDecayKine->SetBranchAddress("pdgD",&pdgD);
   fDecayKine->SetBranchAddress("pxD",&pxD);
   fDecayKine->SetBranchAddress("pyD",&pyD);
   fDecayKine->SetBranchAddress("pzD",&pzD);
   fDecayKine->SetBranchAddress("yD",&yD);
   fDecayKine->SetBranchAddress("etaD",&etaD);
   Double_t ptD,ptM;

   Int_t nentries = (Int_t)fDecayKine->GetEntries();
   Int_t fNcurrent=0;
   Int_t nb = (Int_t)fDecayKine->GetEvent(fNcurrent);

   for (Int_t iev=0; iev<nentries; iev++){
    pdgM = TMath::Abs(pdgM);
    if((TMath::Abs(pdgD))==fDaughter && IsMothers((Int_t)pdgM))
     {
     //check if rapidity or pseudorapidity range is set
     if(fUseEta == kFALSE) {cutVarD = yD; cutVarM = yM;}
     else {cutVarD = etaD; cutVarM = etaM;}
     ptD=TMath::Sqrt(pxD*pxD+pyD*pyD);
     ptM=TMath::Sqrt(pxM*pxM+pyM*pyM);
     i=GiveBinIndex(ptM,fHistoPtMothers);
     iMin=GiveBinIndex(ptM,fHistoPtMinMothers);
        if(cutVarM<fyMothMin || cutVarM>fyMothMax){ 
        fNcurrent++; nb = (Int_t)fDecayKine->GetEvent(fNcurrent); continue;}
     entries[i]++;
     for(Int_t k=0; k<iMin+1;k++) {entries1[k]+=1;}
     if(!CutDaugh(cutVarD,ptD)) {fNcurrent++; nb = (Int_t)fDecayKine->GetEvent(fNcurrent); continue;}
     fFi[i]+=1.;
     for(Int_t k=0; k<iMin+1;k++) {fFiMin[k]+=1.;}
     }
    fNcurrent++;
    nb = (Int_t) fDecayKine->GetEvent(fNcurrent);
    }

  for(Int_t ii=0;ii<nbinsM;ii++){
       if(entries[ii]>0){
         fFi[ii]/=entries[ii];
         fFi[ii+nbinsM]=(TMath::Sqrt(fFi[ii]*(1-(fFi[ii]/entries[ii]))))/entries[ii];
        }
      else{
         fFi[ii]=-1;
         fFi[ii+nbinsM]=999;
        }
      }
  for(Int_t ii=0;ii<nbinsMmin;ii++){
     if(entries1[ii]>0){
        fFiMin[ii]/=entries1[ii];
        fFiMin[ii+nbinsMmin]=(TMath::Sqrt(fFiMin[ii]*(1-(fFiMin[ii]/entries1[ii]))))/entries1[ii];
       }
      else {fFiMin[ii]=-1; fFiMin[ii+nbinsMmin]=999;}
     }
   delete [] entries;
   delete [] entries1;
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluatePtMothRaw(TH1F *histoPt, TH1F *histoPtMin)
  {
   //Apply the fWij and fWijMin on the daughter distribution
   //in order to evaluate the pt and ptMin raw distributions for mothers

   if(!fHistoPtMothers || !fHistoPtDaughter || !fHistoPtMinMothers)
   {AliError("Control mother and daughter histograms!\n"); return kFALSE;}

   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Double_t lowedge=0.;
   Double_t wfill=0.;
   Double_t wMinfill=0.;
   for(Int_t j=0;j<nbinsD;j++){
     for(Int_t i=0;i<nbinsM;i++){
         lowedge=fHistoPtMothers->GetBinCenter(i);
         if(fWij[i][j]>=0) wfill=fWij[i][j]*(fHistoPtDaughter->GetBinContent(j));
         histoPt->Fill(lowedge,wfill);
         }
     for(Int_t i=0;i<nbinsMmin;i++){
         lowedge=fHistoPtMinMothers->GetBinLowEdge(i);
         if(fWijMin[i][j]>=0) wMinfill=fWijMin[i][j]*(fHistoPtDaughter->GetBinContent(j));
         histoPtMin->Fill(lowedge,wMinfill);
         }
      }
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluateErrPt(Double_t *erStat)
  {
   // Evaluate the statistical error on the pt-mothers distribution.
   // sigmaX: contribution that came from the measured distibution
   // sigmaWij: contribution that came from the fWij factors
   // sigmaFi: contribution that came from the fFi factors 

   if(!fHistoPtMothers || !fHistoPtDaughter)
   {AliError("Control mother(pt) and daughter histograms!\n"); return kFALSE;}
 
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   Double_t m=0;
   Double_t  sigmaX, sigmaWij, sigmaFi;
   for(Int_t i=0;i<nbinsM;i++){
       sigmaX=0.;
       sigmaWij=0.;
       sigmaFi=0;
       for(Int_t j=0;j<nbinsD;j++){
          if(fWij[i][j]>=0){
          sigmaX+=(((fWij[i][j])*(fWij[i][j]))*((fHistoPtDaughter->GetBinError(j))*(fHistoPtDaughter->GetBinError(j))));
          sigmaWij+=((fHistoPtDaughter->GetBinContent(j))*(fHistoPtDaughter->GetBinContent(j)))*(fWij[i+nbinsM][j]*fWij[i+nbinsM][j]);
          sigmaFi+=(fWij[i][j])*(fHistoPtDaughter->GetBinContent(j));
          }
       } 
      if(fFi[i]>0) sigmaFi=((sigmaFi*sigmaFi)*(fFi[i+nbinsM]*fFi[i+nbinsM]))/(fFi[i]*fFi[i]);
      m=TMath::Sqrt(sigmaX+sigmaWij+sigmaFi);
      if(fFi[i]>0) erStat[i]=(1/fFi[i])*m;
      else erStat[i]=999;  
    }
   return kTRUE;
  }

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluateErrPtMin(Double_t *erStat)
  {
   // Evaluate statistical error on ptMin mothers distribution
   // sigmaMinX: contribution that came from the measured distibution
   // sigmaMinWij: contribution that came from the fWijMin factors
   // sigmaMinFi: contribution that came from the fFiMin factors  
   
   if(!fHistoPtDaughter || !fHistoPtMinMothers)
   {AliError("Control mother(ptMin) and daughter histograms!\n"); return kFALSE;}

   Int_t nbinsD=fHistoPtDaughter->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Double_t m1=0;
   Double_t sigmaMinX;
   Double_t sigmaMinWij;
   Double_t sigmaMinFi;
      for(Int_t i=0;i<nbinsMmin;i++){
        sigmaMinX=0.;
        sigmaMinWij=0.;
        sigmaMinFi=0.;
           for(Int_t j=0;j<nbinsD;j++){
            if(fWijMin[i][j]>=0){
            sigmaMinX+=(((fWijMin[i][j])*(fWijMin[i][j]))*((fHistoPtDaughter->GetBinError(j))*(fHistoPtDaughter->GetBinError(j))));
            sigmaMinWij+=((fHistoPtDaughter->GetBinContent(j))*(fHistoPtDaughter->GetBinContent(j)))*(fWijMin[i+nbinsMmin][j]*fWijMin[i+nbinsMmin][j]);
            sigmaMinFi+=(fWijMin[i][j])*(fHistoPtDaughter->GetBinContent(j));
           }
         }
    if(fFiMin[i]>0) sigmaMinFi=((sigmaMinFi*sigmaMinFi)*(fFiMin[i+nbinsMmin]*fFiMin[i+nbinsMmin]))/(fFiMin[i]*fFiMin[i]);
    m1=TMath::Sqrt(sigmaMinX+sigmaMinWij+sigmaMinFi);
    if(fFiMin[i]>0) erStat[i]=(1/fFiMin[i])*m1;
    else erStat[i]=999;
    }

   return kTRUE;
  } 

//______________________________________________________________________________________
Bool_t AliPtMothFromPtDaugh::EvaluatePtMoth()
  {
   // Evaluate pt and ptMin distribution for mothers
   // First evaluate the sigma raw distribution by calling EvaluatePtMothRaw
   // then evaluate the pt and ptMin mothers distribution.  
   // Statistical errors on those distributions are evaluated too.  

   if(!EvaluateWij()) return kFALSE;
   if(!EvaluateFi())  return kFALSE;

   // reset pt and ptMin mothers histograms
   fHistoPtMothers->Reset();
   fHistoPtMinMothers->Reset();

   TH1F *histoPt=(TH1F*)fHistoPtMothers->Clone();
   TH1F *histoPtMin=(TH1F*)fHistoPtMinMothers->Clone();
   EvaluatePtMothRaw(histoPt,histoPtMin);
   Int_t nbinsM=fHistoPtMothers->GetNbinsX()+2;
   Int_t nbinsMmin=fHistoPtMinMothers->GetNbinsX()+2;
   Double_t *erPtStat=new Double_t[nbinsM+2];
   EvaluateErrPt(erPtStat);
   Double_t *erPtMinStat=new Double_t[nbinsMmin+2];
   EvaluateErrPtMin(erPtMinStat); 
   Double_t lowedge=0;
   Double_t fwfill;
   Double_t fMinfill;
   for(Int_t i=0;i<nbinsM;i++){
      fwfill=0.;
      lowedge=fHistoPtMothers->GetBinCenter(i);
       if(fFi[i]>0){
        fwfill=(histoPt->GetBinContent(i))/fFi[i];
        fHistoPtMothers->Fill(lowedge,fwfill);
        fHistoPtMothers->SetBinError(i,erPtStat[i]);
       }
      }
   for(Int_t i=0;i<nbinsMmin;i++){
      fMinfill=0.;
      lowedge=fHistoPtMinMothers->GetBinCenter(i);
      if(fFiMin[i]>0){
         fMinfill=(histoPtMin->GetBinContent(i))/fFiMin[i];
         fHistoPtMinMothers->Fill(lowedge,fMinfill);
         fHistoPtMinMothers->SetBinError(i,erPtMinStat[i]);
        }
      }
   delete [] erPtStat;
   delete [] erPtMinStat;  
   return kTRUE;
  }

//______________________________________________________________________________________
void AliPtMothFromPtDaugh::WritePtMothHistoToFile(TString fileOutName)
  {
   // Write pt and ptMin histograms of mothers in a file 
   // with name fileOutName. Default name is "Mothers.root".
    AliError(Form("Write mothers histograms in the file %s \n",fileOutName.Data()));
   if(!fHistoPtMothers) {AliError("Cannot write pt-mothers histogram! It doesn't exists!"); return;}
   if(!fHistoPtMinMothers)  { AliError("Cannot write ptMin-mothers histogram! It doesn't exists!"); return;} 
   TFile *outFile = TFile::Open(fileOutName.Data(),"RECREATE");
   outFile->cd();
   fHistoPtMothers->Write();
   fHistoPtMinMothers->Write();
   outFile->Close();
   return;  
  }
 AliPtMothFromPtDaugh.cxx:1
 AliPtMothFromPtDaugh.cxx:2
 AliPtMothFromPtDaugh.cxx:3
 AliPtMothFromPtDaugh.cxx:4
 AliPtMothFromPtDaugh.cxx:5
 AliPtMothFromPtDaugh.cxx:6
 AliPtMothFromPtDaugh.cxx:7
 AliPtMothFromPtDaugh.cxx:8
 AliPtMothFromPtDaugh.cxx:9
 AliPtMothFromPtDaugh.cxx:10
 AliPtMothFromPtDaugh.cxx:11
 AliPtMothFromPtDaugh.cxx:12
 AliPtMothFromPtDaugh.cxx:13
 AliPtMothFromPtDaugh.cxx:14
 AliPtMothFromPtDaugh.cxx:15
 AliPtMothFromPtDaugh.cxx:16
 AliPtMothFromPtDaugh.cxx:17
 AliPtMothFromPtDaugh.cxx:18
 AliPtMothFromPtDaugh.cxx:19
 AliPtMothFromPtDaugh.cxx:20
 AliPtMothFromPtDaugh.cxx:21
 AliPtMothFromPtDaugh.cxx:22
 AliPtMothFromPtDaugh.cxx:23
 AliPtMothFromPtDaugh.cxx:24
 AliPtMothFromPtDaugh.cxx:25
 AliPtMothFromPtDaugh.cxx:26
 AliPtMothFromPtDaugh.cxx:27
 AliPtMothFromPtDaugh.cxx:28
 AliPtMothFromPtDaugh.cxx:29
 AliPtMothFromPtDaugh.cxx:30
 AliPtMothFromPtDaugh.cxx:31
 AliPtMothFromPtDaugh.cxx:32
 AliPtMothFromPtDaugh.cxx:33
 AliPtMothFromPtDaugh.cxx:34
 AliPtMothFromPtDaugh.cxx:35
 AliPtMothFromPtDaugh.cxx:36
 AliPtMothFromPtDaugh.cxx:37
 AliPtMothFromPtDaugh.cxx:38
 AliPtMothFromPtDaugh.cxx:39
 AliPtMothFromPtDaugh.cxx:40
 AliPtMothFromPtDaugh.cxx:41
 AliPtMothFromPtDaugh.cxx:42
 AliPtMothFromPtDaugh.cxx:43
 AliPtMothFromPtDaugh.cxx:44
 AliPtMothFromPtDaugh.cxx:45
 AliPtMothFromPtDaugh.cxx:46
 AliPtMothFromPtDaugh.cxx:47
 AliPtMothFromPtDaugh.cxx:48
 AliPtMothFromPtDaugh.cxx:49
 AliPtMothFromPtDaugh.cxx:50
 AliPtMothFromPtDaugh.cxx:51
 AliPtMothFromPtDaugh.cxx:52
 AliPtMothFromPtDaugh.cxx:53
 AliPtMothFromPtDaugh.cxx:54
 AliPtMothFromPtDaugh.cxx:55
 AliPtMothFromPtDaugh.cxx:56
 AliPtMothFromPtDaugh.cxx:57
 AliPtMothFromPtDaugh.cxx:58
 AliPtMothFromPtDaugh.cxx:59
 AliPtMothFromPtDaugh.cxx:60
 AliPtMothFromPtDaugh.cxx:61
 AliPtMothFromPtDaugh.cxx:62
 AliPtMothFromPtDaugh.cxx:63
 AliPtMothFromPtDaugh.cxx:64
 AliPtMothFromPtDaugh.cxx:65
 AliPtMothFromPtDaugh.cxx:66
 AliPtMothFromPtDaugh.cxx:67
 AliPtMothFromPtDaugh.cxx:68
 AliPtMothFromPtDaugh.cxx:69
 AliPtMothFromPtDaugh.cxx:70
 AliPtMothFromPtDaugh.cxx:71
 AliPtMothFromPtDaugh.cxx:72
 AliPtMothFromPtDaugh.cxx:73
 AliPtMothFromPtDaugh.cxx:74
 AliPtMothFromPtDaugh.cxx:75
 AliPtMothFromPtDaugh.cxx:76
 AliPtMothFromPtDaugh.cxx:77
 AliPtMothFromPtDaugh.cxx:78
 AliPtMothFromPtDaugh.cxx:79
 AliPtMothFromPtDaugh.cxx:80
 AliPtMothFromPtDaugh.cxx:81
 AliPtMothFromPtDaugh.cxx:82
 AliPtMothFromPtDaugh.cxx:83
 AliPtMothFromPtDaugh.cxx:84
 AliPtMothFromPtDaugh.cxx:85
 AliPtMothFromPtDaugh.cxx:86
 AliPtMothFromPtDaugh.cxx:87
 AliPtMothFromPtDaugh.cxx:88
 AliPtMothFromPtDaugh.cxx:89
 AliPtMothFromPtDaugh.cxx:90
 AliPtMothFromPtDaugh.cxx:91
 AliPtMothFromPtDaugh.cxx:92
 AliPtMothFromPtDaugh.cxx:93
 AliPtMothFromPtDaugh.cxx:94
 AliPtMothFromPtDaugh.cxx:95
 AliPtMothFromPtDaugh.cxx:96
 AliPtMothFromPtDaugh.cxx:97
 AliPtMothFromPtDaugh.cxx:98
 AliPtMothFromPtDaugh.cxx:99
 AliPtMothFromPtDaugh.cxx:100
 AliPtMothFromPtDaugh.cxx:101
 AliPtMothFromPtDaugh.cxx:102
 AliPtMothFromPtDaugh.cxx:103
 AliPtMothFromPtDaugh.cxx:104
 AliPtMothFromPtDaugh.cxx:105
 AliPtMothFromPtDaugh.cxx:106
 AliPtMothFromPtDaugh.cxx:107
 AliPtMothFromPtDaugh.cxx:108
 AliPtMothFromPtDaugh.cxx:109
 AliPtMothFromPtDaugh.cxx:110
 AliPtMothFromPtDaugh.cxx:111
 AliPtMothFromPtDaugh.cxx:112
 AliPtMothFromPtDaugh.cxx:113
 AliPtMothFromPtDaugh.cxx:114
 AliPtMothFromPtDaugh.cxx:115
 AliPtMothFromPtDaugh.cxx:116
 AliPtMothFromPtDaugh.cxx:117
 AliPtMothFromPtDaugh.cxx:118
 AliPtMothFromPtDaugh.cxx:119
 AliPtMothFromPtDaugh.cxx:120
 AliPtMothFromPtDaugh.cxx:121
 AliPtMothFromPtDaugh.cxx:122
 AliPtMothFromPtDaugh.cxx:123
 AliPtMothFromPtDaugh.cxx:124
 AliPtMothFromPtDaugh.cxx:125
 AliPtMothFromPtDaugh.cxx:126
 AliPtMothFromPtDaugh.cxx:127
 AliPtMothFromPtDaugh.cxx:128
 AliPtMothFromPtDaugh.cxx:129
 AliPtMothFromPtDaugh.cxx:130
 AliPtMothFromPtDaugh.cxx:131
 AliPtMothFromPtDaugh.cxx:132
 AliPtMothFromPtDaugh.cxx:133
 AliPtMothFromPtDaugh.cxx:134
 AliPtMothFromPtDaugh.cxx:135
 AliPtMothFromPtDaugh.cxx:136
 AliPtMothFromPtDaugh.cxx:137
 AliPtMothFromPtDaugh.cxx:138
 AliPtMothFromPtDaugh.cxx:139
 AliPtMothFromPtDaugh.cxx:140
 AliPtMothFromPtDaugh.cxx:141
 AliPtMothFromPtDaugh.cxx:142
 AliPtMothFromPtDaugh.cxx:143
 AliPtMothFromPtDaugh.cxx:144
 AliPtMothFromPtDaugh.cxx:145
 AliPtMothFromPtDaugh.cxx:146
 AliPtMothFromPtDaugh.cxx:147
 AliPtMothFromPtDaugh.cxx:148
 AliPtMothFromPtDaugh.cxx:149
 AliPtMothFromPtDaugh.cxx:150
 AliPtMothFromPtDaugh.cxx:151
 AliPtMothFromPtDaugh.cxx:152
 AliPtMothFromPtDaugh.cxx:153
 AliPtMothFromPtDaugh.cxx:154
 AliPtMothFromPtDaugh.cxx:155
 AliPtMothFromPtDaugh.cxx:156
 AliPtMothFromPtDaugh.cxx:157
 AliPtMothFromPtDaugh.cxx:158
 AliPtMothFromPtDaugh.cxx:159
 AliPtMothFromPtDaugh.cxx:160
 AliPtMothFromPtDaugh.cxx:161
 AliPtMothFromPtDaugh.cxx:162
 AliPtMothFromPtDaugh.cxx:163
 AliPtMothFromPtDaugh.cxx:164
 AliPtMothFromPtDaugh.cxx:165
 AliPtMothFromPtDaugh.cxx:166
 AliPtMothFromPtDaugh.cxx:167
 AliPtMothFromPtDaugh.cxx:168
 AliPtMothFromPtDaugh.cxx:169
 AliPtMothFromPtDaugh.cxx:170
 AliPtMothFromPtDaugh.cxx:171
 AliPtMothFromPtDaugh.cxx:172
 AliPtMothFromPtDaugh.cxx:173
 AliPtMothFromPtDaugh.cxx:174
 AliPtMothFromPtDaugh.cxx:175
 AliPtMothFromPtDaugh.cxx:176
 AliPtMothFromPtDaugh.cxx:177
 AliPtMothFromPtDaugh.cxx:178
 AliPtMothFromPtDaugh.cxx:179
 AliPtMothFromPtDaugh.cxx:180
 AliPtMothFromPtDaugh.cxx:181
 AliPtMothFromPtDaugh.cxx:182
 AliPtMothFromPtDaugh.cxx:183
 AliPtMothFromPtDaugh.cxx:184
 AliPtMothFromPtDaugh.cxx:185
 AliPtMothFromPtDaugh.cxx:186
 AliPtMothFromPtDaugh.cxx:187
 AliPtMothFromPtDaugh.cxx:188
 AliPtMothFromPtDaugh.cxx:189
 AliPtMothFromPtDaugh.cxx:190
 AliPtMothFromPtDaugh.cxx:191
 AliPtMothFromPtDaugh.cxx:192
 AliPtMothFromPtDaugh.cxx:193
 AliPtMothFromPtDaugh.cxx:194
 AliPtMothFromPtDaugh.cxx:195
 AliPtMothFromPtDaugh.cxx:196
 AliPtMothFromPtDaugh.cxx:197
 AliPtMothFromPtDaugh.cxx:198
 AliPtMothFromPtDaugh.cxx:199
 AliPtMothFromPtDaugh.cxx:200
 AliPtMothFromPtDaugh.cxx:201
 AliPtMothFromPtDaugh.cxx:202
 AliPtMothFromPtDaugh.cxx:203
 AliPtMothFromPtDaugh.cxx:204
 AliPtMothFromPtDaugh.cxx:205
 AliPtMothFromPtDaugh.cxx:206
 AliPtMothFromPtDaugh.cxx:207
 AliPtMothFromPtDaugh.cxx:208
 AliPtMothFromPtDaugh.cxx:209
 AliPtMothFromPtDaugh.cxx:210
 AliPtMothFromPtDaugh.cxx:211
 AliPtMothFromPtDaugh.cxx:212
 AliPtMothFromPtDaugh.cxx:213
 AliPtMothFromPtDaugh.cxx:214
 AliPtMothFromPtDaugh.cxx:215
 AliPtMothFromPtDaugh.cxx:216
 AliPtMothFromPtDaugh.cxx:217
 AliPtMothFromPtDaugh.cxx:218
 AliPtMothFromPtDaugh.cxx:219
 AliPtMothFromPtDaugh.cxx:220
 AliPtMothFromPtDaugh.cxx:221
 AliPtMothFromPtDaugh.cxx:222
 AliPtMothFromPtDaugh.cxx:223
 AliPtMothFromPtDaugh.cxx:224
 AliPtMothFromPtDaugh.cxx:225
 AliPtMothFromPtDaugh.cxx:226
 AliPtMothFromPtDaugh.cxx:227
 AliPtMothFromPtDaugh.cxx:228
 AliPtMothFromPtDaugh.cxx:229
 AliPtMothFromPtDaugh.cxx:230
 AliPtMothFromPtDaugh.cxx:231
 AliPtMothFromPtDaugh.cxx:232
 AliPtMothFromPtDaugh.cxx:233
 AliPtMothFromPtDaugh.cxx:234
 AliPtMothFromPtDaugh.cxx:235
 AliPtMothFromPtDaugh.cxx:236
 AliPtMothFromPtDaugh.cxx:237
 AliPtMothFromPtDaugh.cxx:238
 AliPtMothFromPtDaugh.cxx:239
 AliPtMothFromPtDaugh.cxx:240
 AliPtMothFromPtDaugh.cxx:241
 AliPtMothFromPtDaugh.cxx:242
 AliPtMothFromPtDaugh.cxx:243
 AliPtMothFromPtDaugh.cxx:244
 AliPtMothFromPtDaugh.cxx:245
 AliPtMothFromPtDaugh.cxx:246
 AliPtMothFromPtDaugh.cxx:247
 AliPtMothFromPtDaugh.cxx:248
 AliPtMothFromPtDaugh.cxx:249
 AliPtMothFromPtDaugh.cxx:250
 AliPtMothFromPtDaugh.cxx:251
 AliPtMothFromPtDaugh.cxx:252
 AliPtMothFromPtDaugh.cxx:253
 AliPtMothFromPtDaugh.cxx:254
 AliPtMothFromPtDaugh.cxx:255
 AliPtMothFromPtDaugh.cxx:256
 AliPtMothFromPtDaugh.cxx:257
 AliPtMothFromPtDaugh.cxx:258
 AliPtMothFromPtDaugh.cxx:259
 AliPtMothFromPtDaugh.cxx:260
 AliPtMothFromPtDaugh.cxx:261
 AliPtMothFromPtDaugh.cxx:262
 AliPtMothFromPtDaugh.cxx:263
 AliPtMothFromPtDaugh.cxx:264
 AliPtMothFromPtDaugh.cxx:265
 AliPtMothFromPtDaugh.cxx:266
 AliPtMothFromPtDaugh.cxx:267
 AliPtMothFromPtDaugh.cxx:268
 AliPtMothFromPtDaugh.cxx:269
 AliPtMothFromPtDaugh.cxx:270
 AliPtMothFromPtDaugh.cxx:271
 AliPtMothFromPtDaugh.cxx:272
 AliPtMothFromPtDaugh.cxx:273
 AliPtMothFromPtDaugh.cxx:274
 AliPtMothFromPtDaugh.cxx:275
 AliPtMothFromPtDaugh.cxx:276
 AliPtMothFromPtDaugh.cxx:277
 AliPtMothFromPtDaugh.cxx:278
 AliPtMothFromPtDaugh.cxx:279
 AliPtMothFromPtDaugh.cxx:280
 AliPtMothFromPtDaugh.cxx:281
 AliPtMothFromPtDaugh.cxx:282
 AliPtMothFromPtDaugh.cxx:283
 AliPtMothFromPtDaugh.cxx:284
 AliPtMothFromPtDaugh.cxx:285
 AliPtMothFromPtDaugh.cxx:286
 AliPtMothFromPtDaugh.cxx:287
 AliPtMothFromPtDaugh.cxx:288
 AliPtMothFromPtDaugh.cxx:289
 AliPtMothFromPtDaugh.cxx:290
 AliPtMothFromPtDaugh.cxx:291
 AliPtMothFromPtDaugh.cxx:292
 AliPtMothFromPtDaugh.cxx:293
 AliPtMothFromPtDaugh.cxx:294
 AliPtMothFromPtDaugh.cxx:295
 AliPtMothFromPtDaugh.cxx:296
 AliPtMothFromPtDaugh.cxx:297
 AliPtMothFromPtDaugh.cxx:298
 AliPtMothFromPtDaugh.cxx:299
 AliPtMothFromPtDaugh.cxx:300
 AliPtMothFromPtDaugh.cxx:301
 AliPtMothFromPtDaugh.cxx:302
 AliPtMothFromPtDaugh.cxx:303
 AliPtMothFromPtDaugh.cxx:304
 AliPtMothFromPtDaugh.cxx:305
 AliPtMothFromPtDaugh.cxx:306
 AliPtMothFromPtDaugh.cxx:307
 AliPtMothFromPtDaugh.cxx:308
 AliPtMothFromPtDaugh.cxx:309
 AliPtMothFromPtDaugh.cxx:310
 AliPtMothFromPtDaugh.cxx:311
 AliPtMothFromPtDaugh.cxx:312
 AliPtMothFromPtDaugh.cxx:313
 AliPtMothFromPtDaugh.cxx:314
 AliPtMothFromPtDaugh.cxx:315
 AliPtMothFromPtDaugh.cxx:316
 AliPtMothFromPtDaugh.cxx:317
 AliPtMothFromPtDaugh.cxx:318
 AliPtMothFromPtDaugh.cxx:319
 AliPtMothFromPtDaugh.cxx:320
 AliPtMothFromPtDaugh.cxx:321
 AliPtMothFromPtDaugh.cxx:322
 AliPtMothFromPtDaugh.cxx:323
 AliPtMothFromPtDaugh.cxx:324
 AliPtMothFromPtDaugh.cxx:325
 AliPtMothFromPtDaugh.cxx:326
 AliPtMothFromPtDaugh.cxx:327
 AliPtMothFromPtDaugh.cxx:328
 AliPtMothFromPtDaugh.cxx:329
 AliPtMothFromPtDaugh.cxx:330
 AliPtMothFromPtDaugh.cxx:331
 AliPtMothFromPtDaugh.cxx:332
 AliPtMothFromPtDaugh.cxx:333
 AliPtMothFromPtDaugh.cxx:334
 AliPtMothFromPtDaugh.cxx:335
 AliPtMothFromPtDaugh.cxx:336
 AliPtMothFromPtDaugh.cxx:337
 AliPtMothFromPtDaugh.cxx:338
 AliPtMothFromPtDaugh.cxx:339
 AliPtMothFromPtDaugh.cxx:340
 AliPtMothFromPtDaugh.cxx:341
 AliPtMothFromPtDaugh.cxx:342
 AliPtMothFromPtDaugh.cxx:343
 AliPtMothFromPtDaugh.cxx:344
 AliPtMothFromPtDaugh.cxx:345
 AliPtMothFromPtDaugh.cxx:346
 AliPtMothFromPtDaugh.cxx:347
 AliPtMothFromPtDaugh.cxx:348
 AliPtMothFromPtDaugh.cxx:349
 AliPtMothFromPtDaugh.cxx:350
 AliPtMothFromPtDaugh.cxx:351
 AliPtMothFromPtDaugh.cxx:352
 AliPtMothFromPtDaugh.cxx:353
 AliPtMothFromPtDaugh.cxx:354
 AliPtMothFromPtDaugh.cxx:355
 AliPtMothFromPtDaugh.cxx:356
 AliPtMothFromPtDaugh.cxx:357
 AliPtMothFromPtDaugh.cxx:358
 AliPtMothFromPtDaugh.cxx:359
 AliPtMothFromPtDaugh.cxx:360
 AliPtMothFromPtDaugh.cxx:361
 AliPtMothFromPtDaugh.cxx:362
 AliPtMothFromPtDaugh.cxx:363
 AliPtMothFromPtDaugh.cxx:364
 AliPtMothFromPtDaugh.cxx:365
 AliPtMothFromPtDaugh.cxx:366
 AliPtMothFromPtDaugh.cxx:367
 AliPtMothFromPtDaugh.cxx:368
 AliPtMothFromPtDaugh.cxx:369
 AliPtMothFromPtDaugh.cxx:370
 AliPtMothFromPtDaugh.cxx:371
 AliPtMothFromPtDaugh.cxx:372
 AliPtMothFromPtDaugh.cxx:373
 AliPtMothFromPtDaugh.cxx:374
 AliPtMothFromPtDaugh.cxx:375
 AliPtMothFromPtDaugh.cxx:376
 AliPtMothFromPtDaugh.cxx:377
 AliPtMothFromPtDaugh.cxx:378
 AliPtMothFromPtDaugh.cxx:379
 AliPtMothFromPtDaugh.cxx:380
 AliPtMothFromPtDaugh.cxx:381
 AliPtMothFromPtDaugh.cxx:382
 AliPtMothFromPtDaugh.cxx:383
 AliPtMothFromPtDaugh.cxx:384
 AliPtMothFromPtDaugh.cxx:385
 AliPtMothFromPtDaugh.cxx:386
 AliPtMothFromPtDaugh.cxx:387
 AliPtMothFromPtDaugh.cxx:388
 AliPtMothFromPtDaugh.cxx:389
 AliPtMothFromPtDaugh.cxx:390
 AliPtMothFromPtDaugh.cxx:391
 AliPtMothFromPtDaugh.cxx:392
 AliPtMothFromPtDaugh.cxx:393
 AliPtMothFromPtDaugh.cxx:394
 AliPtMothFromPtDaugh.cxx:395
 AliPtMothFromPtDaugh.cxx:396
 AliPtMothFromPtDaugh.cxx:397
 AliPtMothFromPtDaugh.cxx:398
 AliPtMothFromPtDaugh.cxx:399
 AliPtMothFromPtDaugh.cxx:400
 AliPtMothFromPtDaugh.cxx:401
 AliPtMothFromPtDaugh.cxx:402
 AliPtMothFromPtDaugh.cxx:403
 AliPtMothFromPtDaugh.cxx:404
 AliPtMothFromPtDaugh.cxx:405
 AliPtMothFromPtDaugh.cxx:406
 AliPtMothFromPtDaugh.cxx:407
 AliPtMothFromPtDaugh.cxx:408
 AliPtMothFromPtDaugh.cxx:409
 AliPtMothFromPtDaugh.cxx:410
 AliPtMothFromPtDaugh.cxx:411
 AliPtMothFromPtDaugh.cxx:412
 AliPtMothFromPtDaugh.cxx:413
 AliPtMothFromPtDaugh.cxx:414
 AliPtMothFromPtDaugh.cxx:415
 AliPtMothFromPtDaugh.cxx:416
 AliPtMothFromPtDaugh.cxx:417
 AliPtMothFromPtDaugh.cxx:418
 AliPtMothFromPtDaugh.cxx:419
 AliPtMothFromPtDaugh.cxx:420
 AliPtMothFromPtDaugh.cxx:421
 AliPtMothFromPtDaugh.cxx:422
 AliPtMothFromPtDaugh.cxx:423
 AliPtMothFromPtDaugh.cxx:424
 AliPtMothFromPtDaugh.cxx:425
 AliPtMothFromPtDaugh.cxx:426
 AliPtMothFromPtDaugh.cxx:427
 AliPtMothFromPtDaugh.cxx:428
 AliPtMothFromPtDaugh.cxx:429
 AliPtMothFromPtDaugh.cxx:430
 AliPtMothFromPtDaugh.cxx:431
 AliPtMothFromPtDaugh.cxx:432
 AliPtMothFromPtDaugh.cxx:433
 AliPtMothFromPtDaugh.cxx:434
 AliPtMothFromPtDaugh.cxx:435
 AliPtMothFromPtDaugh.cxx:436
 AliPtMothFromPtDaugh.cxx:437
 AliPtMothFromPtDaugh.cxx:438
 AliPtMothFromPtDaugh.cxx:439
 AliPtMothFromPtDaugh.cxx:440
 AliPtMothFromPtDaugh.cxx:441
 AliPtMothFromPtDaugh.cxx:442
 AliPtMothFromPtDaugh.cxx:443
 AliPtMothFromPtDaugh.cxx:444
 AliPtMothFromPtDaugh.cxx:445
 AliPtMothFromPtDaugh.cxx:446
 AliPtMothFromPtDaugh.cxx:447
 AliPtMothFromPtDaugh.cxx:448
 AliPtMothFromPtDaugh.cxx:449
 AliPtMothFromPtDaugh.cxx:450
 AliPtMothFromPtDaugh.cxx:451
 AliPtMothFromPtDaugh.cxx:452
 AliPtMothFromPtDaugh.cxx:453
 AliPtMothFromPtDaugh.cxx:454
 AliPtMothFromPtDaugh.cxx:455
 AliPtMothFromPtDaugh.cxx:456
 AliPtMothFromPtDaugh.cxx:457
 AliPtMothFromPtDaugh.cxx:458
 AliPtMothFromPtDaugh.cxx:459
 AliPtMothFromPtDaugh.cxx:460
 AliPtMothFromPtDaugh.cxx:461
 AliPtMothFromPtDaugh.cxx:462
 AliPtMothFromPtDaugh.cxx:463
 AliPtMothFromPtDaugh.cxx:464
 AliPtMothFromPtDaugh.cxx:465
 AliPtMothFromPtDaugh.cxx:466
 AliPtMothFromPtDaugh.cxx:467
 AliPtMothFromPtDaugh.cxx:468
 AliPtMothFromPtDaugh.cxx:469
 AliPtMothFromPtDaugh.cxx:470
 AliPtMothFromPtDaugh.cxx:471
 AliPtMothFromPtDaugh.cxx:472
 AliPtMothFromPtDaugh.cxx:473
 AliPtMothFromPtDaugh.cxx:474
 AliPtMothFromPtDaugh.cxx:475
 AliPtMothFromPtDaugh.cxx:476
 AliPtMothFromPtDaugh.cxx:477
 AliPtMothFromPtDaugh.cxx:478
 AliPtMothFromPtDaugh.cxx:479
 AliPtMothFromPtDaugh.cxx:480
 AliPtMothFromPtDaugh.cxx:481
 AliPtMothFromPtDaugh.cxx:482
 AliPtMothFromPtDaugh.cxx:483
 AliPtMothFromPtDaugh.cxx:484
 AliPtMothFromPtDaugh.cxx:485
 AliPtMothFromPtDaugh.cxx:486
 AliPtMothFromPtDaugh.cxx:487
 AliPtMothFromPtDaugh.cxx:488
 AliPtMothFromPtDaugh.cxx:489
 AliPtMothFromPtDaugh.cxx:490
 AliPtMothFromPtDaugh.cxx:491
 AliPtMothFromPtDaugh.cxx:492
 AliPtMothFromPtDaugh.cxx:493
 AliPtMothFromPtDaugh.cxx:494
 AliPtMothFromPtDaugh.cxx:495
 AliPtMothFromPtDaugh.cxx:496
 AliPtMothFromPtDaugh.cxx:497
 AliPtMothFromPtDaugh.cxx:498
 AliPtMothFromPtDaugh.cxx:499
 AliPtMothFromPtDaugh.cxx:500
 AliPtMothFromPtDaugh.cxx:501
 AliPtMothFromPtDaugh.cxx:502
 AliPtMothFromPtDaugh.cxx:503
 AliPtMothFromPtDaugh.cxx:504
 AliPtMothFromPtDaugh.cxx:505
 AliPtMothFromPtDaugh.cxx:506
 AliPtMothFromPtDaugh.cxx:507
 AliPtMothFromPtDaugh.cxx:508
 AliPtMothFromPtDaugh.cxx:509
 AliPtMothFromPtDaugh.cxx:510
 AliPtMothFromPtDaugh.cxx:511
 AliPtMothFromPtDaugh.cxx:512
 AliPtMothFromPtDaugh.cxx:513
 AliPtMothFromPtDaugh.cxx:514
 AliPtMothFromPtDaugh.cxx:515
 AliPtMothFromPtDaugh.cxx:516
 AliPtMothFromPtDaugh.cxx:517
 AliPtMothFromPtDaugh.cxx:518
 AliPtMothFromPtDaugh.cxx:519
 AliPtMothFromPtDaugh.cxx:520
 AliPtMothFromPtDaugh.cxx:521
 AliPtMothFromPtDaugh.cxx:522
 AliPtMothFromPtDaugh.cxx:523
 AliPtMothFromPtDaugh.cxx:524
 AliPtMothFromPtDaugh.cxx:525
 AliPtMothFromPtDaugh.cxx:526
 AliPtMothFromPtDaugh.cxx:527
 AliPtMothFromPtDaugh.cxx:528
 AliPtMothFromPtDaugh.cxx:529
 AliPtMothFromPtDaugh.cxx:530
 AliPtMothFromPtDaugh.cxx:531
 AliPtMothFromPtDaugh.cxx:532
 AliPtMothFromPtDaugh.cxx:533
 AliPtMothFromPtDaugh.cxx:534
 AliPtMothFromPtDaugh.cxx:535
 AliPtMothFromPtDaugh.cxx:536
 AliPtMothFromPtDaugh.cxx:537
 AliPtMothFromPtDaugh.cxx:538
 AliPtMothFromPtDaugh.cxx:539
 AliPtMothFromPtDaugh.cxx:540
 AliPtMothFromPtDaugh.cxx:541
 AliPtMothFromPtDaugh.cxx:542
 AliPtMothFromPtDaugh.cxx:543
 AliPtMothFromPtDaugh.cxx:544
 AliPtMothFromPtDaugh.cxx:545
 AliPtMothFromPtDaugh.cxx:546
 AliPtMothFromPtDaugh.cxx:547
 AliPtMothFromPtDaugh.cxx:548
 AliPtMothFromPtDaugh.cxx:549
 AliPtMothFromPtDaugh.cxx:550
 AliPtMothFromPtDaugh.cxx:551
 AliPtMothFromPtDaugh.cxx:552
 AliPtMothFromPtDaugh.cxx:553
 AliPtMothFromPtDaugh.cxx:554
 AliPtMothFromPtDaugh.cxx:555
 AliPtMothFromPtDaugh.cxx:556
 AliPtMothFromPtDaugh.cxx:557
 AliPtMothFromPtDaugh.cxx:558
 AliPtMothFromPtDaugh.cxx:559
 AliPtMothFromPtDaugh.cxx:560
 AliPtMothFromPtDaugh.cxx:561
 AliPtMothFromPtDaugh.cxx:562
 AliPtMothFromPtDaugh.cxx:563
 AliPtMothFromPtDaugh.cxx:564
 AliPtMothFromPtDaugh.cxx:565
 AliPtMothFromPtDaugh.cxx:566
 AliPtMothFromPtDaugh.cxx:567
 AliPtMothFromPtDaugh.cxx:568
 AliPtMothFromPtDaugh.cxx:569
 AliPtMothFromPtDaugh.cxx:570
 AliPtMothFromPtDaugh.cxx:571
 AliPtMothFromPtDaugh.cxx:572
 AliPtMothFromPtDaugh.cxx:573
 AliPtMothFromPtDaugh.cxx:574
 AliPtMothFromPtDaugh.cxx:575
 AliPtMothFromPtDaugh.cxx:576
 AliPtMothFromPtDaugh.cxx:577
 AliPtMothFromPtDaugh.cxx:578
 AliPtMothFromPtDaugh.cxx:579
 AliPtMothFromPtDaugh.cxx:580
 AliPtMothFromPtDaugh.cxx:581
 AliPtMothFromPtDaugh.cxx:582
 AliPtMothFromPtDaugh.cxx:583
 AliPtMothFromPtDaugh.cxx:584
 AliPtMothFromPtDaugh.cxx:585
 AliPtMothFromPtDaugh.cxx:586
 AliPtMothFromPtDaugh.cxx:587
 AliPtMothFromPtDaugh.cxx:588
 AliPtMothFromPtDaugh.cxx:589
 AliPtMothFromPtDaugh.cxx:590
 AliPtMothFromPtDaugh.cxx:591
 AliPtMothFromPtDaugh.cxx:592
 AliPtMothFromPtDaugh.cxx:593
 AliPtMothFromPtDaugh.cxx:594
 AliPtMothFromPtDaugh.cxx:595
 AliPtMothFromPtDaugh.cxx:596
 AliPtMothFromPtDaugh.cxx:597
 AliPtMothFromPtDaugh.cxx:598
 AliPtMothFromPtDaugh.cxx:599
 AliPtMothFromPtDaugh.cxx:600
 AliPtMothFromPtDaugh.cxx:601
 AliPtMothFromPtDaugh.cxx:602
 AliPtMothFromPtDaugh.cxx:603
 AliPtMothFromPtDaugh.cxx:604
 AliPtMothFromPtDaugh.cxx:605
 AliPtMothFromPtDaugh.cxx:606
 AliPtMothFromPtDaugh.cxx:607
 AliPtMothFromPtDaugh.cxx:608
 AliPtMothFromPtDaugh.cxx:609
 AliPtMothFromPtDaugh.cxx:610
 AliPtMothFromPtDaugh.cxx:611
 AliPtMothFromPtDaugh.cxx:612
 AliPtMothFromPtDaugh.cxx:613
 AliPtMothFromPtDaugh.cxx:614
 AliPtMothFromPtDaugh.cxx:615
 AliPtMothFromPtDaugh.cxx:616
 AliPtMothFromPtDaugh.cxx:617
 AliPtMothFromPtDaugh.cxx:618
 AliPtMothFromPtDaugh.cxx:619
 AliPtMothFromPtDaugh.cxx:620
 AliPtMothFromPtDaugh.cxx:621
 AliPtMothFromPtDaugh.cxx:622
 AliPtMothFromPtDaugh.cxx:623
 AliPtMothFromPtDaugh.cxx:624
 AliPtMothFromPtDaugh.cxx:625
 AliPtMothFromPtDaugh.cxx:626
 AliPtMothFromPtDaugh.cxx:627
 AliPtMothFromPtDaugh.cxx:628
 AliPtMothFromPtDaugh.cxx:629
 AliPtMothFromPtDaugh.cxx:630
 AliPtMothFromPtDaugh.cxx:631
 AliPtMothFromPtDaugh.cxx:632
 AliPtMothFromPtDaugh.cxx:633
 AliPtMothFromPtDaugh.cxx:634
 AliPtMothFromPtDaugh.cxx:635
 AliPtMothFromPtDaugh.cxx:636
 AliPtMothFromPtDaugh.cxx:637
 AliPtMothFromPtDaugh.cxx:638
 AliPtMothFromPtDaugh.cxx:639
 AliPtMothFromPtDaugh.cxx:640
 AliPtMothFromPtDaugh.cxx:641
 AliPtMothFromPtDaugh.cxx:642
 AliPtMothFromPtDaugh.cxx:643
 AliPtMothFromPtDaugh.cxx:644
 AliPtMothFromPtDaugh.cxx:645
 AliPtMothFromPtDaugh.cxx:646
 AliPtMothFromPtDaugh.cxx:647
 AliPtMothFromPtDaugh.cxx:648
 AliPtMothFromPtDaugh.cxx:649
 AliPtMothFromPtDaugh.cxx:650
 AliPtMothFromPtDaugh.cxx:651
 AliPtMothFromPtDaugh.cxx:652
 AliPtMothFromPtDaugh.cxx:653
 AliPtMothFromPtDaugh.cxx:654
 AliPtMothFromPtDaugh.cxx:655
 AliPtMothFromPtDaugh.cxx:656
 AliPtMothFromPtDaugh.cxx:657
 AliPtMothFromPtDaugh.cxx:658
 AliPtMothFromPtDaugh.cxx:659
 AliPtMothFromPtDaugh.cxx:660
 AliPtMothFromPtDaugh.cxx:661
 AliPtMothFromPtDaugh.cxx:662
 AliPtMothFromPtDaugh.cxx:663
 AliPtMothFromPtDaugh.cxx:664
 AliPtMothFromPtDaugh.cxx:665
 AliPtMothFromPtDaugh.cxx:666
 AliPtMothFromPtDaugh.cxx:667
 AliPtMothFromPtDaugh.cxx:668
 AliPtMothFromPtDaugh.cxx:669
 AliPtMothFromPtDaugh.cxx:670
 AliPtMothFromPtDaugh.cxx:671
 AliPtMothFromPtDaugh.cxx:672
 AliPtMothFromPtDaugh.cxx:673
 AliPtMothFromPtDaugh.cxx:674
 AliPtMothFromPtDaugh.cxx:675
 AliPtMothFromPtDaugh.cxx:676
 AliPtMothFromPtDaugh.cxx:677
 AliPtMothFromPtDaugh.cxx:678
 AliPtMothFromPtDaugh.cxx:679
 AliPtMothFromPtDaugh.cxx:680
 AliPtMothFromPtDaugh.cxx:681
 AliPtMothFromPtDaugh.cxx:682
 AliPtMothFromPtDaugh.cxx:683
 AliPtMothFromPtDaugh.cxx:684
 AliPtMothFromPtDaugh.cxx:685
 AliPtMothFromPtDaugh.cxx:686
 AliPtMothFromPtDaugh.cxx:687
 AliPtMothFromPtDaugh.cxx:688
 AliPtMothFromPtDaugh.cxx:689
 AliPtMothFromPtDaugh.cxx:690
 AliPtMothFromPtDaugh.cxx:691
 AliPtMothFromPtDaugh.cxx:692
 AliPtMothFromPtDaugh.cxx:693
 AliPtMothFromPtDaugh.cxx:694
 AliPtMothFromPtDaugh.cxx:695
 AliPtMothFromPtDaugh.cxx:696
 AliPtMothFromPtDaugh.cxx:697
 AliPtMothFromPtDaugh.cxx:698
 AliPtMothFromPtDaugh.cxx:699
 AliPtMothFromPtDaugh.cxx:700
 AliPtMothFromPtDaugh.cxx:701
 AliPtMothFromPtDaugh.cxx:702
 AliPtMothFromPtDaugh.cxx:703
 AliPtMothFromPtDaugh.cxx:704
 AliPtMothFromPtDaugh.cxx:705
 AliPtMothFromPtDaugh.cxx:706
 AliPtMothFromPtDaugh.cxx:707
 AliPtMothFromPtDaugh.cxx:708
 AliPtMothFromPtDaugh.cxx:709
 AliPtMothFromPtDaugh.cxx:710
 AliPtMothFromPtDaugh.cxx:711
 AliPtMothFromPtDaugh.cxx:712
 AliPtMothFromPtDaugh.cxx:713
 AliPtMothFromPtDaugh.cxx:714
 AliPtMothFromPtDaugh.cxx:715
 AliPtMothFromPtDaugh.cxx:716
 AliPtMothFromPtDaugh.cxx:717
 AliPtMothFromPtDaugh.cxx:718
 AliPtMothFromPtDaugh.cxx:719
 AliPtMothFromPtDaugh.cxx:720
 AliPtMothFromPtDaugh.cxx:721
 AliPtMothFromPtDaugh.cxx:722
 AliPtMothFromPtDaugh.cxx:723
 AliPtMothFromPtDaugh.cxx:724
 AliPtMothFromPtDaugh.cxx:725
 AliPtMothFromPtDaugh.cxx:726
 AliPtMothFromPtDaugh.cxx:727
 AliPtMothFromPtDaugh.cxx:728
 AliPtMothFromPtDaugh.cxx:729
 AliPtMothFromPtDaugh.cxx:730
 AliPtMothFromPtDaugh.cxx:731
 AliPtMothFromPtDaugh.cxx:732
 AliPtMothFromPtDaugh.cxx:733
 AliPtMothFromPtDaugh.cxx:734
 AliPtMothFromPtDaugh.cxx:735
 AliPtMothFromPtDaugh.cxx:736
 AliPtMothFromPtDaugh.cxx:737
 AliPtMothFromPtDaugh.cxx:738
 AliPtMothFromPtDaugh.cxx:739
 AliPtMothFromPtDaugh.cxx:740
 AliPtMothFromPtDaugh.cxx:741
 AliPtMothFromPtDaugh.cxx:742
 AliPtMothFromPtDaugh.cxx:743
 AliPtMothFromPtDaugh.cxx:744
 AliPtMothFromPtDaugh.cxx:745
 AliPtMothFromPtDaugh.cxx:746
 AliPtMothFromPtDaugh.cxx:747
 AliPtMothFromPtDaugh.cxx:748
 AliPtMothFromPtDaugh.cxx:749
 AliPtMothFromPtDaugh.cxx:750
 AliPtMothFromPtDaugh.cxx:751
 AliPtMothFromPtDaugh.cxx:752
 AliPtMothFromPtDaugh.cxx:753
 AliPtMothFromPtDaugh.cxx:754
 AliPtMothFromPtDaugh.cxx:755
 AliPtMothFromPtDaugh.cxx:756
 AliPtMothFromPtDaugh.cxx:757
 AliPtMothFromPtDaugh.cxx:758
 AliPtMothFromPtDaugh.cxx:759
 AliPtMothFromPtDaugh.cxx:760
 AliPtMothFromPtDaugh.cxx:761
 AliPtMothFromPtDaugh.cxx:762
 AliPtMothFromPtDaugh.cxx:763
 AliPtMothFromPtDaugh.cxx:764
 AliPtMothFromPtDaugh.cxx:765
 AliPtMothFromPtDaugh.cxx:766
 AliPtMothFromPtDaugh.cxx:767
 AliPtMothFromPtDaugh.cxx:768
 AliPtMothFromPtDaugh.cxx:769
 AliPtMothFromPtDaugh.cxx:770
 AliPtMothFromPtDaugh.cxx:771
 AliPtMothFromPtDaugh.cxx:772
 AliPtMothFromPtDaugh.cxx:773
 AliPtMothFromPtDaugh.cxx:774
 AliPtMothFromPtDaugh.cxx:775
 AliPtMothFromPtDaugh.cxx:776
 AliPtMothFromPtDaugh.cxx:777
 AliPtMothFromPtDaugh.cxx:778
 AliPtMothFromPtDaugh.cxx:779
 AliPtMothFromPtDaugh.cxx:780
 AliPtMothFromPtDaugh.cxx:781
 AliPtMothFromPtDaugh.cxx:782
 AliPtMothFromPtDaugh.cxx:783
 AliPtMothFromPtDaugh.cxx:784
 AliPtMothFromPtDaugh.cxx:785
 AliPtMothFromPtDaugh.cxx:786
 AliPtMothFromPtDaugh.cxx:787
 AliPtMothFromPtDaugh.cxx:788
 AliPtMothFromPtDaugh.cxx:789
 AliPtMothFromPtDaugh.cxx:790
 AliPtMothFromPtDaugh.cxx:791
 AliPtMothFromPtDaugh.cxx:792
 AliPtMothFromPtDaugh.cxx:793
 AliPtMothFromPtDaugh.cxx:794
 AliPtMothFromPtDaugh.cxx:795
 AliPtMothFromPtDaugh.cxx:796
 AliPtMothFromPtDaugh.cxx:797
 AliPtMothFromPtDaugh.cxx:798
 AliPtMothFromPtDaugh.cxx:799
 AliPtMothFromPtDaugh.cxx:800
 AliPtMothFromPtDaugh.cxx:801
 AliPtMothFromPtDaugh.cxx:802
 AliPtMothFromPtDaugh.cxx:803
 AliPtMothFromPtDaugh.cxx:804
 AliPtMothFromPtDaugh.cxx:805
 AliPtMothFromPtDaugh.cxx:806
 AliPtMothFromPtDaugh.cxx:807
 AliPtMothFromPtDaugh.cxx:808
 AliPtMothFromPtDaugh.cxx:809
 AliPtMothFromPtDaugh.cxx:810
 AliPtMothFromPtDaugh.cxx:811
 AliPtMothFromPtDaugh.cxx:812
 AliPtMothFromPtDaugh.cxx:813
 AliPtMothFromPtDaugh.cxx:814
 AliPtMothFromPtDaugh.cxx:815
 AliPtMothFromPtDaugh.cxx:816
 AliPtMothFromPtDaugh.cxx:817
 AliPtMothFromPtDaugh.cxx:818
 AliPtMothFromPtDaugh.cxx:819
 AliPtMothFromPtDaugh.cxx:820
 AliPtMothFromPtDaugh.cxx:821
 AliPtMothFromPtDaugh.cxx:822
 AliPtMothFromPtDaugh.cxx:823
 AliPtMothFromPtDaugh.cxx:824
 AliPtMothFromPtDaugh.cxx:825
 AliPtMothFromPtDaugh.cxx:826
 AliPtMothFromPtDaugh.cxx:827
 AliPtMothFromPtDaugh.cxx:828
 AliPtMothFromPtDaugh.cxx:829
 AliPtMothFromPtDaugh.cxx:830
 AliPtMothFromPtDaugh.cxx:831
 AliPtMothFromPtDaugh.cxx:832
 AliPtMothFromPtDaugh.cxx:833
 AliPtMothFromPtDaugh.cxx:834
 AliPtMothFromPtDaugh.cxx:835
 AliPtMothFromPtDaugh.cxx:836
 AliPtMothFromPtDaugh.cxx:837
 AliPtMothFromPtDaugh.cxx:838
 AliPtMothFromPtDaugh.cxx:839
 AliPtMothFromPtDaugh.cxx:840
 AliPtMothFromPtDaugh.cxx:841
 AliPtMothFromPtDaugh.cxx:842
 AliPtMothFromPtDaugh.cxx:843
 AliPtMothFromPtDaugh.cxx:844
 AliPtMothFromPtDaugh.cxx:845
 AliPtMothFromPtDaugh.cxx:846
 AliPtMothFromPtDaugh.cxx:847
 AliPtMothFromPtDaugh.cxx:848
 AliPtMothFromPtDaugh.cxx:849
 AliPtMothFromPtDaugh.cxx:850
 AliPtMothFromPtDaugh.cxx:851
 AliPtMothFromPtDaugh.cxx:852
 AliPtMothFromPtDaugh.cxx:853
 AliPtMothFromPtDaugh.cxx:854
 AliPtMothFromPtDaugh.cxx:855
 AliPtMothFromPtDaugh.cxx:856
 AliPtMothFromPtDaugh.cxx:857
 AliPtMothFromPtDaugh.cxx:858
 AliPtMothFromPtDaugh.cxx:859
 AliPtMothFromPtDaugh.cxx:860
 AliPtMothFromPtDaugh.cxx:861
 AliPtMothFromPtDaugh.cxx:862
 AliPtMothFromPtDaugh.cxx:863
 AliPtMothFromPtDaugh.cxx:864
 AliPtMothFromPtDaugh.cxx:865
 AliPtMothFromPtDaugh.cxx:866
 AliPtMothFromPtDaugh.cxx:867
 AliPtMothFromPtDaugh.cxx:868
 AliPtMothFromPtDaugh.cxx:869
 AliPtMothFromPtDaugh.cxx:870
 AliPtMothFromPtDaugh.cxx:871
 AliPtMothFromPtDaugh.cxx:872
 AliPtMothFromPtDaugh.cxx:873
 AliPtMothFromPtDaugh.cxx:874
 AliPtMothFromPtDaugh.cxx:875
 AliPtMothFromPtDaugh.cxx:876
 AliPtMothFromPtDaugh.cxx:877
 AliPtMothFromPtDaugh.cxx:878
 AliPtMothFromPtDaugh.cxx:879
 AliPtMothFromPtDaugh.cxx:880
 AliPtMothFromPtDaugh.cxx:881
 AliPtMothFromPtDaugh.cxx:882
 AliPtMothFromPtDaugh.cxx:883
 AliPtMothFromPtDaugh.cxx:884
 AliPtMothFromPtDaugh.cxx:885
 AliPtMothFromPtDaugh.cxx:886
 AliPtMothFromPtDaugh.cxx:887
 AliPtMothFromPtDaugh.cxx:888
 AliPtMothFromPtDaugh.cxx:889
 AliPtMothFromPtDaugh.cxx:890
 AliPtMothFromPtDaugh.cxx:891
 AliPtMothFromPtDaugh.cxx:892
 AliPtMothFromPtDaugh.cxx:893
 AliPtMothFromPtDaugh.cxx:894
 AliPtMothFromPtDaugh.cxx:895
 AliPtMothFromPtDaugh.cxx:896
 AliPtMothFromPtDaugh.cxx:897
 AliPtMothFromPtDaugh.cxx:898
 AliPtMothFromPtDaugh.cxx:899
 AliPtMothFromPtDaugh.cxx:900
 AliPtMothFromPtDaugh.cxx:901
 AliPtMothFromPtDaugh.cxx:902
 AliPtMothFromPtDaugh.cxx:903
 AliPtMothFromPtDaugh.cxx:904
 AliPtMothFromPtDaugh.cxx:905
 AliPtMothFromPtDaugh.cxx:906
 AliPtMothFromPtDaugh.cxx:907
 AliPtMothFromPtDaugh.cxx:908
 AliPtMothFromPtDaugh.cxx:909
 AliPtMothFromPtDaugh.cxx:910
 AliPtMothFromPtDaugh.cxx:911
 AliPtMothFromPtDaugh.cxx:912
 AliPtMothFromPtDaugh.cxx:913
 AliPtMothFromPtDaugh.cxx:914
 AliPtMothFromPtDaugh.cxx:915
 AliPtMothFromPtDaugh.cxx:916
 AliPtMothFromPtDaugh.cxx:917
 AliPtMothFromPtDaugh.cxx:918
 AliPtMothFromPtDaugh.cxx:919
 AliPtMothFromPtDaugh.cxx:920
 AliPtMothFromPtDaugh.cxx:921
 AliPtMothFromPtDaugh.cxx:922
 AliPtMothFromPtDaugh.cxx:923
 AliPtMothFromPtDaugh.cxx:924
 AliPtMothFromPtDaugh.cxx:925
 AliPtMothFromPtDaugh.cxx:926
 AliPtMothFromPtDaugh.cxx:927
 AliPtMothFromPtDaugh.cxx:928
 AliPtMothFromPtDaugh.cxx:929
 AliPtMothFromPtDaugh.cxx:930
 AliPtMothFromPtDaugh.cxx:931
 AliPtMothFromPtDaugh.cxx:932
 AliPtMothFromPtDaugh.cxx:933
 AliPtMothFromPtDaugh.cxx:934
 AliPtMothFromPtDaugh.cxx:935
 AliPtMothFromPtDaugh.cxx:936
 AliPtMothFromPtDaugh.cxx:937
 AliPtMothFromPtDaugh.cxx:938
 AliPtMothFromPtDaugh.cxx:939
 AliPtMothFromPtDaugh.cxx:940
 AliPtMothFromPtDaugh.cxx:941
 AliPtMothFromPtDaugh.cxx:942
 AliPtMothFromPtDaugh.cxx:943
 AliPtMothFromPtDaugh.cxx:944
 AliPtMothFromPtDaugh.cxx:945
 AliPtMothFromPtDaugh.cxx:946
 AliPtMothFromPtDaugh.cxx:947
 AliPtMothFromPtDaugh.cxx:948
 AliPtMothFromPtDaugh.cxx:949
 AliPtMothFromPtDaugh.cxx:950
 AliPtMothFromPtDaugh.cxx:951
 AliPtMothFromPtDaugh.cxx:952
 AliPtMothFromPtDaugh.cxx:953
 AliPtMothFromPtDaugh.cxx:954
 AliPtMothFromPtDaugh.cxx:955
 AliPtMothFromPtDaugh.cxx:956
 AliPtMothFromPtDaugh.cxx:957
 AliPtMothFromPtDaugh.cxx:958
 AliPtMothFromPtDaugh.cxx:959
 AliPtMothFromPtDaugh.cxx:960
 AliPtMothFromPtDaugh.cxx:961
 AliPtMothFromPtDaugh.cxx:962
 AliPtMothFromPtDaugh.cxx:963
 AliPtMothFromPtDaugh.cxx:964
 AliPtMothFromPtDaugh.cxx:965
 AliPtMothFromPtDaugh.cxx:966
 AliPtMothFromPtDaugh.cxx:967
 AliPtMothFromPtDaugh.cxx:968
 AliPtMothFromPtDaugh.cxx:969
 AliPtMothFromPtDaugh.cxx:970
 AliPtMothFromPtDaugh.cxx:971
 AliPtMothFromPtDaugh.cxx:972
 AliPtMothFromPtDaugh.cxx:973
 AliPtMothFromPtDaugh.cxx:974
 AliPtMothFromPtDaugh.cxx:975
 AliPtMothFromPtDaugh.cxx:976
 AliPtMothFromPtDaugh.cxx:977
 AliPtMothFromPtDaugh.cxx:978
 AliPtMothFromPtDaugh.cxx:979
 AliPtMothFromPtDaugh.cxx:980
 AliPtMothFromPtDaugh.cxx:981
 AliPtMothFromPtDaugh.cxx:982
 AliPtMothFromPtDaugh.cxx:983
 AliPtMothFromPtDaugh.cxx:984
 AliPtMothFromPtDaugh.cxx:985
 AliPtMothFromPtDaugh.cxx:986
 AliPtMothFromPtDaugh.cxx:987
 AliPtMothFromPtDaugh.cxx:988
 AliPtMothFromPtDaugh.cxx:989
 AliPtMothFromPtDaugh.cxx:990
 AliPtMothFromPtDaugh.cxx:991
 AliPtMothFromPtDaugh.cxx:992
 AliPtMothFromPtDaugh.cxx:993
 AliPtMothFromPtDaugh.cxx:994
 AliPtMothFromPtDaugh.cxx:995
 AliPtMothFromPtDaugh.cxx:996
 AliPtMothFromPtDaugh.cxx:997
 AliPtMothFromPtDaugh.cxx:998
 AliPtMothFromPtDaugh.cxx:999
 AliPtMothFromPtDaugh.cxx:1000
 AliPtMothFromPtDaugh.cxx:1001
 AliPtMothFromPtDaugh.cxx:1002
 AliPtMothFromPtDaugh.cxx:1003
 AliPtMothFromPtDaugh.cxx:1004
 AliPtMothFromPtDaugh.cxx:1005
 AliPtMothFromPtDaugh.cxx:1006
 AliPtMothFromPtDaugh.cxx:1007
 AliPtMothFromPtDaugh.cxx:1008
 AliPtMothFromPtDaugh.cxx:1009
 AliPtMothFromPtDaugh.cxx:1010
 AliPtMothFromPtDaugh.cxx:1011
 AliPtMothFromPtDaugh.cxx:1012
 AliPtMothFromPtDaugh.cxx:1013
 AliPtMothFromPtDaugh.cxx:1014
 AliPtMothFromPtDaugh.cxx:1015
 AliPtMothFromPtDaugh.cxx:1016
 AliPtMothFromPtDaugh.cxx:1017
 AliPtMothFromPtDaugh.cxx:1018
 AliPtMothFromPtDaugh.cxx:1019
 AliPtMothFromPtDaugh.cxx:1020
 AliPtMothFromPtDaugh.cxx:1021
 AliPtMothFromPtDaugh.cxx:1022
 AliPtMothFromPtDaugh.cxx:1023
 AliPtMothFromPtDaugh.cxx:1024
 AliPtMothFromPtDaugh.cxx:1025
 AliPtMothFromPtDaugh.cxx:1026
 AliPtMothFromPtDaugh.cxx:1027
 AliPtMothFromPtDaugh.cxx:1028
 AliPtMothFromPtDaugh.cxx:1029
 AliPtMothFromPtDaugh.cxx:1030
 AliPtMothFromPtDaugh.cxx:1031
 AliPtMothFromPtDaugh.cxx:1032
 AliPtMothFromPtDaugh.cxx:1033
 AliPtMothFromPtDaugh.cxx:1034
 AliPtMothFromPtDaugh.cxx:1035
 AliPtMothFromPtDaugh.cxx:1036
 AliPtMothFromPtDaugh.cxx:1037
 AliPtMothFromPtDaugh.cxx:1038
 AliPtMothFromPtDaugh.cxx:1039
 AliPtMothFromPtDaugh.cxx:1040
 AliPtMothFromPtDaugh.cxx:1041
 AliPtMothFromPtDaugh.cxx:1042
 AliPtMothFromPtDaugh.cxx:1043
 AliPtMothFromPtDaugh.cxx:1044
 AliPtMothFromPtDaugh.cxx:1045
 AliPtMothFromPtDaugh.cxx:1046
 AliPtMothFromPtDaugh.cxx:1047
 AliPtMothFromPtDaugh.cxx:1048
 AliPtMothFromPtDaugh.cxx:1049
 AliPtMothFromPtDaugh.cxx:1050
 AliPtMothFromPtDaugh.cxx:1051
 AliPtMothFromPtDaugh.cxx:1052
 AliPtMothFromPtDaugh.cxx:1053
 AliPtMothFromPtDaugh.cxx:1054
 AliPtMothFromPtDaugh.cxx:1055
 AliPtMothFromPtDaugh.cxx:1056
 AliPtMothFromPtDaugh.cxx:1057
 AliPtMothFromPtDaugh.cxx:1058
 AliPtMothFromPtDaugh.cxx:1059
 AliPtMothFromPtDaugh.cxx:1060
 AliPtMothFromPtDaugh.cxx:1061
 AliPtMothFromPtDaugh.cxx:1062
 AliPtMothFromPtDaugh.cxx:1063
 AliPtMothFromPtDaugh.cxx:1064
 AliPtMothFromPtDaugh.cxx:1065
 AliPtMothFromPtDaugh.cxx:1066
 AliPtMothFromPtDaugh.cxx:1067
 AliPtMothFromPtDaugh.cxx:1068
 AliPtMothFromPtDaugh.cxx:1069
 AliPtMothFromPtDaugh.cxx:1070
 AliPtMothFromPtDaugh.cxx:1071
 AliPtMothFromPtDaugh.cxx:1072
 AliPtMothFromPtDaugh.cxx:1073
 AliPtMothFromPtDaugh.cxx:1074
 AliPtMothFromPtDaugh.cxx:1075
 AliPtMothFromPtDaugh.cxx:1076
 AliPtMothFromPtDaugh.cxx:1077
 AliPtMothFromPtDaugh.cxx:1078
 AliPtMothFromPtDaugh.cxx:1079
 AliPtMothFromPtDaugh.cxx:1080
 AliPtMothFromPtDaugh.cxx:1081
 AliPtMothFromPtDaugh.cxx:1082
 AliPtMothFromPtDaugh.cxx:1083
 AliPtMothFromPtDaugh.cxx:1084
 AliPtMothFromPtDaugh.cxx:1085
 AliPtMothFromPtDaugh.cxx:1086
 AliPtMothFromPtDaugh.cxx:1087
 AliPtMothFromPtDaugh.cxx:1088
 AliPtMothFromPtDaugh.cxx:1089
 AliPtMothFromPtDaugh.cxx:1090
 AliPtMothFromPtDaugh.cxx:1091
 AliPtMothFromPtDaugh.cxx:1092
 AliPtMothFromPtDaugh.cxx:1093
 AliPtMothFromPtDaugh.cxx:1094
 AliPtMothFromPtDaugh.cxx:1095
 AliPtMothFromPtDaugh.cxx:1096
 AliPtMothFromPtDaugh.cxx:1097
 AliPtMothFromPtDaugh.cxx:1098
 AliPtMothFromPtDaugh.cxx:1099
 AliPtMothFromPtDaugh.cxx:1100
 AliPtMothFromPtDaugh.cxx:1101
 AliPtMothFromPtDaugh.cxx:1102
 AliPtMothFromPtDaugh.cxx:1103
 AliPtMothFromPtDaugh.cxx:1104
 AliPtMothFromPtDaugh.cxx:1105
 AliPtMothFromPtDaugh.cxx:1106
 AliPtMothFromPtDaugh.cxx:1107
 AliPtMothFromPtDaugh.cxx:1108
 AliPtMothFromPtDaugh.cxx:1109
 AliPtMothFromPtDaugh.cxx:1110
 AliPtMothFromPtDaugh.cxx:1111
 AliPtMothFromPtDaugh.cxx:1112
 AliPtMothFromPtDaugh.cxx:1113
 AliPtMothFromPtDaugh.cxx:1114
 AliPtMothFromPtDaugh.cxx:1115
 AliPtMothFromPtDaugh.cxx:1116
 AliPtMothFromPtDaugh.cxx:1117
 AliPtMothFromPtDaugh.cxx:1118
 AliPtMothFromPtDaugh.cxx:1119
 AliPtMothFromPtDaugh.cxx:1120
 AliPtMothFromPtDaugh.cxx:1121
 AliPtMothFromPtDaugh.cxx:1122
 AliPtMothFromPtDaugh.cxx:1123
 AliPtMothFromPtDaugh.cxx:1124
 AliPtMothFromPtDaugh.cxx:1125
 AliPtMothFromPtDaugh.cxx:1126
 AliPtMothFromPtDaugh.cxx:1127
 AliPtMothFromPtDaugh.cxx:1128
 AliPtMothFromPtDaugh.cxx:1129
 AliPtMothFromPtDaugh.cxx:1130
 AliPtMothFromPtDaugh.cxx:1131
 AliPtMothFromPtDaugh.cxx:1132
 AliPtMothFromPtDaugh.cxx:1133
 AliPtMothFromPtDaugh.cxx:1134
 AliPtMothFromPtDaugh.cxx:1135
 AliPtMothFromPtDaugh.cxx:1136