ROOT logo
/**************************************************************************
 * Copyright(c) 1998-1999, 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$ */
//
// Utility class to make simple Glauber type calculations 
//           for SYMMETRIC collision geometries (AA):
// Impact parameter, production points, reaction plane dependence
//
// The SimulateTrigger method can be used for simple MB and hard-process
// (binary scaling) trigger studies.
// 
// Some basic quantities can be visualized directly.
//
// The default set-up for PbPb or AUAu collisions can be read from a file 
// calling Init(1) or Init(2) if you want to read Almonds too.
//
// ***** If you change settings dont forget to call init afterwards, *****
// ***** in order to update the formulas with the new parameters.    *****
// 
// Author: andreas.morsch@cern.ch
//=================== Added by A. Dainese 11/02/04 ===========================
// Calculate path length for a parton with production point (x0,y0)
// and propagation direction (ux=cos(phi0),uy=sin(phi0)) 
// in a collision with impact parameter b and functions that make use
// of it.
//=================== Added by A. Dainese 05/03/04 ===========================
// Calculation of line integrals I0 and I1
//  integral0 =    \int_0^ellCut dl*(T_A*T_B)(x0+l*ux,y0+l*uy)
//  integral1 =    \int_0^ellCut dl*l*(T_A*T_B)(x0+l*ux,y0+l*uy)
// mostly for use in the Quenching class
//=================== Added by C. Loizdes 27/03/04 ===========================
// Handling of AuAu collisions
// More get/set functions
// Comments, units and clearing of code
//

// from AliRoot
#include "AliFastGlauber.h"
// from root
#include <TCanvas.h>
#include <TF1.h>
#include <TF2.h>
#include <TFile.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TLegend.h>
#include <TMath.h>
#include <TRandom.h>
#include <TStyle.h>

ClassImp(AliFastGlauber)

Float_t AliFastGlauber::fgBMax           = 0.;
TF1*    AliFastGlauber::fgWSb            = NULL;     
TF1*    AliFastGlauber::fgRWSb           = NULL;     
TF2*    AliFastGlauber::fgWSbz           = NULL;    
TF1*    AliFastGlauber::fgWSz            = NULL;     
TF1*    AliFastGlauber::fgWSta           = NULL;    
TF2*    AliFastGlauber::fgWStarfi        = NULL; 
TF2*    AliFastGlauber::fgWAlmond        = NULL; 
TF1*    AliFastGlauber::fgWStaa          = NULL;   
TF1*    AliFastGlauber::fgWSgeo          = NULL;   
TF1*    AliFastGlauber::fgWSbinary       = NULL;   
TF1*    AliFastGlauber::fgWSN            = NULL;   
TF1*    AliFastGlauber::fgWPathLength0   = NULL;   
TF1*    AliFastGlauber::fgWPathLength    = NULL;
TF1*    AliFastGlauber::fgWEnergyDensity = NULL;   
TF1*    AliFastGlauber::fgWIntRadius     = NULL;   
TF2*    AliFastGlauber::fgWKParticipants = NULL; 
TF1*    AliFastGlauber::fgWParticipants  = NULL; 
TF2*    AliFastGlauber::fgWAlmondCurrent = NULL;    
TF2*    AliFastGlauber::fgWAlmondFixedB[40]; 
const Int_t AliFastGlauber::fgkMCInts = 100000;
AliFastGlauber* AliFastGlauber::fgGlauber = NULL;


AliFastGlauber::AliFastGlauber(): 
    fWSr0(0.),
    fWSd(0.), 
    fWSw(0.), 
    fWSn(0.), 
    fSigmaHard(0.),
    fSigmaNN(0.),  
    fA(0),         
    fBmin(0.),     
    fBmax(0.),     
    fEllDef(0),    
    fName()     
{
  //  Default Constructor 
  //  Defaults for Pb
  SetMaxImpact();
  SetLengthDefinition();
  SetPbPbLHC();
  fXY[0] = fXY[1] = 0;
  fI0I1[0] = fI0I1[1] = 0;
}

AliFastGlauber::AliFastGlauber(const AliFastGlauber & gl)
    :TObject(gl),
     fWSr0(0.),
     fWSd(0.), 
     fWSw(0.), 
     fWSn(0.), 
     fSigmaHard(0.),
     fSigmaNN(0.),  
     fA(0),         
     fBmin(0.),     
     fBmax(0.),     
     fEllDef(0),    
     fName()     
{
// Copy constructor
    gl.Copy(*this);
    fXY[0] = fXY[1] = 0;
    fI0I1[0] = fI0I1[1] = 0;
}

AliFastGlauber* AliFastGlauber::Instance()
{ 
// Set random number generator 
    if (fgGlauber) {
	return fgGlauber;
    } else {
	fgGlauber = new AliFastGlauber();
	return fgGlauber;
    }
}

AliFastGlauber::~AliFastGlauber()
{
// Destructor
  for(Int_t k=0; k<40; k++) delete fgWAlmondFixedB[k];
}

void AliFastGlauber::SetAuAuRhic()
{
  //Set all parameters for RHIC
  SetWoodSaxonParametersAu();
  SetHardCrossSection();
  SetNNCrossSection(42);
  SetNucleus(197);
  SetFileName("$(ALICE_ROOT)/FASTSIM/data/glauberAuAu.root");
}

void AliFastGlauber::SetPbPbLHC()
{
  //Set all parameters for LHC
  SetWoodSaxonParametersPb();
  SetHardCrossSection();
  SetNNCrossSection();
  SetNucleus();
  SetFileName();
}

void AliFastGlauber::Init(Int_t mode)
{
  // Initialisation
  //    mode = 0; all functions are calculated 
  //    mode = 1; overlap function is read from file (for Pb-Pb only)
  //    mode = 2; interaction almond functions are read from file 
  //              USE THIS FOR PATH LENGTH CALC.!
  //

  // 
  Reset(); 
  //

  //
  //  Wood-Saxon
  //
  fgWSb = new TF1("WSb", WSb, 0, fgBMax, 4);
  fgWSb->SetParameter(0, fWSr0);
  fgWSb->SetParameter(1, fWSd);
  fgWSb->SetParameter(2, fWSw);
  fgWSb->SetParameter(3, fWSn);

  fgRWSb = new TF1("RWSb", RWSb, 0, fgBMax, 4);
  fgRWSb->SetParameter(0, fWSr0);
  fgRWSb->SetParameter(1, fWSd);
  fgRWSb->SetParameter(2, fWSw);
  fgRWSb->SetParameter(3, fWSn);

  fgWSbz = new TF2("WSbz", WSbz, 0, fgBMax, 0, fgBMax, 4);
  fgWSbz->SetParameter(0, fWSr0);
  fgWSbz->SetParameter(1, fWSd);
  fgWSbz->SetParameter(2, fWSw);
  fgWSbz->SetParameter(3, fWSn);

  fgWSz = new TF1("WSz", WSz, 0, fgBMax, 5);
  fgWSz->SetParameter(0, fWSr0);
  fgWSz->SetParameter(1, fWSd);
  fgWSz->SetParameter(2, fWSw);
  fgWSz->SetParameter(3, fWSn);

  //
  //  Thickness
  //
  fgWSta = new TF1("WSta", WSta, 0., fgBMax, 0);
    
  //
  //  Overlap Kernel
  //
  fgWStarfi = new TF2("WStarfi", WStarfi, 0., fgBMax, 0., TMath::Pi(), 1);
  fgWStarfi->SetParameter(0, 0.);     
  fgWStarfi->SetNpx(200);     
  fgWStarfi->SetNpy(20);    

  //
  //  Participants Kernel
  //
  fgWKParticipants = new TF2("WKParticipants", WKParticipants, 0., fgBMax, 0., TMath::Pi(), 3);
  fgWKParticipants->SetParameter(0, 0.);     
  fgWKParticipants->SetParameter(1, fSigmaNN);     
  fgWKParticipants->SetParameter(2, fA);     
  fgWKParticipants->SetNpx(200);     
  fgWKParticipants->SetNpy(20);      

  //
  //  Overlap and Participants
  //
  if (!mode) {
    fgWStaa = new TF1("WStaa", WStaa, 0., fgBMax, 1);
    fgWStaa->SetNpx(100);
    fgWStaa->SetParameter(0,fA);
    fgWStaa->SetNpx(100);
    fgWParticipants = new TF1("WParticipants", WParticipants, 0., fgBMax, 2);
    fgWParticipants->SetParameter(0, fSigmaNN);     
    fgWParticipants->SetParameter(1, fA);     
    fgWParticipants->SetNpx(100);
  } else {
    Info("Init","Reading overlap function from file %s",fName.Data()); 
    TFile* f = new TFile(fName.Data());
    if(!f->IsOpen()){
      Fatal("Init", "Could not open file %s",fName.Data());
    }
    fgWStaa = (TF1*) f->Get("WStaa");
    fgWParticipants = (TF1*) f->Get("WParticipants");
    delete f;
  }

  //
  //  Energy Density
  //
  fgWEnergyDensity = new TF1("WEnergyDensity", WEnergyDensity, 0., 2. * fWSr0, 1);
  fgWEnergyDensity->SetParameter(0, fWSr0 + 1.);
    
  //
  //  Geometrical Cross-Section
  //
  fgWSgeo = new TF1("WSgeo", WSgeo, 0., fgBMax, 1);
  fgWSgeo->SetParameter(0,fSigmaNN); //mbarn
  fgWSgeo->SetNpx(100);

  //
  //  Hard cross section (binary collisions)
  //
  fgWSbinary = new TF1("WSbinary", WSbinary, 0., fgBMax, 1);
  fgWSbinary->SetParameter(0, fSigmaHard); //mbarn
  fgWSbinary->SetNpx(100);

  //
  // Hard collisions per event
  //
  fgWSN = new TF1("WSN", WSN, 0.01, fgBMax, 1);
  fgWSN->SetNpx(100);

  //
  //  Almond shaped interaction region
  //
  fgWAlmond = new TF2("WAlmond", WAlmond, -fgBMax, fgBMax, -fgBMax, fgBMax, 1);
  fgWAlmond->SetParameter(0, 0.);     
  fgWAlmond->SetNpx(200);     
  fgWAlmond->SetNpy(200);  
  
  if(mode==2) {
    Info("Init","Reading interaction almonds from file: %s",fName.Data());
    Char_t almondName[100];
    TFile* ff = new TFile(fName.Data());
    for(Int_t k=0; k<40; k++) {
      snprintf(almondName,100, "WAlmondFixedB%d",k);
      fgWAlmondCurrent = (TF2*)ff->Get(almondName);
      fgWAlmondFixedB[k] = fgWAlmondCurrent;
    }
    delete ff;
  }

  fgWIntRadius = new TF1("WIntRadius", WIntRadius, 0., fgBMax, 1);
  fgWIntRadius->SetParameter(0, 0.);    

  //
  //  Path Length as a function of Phi
  //    
  fgWPathLength0 = new TF1("WPathLength0", WPathLength0, -TMath::Pi(), TMath::Pi(), 2);
  fgWPathLength0->SetParameter(0, 0.);
  fgWPathLength0->SetParameter(1, 0.); //Pathlength definition     

  fgWPathLength = new TF1("WPathLength", WPathLength, -TMath::Pi(), TMath::Pi(), 3);
  fgWPathLength->SetParameter(0, 0.);    //Impact Parameter
  fgWPathLength->SetParameter(1, 1000.); //Number of interactions used for average
  fgWPathLength->SetParameter(2, 0);     //Pathlength definition
}

void AliFastGlauber::Reset() const
{
  //
  // Reset dynamic allocated formulas
  // in case init is called twice

  if(fgWSb)            delete fgWSb;     
  if(fgRWSb)           delete fgRWSb;     
  if(fgWSbz)           delete fgWSbz;
  if(fgWSz)            delete fgWSz;
  if(fgWSta)           delete fgWSta;
  if(fgWStarfi)        delete fgWStarfi;
  if(fgWAlmond)        delete fgWAlmond;
  if(fgWStaa)          delete fgWStaa;
  if(fgWSgeo)          delete fgWSgeo;
  if(fgWSbinary)       delete fgWSbinary;
  if(fgWSN)            delete fgWSN;
  if(fgWPathLength0)   delete fgWPathLength0;
  if(fgWPathLength)    delete fgWPathLength;
  if(fgWEnergyDensity) delete fgWEnergyDensity;
  if(fgWIntRadius)     delete fgWIntRadius; 
  if(fgWKParticipants) delete fgWKParticipants; 
  if(fgWParticipants)  delete fgWParticipants;
}

void AliFastGlauber::DrawWSb() const
{
  //
  //  Draw Wood-Saxon Nuclear Density Function
  //
  TCanvas *c1 = new TCanvas("c1","Wood Saxon",400,10,600,700);
  c1->cd();
  Double_t max=fgWSb->GetMaximum(0.,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fwsb","Wood Saxon: #rho(r) = n (1-#omega(r/r_{0})^2)/(1+exp((r-r_{0})/d)) [fm^{-3}]",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("r [fm]");
  h2f->GetYaxis()->SetNoExponent(kTRUE);
  h2f->GetYaxis()->SetTitle("#rho [fm^{-3}]");
  h2f->Draw(); 
  fgWSb->Draw("same");
  TLegend *l1a = new TLegend(0.45,0.6,.90,0.8);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label,100, "r_{0} = %.2f fm",fWSr0);
  l1a->AddEntry(fgWSb,label,"");
  snprintf(label,100, "d = %.2f fm",fWSd);
  l1a->AddEntry(fgWSb,label,"");
  snprintf(label,100, "n = %.2e fm^{-3}",fWSn);
  l1a->AddEntry(fgWSb,label,"");
  snprintf(label,100, "#omega = %.2f",fWSw);
  l1a->AddEntry(fgWSb,label,"");
  l1a->Draw();
  c1->Update();
}

void AliFastGlauber::DrawOverlap() const
{
  //
  //  Draw Overlap Function
  //
  TCanvas *c2 = new TCanvas("c2","Overlap",400,10,600,700);
  c2->cd();
  Double_t max=fgWStaa->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2ftaa","Overlap function: T_{AB} [mbarn^{-1}]",2,0,fgBMax,2,0, max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("T_{AB} [mbarn^{-1}]");
  h2f->Draw(); 
  fgWStaa->Draw("same");
}

void AliFastGlauber::DrawParticipants() const
{
  //
  //  Draw Number of Participants Npart
  //
  TCanvas *c3 = new TCanvas("c3","Participants",400,10,600,700);
  c3->cd();
  Double_t max=fgWParticipants->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fpart","Number of Participants",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("N_{part}");
  h2f->Draw(); 
  fgWParticipants->Draw("same");
  TLegend *l1a = new TLegend(0.50,0.75,.90,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label,100, "#sigma^{inel.}_{NN} = %.1f mbarn",fSigmaNN);
  l1a->AddEntry(fgWParticipants,label,"");
  l1a->Draw();
  c3->Update();
}

void AliFastGlauber::DrawThickness() const
{
  //
  //  Draw Thickness Function
  //
  TCanvas *c4 = new TCanvas("c4","Thickness",400,10,600,700);
  c4->cd();
  Double_t max=fgWSta->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fta","Thickness function: T_{A} [fm^{-2}]",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("T_{A} [fm^{-2}]");
  h2f->Draw(); 
  fgWSta->Draw("same");
}

void AliFastGlauber::DrawGeo() const
{
  //
  //  Draw Geometrical Cross-Section
  //
  TCanvas *c5 = new TCanvas("c5","Geometrical Cross-Section",400,10,600,700);
  c5->cd();
  Double_t max=fgWSgeo->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fgeo","Differential Geometrical Cross-Section: d#sigma^{geo}_{AB}/db [fm]",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("d#sigma^{geo}_{AB}/db [fm]");
  h2f->Draw(); 
  fgWSgeo->Draw("same");
  TLegend *l1a = new TLegend(0.10,0.8,.40,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label,100, "#sigma_{NN}^{inel.} = %.1f mbarn",fSigmaNN);
  l1a->AddEntry(fgWSgeo,label,"");
  l1a->Draw();
  c5->Update();
}

void AliFastGlauber::DrawBinary() const
{
  //
  //  Draw Binary Cross-Section
  //
  TCanvas *c6 = new TCanvas("c6","Binary Cross-Section",400,10,600,700);
  c6->cd();
  Double_t max=fgWSbinary->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fbinary","Differential Binary Cross-Section: #sigma^{hard}_{NN} dT_{AB}/db [fm]",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("d#sigma^{hard}_{AB}/db [fm]");
  h2f->Draw(); 
  fgWSbinary->Draw("same");
  TLegend *l1a = new TLegend(0.50,0.8,.90,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label,100, "#sigma_{NN}^{hard} = %.1f mbarn",fSigmaHard);
  l1a->AddEntry(fgWSb,label,"");
  l1a->Draw();
  c6->Update();
}

void AliFastGlauber::DrawN() const
{
  //
  //  Draw Binaries per event (Ncoll)
  //
  TCanvas *c7 = new TCanvas("c7","Binaries per event",400,10,600,700);
  c7->cd();
  Double_t max=fgWSN->GetMaximum(0.01,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fhardcols","Number of hard collisions: T_{AB} #sigma^{hard}_{NN}/#sigma_{AB}^{geo}",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("N_{coll}");
  h2f->Draw(); 
  fgWSN->Draw("same");
  TLegend *l1a = new TLegend(0.50,0.75,.90,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label,100, "#sigma^{hard}_{NN} = %.1f mbarn",fSigmaHard);
  l1a->AddEntry(fgWSN,label,"");
  snprintf(label,100, "#sigma^{inel.}_{NN} = %.1f mbarn",fSigmaNN);
  l1a->AddEntry(fgWSN,label,"");
  l1a->Draw();
  c7->Update();
}

void AliFastGlauber::DrawKernel(Double_t b) const
{
  //
  //  Draw Kernel
  //
  TCanvas *c8 = new TCanvas("c8","Kernel",400,10,600,700);
  c8->cd();
  fgWStarfi->SetParameter(0, b);
  TH2F *h2f=new TH2F("h2fkernel","Kernel of Overlap function: d^{2}T_{AB}/dr/d#phi [fm^{-3}]",2,0,fgBMax,2,0,TMath::Pi());
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("r [fm]");
  h2f->GetYaxis()->SetTitle("#phi [rad]");
  h2f->Draw(); 
  fgWStarfi->Draw("same");
  TLegend *l1a = new TLegend(0.65,0.8,.90,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label, 100, "b = %.1f fm",b);
  l1a->AddEntry(fgWStarfi,label,"");
  l1a->Draw();
  c8->Update();
}

void AliFastGlauber::DrawAlmond(Double_t b) const
{
  //
  //  Draw Interaction Almond
  //
  TCanvas *c9 = new TCanvas("c9","Almond",400,10,600,700);
  c9->cd();
  fgWAlmond->SetParameter(0, b);
  TH2F *h2f=new TH2F("h2falmond","Interaction Almond [fm^{-4}]",2,-fgBMax, fgBMax, 2, -fgBMax, fgBMax);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("x [fm]");
  h2f->GetYaxis()->SetTitle("y [fm]");
  h2f->Draw("");
  gStyle->SetPalette(1);
  fgWAlmond->Draw("colzsame");
  TLegend *l1a = new TLegend(0.65,0.8,.90,0.9);
  l1a->SetFillStyle(0);
  l1a->SetBorderSize(0);
  Char_t label[100];
  snprintf(label, 100, "b = %.1f fm",b);
  l1a->AddEntry(fgWAlmond,label,"");
  l1a->Draw();
  c9->Update();
}

void AliFastGlauber::DrawEnergyDensity() const
{
  //
  //  Draw energy density
  //
  TCanvas *c10 = new TCanvas("c10","Energy Density",400, 10, 600, 700);
  c10->cd();
  fgWEnergyDensity->SetMinimum(0.);
  Double_t max=fgWEnergyDensity->GetMaximum(0,fgWEnergyDensity->GetParameter(0))*1.01;
  TH2F *h2f=new TH2F("h2fenergydens","Energy density",2,0,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("b [fm]");
  h2f->GetYaxis()->SetTitle("fm^{-4}");
  h2f->Draw(); 
  fgWEnergyDensity->Draw("same");
  c10->Update();
}

void AliFastGlauber::DrawPathLength0(Double_t b, Int_t iopt) const
{
  //
  //  Draw Path Length
  //
  TCanvas *c11 = new TCanvas("c11","Path Length",400,10,600,700);
  c11->cd();
  fgWPathLength0->SetParameter(0, b);
  fgWPathLength0->SetParameter(1, Double_t(iopt));
  fgWPathLength0->SetMinimum(0.); 
  fgWPathLength0->SetMaximum(10.); 
  TH2F *h2f=new TH2F("h2fpathlength0","Path length",2,-TMath::Pi(), TMath::Pi(),2,0,10.);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("#phi [rad]");
  h2f->GetYaxis()->SetTitle("l [fm]");
  h2f->Draw(); 
  fgWPathLength0->Draw("same");
}

void AliFastGlauber::DrawPathLength(Double_t b , Int_t ni, Int_t iopt) const
{
  //
  //  Draw Path Length
  //
  TCanvas *c12 = new TCanvas("c12","Path Length",400,10,600,700);
  c12->cd();
  fgWAlmond->SetParameter(0, b);
  fgWPathLength->SetParameter(0, b);
  fgWPathLength->SetParameter(1, Double_t (ni));
  fgWPathLength->SetParameter(2, Double_t (iopt));
  fgWPathLength->SetMinimum(0.); 
  fgWPathLength->SetMaximum(10.); 
  TH2F *h2f=new TH2F("h2fpathlength","Path length",2,-TMath::Pi(), TMath::Pi(),2,0,10.);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("#phi [rad]");
  h2f->GetYaxis()->SetTitle("l [fm]");
  h2f->Draw(); 
  fgWPathLength->Draw("same");
}

void AliFastGlauber::DrawIntRadius(Double_t b) const
{
  //
  //  Draw Interaction Radius
  //
  TCanvas *c13 = new TCanvas("c13","Interaction Radius",400,10,600,700);
  c13->cd();
  fgWIntRadius->SetParameter(0, b);
  fgWIntRadius->SetMinimum(0);
  Double_t max=fgWIntRadius->GetMaximum(0,fgBMax)*1.01;
  TH2F *h2f=new TH2F("h2fintradius","Interaction Density",2,0.,fgBMax,2,0,max);
  h2f->SetStats(0);
  h2f->GetXaxis()->SetTitle("r [fm]");
  h2f->GetYaxis()->SetTitle("[fm^{-3}]");
  h2f->Draw(); 
  fgWIntRadius->Draw("same");
}

Double_t AliFastGlauber::WSb(const Double_t* x, const Double_t* par)
{
  //
  //  Woods-Saxon Parameterisation
  //  as a function of radius (xx)
  //
  const Double_t kxx  = x[0];   //fm
  const Double_t kr0  = par[0]; //fm
  const Double_t kd   = par[1]; //fm   
  const Double_t kw   = par[2]; //no units
  const Double_t kn   = par[3]; //fm^-3 (used to normalize integral to one)
  Double_t y   = kn * (1.+kw*(kxx/kr0)*(kxx/kr0))/(1.+TMath::Exp((kxx-kr0)/kd));
  return y; //fm^-3
}

Double_t AliFastGlauber::RWSb(const Double_t* x, const Double_t* par)
{
  //
  //  Woods-Saxon Parameterisation
  //  as a function of radius (xx)
  //  times r**2
  const Double_t kxx  = x[0];   //fm
  const Double_t kr0  = par[0]; //fm
  const Double_t kd   = par[1]; //fm   
  const Double_t kw   = par[2]; //no units
  const Double_t kn   = par[3]; //fm^-3 (used to normalize integral to one)
  Double_t y   = kxx * kxx * kn * (1.+kw*(kxx/kr0)*(kxx/kr0))/(1.+TMath::Exp((kxx-kr0)/kd));

  return y; //fm^-1
}

Double_t AliFastGlauber::WSbz(const Double_t* x, const Double_t* par)
{
  //
  //  Wood Saxon Parameterisation
  //  as a function of z and  b
  //
  const Double_t kbb  = x[0];   //fm
  const Double_t kzz  = x[1];   //fm
  const Double_t kr0  = par[0]; //fm
  const Double_t kd   = par[1]; //fm
  const Double_t kw   = par[2]; //no units
  const Double_t kn   = par[3]; //fm^-3 (used to normalize integral to one)
  const Double_t kxx  = TMath::Sqrt(kbb*kbb+kzz*kzz);
  Double_t y  = kn * (1.+kw*(kxx/kr0)*(kxx/kr0))/(1.+TMath::Exp((kxx-kr0)/kd));
  return y; //fm^-3
}

Double_t AliFastGlauber::WSz(const Double_t* x, const Double_t* par)
{
  //
  //  Wood Saxon Parameterisation
  //  as a function of z for fixed b
  //
  const Double_t kzz  = x[0];   //fm
  const Double_t kr0  = par[0]; //fm
  const Double_t kd   = par[1]; //fm
  const Double_t kw   = par[2]; //no units
  const Double_t kn   = par[3]; //fm^-3 (used to normalize integral to one)
  const Double_t kbb  = par[4]; //fm
  const Double_t kxx  = TMath::Sqrt(kbb*kbb+kzz*kzz);
  Double_t y  = kn * (1.+kw*(kxx/kr0)*(kxx/kr0))/(1.+TMath::Exp((kxx-kr0)/kd));
  return y; //fm^-3
}

Double_t AliFastGlauber::WSta(const Double_t* x, const Double_t* /*par*/)
{
  //
  //  Thickness function T_A
  //  as a function of b
  //
  const Double_t kb  = x[0];
  fgWSz->SetParameter(4, kb);
  Double_t y  = 2. * fgWSz->Integral(0., fgBMax);
  return y; //fm^-2
}

Double_t AliFastGlauber::WStarfi(const Double_t* x, const Double_t* par)
{
  //
  //  Kernel for overlap function: T_A(s)*T_A(s-b)
  //  as a function of r and phi
  const Double_t kr1  = x[0];
  const Double_t kphi = x[1];
  const Double_t kb   = par[0];
  const Double_t kr2  = TMath::Sqrt(kr1*kr1 + kb*kb - 2.*kr1*kb*TMath::Cos(kphi)); 
  Double_t y = kr1 * fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
  return y; //fm^-3
}

Double_t AliFastGlauber::WStaa(const Double_t* x, const Double_t* par)
{
  //
  //  Overlap function 
  //  T_{AB}=Int d2s T_A(s)*T_B(s-b)
  //  as a function of b
  // (normalized to fA*fB)
  //
  const Double_t kb  = x[0];
  const Double_t ka = par[0];
  fgWStarfi->SetParameter(0, kb);

  // root integration seems to fail
  /* 
     Double_t al[2];
     Double_t bl[2];
     al[0] = 1e-6;
     al[1] = fgBMax;
     bl[0] = 0.;
     bl[1] = TMath::Pi();
     Double_t err;
     
     Double_t y =  2. *  208. * 208. * fgWStarfi->IntegralMultiple(2, al, bl, 0.001, err);
     printf("WStaa: %.5e %.5e %.5e\n", b, y, err);
  */

  //
  //  MC Integration
  //
  Double_t y = 0;
  

  for (Int_t i = 0; i < fgkMCInts; i++)
    {
	
      const Double_t kphi = TMath::Pi() * gRandom->Rndm();
      const Double_t kb1  = fgBMax      * gRandom->Rndm();	
      y += fgWStarfi->Eval(kb1, kphi);
    }
  y *= 2. * TMath::Pi() * fgBMax / fgkMCInts; //fm^-2
  y *= ka * ka * 0.1; //mbarn^-1
  return y;
}

Double_t AliFastGlauber::WKParticipants(const Double_t* x, const Double_t* par)
{
  //
  //  Kernel for number of participants
  //  as a function of r and phi
  //
  const Double_t kr1   = x[0];
  const Double_t kphi  = x[1];
  const Double_t kb    = par[0]; //fm
  const Double_t ksig  = par[1]; //mbarn
  const Double_t ka    = par[2]; //mass number
  const Double_t kr2   = TMath::Sqrt(kr1*kr1 +kb*kb - 2.*kr1*kb*TMath::Cos(kphi)); 
  const Double_t kxsi  = fgWSta->Eval(kr2) * ksig * 0.1; //no units
  /*
    Double_t y=(1-TMath::Power((1-xsi),aa))
   */
  Double_t a   = ka;
  Double_t sum = ka * kxsi;
  Double_t y   = sum;
  for (Int_t i = 1; i <= ka; i++)
    {
      a--;
      sum *= (-kxsi) * a / Float_t(i+1);
      y  += sum;
    }
  y *= kr1 * fgWSta->Eval(kr1);
  return y; //fm^-1
}

Double_t AliFastGlauber::WParticipants(const Double_t* x, const Double_t* par)
{
  //
  //  Number of Participants as 
  //  a function of b
  //
  const Double_t kb = x[0];
  const Double_t ksig  = par[0]; //mbarn
  const Double_t ka   = par[1];  //mass number
  fgWKParticipants->SetParameter(0, kb);
  fgWKParticipants->SetParameter(1, ksig);
  fgWKParticipants->SetParameter(2, ka);

  //
  //  MC Integration
  //
  Double_t y = 0;
  for (Int_t i = 0; i < fgkMCInts; i++)
    {
      const Double_t kphi = TMath::Pi() * gRandom->Rndm();
      const Double_t kb1  = fgBMax      * gRandom->Rndm();	
      y += fgWKParticipants->Eval(kb1, kphi);
    }
  y *= 2. *  ka * 2. * TMath::Pi() * fgBMax / fgkMCInts;
  return y; //no units
}

Double_t AliFastGlauber::WSgeo(const Double_t* x, const Double_t* par)
{
  //
  //  Geometrical Cross-Section
  //  as a function of b
  //
  const Double_t kb     = x[0];              //fm
  const Double_t ksigNN = par[0];            //mbarn
  const Double_t ktaa   = fgWStaa->Eval(kb); //mbarn^-1
  Double_t y     = 2. * TMath::Pi() * kb * (1. - TMath::Exp(- ksigNN * ktaa)); 
  return y; //fm
}

Double_t AliFastGlauber::WSbinary(const Double_t* x, const Double_t* par)
{
  //
  //  Number of binary hard collisions
  //  as a function of b
  //
  const Double_t kb     = x[0];              //fm
  const Double_t ksig   = par[0];            //mbarn
  const Double_t ktaa   = fgWStaa->Eval(kb); //mbarn^-1
  Double_t y = 2. * TMath::Pi() * kb * ksig * ktaa; 
  return y; //fm
}

Double_t AliFastGlauber::WSN(const Double_t* x, const Double_t* /*par*/)
{
  //
  //  Number of hard processes per event
  //  as a function of b
  const Double_t kb = x[0];
  Double_t y = fgWSbinary->Eval(kb)/fgWSgeo->Eval(kb);
  return y; //no units
}

Double_t AliFastGlauber::WEnergyDensity(const Double_t* x, const Double_t* par)
{
  //
  //  Initial energy density 
  //  as a function of the impact parameter
  //
  const Double_t kb     = x[0];
  const Double_t krA    = par[0];
  //
  //  Attention: area of transverse reaction zone in hard-sphere approximation !     
  const Double_t krA2=krA*krA;
  const Double_t kb2=kb*kb;  
  const Double_t ksaa = (TMath::Pi() - 2. * TMath::ASin(kb/ 2./ krA)) * krA2 
                      - kb * TMath::Sqrt(krA2 - kb2/ 4.); //fm^2
  const Double_t ktaa = fgWStaa->Eval(kb); //mbarn^-1
  Double_t y=ktaa/ksaa*10;
  return y; //fm^-4
}

Double_t AliFastGlauber::WAlmond(const Double_t* x, const Double_t* par)
{
  //
  //  Almond shaped interaction region
  //  as a function of cartesian x,y.
  //
  const Double_t kb    = par[0];
  const Double_t kxx   = x[0] + kb/2.;
  const Double_t kyy   = x[1];
  const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
  const Double_t kphi  = TMath::ATan2(kyy,kxx);
  const Double_t kr2   = TMath::Sqrt(kr1*kr1 + kb*kb - 2.*kr1*kb*TMath::Cos(kphi)); 
  //
  //  Interaction probability calculated as product of thicknesses
  //
  Double_t y    = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
  return y; //fm^-4
}

Double_t AliFastGlauber::WIntRadius(const Double_t* x, const Double_t* par)
{
  //
  //  Average interaction density over radius 
  //  at which interaction takes place
  //  as a function of radius
  //
  const Double_t kr    = x[0];
  const Double_t kb    = par[0];
  fgWAlmond->SetParameter(0, kb);
  //  Average over phi in small steps   
  const Double_t kdphi = 2. * TMath::Pi() / 100.;
  Double_t phi  = 0.;
  Double_t y    = 0.;
  for (Int_t i = 0; i < 100; i++) {
    const Double_t kxx = kr * TMath::Cos(phi);
    const Double_t kyy = kr * TMath::Sin(phi);
    y   += fgWAlmond->Eval(kxx,kyy);
    phi += kdphi;
  } // phi loop
  // Result multiplied by Jacobian (2 pi r)     
  y *= 2. * TMath::Pi() * kr / 100.;
  return y; //fm^-3
}

Double_t AliFastGlauber::WPathLength0(const Double_t* x, const Double_t* par)
{
  //
  //  Path Length as a function of phi 
  //  for interaction point fixed at (0,0)
  //  as a function of phi-direction
  //
  //  Phi direction in Almond
  const Double_t kphi0   = x[0];
  const Double_t kb      = par[0];
  //  Path Length definition
  const Int_t    kiopt   = Int_t(par[1]);

  //  Step along radial direction phi   
  const Int_t    kNp  = 100; // Steps in r 
  const Double_t kDr  = fgBMax/kNp;
  Double_t r  = 0.;
  Double_t rw = 0.;
  Double_t w  = 0.;
  for (Int_t i = 0; i < kNp; i++) {
    //
    //  Transform into target frame
    //
    const Double_t kxx   = r * TMath::Cos(kphi0) + kb / 2.;
    const Double_t kyy   = r * TMath::Sin(kphi0);
    const Double_t kphi  = TMath::ATan2(kyy, kxx);
    const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
    // Radius in projectile frame
    const Double_t kr2   = TMath::Sqrt(kr1*kr1 + kb*kb - 2.*kr1*kb*TMath::Cos(kphi)); 
    const Double_t ky    = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);

    rw += ky * r;
    w  += ky;
    r  += kDr;
  } // radial steps

  Double_t y=0.;
  if (!kiopt) { // My length definition (is exact for hard disk)
      if(w) y= 2. * rw / w; 
  } else {
      const Double_t knorm=fgWSta->Eval(1e-4);
      if(knorm) y =  TMath::Sqrt(2. * rw * kDr / knorm / knorm);
  }
  return y; //fm
}

Double_t AliFastGlauber::WPathLength(const Double_t* x, const Double_t* par)
{
  //
  //  Path Length as a function of phi 
  //  Interaction point from random distribution
  //  as a function of the phi-direction
  const Double_t kphi0   = x[0];
  const Double_t kb      = par[0];
  fgWAlmond->SetParameter(0, kb); 
  const Int_t    kNpi    = Int_t (par[1]); //Number of interactions
  const Int_t    kiopt   = Int_t(par[2]);  //Path Length definition 

  //
  //  r-steps
  // 
  const Int_t    kNp   = 100;
  const Double_t kDr  = fgBMax/Double_t(kNp);
  Double_t l = 0.;  //  Path length 
  for (Int_t in = 0; in < kNpi; in ++) {
    Double_t rw = 0.;
    Double_t w  = 0.;
    // Interaction point
    Double_t x0, y0;
    fgWAlmond->GetRandom2(x0, y0);
    // Initial radius
    const Double_t kr0  = TMath::Sqrt(x0*x0 + y0*y0);
    const Int_t    knps = Int_t ((fgBMax - kr0)/kDr) - 1;
	
    // Radial steps
    Double_t r  = 0.;
    for (Int_t i = 0; (i < knps ); i++) {
      // Transform into target frame
      const Double_t kxx   = x0 + r * TMath::Cos(kphi0) + kb / 2.;
      const Double_t kyy   = y0 + r * TMath::Sin(kphi0);
      const Double_t kphi  = TMath::ATan2(kyy, kxx);
      const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
      // Radius in projectile frame
      const Double_t kr2   = TMath::Sqrt(kr1*kr1 + kb*kb - 2.*kr1*kb*TMath::Cos(kphi)); 
      const Double_t ky    = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
	    
      rw += ky * r;
      w  += ky;
      r  += kDr;
    } // steps
    // Average over interactions
    if (!kiopt) {
      if(w) l += (2. * rw / w);
    } else {
      const Double_t knorm=fgWSta->Eval(1e-4);
      if(knorm) l+= 2. * rw * kDr / knorm / knorm;
    }
  } // interactions
  Double_t ret=0;
  if (!kiopt) 
    ret= l / kNpi;
  else 
    ret=TMath::Sqrt( l / kNpi);
  return ret; //fm
}

Double_t AliFastGlauber::CrossSection(Double_t b1, Double_t b2) const
{
  //
  // Return the geometrical cross-section integrated from b1 to b2 
  //
  return fgWSgeo->Integral(b1, b2)*10.; //mbarn
}

Double_t AliFastGlauber::HardCrossSection(Double_t b1, Double_t b2) const
{
  //
  // Return the hard cross-section integrated from b1 to b2 
  //
  return fgWSbinary->Integral(b1, b2)*10.; //mbarn
}

Double_t AliFastGlauber::FractionOfHardCrossSection(Double_t b1, Double_t b2) const
{
  //
  // Return fraction of hard cross-section integrated from b1 to b2 
  //
  return fgWSbinary->Integral(b1, b2)/fgWSbinary->Integral(0., 100.);
}

Double_t AliFastGlauber::NHard(const Double_t b1, const Double_t b2) const
{
  //
  //  Number of binary hard collisions 
  //  as a function of b (nucl/ex/0302016 eq. 19)
  //
  const Double_t kshard=HardCrossSection(b1,b2);
  const Double_t ksgeo=CrossSection(b1,b2); 
  if(ksgeo>0)
    return kshard/ksgeo;
  else return -1; 
}

Double_t AliFastGlauber::Binaries(Double_t b) const
{
  //
  // Return number of binary hard collisions normalized to 1 at b=0
  //
  if(b < 1.e-4) b = 1e-4;
  return fgWSN->Eval(b)/fgWSN->Eval(1e-4);
}

Double_t AliFastGlauber::MeanOverlap(Double_t b1, Double_t b2)
{
//
// Calculate the mean overlap for impact parameter range b1 .. b2
//
    Double_t sum  = 0.;
    Double_t sumc = 0.;
    Double_t b    = b1;
    
    while (b < b2-0.005) {
	Double_t  nc = GetNumberOfCollisions(b);
	sum  += 10. * fgWStaa->Eval(b) * fgWSgeo->Eval(b) * 0.01 / (1. - TMath::Exp(-nc));
	sumc += 10. * fgWSgeo->Eval(b) * 0.01;
	b += 0.01;
    }
    return (sum / CrossSection(b1, b2));
}


Double_t AliFastGlauber::MeanNumberOfCollisionsPerEvent(Double_t b1, Double_t b2)
{
//
// Calculate the mean number of collisions per event for impact parameter range b1 .. b2
//
    Double_t sum  = 0.;
    Double_t sumc = 0.;
    Double_t b    = b1;
    
    while (b < b2-0.005) {
	Double_t  nc = GetNumberOfCollisions(b);
	sum  += nc / (1. - TMath::Exp(-nc)) * 10. * fgWSgeo->Eval(b) * 0.01;
	sumc += 10. * fgWSgeo->Eval(b) * 0.01;
	b += 0.01;
    }
    return (sum / CrossSection(b1, b2));
}


Double_t AliFastGlauber::GetNumberOfBinaries(Double_t b) const
{
  //
  // Return number of binary hard collisions at b
  //
  if(b<1.e-4) b=1e-4;
  return fgWSN->Eval(b);
}

Double_t AliFastGlauber::Participants(Double_t  b) const
{
  //
  // Return the number of participants normalized to 1 at b=0
  //
  if(b<1.e-4) b=1e-4;
  return (fgWParticipants->Eval(b)/fgWParticipants->Eval(1e-4));
}

Double_t AliFastGlauber::GetNumberOfParticipants(Double_t  b) const
{
  //
  // Return the number of participants for impact parameter b
  //
  if(b<1.e-4) b=1e-4;
  return (fgWParticipants->Eval(b));
}

Double_t AliFastGlauber::GetNumberOfCollisions(Double_t  b) const
{
  //
  // Return the number of collisions for impact parameter b
  //
  if(b<1.e-4) b=1e-4;
  return (fgWStaa->Eval(b)*fSigmaNN);
}

Double_t AliFastGlauber::GetNumberOfCollisionsPerEvent(Double_t  b) const
{
  //
  // Return the number of collisions per event (at least one collision)
  // for impact parameter b
  //
    Double_t n = GetNumberOfCollisions(b);
    if (n > 0.) {
	return (n / (1. - TMath::Exp(- n)));
    } else {
	return (0.);
    }
}

void AliFastGlauber::SimulateTrigger(Int_t n)
{
  //
  //  Simulates Trigger
  //
  TH1F* mbtH = new TH1F("mbtH", "MB Trigger b-Distribution",   100, 0., 20.);
  TH1F* hdtH = new TH1F("hdtH", "Hard Trigger b-Distribution", 100, 0., 20.);   
  TH1F* mbmH = new TH1F("mbmH", "MB Trigger Multiplicity Distribution",   100, 0., 8000.);
  TH1F* hdmH = new TH1F("hdmH", "Hard Trigger Multiplicity Distribution", 100, 0., 8000.);   

  mbtH->SetXTitle("b [fm]");
  hdtH->SetXTitle("b [fm]");    
  mbmH->SetXTitle("Multiplicity");
  hdmH->SetXTitle("Multiplicity");    

  TCanvas *c0 = new TCanvas("c0","Trigger Simulation",400,10,600,700);    
  c0->Divide(2,1);
  TCanvas *c1 = new TCanvas("c1","Trigger Simulation",400,10,600,700);    
  c1->Divide(1,2);

  //
  //
  Init(1);
  for (Int_t iev = 0; iev < n; iev++)
    {
      Float_t b, p, mult;
      GetRandom(b, p, mult);
      mbtH->Fill(b,1.);
      hdtH->Fill(b, p);
      mbmH->Fill(mult, 1.);
      hdmH->Fill(mult, p);

      c0->cd(1);
      mbtH->Draw();
      c0->cd(2);
      hdtH->Draw();	
      c0->Update();

      c1->cd(1);
      mbmH->Draw();
      c1->cd(2);
      hdmH->Draw();	
      c1->Update();
    }
}

void AliFastGlauber::GetRandom(Float_t& b, Float_t& p, Float_t& mult)
{
  //
  // Gives back a random impact parameter, hard trigger probability and multiplicity
  //
  b = fgWSgeo->GetRandom();
  const Float_t kmu = fgWSN->Eval(b);
  p = 1.-TMath::Exp(-kmu);
  mult = 6000./fgWSN->Eval(1.) * kmu;
}

void AliFastGlauber::GetRandom(Int_t& bin, Bool_t& hard)
{
  //
  // Gives back a random impact parameter bin, and hard trigger decission
  //
  const Float_t kb  = fgWSgeo->GetRandom();
  const Float_t kmu = fgWSN->Eval(kb) * fSigmaHard;
  const Float_t kp  = 1.-TMath::Exp(-kmu);
  if (kb < 5.) {
    bin = 1;
  } else if (kb <  8.6) {
    bin = 2;
  } else if (kb < 11.2) {
    bin = 3;
  } else if (kb < 13.2) {
    bin = 4;
  } else if (kb < 15.0) {
    bin = 5;
  } else {
    bin = 6;
  }
  hard = kFALSE;
  const Float_t kr = gRandom->Rndm();
  if (kr < kp) hard = kTRUE;
}

Double_t  AliFastGlauber::GetRandomImpactParameter(Double_t bmin, Double_t bmax)
{
  //
  // Gives back a random impact parameter in the range bmin .. bmax
  //
  Float_t b = -1.;
  while(b < bmin || b > bmax)
    b = fgWSgeo->GetRandom();
  return b;
}

void AliFastGlauber::StoreFunctions() const
{
  //
  // Store in file functions
  //
  TFile* ff = new TFile(fName.Data(),"recreate");
  fgWStaa->Write("WStaa");
  fgWParticipants->Write("WParticipants");
  ff->Close();
  return;
}

//=================== Added by A. Dainese 11/02/04 ===========================

void AliFastGlauber::StoreAlmonds() const
{
  //
  // Store in file 
  // 40 almonds for b = (0.25+k*0.5) fm (k=0->39)
  //
  Char_t almondName[100];
  TFile* ff = new TFile(fName.Data(),"update");
  for(Int_t k=0; k<40; k++) {
    snprintf(almondName, 100, "WAlmondFixedB%d",k);
    Double_t b = 0.25+k*0.5;
    Info("StoreAlmonds"," b = %f\n",b); 
    fgWAlmond->SetParameter(0,b);
    fgWAlmond->Write(almondName);
  }
  ff->Close();
  return;
}

void AliFastGlauber::SetCentralityClass(Double_t xsecFrLow,Double_t xsecFrUp)
{
  //
  // Set limits of centrality class as fractions 
  // of the geomtrical cross section  
  //
  if(xsecFrLow>1. || xsecFrUp>1. || xsecFrLow>xsecFrUp) {
    Error("SetCentralityClass", "Please set 0 <= xsecFrLow <= xsecFrUp <= 1\n");
    return;
  }

  Double_t bLow=0.,bUp=0.;
  Double_t xsecFr=0.;
  const Double_t knorm=fgWSgeo->Integral(0.,100.);
  while(xsecFr<xsecFrLow) {
    xsecFr = fgWSgeo->Integral(0.,bLow)/knorm;
    bLow += 0.1;
  }
  bUp = bLow;
  while(xsecFr<xsecFrUp) {
    xsecFr = fgWSgeo->Integral(0.,bUp)/knorm;
    bUp += 0.1;
  }

  Info("SetCentralityClass", "Centrality class: %4.2f-%4.2f; %4.1f < b < %4.1f fm",
	 xsecFrLow,xsecFrUp,bLow,bUp);
  fgWSbinary->SetRange(bLow,bUp);
  fBmin=bLow;
  fBmax=bUp;
  return;
}

void AliFastGlauber::GetRandomBHard(Double_t& b)
{
  //
  // Get random impact parameter according to distribution of 
  // hard (binary) cross-section, in the range defined by the centrality class
  //
  b = fgWSbinary->GetRandom();
  Int_t bin = 2*(Int_t)b;
  if( (b-(Int_t)b) > 0.5) bin++;
  fgWAlmondCurrent = fgWAlmondFixedB[bin]; 
  return;
}

void AliFastGlauber::GetRandomXY(Double_t& x,Double_t& y)
{
  //
  // Get random position of parton production point according to 
  // product of thickness functions
  //
  fgWAlmondCurrent->GetRandom2(x,y);
  return;
}

void AliFastGlauber::GetRandomPhi(Double_t& phi) 
{
  //
  // Get random parton azimuthal propagation direction
  //
  phi = 2.*TMath::Pi()*gRandom->Rndm();
  return;
}

Double_t AliFastGlauber::CalculateLength(Double_t b,Double_t x0,Double_t y0,Double_t phi0)
{
  // 
  // Calculate path length for a parton with production point (x0,y0)
  // and propagation direction (ux=cos(phi0),uy=sin(phi0)) 
  // in a collision with impact parameter b
  //

  // number of steps in l
  const Int_t    kNp  = 100;
  const Double_t kDl  = fgBMax/Double_t(kNp);

  if(fEllDef==1) {
    //
    // Definition 1:
    // 
    // ell = 2 * \int_0^\infty dl*l*(T_A*T_B)(x0+l*ux,y0+l*uy) /
    //           \int_0^\infty dl*(T_A*T_B)(x0+l*ux,y0+l*uy)
    //

    // Initial radius
    const Double_t kr0 = TMath::Sqrt(x0*x0 + y0*y0);
    const Int_t knps = Int_t ((fgBMax - kr0)/kDl) - 1;
    Double_t l  = 0.;
    Double_t integral1 = 0.;
    Double_t integral2 = 0.;
    // Radial steps
    for (Int_t i = 0; i < knps; i++) {
      
      // Transform into target frame
      const Double_t kxx   = x0 + l * TMath::Cos(phi0) + b / 2.;
      const Double_t kyy   = y0 + l * TMath::Sin(phi0);
      const Double_t kphi  = TMath::ATan2(kyy, kxx);
      const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
      // Radius in projectile frame
      const Double_t kr2   = TMath::Sqrt(kr1*kr1 + b*b - 2.*kr1*b*TMath::Cos(kphi)); 
      const Double_t kprodTATB = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
      
      integral1 += kprodTATB * l * kDl;
      integral2 += kprodTATB * kDl;
      l  += kDl;
    } // steps
    
    Double_t ell=0.;
    if(integral2)
      ell = (2. * integral1 / integral2);
    return ell;
  } else if(fEllDef==2) {
    //
    // Definition 2:
    // 
    // ell = \int_0^\infty dl*
    //       \Theta((T_A*T_B)(x0+l*ux,y0+l*uy)-0.5*(T_A*T_B)(0,0))
    //

    // Initial radius
    const Double_t kr0  = TMath::Sqrt(x0*x0 + y0*y0);
    const Int_t knps = Int_t ((fgBMax - kr0)/kDl) - 1;
    const Double_t kprodTATBHalfMax = 0.5*fgWAlmondCurrent->Eval(0.,0.);
    // Radial steps
    Double_t l  = 0.;
    Double_t integral = 0.;
    for (Int_t i = 0; i < knps; i++) {
      // Transform into target frame
      const Double_t kxx   = x0 + l * TMath::Cos(phi0) + b / 2.;
      const Double_t kyy   = y0 + l * TMath::Sin(phi0);
      const Double_t kphi  = TMath::ATan2(kyy, kxx);
      const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
      // Radius in projectile frame
      const Double_t kr2   = TMath::Sqrt(kr1*kr1 + b*b - 2.*kr1*b*TMath::Cos(kphi)); 
      const Double_t kprodTATB = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
      if(kprodTATB>kprodTATBHalfMax) integral += kDl;
      l  += kDl;
    } // steps
    Double_t ell = integral;
    return ell;
  } else {
    Error("CalculateLength","Wrong length definition setting: %d !\n",fEllDef);
    return -1.;
  }
}

void AliFastGlauber::GetLengthAndPhi(Double_t& ell,Double_t& phi,Double_t b)
{
  //
  // Return length from random b, x0, y0, phi0 
  // Return also phi0
  //
  Double_t x0,y0,phi0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  GetRandomPhi(phi0);
  phi = phi0;
  ell = CalculateLength(b,x0,y0,phi0);
  return;
}

void AliFastGlauber::GetLength(Double_t& ell,Double_t b)
{
  //
  // Return length from random b, x0, y0, phi0 
  //
  Double_t phi;
  GetLengthAndPhi(ell,phi,b);
  return;
}

void AliFastGlauber::GetLengthsBackToBackAndPhi(Double_t& ell1,Double_t& ell2,Double_t &phi,Double_t b)
{
  //
  // Return 2 lengths back to back from random b, x0, y0, phi0 
  // Return also phi0 
 // 
  Double_t x0,y0,phi0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  GetRandomPhi(phi0);
  const Double_t kphi0plusPi = phi0+TMath::Pi();
  phi = phi0;
  ell1 = CalculateLength(b,x0,y0,phi0);
  ell2 = CalculateLength(b,x0,y0,kphi0plusPi);
  return;
}

void AliFastGlauber::GetLengthsBackToBack(Double_t& ell1,Double_t& ell2,
					  Double_t b)
{
  //
  // Return 2 lengths back to back from random b, x0, y0, phi0 
  // 
  Double_t phi;
  GetLengthsBackToBackAndPhi(ell1,ell2,phi,b);
  return;
}

void AliFastGlauber::GetLengthsForPythia(Int_t n,Double_t* const phi,Double_t* ell, Double_t b)
{
  //
  // Returns lenghts for n partons with azimuthal angles phi[n] 
  // from random b, x0, y0
  //
  Double_t x0, y0;
  if(b < 0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  for(Int_t i = 0; i< n; i++) ell[i] = CalculateLength(b,x0,y0,phi[i]);
  return;
}

void AliFastGlauber::PlotBDistr(Int_t n)
{  
  // 
  // Plot distribution of n impact parameters
  //
  Double_t b;
  TH1F *hB = new TH1F("hB","dN/db",100,0,fgBMax); 
  hB->SetXTitle("b [fm]");
  hB->SetYTitle("dN/db [a.u.]");
  hB->SetFillColor(3);
  for(Int_t i=0; i<n; i++) {
    GetRandomBHard(b);
    hB->Fill(b);
  }
  TCanvas *cB = new TCanvas("cB","Impact parameter distribution",0,0,500,500);
  cB->cd();
  hB->Draw();
  return;
}

void AliFastGlauber::PlotLengthDistr(Int_t n,Bool_t save,const char *fname)
{
  //
  // Plot length distribution
  //
  Double_t ell;
  TH1F *hEll = new TH1F("hEll","Length distribution",64,-0.5,15); 
  hEll->SetXTitle("Transverse path length, L [fm]");
  hEll->SetYTitle("Probability");
  hEll->SetFillColor(2);
  for(Int_t i=0; i<n; i++) {
    GetLength(ell);
    hEll->Fill(ell);
  }
  hEll->Scale(1/(Double_t)n);
  TCanvas *cL = new TCanvas("cL","Length distribution",0,0,500,500);
  cL->cd();
  hEll->Draw();

  if(save) {
    TFile *f = new TFile(fname,"recreate");
    hEll->Write();
    f->Close();
  }
  return;
}

void AliFastGlauber::PlotLengthB2BDistr(Int_t n,Bool_t save,const char *fname)
{
  //
  // Plot lengths back-to-back distributions
  //
  Double_t ell1,ell2;
  TH2F *hElls = new TH2F("hElls","Lengths back-to-back",100,0,15,100,0,15); 
  hElls->SetXTitle("Transverse path length, L [fm]");
  hElls->SetYTitle("Transverse path length, L [fm]");
  for(Int_t i=0; i<n; i++) {
    GetLengthsBackToBack(ell1,ell2);
    hElls->Fill(ell1,ell2);
  }
  hElls->Scale(1/(Double_t)n);
  TCanvas *cLs = new TCanvas("cLs","Length back-to-back distribution",0,0,500,500);
  gStyle->SetPalette(1,0);
  cLs->cd();
  hElls->Draw("col,Z");
  if(save) {
    TFile *f = new TFile(fname,"recreate");
    hElls->Write();
    f->Close();
  }
  return;
}

void AliFastGlauber::PlotAlmonds() const
{
  //
  // Plot almonds for some impact parameters
  //
  TCanvas *c = new TCanvas("c","Almonds",0,0,500,500);
  gStyle->SetPalette(1,0);
  c->Divide(2,2);
  c->cd(1);
  fgWAlmondFixedB[0]->Draw("cont1");
  c->cd(2);
  fgWAlmondFixedB[10]->Draw("cont1");
  c->cd(3);
  fgWAlmondFixedB[20]->Draw("cont1");
  c->cd(4);
  fgWAlmondFixedB[30]->Draw("cont1");
  return;
}

//=================== Added by A. Dainese 05/03/04 ===========================

void AliFastGlauber::CalculateI0I1(Double_t& integral0,Double_t& integral1,
				   Double_t b,Double_t x0,Double_t y0,
                                   Double_t phi0,Double_t ellCut) const
{
  // 
  // Calculate integrals: 
  //  integral0 = \int_0^ellCut dl*(T_A*T_B)(x0+l*ux,y0+l*uy)
  //  integral1 = \int_0^ellCut dl*l*(T_A*T_B)(x0+l*ux,y0+l*uy)
  //
  // for a parton with production point (x0,y0)
  // and propagation direction (ux=cos(phi0),uy=sin(phi0)) 
  // in a collision with impact parameter b
  // 

  // number of steps in l
  const Int_t    kNp  = 100;
  const Double_t kDl  = fgBMax/Double_t(kNp);
    
  // Initial radius
  const Double_t kr0  = TMath::Sqrt(x0 * x0 + y0 * y0);
  const Int_t knps = Int_t ((fgBMax - kr0)/kDl) - 1;
    
  // Radial steps
  Double_t l  = 0.;
  integral0 = 0.;
  integral1 = 0.;
  Int_t i = 0;
  while((i < knps) && (l < ellCut)) {
    // Transform into target frame
    const Double_t kxx   = x0 + l * TMath::Cos(phi0) + b / 2.;
    const Double_t kyy   = y0 + l * TMath::Sin(phi0);
    const Double_t kphi  = TMath::ATan2(kyy, kxx);
    const Double_t kr1   = TMath::Sqrt(kxx*kxx + kyy*kyy);
    // Radius in projectile frame
    const Double_t kr2   = TMath::Sqrt(kr1*kr1 + b*b - 2.*kr1*b*TMath::Cos(kphi)); 
    const Double_t kprodTATB = fgWSta->Eval(kr1) * fgWSta->Eval(kr2);
    integral0 += kprodTATB * kDl;
    integral1 += kprodTATB * l * kDl;
    l  += kDl;
    i++;
  } // steps
  return;  
}

void AliFastGlauber::GetI0I1AndPhi(Double_t& integral0,Double_t& integral1,
				   Double_t& phi,
				   Double_t ellCut,Double_t b)
{
  //
  // Return I0 and I1 from random b, x0, y0, phi0 
  // Return also phi
  //
  Double_t x0,y0,phi0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  GetRandomPhi(phi0);
  phi = phi0;
  CalculateI0I1(integral0,integral1,b,x0,y0,phi0,ellCut);
  return;
}

void AliFastGlauber::GetI0I1(Double_t& integral0,Double_t& integral1,
			     Double_t ellCut,Double_t b)
{
  //
  // Return I0 and I1 from random b, x0, y0, phi0 
  //
  Double_t phi;
  GetI0I1AndPhi(integral0,integral1,phi,ellCut,b);
  return;
}

void AliFastGlauber::GetI0I1BackToBackAndPhi(Double_t& integral01,Double_t& integral11,
					     Double_t& integral02,Double_t& integral12,
					     Double_t& phi,
					     Double_t ellCut,Double_t b)
{
  //
  // Return 2 pairs of I0 and I1 back to back from random b, x0, y0, phi0 
  // Return also phi0
  //
  Double_t x0,y0,phi0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  GetRandomPhi(phi0);
  phi = phi0;
  const Double_t kphi0plusPi = phi0+TMath::Pi();
  CalculateI0I1(integral01,integral11,b,x0,y0,phi0,ellCut);
  CalculateI0I1(integral02,integral12,b,x0,y0,kphi0plusPi,ellCut);
  return;
}

void AliFastGlauber::GetI0I1BackToBackAndPhiAndXY(Double_t& integral01,Double_t& integral11,
					          Double_t& integral02,Double_t& integral12,
						  Double_t& phi,Double_t &x,Double_t &y,
						  Double_t ellCut,Double_t b)
{
  //
  // Return 2 pairs of I0 and I1 back to back from random b, x0, y0, phi0 
  // Return also phi0
  //
  Double_t x0,y0,phi0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  GetRandomPhi(phi0);
  phi = phi0; x=x0; y=y0;
  const Double_t kphi0plusPi = phi0+TMath::Pi();
  CalculateI0I1(integral01,integral11,b,x0,y0,phi0,ellCut);
  CalculateI0I1(integral02,integral12,b,x0,y0,kphi0plusPi,ellCut);
  return;
}

void AliFastGlauber::GetI0I1BackToBack(Double_t& integral01,Double_t& integral11,
				       Double_t& integral02,Double_t& integral12,
				       Double_t ellCut,Double_t b)
{
  //
  // Return 2 pairs of I0 and I1 back to back from random b, x0, y0, phi0 
  //
  Double_t phi;
  GetI0I1BackToBackAndPhi(integral01,integral11,integral02,integral12,
			  phi,ellCut,b);
  return;
}

void AliFastGlauber::GetI0I1ForPythia(Int_t n,Double_t* phi,
				      Double_t* integral0,Double_t* integral1,
				      Double_t ellCut,Double_t b)
{
  //
  // Returns I0 and I1 pairs for n partons with azimuthal angles phi[n] 
  // from random b, x0, y0
  //
  Double_t x0,y0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  for(Int_t i=0; i<n; i++) 
    CalculateI0I1(integral0[i],integral1[i],b,x0,y0,phi[i],ellCut);
  return;
}

void AliFastGlauber::GetI0I1ForPythiaAndXY(Int_t n,Double_t* phi,
				      Double_t* integral0,Double_t* integral1,
				      Double_t &x,Double_t& y,
				      Double_t ellCut,Double_t b)
{
  //
  // Returns I0 and I1 pairs for n partons with azimuthal angles phi[n] 
  // from random b, x0, y0 and return x0,y0
  //
  Double_t x0,y0;
  if(b<0.) GetRandomBHard(b);
  GetRandomXY(x0,y0);
  for(Int_t i=0; i<n; i++) 
    CalculateI0I1(integral0[i],integral1[i],b,x0,y0,phi[i],ellCut);
  x=x0;
  y=y0;
  return;
}

void AliFastGlauber::PlotI0I1Distr(Int_t n,Double_t ellCut,
				   Bool_t save,const char *fname)
{
  //
  // Plot I0-I1 distribution
  //
  Double_t i0,i1;
  TH2F *hI0I1s = new TH2F("hI0I1s","I_{0} versus I_{1}",1000,0,0.001,1000,0,0.01); 
  hI0I1s->SetXTitle("I_{0} [fm^{-3}]");
  hI0I1s->SetYTitle("I_{1} [fm^{-2}]");

  TH1F *hI0 = new TH1F("hI0","I_{0} = #hat{q}L / k",
		       1000,0,0.001); 
  hI0->SetXTitle("I_{0} [fm^{-3}]");
  hI0->SetYTitle("Probability");
  hI0->SetFillColor(3);
  TH1F *hI1 = new TH1F("hI1","I_{1} = #omega_{c} / k",
		       1000,0,0.01); 
  hI1->SetXTitle("I_{1} [fm^{-2}]");
  hI1->SetYTitle("Probability");
  hI1->SetFillColor(4);
  TH1F *h2 = new TH1F("h2","2 I_{1}^{2}/I_{0} = R / k",
		      100,0,0.02); 
  h2->SetXTitle("2 I_{1}^{2}/I_{0} [fm^{-1}]");
  h2->SetYTitle("Probability");
  h2->SetFillColor(2);
  TH1F *h3 = new TH1F("h3","2 I_{1}/I_{0} = L",
		      100,0,15); 
  h3->SetXTitle("2 I_{1}/I_{0} [fm]");
  h3->SetYTitle("Probability");
  h3->SetFillColor(5);
  TH1F *h4 = new TH1F("h4","I_{0}^{2}/(2 I_{1}) = #hat{q} / k",
		      100,0,0.00015); 
  h4->SetXTitle("I_{0}^{2}/(2 I_{1}) [fm^{-4}]");
  h4->SetYTitle("Probability");
  h4->SetFillColor(7);

  for(Int_t i=0; i<n; i++) {
    GetI0I1(i0,i1,ellCut);
    hI0I1s->Fill(i0,i1);
    hI0->Fill(i0);
    hI1->Fill(i1);
    h2->Fill(2.*i1*i1/i0);
    h3->Fill(2.*i1/i0);
    h4->Fill(i0*i0/2./i1);
  }
  hI0->Scale(1/(Double_t)n);
  hI1->Scale(1/(Double_t)n);
  h2->Scale(1/(Double_t)n);
  h3->Scale(1/(Double_t)n);
  h4->Scale(1/(Double_t)n);
  hI0I1s->Scale(1/(Double_t)n);

  TCanvas *cI0I1 = new TCanvas("cI0I1","I0 and I1",0,0,900,700);
  cI0I1->Divide(3,2);
  cI0I1->cd(1);
  hI0->Draw();
  cI0I1->cd(2);
  hI1->Draw();
  cI0I1->cd(3);
  h2->Draw();
  cI0I1->cd(4);
  h3->Draw();
  cI0I1->cd(5);
  h4->Draw();
  cI0I1->cd(6);
  gStyle->SetPalette(1,0);
  hI0I1s->Draw("col,Z");

  if(save) {
    TFile *f = new TFile(fname,"recreate");
    hI0I1s->Write();
    hI0->Write();
    hI1->Write();
    h2->Write();
    h3->Write();
    h4->Write();
    f->Close();
  }
  return;
}

void AliFastGlauber::PlotI0I1B2BDistr(Int_t n,Double_t ellCut,
				      Bool_t save,const char *fname)
{
  //
  // Plot I0-I1 back-to-back distributions
  //
  Double_t i01,i11,i02,i12;
  TH2F *hI0s = new TH2F("hI0s","I_{0}'s back-to-back",100,0,100,100,0,100); 
  hI0s->SetXTitle("I_{0} [fm^{-3}]");
  hI0s->SetYTitle("I_{0} [fm^{-3}]");
  TH2F *hI1s = new TH2F("hI1s","I_{1}'s back-to-back",100,0,100,100,0,100); 
  hI1s->SetXTitle("I_{1} [fm^{-2}]");
  hI1s->SetYTitle("I_{1} [fm^{-2}]");

  for(Int_t i=0; i<n; i++) {
    GetI0I1BackToBack(i01,i11,i02,i12,ellCut);
    hI0s->Fill(i01,i02);
    hI1s->Fill(i11,i12);
  }
  hI0s->Scale(1/(Double_t)n);
  hI1s->Scale(1/(Double_t)n);

  TCanvas *cI0I1s = new TCanvas("cI0I1s","I0 and I1 back-to-back distributions",0,0,800,400);
  gStyle->SetPalette(1,0);
  cI0I1s->Divide(2,1);
  cI0I1s->cd(1);
  hI0s->Draw("col,Z");
  cI0I1s->cd(2);
  hI1s->Draw("col,Z");

  if(save) {
    TFile *f = new TFile(fname,"recreate");
    hI0s->Write();
    hI1s->Write();
    f->Close();
  }
  return;
}

AliFastGlauber& AliFastGlauber::operator=(const  AliFastGlauber& rhs)
{
// Assignment operator
    rhs.Copy(*this);
    return *this;
}

void AliFastGlauber::Copy(TObject&) const
{
    //
    // Copy 
    //
    Fatal("Copy","Not implemented!\n");
}

 AliFastGlauber.cxx:1
 AliFastGlauber.cxx:2
 AliFastGlauber.cxx:3
 AliFastGlauber.cxx:4
 AliFastGlauber.cxx:5
 AliFastGlauber.cxx:6
 AliFastGlauber.cxx:7
 AliFastGlauber.cxx:8
 AliFastGlauber.cxx:9
 AliFastGlauber.cxx:10
 AliFastGlauber.cxx:11
 AliFastGlauber.cxx:12
 AliFastGlauber.cxx:13
 AliFastGlauber.cxx:14
 AliFastGlauber.cxx:15
 AliFastGlauber.cxx:16
 AliFastGlauber.cxx:17
 AliFastGlauber.cxx:18
 AliFastGlauber.cxx:19
 AliFastGlauber.cxx:20
 AliFastGlauber.cxx:21
 AliFastGlauber.cxx:22
 AliFastGlauber.cxx:23
 AliFastGlauber.cxx:24
 AliFastGlauber.cxx:25
 AliFastGlauber.cxx:26
 AliFastGlauber.cxx:27
 AliFastGlauber.cxx:28
 AliFastGlauber.cxx:29
 AliFastGlauber.cxx:30
 AliFastGlauber.cxx:31
 AliFastGlauber.cxx:32
 AliFastGlauber.cxx:33
 AliFastGlauber.cxx:34
 AliFastGlauber.cxx:35
 AliFastGlauber.cxx:36
 AliFastGlauber.cxx:37
 AliFastGlauber.cxx:38
 AliFastGlauber.cxx:39
 AliFastGlauber.cxx:40
 AliFastGlauber.cxx:41
 AliFastGlauber.cxx:42
 AliFastGlauber.cxx:43
 AliFastGlauber.cxx:44
 AliFastGlauber.cxx:45
 AliFastGlauber.cxx:46
 AliFastGlauber.cxx:47
 AliFastGlauber.cxx:48
 AliFastGlauber.cxx:49
 AliFastGlauber.cxx:50
 AliFastGlauber.cxx:51
 AliFastGlauber.cxx:52
 AliFastGlauber.cxx:53
 AliFastGlauber.cxx:54
 AliFastGlauber.cxx:55
 AliFastGlauber.cxx:56
 AliFastGlauber.cxx:57
 AliFastGlauber.cxx:58
 AliFastGlauber.cxx:59
 AliFastGlauber.cxx:60
 AliFastGlauber.cxx:61
 AliFastGlauber.cxx:62
 AliFastGlauber.cxx:63
 AliFastGlauber.cxx:64
 AliFastGlauber.cxx:65
 AliFastGlauber.cxx:66
 AliFastGlauber.cxx:67
 AliFastGlauber.cxx:68
 AliFastGlauber.cxx:69
 AliFastGlauber.cxx:70
 AliFastGlauber.cxx:71
 AliFastGlauber.cxx:72
 AliFastGlauber.cxx:73
 AliFastGlauber.cxx:74
 AliFastGlauber.cxx:75
 AliFastGlauber.cxx:76
 AliFastGlauber.cxx:77
 AliFastGlauber.cxx:78
 AliFastGlauber.cxx:79
 AliFastGlauber.cxx:80
 AliFastGlauber.cxx:81
 AliFastGlauber.cxx:82
 AliFastGlauber.cxx:83
 AliFastGlauber.cxx:84
 AliFastGlauber.cxx:85
 AliFastGlauber.cxx:86
 AliFastGlauber.cxx:87
 AliFastGlauber.cxx:88
 AliFastGlauber.cxx:89
 AliFastGlauber.cxx:90
 AliFastGlauber.cxx:91
 AliFastGlauber.cxx:92
 AliFastGlauber.cxx:93
 AliFastGlauber.cxx:94
 AliFastGlauber.cxx:95
 AliFastGlauber.cxx:96
 AliFastGlauber.cxx:97
 AliFastGlauber.cxx:98
 AliFastGlauber.cxx:99
 AliFastGlauber.cxx:100
 AliFastGlauber.cxx:101
 AliFastGlauber.cxx:102
 AliFastGlauber.cxx:103
 AliFastGlauber.cxx:104
 AliFastGlauber.cxx:105
 AliFastGlauber.cxx:106
 AliFastGlauber.cxx:107
 AliFastGlauber.cxx:108
 AliFastGlauber.cxx:109
 AliFastGlauber.cxx:110
 AliFastGlauber.cxx:111
 AliFastGlauber.cxx:112
 AliFastGlauber.cxx:113
 AliFastGlauber.cxx:114
 AliFastGlauber.cxx:115
 AliFastGlauber.cxx:116
 AliFastGlauber.cxx:117
 AliFastGlauber.cxx:118
 AliFastGlauber.cxx:119
 AliFastGlauber.cxx:120
 AliFastGlauber.cxx:121
 AliFastGlauber.cxx:122
 AliFastGlauber.cxx:123
 AliFastGlauber.cxx:124
 AliFastGlauber.cxx:125
 AliFastGlauber.cxx:126
 AliFastGlauber.cxx:127
 AliFastGlauber.cxx:128
 AliFastGlauber.cxx:129
 AliFastGlauber.cxx:130
 AliFastGlauber.cxx:131
 AliFastGlauber.cxx:132
 AliFastGlauber.cxx:133
 AliFastGlauber.cxx:134
 AliFastGlauber.cxx:135
 AliFastGlauber.cxx:136
 AliFastGlauber.cxx:137
 AliFastGlauber.cxx:138
 AliFastGlauber.cxx:139
 AliFastGlauber.cxx:140
 AliFastGlauber.cxx:141
 AliFastGlauber.cxx:142
 AliFastGlauber.cxx:143
 AliFastGlauber.cxx:144
 AliFastGlauber.cxx:145
 AliFastGlauber.cxx:146
 AliFastGlauber.cxx:147
 AliFastGlauber.cxx:148
 AliFastGlauber.cxx:149
 AliFastGlauber.cxx:150
 AliFastGlauber.cxx:151
 AliFastGlauber.cxx:152
 AliFastGlauber.cxx:153
 AliFastGlauber.cxx:154
 AliFastGlauber.cxx:155
 AliFastGlauber.cxx:156
 AliFastGlauber.cxx:157
 AliFastGlauber.cxx:158
 AliFastGlauber.cxx:159
 AliFastGlauber.cxx:160
 AliFastGlauber.cxx:161
 AliFastGlauber.cxx:162
 AliFastGlauber.cxx:163
 AliFastGlauber.cxx:164
 AliFastGlauber.cxx:165
 AliFastGlauber.cxx:166
 AliFastGlauber.cxx:167
 AliFastGlauber.cxx:168
 AliFastGlauber.cxx:169
 AliFastGlauber.cxx:170
 AliFastGlauber.cxx:171
 AliFastGlauber.cxx:172
 AliFastGlauber.cxx:173
 AliFastGlauber.cxx:174
 AliFastGlauber.cxx:175
 AliFastGlauber.cxx:176
 AliFastGlauber.cxx:177
 AliFastGlauber.cxx:178
 AliFastGlauber.cxx:179
 AliFastGlauber.cxx:180
 AliFastGlauber.cxx:181
 AliFastGlauber.cxx:182
 AliFastGlauber.cxx:183
 AliFastGlauber.cxx:184
 AliFastGlauber.cxx:185
 AliFastGlauber.cxx:186
 AliFastGlauber.cxx:187
 AliFastGlauber.cxx:188
 AliFastGlauber.cxx:189
 AliFastGlauber.cxx:190
 AliFastGlauber.cxx:191
 AliFastGlauber.cxx:192
 AliFastGlauber.cxx:193
 AliFastGlauber.cxx:194
 AliFastGlauber.cxx:195
 AliFastGlauber.cxx:196
 AliFastGlauber.cxx:197
 AliFastGlauber.cxx:198
 AliFastGlauber.cxx:199
 AliFastGlauber.cxx:200
 AliFastGlauber.cxx:201
 AliFastGlauber.cxx:202
 AliFastGlauber.cxx:203
 AliFastGlauber.cxx:204
 AliFastGlauber.cxx:205
 AliFastGlauber.cxx:206
 AliFastGlauber.cxx:207
 AliFastGlauber.cxx:208
 AliFastGlauber.cxx:209
 AliFastGlauber.cxx:210
 AliFastGlauber.cxx:211
 AliFastGlauber.cxx:212
 AliFastGlauber.cxx:213
 AliFastGlauber.cxx:214
 AliFastGlauber.cxx:215
 AliFastGlauber.cxx:216
 AliFastGlauber.cxx:217
 AliFastGlauber.cxx:218
 AliFastGlauber.cxx:219
 AliFastGlauber.cxx:220
 AliFastGlauber.cxx:221
 AliFastGlauber.cxx:222
 AliFastGlauber.cxx:223
 AliFastGlauber.cxx:224
 AliFastGlauber.cxx:225
 AliFastGlauber.cxx:226
 AliFastGlauber.cxx:227
 AliFastGlauber.cxx:228
 AliFastGlauber.cxx:229
 AliFastGlauber.cxx:230
 AliFastGlauber.cxx:231
 AliFastGlauber.cxx:232
 AliFastGlauber.cxx:233
 AliFastGlauber.cxx:234
 AliFastGlauber.cxx:235
 AliFastGlauber.cxx:236
 AliFastGlauber.cxx:237
 AliFastGlauber.cxx:238
 AliFastGlauber.cxx:239
 AliFastGlauber.cxx:240
 AliFastGlauber.cxx:241
 AliFastGlauber.cxx:242
 AliFastGlauber.cxx:243
 AliFastGlauber.cxx:244
 AliFastGlauber.cxx:245
 AliFastGlauber.cxx:246
 AliFastGlauber.cxx:247
 AliFastGlauber.cxx:248
 AliFastGlauber.cxx:249
 AliFastGlauber.cxx:250
 AliFastGlauber.cxx:251
 AliFastGlauber.cxx:252
 AliFastGlauber.cxx:253
 AliFastGlauber.cxx:254
 AliFastGlauber.cxx:255
 AliFastGlauber.cxx:256
 AliFastGlauber.cxx:257
 AliFastGlauber.cxx:258
 AliFastGlauber.cxx:259
 AliFastGlauber.cxx:260
 AliFastGlauber.cxx:261
 AliFastGlauber.cxx:262
 AliFastGlauber.cxx:263
 AliFastGlauber.cxx:264
 AliFastGlauber.cxx:265
 AliFastGlauber.cxx:266
 AliFastGlauber.cxx:267
 AliFastGlauber.cxx:268
 AliFastGlauber.cxx:269
 AliFastGlauber.cxx:270
 AliFastGlauber.cxx:271
 AliFastGlauber.cxx:272
 AliFastGlauber.cxx:273
 AliFastGlauber.cxx:274
 AliFastGlauber.cxx:275
 AliFastGlauber.cxx:276
 AliFastGlauber.cxx:277
 AliFastGlauber.cxx:278
 AliFastGlauber.cxx:279
 AliFastGlauber.cxx:280
 AliFastGlauber.cxx:281
 AliFastGlauber.cxx:282
 AliFastGlauber.cxx:283
 AliFastGlauber.cxx:284
 AliFastGlauber.cxx:285
 AliFastGlauber.cxx:286
 AliFastGlauber.cxx:287
 AliFastGlauber.cxx:288
 AliFastGlauber.cxx:289
 AliFastGlauber.cxx:290
 AliFastGlauber.cxx:291
 AliFastGlauber.cxx:292
 AliFastGlauber.cxx:293
 AliFastGlauber.cxx:294
 AliFastGlauber.cxx:295
 AliFastGlauber.cxx:296
 AliFastGlauber.cxx:297
 AliFastGlauber.cxx:298
 AliFastGlauber.cxx:299
 AliFastGlauber.cxx:300
 AliFastGlauber.cxx:301
 AliFastGlauber.cxx:302
 AliFastGlauber.cxx:303
 AliFastGlauber.cxx:304
 AliFastGlauber.cxx:305
 AliFastGlauber.cxx:306
 AliFastGlauber.cxx:307
 AliFastGlauber.cxx:308
 AliFastGlauber.cxx:309
 AliFastGlauber.cxx:310
 AliFastGlauber.cxx:311
 AliFastGlauber.cxx:312
 AliFastGlauber.cxx:313
 AliFastGlauber.cxx:314
 AliFastGlauber.cxx:315
 AliFastGlauber.cxx:316
 AliFastGlauber.cxx:317
 AliFastGlauber.cxx:318
 AliFastGlauber.cxx:319
 AliFastGlauber.cxx:320
 AliFastGlauber.cxx:321
 AliFastGlauber.cxx:322
 AliFastGlauber.cxx:323
 AliFastGlauber.cxx:324
 AliFastGlauber.cxx:325
 AliFastGlauber.cxx:326
 AliFastGlauber.cxx:327
 AliFastGlauber.cxx:328
 AliFastGlauber.cxx:329
 AliFastGlauber.cxx:330
 AliFastGlauber.cxx:331
 AliFastGlauber.cxx:332
 AliFastGlauber.cxx:333
 AliFastGlauber.cxx:334
 AliFastGlauber.cxx:335
 AliFastGlauber.cxx:336
 AliFastGlauber.cxx:337
 AliFastGlauber.cxx:338
 AliFastGlauber.cxx:339
 AliFastGlauber.cxx:340
 AliFastGlauber.cxx:341
 AliFastGlauber.cxx:342
 AliFastGlauber.cxx:343
 AliFastGlauber.cxx:344
 AliFastGlauber.cxx:345
 AliFastGlauber.cxx:346
 AliFastGlauber.cxx:347
 AliFastGlauber.cxx:348
 AliFastGlauber.cxx:349
 AliFastGlauber.cxx:350
 AliFastGlauber.cxx:351
 AliFastGlauber.cxx:352
 AliFastGlauber.cxx:353
 AliFastGlauber.cxx:354
 AliFastGlauber.cxx:355
 AliFastGlauber.cxx:356
 AliFastGlauber.cxx:357
 AliFastGlauber.cxx:358
 AliFastGlauber.cxx:359
 AliFastGlauber.cxx:360
 AliFastGlauber.cxx:361
 AliFastGlauber.cxx:362
 AliFastGlauber.cxx:363
 AliFastGlauber.cxx:364
 AliFastGlauber.cxx:365
 AliFastGlauber.cxx:366
 AliFastGlauber.cxx:367
 AliFastGlauber.cxx:368
 AliFastGlauber.cxx:369
 AliFastGlauber.cxx:370
 AliFastGlauber.cxx:371
 AliFastGlauber.cxx:372
 AliFastGlauber.cxx:373
 AliFastGlauber.cxx:374
 AliFastGlauber.cxx:375
 AliFastGlauber.cxx:376
 AliFastGlauber.cxx:377
 AliFastGlauber.cxx:378
 AliFastGlauber.cxx:379
 AliFastGlauber.cxx:380
 AliFastGlauber.cxx:381
 AliFastGlauber.cxx:382
 AliFastGlauber.cxx:383
 AliFastGlauber.cxx:384
 AliFastGlauber.cxx:385
 AliFastGlauber.cxx:386
 AliFastGlauber.cxx:387
 AliFastGlauber.cxx:388
 AliFastGlauber.cxx:389
 AliFastGlauber.cxx:390
 AliFastGlauber.cxx:391
 AliFastGlauber.cxx:392
 AliFastGlauber.cxx:393
 AliFastGlauber.cxx:394
 AliFastGlauber.cxx:395
 AliFastGlauber.cxx:396
 AliFastGlauber.cxx:397
 AliFastGlauber.cxx:398
 AliFastGlauber.cxx:399
 AliFastGlauber.cxx:400
 AliFastGlauber.cxx:401
 AliFastGlauber.cxx:402
 AliFastGlauber.cxx:403
 AliFastGlauber.cxx:404
 AliFastGlauber.cxx:405
 AliFastGlauber.cxx:406
 AliFastGlauber.cxx:407
 AliFastGlauber.cxx:408
 AliFastGlauber.cxx:409
 AliFastGlauber.cxx:410
 AliFastGlauber.cxx:411
 AliFastGlauber.cxx:412
 AliFastGlauber.cxx:413
 AliFastGlauber.cxx:414
 AliFastGlauber.cxx:415
 AliFastGlauber.cxx:416
 AliFastGlauber.cxx:417
 AliFastGlauber.cxx:418
 AliFastGlauber.cxx:419
 AliFastGlauber.cxx:420
 AliFastGlauber.cxx:421
 AliFastGlauber.cxx:422
 AliFastGlauber.cxx:423
 AliFastGlauber.cxx:424
 AliFastGlauber.cxx:425
 AliFastGlauber.cxx:426
 AliFastGlauber.cxx:427
 AliFastGlauber.cxx:428
 AliFastGlauber.cxx:429
 AliFastGlauber.cxx:430
 AliFastGlauber.cxx:431
 AliFastGlauber.cxx:432
 AliFastGlauber.cxx:433
 AliFastGlauber.cxx:434
 AliFastGlauber.cxx:435
 AliFastGlauber.cxx:436
 AliFastGlauber.cxx:437
 AliFastGlauber.cxx:438
 AliFastGlauber.cxx:439
 AliFastGlauber.cxx:440
 AliFastGlauber.cxx:441
 AliFastGlauber.cxx:442
 AliFastGlauber.cxx:443
 AliFastGlauber.cxx:444
 AliFastGlauber.cxx:445
 AliFastGlauber.cxx:446
 AliFastGlauber.cxx:447
 AliFastGlauber.cxx:448
 AliFastGlauber.cxx:449
 AliFastGlauber.cxx:450
 AliFastGlauber.cxx:451
 AliFastGlauber.cxx:452
 AliFastGlauber.cxx:453
 AliFastGlauber.cxx:454
 AliFastGlauber.cxx:455
 AliFastGlauber.cxx:456
 AliFastGlauber.cxx:457
 AliFastGlauber.cxx:458
 AliFastGlauber.cxx:459
 AliFastGlauber.cxx:460
 AliFastGlauber.cxx:461
 AliFastGlauber.cxx:462
 AliFastGlauber.cxx:463
 AliFastGlauber.cxx:464
 AliFastGlauber.cxx:465
 AliFastGlauber.cxx:466
 AliFastGlauber.cxx:467
 AliFastGlauber.cxx:468
 AliFastGlauber.cxx:469
 AliFastGlauber.cxx:470
 AliFastGlauber.cxx:471
 AliFastGlauber.cxx:472
 AliFastGlauber.cxx:473
 AliFastGlauber.cxx:474
 AliFastGlauber.cxx:475
 AliFastGlauber.cxx:476
 AliFastGlauber.cxx:477
 AliFastGlauber.cxx:478
 AliFastGlauber.cxx:479
 AliFastGlauber.cxx:480
 AliFastGlauber.cxx:481
 AliFastGlauber.cxx:482
 AliFastGlauber.cxx:483
 AliFastGlauber.cxx:484
 AliFastGlauber.cxx:485
 AliFastGlauber.cxx:486
 AliFastGlauber.cxx:487
 AliFastGlauber.cxx:488
 AliFastGlauber.cxx:489
 AliFastGlauber.cxx:490
 AliFastGlauber.cxx:491
 AliFastGlauber.cxx:492
 AliFastGlauber.cxx:493
 AliFastGlauber.cxx:494
 AliFastGlauber.cxx:495
 AliFastGlauber.cxx:496
 AliFastGlauber.cxx:497
 AliFastGlauber.cxx:498
 AliFastGlauber.cxx:499
 AliFastGlauber.cxx:500
 AliFastGlauber.cxx:501
 AliFastGlauber.cxx:502
 AliFastGlauber.cxx:503
 AliFastGlauber.cxx:504
 AliFastGlauber.cxx:505
 AliFastGlauber.cxx:506
 AliFastGlauber.cxx:507
 AliFastGlauber.cxx:508
 AliFastGlauber.cxx:509
 AliFastGlauber.cxx:510
 AliFastGlauber.cxx:511
 AliFastGlauber.cxx:512
 AliFastGlauber.cxx:513
 AliFastGlauber.cxx:514
 AliFastGlauber.cxx:515
 AliFastGlauber.cxx:516
 AliFastGlauber.cxx:517
 AliFastGlauber.cxx:518
 AliFastGlauber.cxx:519
 AliFastGlauber.cxx:520
 AliFastGlauber.cxx:521
 AliFastGlauber.cxx:522
 AliFastGlauber.cxx:523
 AliFastGlauber.cxx:524
 AliFastGlauber.cxx:525
 AliFastGlauber.cxx:526
 AliFastGlauber.cxx:527
 AliFastGlauber.cxx:528
 AliFastGlauber.cxx:529
 AliFastGlauber.cxx:530
 AliFastGlauber.cxx:531
 AliFastGlauber.cxx:532
 AliFastGlauber.cxx:533
 AliFastGlauber.cxx:534
 AliFastGlauber.cxx:535
 AliFastGlauber.cxx:536
 AliFastGlauber.cxx:537
 AliFastGlauber.cxx:538
 AliFastGlauber.cxx:539
 AliFastGlauber.cxx:540
 AliFastGlauber.cxx:541
 AliFastGlauber.cxx:542
 AliFastGlauber.cxx:543
 AliFastGlauber.cxx:544
 AliFastGlauber.cxx:545
 AliFastGlauber.cxx:546
 AliFastGlauber.cxx:547
 AliFastGlauber.cxx:548
 AliFastGlauber.cxx:549
 AliFastGlauber.cxx:550
 AliFastGlauber.cxx:551
 AliFastGlauber.cxx:552
 AliFastGlauber.cxx:553
 AliFastGlauber.cxx:554
 AliFastGlauber.cxx:555
 AliFastGlauber.cxx:556
 AliFastGlauber.cxx:557
 AliFastGlauber.cxx:558
 AliFastGlauber.cxx:559
 AliFastGlauber.cxx:560
 AliFastGlauber.cxx:561
 AliFastGlauber.cxx:562
 AliFastGlauber.cxx:563
 AliFastGlauber.cxx:564
 AliFastGlauber.cxx:565
 AliFastGlauber.cxx:566
 AliFastGlauber.cxx:567
 AliFastGlauber.cxx:568
 AliFastGlauber.cxx:569
 AliFastGlauber.cxx:570
 AliFastGlauber.cxx:571
 AliFastGlauber.cxx:572
 AliFastGlauber.cxx:573
 AliFastGlauber.cxx:574
 AliFastGlauber.cxx:575
 AliFastGlauber.cxx:576
 AliFastGlauber.cxx:577
 AliFastGlauber.cxx:578
 AliFastGlauber.cxx:579
 AliFastGlauber.cxx:580
 AliFastGlauber.cxx:581
 AliFastGlauber.cxx:582
 AliFastGlauber.cxx:583
 AliFastGlauber.cxx:584
 AliFastGlauber.cxx:585
 AliFastGlauber.cxx:586
 AliFastGlauber.cxx:587
 AliFastGlauber.cxx:588
 AliFastGlauber.cxx:589
 AliFastGlauber.cxx:590
 AliFastGlauber.cxx:591
 AliFastGlauber.cxx:592
 AliFastGlauber.cxx:593
 AliFastGlauber.cxx:594
 AliFastGlauber.cxx:595
 AliFastGlauber.cxx:596
 AliFastGlauber.cxx:597
 AliFastGlauber.cxx:598
 AliFastGlauber.cxx:599
 AliFastGlauber.cxx:600
 AliFastGlauber.cxx:601
 AliFastGlauber.cxx:602
 AliFastGlauber.cxx:603
 AliFastGlauber.cxx:604
 AliFastGlauber.cxx:605
 AliFastGlauber.cxx:606
 AliFastGlauber.cxx:607
 AliFastGlauber.cxx:608
 AliFastGlauber.cxx:609
 AliFastGlauber.cxx:610
 AliFastGlauber.cxx:611
 AliFastGlauber.cxx:612
 AliFastGlauber.cxx:613
 AliFastGlauber.cxx:614
 AliFastGlauber.cxx:615
 AliFastGlauber.cxx:616
 AliFastGlauber.cxx:617
 AliFastGlauber.cxx:618
 AliFastGlauber.cxx:619
 AliFastGlauber.cxx:620
 AliFastGlauber.cxx:621
 AliFastGlauber.cxx:622
 AliFastGlauber.cxx:623
 AliFastGlauber.cxx:624
 AliFastGlauber.cxx:625
 AliFastGlauber.cxx:626
 AliFastGlauber.cxx:627
 AliFastGlauber.cxx:628
 AliFastGlauber.cxx:629
 AliFastGlauber.cxx:630
 AliFastGlauber.cxx:631
 AliFastGlauber.cxx:632
 AliFastGlauber.cxx:633
 AliFastGlauber.cxx:634
 AliFastGlauber.cxx:635
 AliFastGlauber.cxx:636
 AliFastGlauber.cxx:637
 AliFastGlauber.cxx:638
 AliFastGlauber.cxx:639
 AliFastGlauber.cxx:640
 AliFastGlauber.cxx:641
 AliFastGlauber.cxx:642
 AliFastGlauber.cxx:643
 AliFastGlauber.cxx:644
 AliFastGlauber.cxx:645
 AliFastGlauber.cxx:646
 AliFastGlauber.cxx:647
 AliFastGlauber.cxx:648
 AliFastGlauber.cxx:649
 AliFastGlauber.cxx:650
 AliFastGlauber.cxx:651
 AliFastGlauber.cxx:652
 AliFastGlauber.cxx:653
 AliFastGlauber.cxx:654
 AliFastGlauber.cxx:655
 AliFastGlauber.cxx:656
 AliFastGlauber.cxx:657
 AliFastGlauber.cxx:658
 AliFastGlauber.cxx:659
 AliFastGlauber.cxx:660
 AliFastGlauber.cxx:661
 AliFastGlauber.cxx:662
 AliFastGlauber.cxx:663
 AliFastGlauber.cxx:664
 AliFastGlauber.cxx:665
 AliFastGlauber.cxx:666
 AliFastGlauber.cxx:667
 AliFastGlauber.cxx:668
 AliFastGlauber.cxx:669
 AliFastGlauber.cxx:670
 AliFastGlauber.cxx:671
 AliFastGlauber.cxx:672
 AliFastGlauber.cxx:673
 AliFastGlauber.cxx:674
 AliFastGlauber.cxx:675
 AliFastGlauber.cxx:676
 AliFastGlauber.cxx:677
 AliFastGlauber.cxx:678
 AliFastGlauber.cxx:679
 AliFastGlauber.cxx:680
 AliFastGlauber.cxx:681
 AliFastGlauber.cxx:682
 AliFastGlauber.cxx:683
 AliFastGlauber.cxx:684
 AliFastGlauber.cxx:685
 AliFastGlauber.cxx:686
 AliFastGlauber.cxx:687
 AliFastGlauber.cxx:688
 AliFastGlauber.cxx:689
 AliFastGlauber.cxx:690
 AliFastGlauber.cxx:691
 AliFastGlauber.cxx:692
 AliFastGlauber.cxx:693
 AliFastGlauber.cxx:694
 AliFastGlauber.cxx:695
 AliFastGlauber.cxx:696
 AliFastGlauber.cxx:697
 AliFastGlauber.cxx:698
 AliFastGlauber.cxx:699
 AliFastGlauber.cxx:700
 AliFastGlauber.cxx:701
 AliFastGlauber.cxx:702
 AliFastGlauber.cxx:703
 AliFastGlauber.cxx:704
 AliFastGlauber.cxx:705
 AliFastGlauber.cxx:706
 AliFastGlauber.cxx:707
 AliFastGlauber.cxx:708
 AliFastGlauber.cxx:709
 AliFastGlauber.cxx:710
 AliFastGlauber.cxx:711
 AliFastGlauber.cxx:712
 AliFastGlauber.cxx:713
 AliFastGlauber.cxx:714
 AliFastGlauber.cxx:715
 AliFastGlauber.cxx:716
 AliFastGlauber.cxx:717
 AliFastGlauber.cxx:718
 AliFastGlauber.cxx:719
 AliFastGlauber.cxx:720
 AliFastGlauber.cxx:721
 AliFastGlauber.cxx:722
 AliFastGlauber.cxx:723
 AliFastGlauber.cxx:724
 AliFastGlauber.cxx:725
 AliFastGlauber.cxx:726
 AliFastGlauber.cxx:727
 AliFastGlauber.cxx:728
 AliFastGlauber.cxx:729
 AliFastGlauber.cxx:730
 AliFastGlauber.cxx:731
 AliFastGlauber.cxx:732
 AliFastGlauber.cxx:733
 AliFastGlauber.cxx:734
 AliFastGlauber.cxx:735
 AliFastGlauber.cxx:736
 AliFastGlauber.cxx:737
 AliFastGlauber.cxx:738
 AliFastGlauber.cxx:739
 AliFastGlauber.cxx:740
 AliFastGlauber.cxx:741
 AliFastGlauber.cxx:742
 AliFastGlauber.cxx:743
 AliFastGlauber.cxx:744
 AliFastGlauber.cxx:745
 AliFastGlauber.cxx:746
 AliFastGlauber.cxx:747
 AliFastGlauber.cxx:748
 AliFastGlauber.cxx:749
 AliFastGlauber.cxx:750
 AliFastGlauber.cxx:751
 AliFastGlauber.cxx:752
 AliFastGlauber.cxx:753
 AliFastGlauber.cxx:754
 AliFastGlauber.cxx:755
 AliFastGlauber.cxx:756
 AliFastGlauber.cxx:757
 AliFastGlauber.cxx:758
 AliFastGlauber.cxx:759
 AliFastGlauber.cxx:760
 AliFastGlauber.cxx:761
 AliFastGlauber.cxx:762
 AliFastGlauber.cxx:763
 AliFastGlauber.cxx:764
 AliFastGlauber.cxx:765
 AliFastGlauber.cxx:766
 AliFastGlauber.cxx:767
 AliFastGlauber.cxx:768
 AliFastGlauber.cxx:769
 AliFastGlauber.cxx:770
 AliFastGlauber.cxx:771
 AliFastGlauber.cxx:772
 AliFastGlauber.cxx:773
 AliFastGlauber.cxx:774
 AliFastGlauber.cxx:775
 AliFastGlauber.cxx:776
 AliFastGlauber.cxx:777
 AliFastGlauber.cxx:778
 AliFastGlauber.cxx:779
 AliFastGlauber.cxx:780
 AliFastGlauber.cxx:781
 AliFastGlauber.cxx:782
 AliFastGlauber.cxx:783
 AliFastGlauber.cxx:784
 AliFastGlauber.cxx:785
 AliFastGlauber.cxx:786
 AliFastGlauber.cxx:787
 AliFastGlauber.cxx:788
 AliFastGlauber.cxx:789
 AliFastGlauber.cxx:790
 AliFastGlauber.cxx:791
 AliFastGlauber.cxx:792
 AliFastGlauber.cxx:793
 AliFastGlauber.cxx:794
 AliFastGlauber.cxx:795
 AliFastGlauber.cxx:796
 AliFastGlauber.cxx:797
 AliFastGlauber.cxx:798
 AliFastGlauber.cxx:799
 AliFastGlauber.cxx:800
 AliFastGlauber.cxx:801
 AliFastGlauber.cxx:802
 AliFastGlauber.cxx:803
 AliFastGlauber.cxx:804
 AliFastGlauber.cxx:805
 AliFastGlauber.cxx:806
 AliFastGlauber.cxx:807
 AliFastGlauber.cxx:808
 AliFastGlauber.cxx:809
 AliFastGlauber.cxx:810
 AliFastGlauber.cxx:811
 AliFastGlauber.cxx:812
 AliFastGlauber.cxx:813
 AliFastGlauber.cxx:814
 AliFastGlauber.cxx:815
 AliFastGlauber.cxx:816
 AliFastGlauber.cxx:817
 AliFastGlauber.cxx:818
 AliFastGlauber.cxx:819
 AliFastGlauber.cxx:820
 AliFastGlauber.cxx:821
 AliFastGlauber.cxx:822
 AliFastGlauber.cxx:823
 AliFastGlauber.cxx:824
 AliFastGlauber.cxx:825
 AliFastGlauber.cxx:826
 AliFastGlauber.cxx:827
 AliFastGlauber.cxx:828
 AliFastGlauber.cxx:829
 AliFastGlauber.cxx:830
 AliFastGlauber.cxx:831
 AliFastGlauber.cxx:832
 AliFastGlauber.cxx:833
 AliFastGlauber.cxx:834
 AliFastGlauber.cxx:835
 AliFastGlauber.cxx:836
 AliFastGlauber.cxx:837
 AliFastGlauber.cxx:838
 AliFastGlauber.cxx:839
 AliFastGlauber.cxx:840
 AliFastGlauber.cxx:841
 AliFastGlauber.cxx:842
 AliFastGlauber.cxx:843
 AliFastGlauber.cxx:844
 AliFastGlauber.cxx:845
 AliFastGlauber.cxx:846
 AliFastGlauber.cxx:847
 AliFastGlauber.cxx:848
 AliFastGlauber.cxx:849
 AliFastGlauber.cxx:850
 AliFastGlauber.cxx:851
 AliFastGlauber.cxx:852
 AliFastGlauber.cxx:853
 AliFastGlauber.cxx:854
 AliFastGlauber.cxx:855
 AliFastGlauber.cxx:856
 AliFastGlauber.cxx:857
 AliFastGlauber.cxx:858
 AliFastGlauber.cxx:859
 AliFastGlauber.cxx:860
 AliFastGlauber.cxx:861
 AliFastGlauber.cxx:862
 AliFastGlauber.cxx:863
 AliFastGlauber.cxx:864
 AliFastGlauber.cxx:865
 AliFastGlauber.cxx:866
 AliFastGlauber.cxx:867
 AliFastGlauber.cxx:868
 AliFastGlauber.cxx:869
 AliFastGlauber.cxx:870
 AliFastGlauber.cxx:871
 AliFastGlauber.cxx:872
 AliFastGlauber.cxx:873
 AliFastGlauber.cxx:874
 AliFastGlauber.cxx:875
 AliFastGlauber.cxx:876
 AliFastGlauber.cxx:877
 AliFastGlauber.cxx:878
 AliFastGlauber.cxx:879
 AliFastGlauber.cxx:880
 AliFastGlauber.cxx:881
 AliFastGlauber.cxx:882
 AliFastGlauber.cxx:883
 AliFastGlauber.cxx:884
 AliFastGlauber.cxx:885
 AliFastGlauber.cxx:886
 AliFastGlauber.cxx:887
 AliFastGlauber.cxx:888
 AliFastGlauber.cxx:889
 AliFastGlauber.cxx:890
 AliFastGlauber.cxx:891
 AliFastGlauber.cxx:892
 AliFastGlauber.cxx:893
 AliFastGlauber.cxx:894
 AliFastGlauber.cxx:895
 AliFastGlauber.cxx:896
 AliFastGlauber.cxx:897
 AliFastGlauber.cxx:898
 AliFastGlauber.cxx:899
 AliFastGlauber.cxx:900
 AliFastGlauber.cxx:901
 AliFastGlauber.cxx:902
 AliFastGlauber.cxx:903
 AliFastGlauber.cxx:904
 AliFastGlauber.cxx:905
 AliFastGlauber.cxx:906
 AliFastGlauber.cxx:907
 AliFastGlauber.cxx:908
 AliFastGlauber.cxx:909
 AliFastGlauber.cxx:910
 AliFastGlauber.cxx:911
 AliFastGlauber.cxx:912
 AliFastGlauber.cxx:913
 AliFastGlauber.cxx:914
 AliFastGlauber.cxx:915
 AliFastGlauber.cxx:916
 AliFastGlauber.cxx:917
 AliFastGlauber.cxx:918
 AliFastGlauber.cxx:919
 AliFastGlauber.cxx:920
 AliFastGlauber.cxx:921
 AliFastGlauber.cxx:922
 AliFastGlauber.cxx:923
 AliFastGlauber.cxx:924
 AliFastGlauber.cxx:925
 AliFastGlauber.cxx:926
 AliFastGlauber.cxx:927
 AliFastGlauber.cxx:928
 AliFastGlauber.cxx:929
 AliFastGlauber.cxx:930
 AliFastGlauber.cxx:931
 AliFastGlauber.cxx:932
 AliFastGlauber.cxx:933
 AliFastGlauber.cxx:934
 AliFastGlauber.cxx:935
 AliFastGlauber.cxx:936
 AliFastGlauber.cxx:937
 AliFastGlauber.cxx:938
 AliFastGlauber.cxx:939
 AliFastGlauber.cxx:940
 AliFastGlauber.cxx:941
 AliFastGlauber.cxx:942
 AliFastGlauber.cxx:943
 AliFastGlauber.cxx:944
 AliFastGlauber.cxx:945
 AliFastGlauber.cxx:946
 AliFastGlauber.cxx:947
 AliFastGlauber.cxx:948
 AliFastGlauber.cxx:949
 AliFastGlauber.cxx:950
 AliFastGlauber.cxx:951
 AliFastGlauber.cxx:952
 AliFastGlauber.cxx:953
 AliFastGlauber.cxx:954
 AliFastGlauber.cxx:955
 AliFastGlauber.cxx:956
 AliFastGlauber.cxx:957
 AliFastGlauber.cxx:958
 AliFastGlauber.cxx:959
 AliFastGlauber.cxx:960
 AliFastGlauber.cxx:961
 AliFastGlauber.cxx:962
 AliFastGlauber.cxx:963
 AliFastGlauber.cxx:964
 AliFastGlauber.cxx:965
 AliFastGlauber.cxx:966
 AliFastGlauber.cxx:967
 AliFastGlauber.cxx:968
 AliFastGlauber.cxx:969
 AliFastGlauber.cxx:970
 AliFastGlauber.cxx:971
 AliFastGlauber.cxx:972
 AliFastGlauber.cxx:973
 AliFastGlauber.cxx:974
 AliFastGlauber.cxx:975
 AliFastGlauber.cxx:976
 AliFastGlauber.cxx:977
 AliFastGlauber.cxx:978
 AliFastGlauber.cxx:979
 AliFastGlauber.cxx:980
 AliFastGlauber.cxx:981
 AliFastGlauber.cxx:982
 AliFastGlauber.cxx:983
 AliFastGlauber.cxx:984
 AliFastGlauber.cxx:985
 AliFastGlauber.cxx:986
 AliFastGlauber.cxx:987
 AliFastGlauber.cxx:988
 AliFastGlauber.cxx:989
 AliFastGlauber.cxx:990
 AliFastGlauber.cxx:991
 AliFastGlauber.cxx:992
 AliFastGlauber.cxx:993
 AliFastGlauber.cxx:994
 AliFastGlauber.cxx:995
 AliFastGlauber.cxx:996
 AliFastGlauber.cxx:997
 AliFastGlauber.cxx:998
 AliFastGlauber.cxx:999
 AliFastGlauber.cxx:1000
 AliFastGlauber.cxx:1001
 AliFastGlauber.cxx:1002
 AliFastGlauber.cxx:1003
 AliFastGlauber.cxx:1004
 AliFastGlauber.cxx:1005
 AliFastGlauber.cxx:1006
 AliFastGlauber.cxx:1007
 AliFastGlauber.cxx:1008
 AliFastGlauber.cxx:1009
 AliFastGlauber.cxx:1010
 AliFastGlauber.cxx:1011
 AliFastGlauber.cxx:1012
 AliFastGlauber.cxx:1013
 AliFastGlauber.cxx:1014
 AliFastGlauber.cxx:1015
 AliFastGlauber.cxx:1016
 AliFastGlauber.cxx:1017
 AliFastGlauber.cxx:1018
 AliFastGlauber.cxx:1019
 AliFastGlauber.cxx:1020
 AliFastGlauber.cxx:1021
 AliFastGlauber.cxx:1022
 AliFastGlauber.cxx:1023
 AliFastGlauber.cxx:1024
 AliFastGlauber.cxx:1025
 AliFastGlauber.cxx:1026
 AliFastGlauber.cxx:1027
 AliFastGlauber.cxx:1028
 AliFastGlauber.cxx:1029
 AliFastGlauber.cxx:1030
 AliFastGlauber.cxx:1031
 AliFastGlauber.cxx:1032
 AliFastGlauber.cxx:1033
 AliFastGlauber.cxx:1034
 AliFastGlauber.cxx:1035
 AliFastGlauber.cxx:1036
 AliFastGlauber.cxx:1037
 AliFastGlauber.cxx:1038
 AliFastGlauber.cxx:1039
 AliFastGlauber.cxx:1040
 AliFastGlauber.cxx:1041
 AliFastGlauber.cxx:1042
 AliFastGlauber.cxx:1043
 AliFastGlauber.cxx:1044
 AliFastGlauber.cxx:1045
 AliFastGlauber.cxx:1046
 AliFastGlauber.cxx:1047
 AliFastGlauber.cxx:1048
 AliFastGlauber.cxx:1049
 AliFastGlauber.cxx:1050
 AliFastGlauber.cxx:1051
 AliFastGlauber.cxx:1052
 AliFastGlauber.cxx:1053
 AliFastGlauber.cxx:1054
 AliFastGlauber.cxx:1055
 AliFastGlauber.cxx:1056
 AliFastGlauber.cxx:1057
 AliFastGlauber.cxx:1058
 AliFastGlauber.cxx:1059
 AliFastGlauber.cxx:1060
 AliFastGlauber.cxx:1061
 AliFastGlauber.cxx:1062
 AliFastGlauber.cxx:1063
 AliFastGlauber.cxx:1064
 AliFastGlauber.cxx:1065
 AliFastGlauber.cxx:1066
 AliFastGlauber.cxx:1067
 AliFastGlauber.cxx:1068
 AliFastGlauber.cxx:1069
 AliFastGlauber.cxx:1070
 AliFastGlauber.cxx:1071
 AliFastGlauber.cxx:1072
 AliFastGlauber.cxx:1073
 AliFastGlauber.cxx:1074
 AliFastGlauber.cxx:1075
 AliFastGlauber.cxx:1076
 AliFastGlauber.cxx:1077
 AliFastGlauber.cxx:1078
 AliFastGlauber.cxx:1079
 AliFastGlauber.cxx:1080
 AliFastGlauber.cxx:1081
 AliFastGlauber.cxx:1082
 AliFastGlauber.cxx:1083
 AliFastGlauber.cxx:1084
 AliFastGlauber.cxx:1085
 AliFastGlauber.cxx:1086
 AliFastGlauber.cxx:1087
 AliFastGlauber.cxx:1088
 AliFastGlauber.cxx:1089
 AliFastGlauber.cxx:1090
 AliFastGlauber.cxx:1091
 AliFastGlauber.cxx:1092
 AliFastGlauber.cxx:1093
 AliFastGlauber.cxx:1094
 AliFastGlauber.cxx:1095
 AliFastGlauber.cxx:1096
 AliFastGlauber.cxx:1097
 AliFastGlauber.cxx:1098
 AliFastGlauber.cxx:1099
 AliFastGlauber.cxx:1100
 AliFastGlauber.cxx:1101
 AliFastGlauber.cxx:1102
 AliFastGlauber.cxx:1103
 AliFastGlauber.cxx:1104
 AliFastGlauber.cxx:1105
 AliFastGlauber.cxx:1106
 AliFastGlauber.cxx:1107
 AliFastGlauber.cxx:1108
 AliFastGlauber.cxx:1109
 AliFastGlauber.cxx:1110
 AliFastGlauber.cxx:1111
 AliFastGlauber.cxx:1112
 AliFastGlauber.cxx:1113
 AliFastGlauber.cxx:1114
 AliFastGlauber.cxx:1115
 AliFastGlauber.cxx:1116
 AliFastGlauber.cxx:1117
 AliFastGlauber.cxx:1118
 AliFastGlauber.cxx:1119
 AliFastGlauber.cxx:1120
 AliFastGlauber.cxx:1121
 AliFastGlauber.cxx:1122
 AliFastGlauber.cxx:1123
 AliFastGlauber.cxx:1124
 AliFastGlauber.cxx:1125
 AliFastGlauber.cxx:1126
 AliFastGlauber.cxx:1127
 AliFastGlauber.cxx:1128
 AliFastGlauber.cxx:1129
 AliFastGlauber.cxx:1130
 AliFastGlauber.cxx:1131
 AliFastGlauber.cxx:1132
 AliFastGlauber.cxx:1133
 AliFastGlauber.cxx:1134
 AliFastGlauber.cxx:1135
 AliFastGlauber.cxx:1136
 AliFastGlauber.cxx:1137
 AliFastGlauber.cxx:1138
 AliFastGlauber.cxx:1139
 AliFastGlauber.cxx:1140
 AliFastGlauber.cxx:1141
 AliFastGlauber.cxx:1142
 AliFastGlauber.cxx:1143
 AliFastGlauber.cxx:1144
 AliFastGlauber.cxx:1145
 AliFastGlauber.cxx:1146
 AliFastGlauber.cxx:1147
 AliFastGlauber.cxx:1148
 AliFastGlauber.cxx:1149
 AliFastGlauber.cxx:1150
 AliFastGlauber.cxx:1151
 AliFastGlauber.cxx:1152
 AliFastGlauber.cxx:1153
 AliFastGlauber.cxx:1154
 AliFastGlauber.cxx:1155
 AliFastGlauber.cxx:1156
 AliFastGlauber.cxx:1157
 AliFastGlauber.cxx:1158
 AliFastGlauber.cxx:1159
 AliFastGlauber.cxx:1160
 AliFastGlauber.cxx:1161
 AliFastGlauber.cxx:1162
 AliFastGlauber.cxx:1163
 AliFastGlauber.cxx:1164
 AliFastGlauber.cxx:1165
 AliFastGlauber.cxx:1166
 AliFastGlauber.cxx:1167
 AliFastGlauber.cxx:1168
 AliFastGlauber.cxx:1169
 AliFastGlauber.cxx:1170
 AliFastGlauber.cxx:1171
 AliFastGlauber.cxx:1172
 AliFastGlauber.cxx:1173
 AliFastGlauber.cxx:1174
 AliFastGlauber.cxx:1175
 AliFastGlauber.cxx:1176
 AliFastGlauber.cxx:1177
 AliFastGlauber.cxx:1178
 AliFastGlauber.cxx:1179
 AliFastGlauber.cxx:1180
 AliFastGlauber.cxx:1181
 AliFastGlauber.cxx:1182
 AliFastGlauber.cxx:1183
 AliFastGlauber.cxx:1184
 AliFastGlauber.cxx:1185
 AliFastGlauber.cxx:1186
 AliFastGlauber.cxx:1187
 AliFastGlauber.cxx:1188
 AliFastGlauber.cxx:1189
 AliFastGlauber.cxx:1190
 AliFastGlauber.cxx:1191
 AliFastGlauber.cxx:1192
 AliFastGlauber.cxx:1193
 AliFastGlauber.cxx:1194
 AliFastGlauber.cxx:1195
 AliFastGlauber.cxx:1196
 AliFastGlauber.cxx:1197
 AliFastGlauber.cxx:1198
 AliFastGlauber.cxx:1199
 AliFastGlauber.cxx:1200
 AliFastGlauber.cxx:1201
 AliFastGlauber.cxx:1202
 AliFastGlauber.cxx:1203
 AliFastGlauber.cxx:1204
 AliFastGlauber.cxx:1205
 AliFastGlauber.cxx:1206
 AliFastGlauber.cxx:1207
 AliFastGlauber.cxx:1208
 AliFastGlauber.cxx:1209
 AliFastGlauber.cxx:1210
 AliFastGlauber.cxx:1211
 AliFastGlauber.cxx:1212
 AliFastGlauber.cxx:1213
 AliFastGlauber.cxx:1214
 AliFastGlauber.cxx:1215
 AliFastGlauber.cxx:1216
 AliFastGlauber.cxx:1217
 AliFastGlauber.cxx:1218
 AliFastGlauber.cxx:1219
 AliFastGlauber.cxx:1220
 AliFastGlauber.cxx:1221
 AliFastGlauber.cxx:1222
 AliFastGlauber.cxx:1223
 AliFastGlauber.cxx:1224
 AliFastGlauber.cxx:1225
 AliFastGlauber.cxx:1226
 AliFastGlauber.cxx:1227
 AliFastGlauber.cxx:1228
 AliFastGlauber.cxx:1229
 AliFastGlauber.cxx:1230
 AliFastGlauber.cxx:1231
 AliFastGlauber.cxx:1232
 AliFastGlauber.cxx:1233
 AliFastGlauber.cxx:1234
 AliFastGlauber.cxx:1235
 AliFastGlauber.cxx:1236
 AliFastGlauber.cxx:1237
 AliFastGlauber.cxx:1238
 AliFastGlauber.cxx:1239
 AliFastGlauber.cxx:1240
 AliFastGlauber.cxx:1241
 AliFastGlauber.cxx:1242
 AliFastGlauber.cxx:1243
 AliFastGlauber.cxx:1244
 AliFastGlauber.cxx:1245
 AliFastGlauber.cxx:1246
 AliFastGlauber.cxx:1247
 AliFastGlauber.cxx:1248
 AliFastGlauber.cxx:1249
 AliFastGlauber.cxx:1250
 AliFastGlauber.cxx:1251
 AliFastGlauber.cxx:1252
 AliFastGlauber.cxx:1253
 AliFastGlauber.cxx:1254
 AliFastGlauber.cxx:1255
 AliFastGlauber.cxx:1256
 AliFastGlauber.cxx:1257
 AliFastGlauber.cxx:1258
 AliFastGlauber.cxx:1259
 AliFastGlauber.cxx:1260
 AliFastGlauber.cxx:1261
 AliFastGlauber.cxx:1262
 AliFastGlauber.cxx:1263
 AliFastGlauber.cxx:1264
 AliFastGlauber.cxx:1265
 AliFastGlauber.cxx:1266
 AliFastGlauber.cxx:1267
 AliFastGlauber.cxx:1268
 AliFastGlauber.cxx:1269
 AliFastGlauber.cxx:1270
 AliFastGlauber.cxx:1271
 AliFastGlauber.cxx:1272
 AliFastGlauber.cxx:1273
 AliFastGlauber.cxx:1274
 AliFastGlauber.cxx:1275
 AliFastGlauber.cxx:1276
 AliFastGlauber.cxx:1277
 AliFastGlauber.cxx:1278
 AliFastGlauber.cxx:1279
 AliFastGlauber.cxx:1280
 AliFastGlauber.cxx:1281
 AliFastGlauber.cxx:1282
 AliFastGlauber.cxx:1283
 AliFastGlauber.cxx:1284
 AliFastGlauber.cxx:1285
 AliFastGlauber.cxx:1286
 AliFastGlauber.cxx:1287
 AliFastGlauber.cxx:1288
 AliFastGlauber.cxx:1289
 AliFastGlauber.cxx:1290
 AliFastGlauber.cxx:1291
 AliFastGlauber.cxx:1292
 AliFastGlauber.cxx:1293
 AliFastGlauber.cxx:1294
 AliFastGlauber.cxx:1295
 AliFastGlauber.cxx:1296
 AliFastGlauber.cxx:1297
 AliFastGlauber.cxx:1298
 AliFastGlauber.cxx:1299
 AliFastGlauber.cxx:1300
 AliFastGlauber.cxx:1301
 AliFastGlauber.cxx:1302
 AliFastGlauber.cxx:1303
 AliFastGlauber.cxx:1304
 AliFastGlauber.cxx:1305
 AliFastGlauber.cxx:1306
 AliFastGlauber.cxx:1307
 AliFastGlauber.cxx:1308
 AliFastGlauber.cxx:1309
 AliFastGlauber.cxx:1310
 AliFastGlauber.cxx:1311
 AliFastGlauber.cxx:1312
 AliFastGlauber.cxx:1313
 AliFastGlauber.cxx:1314
 AliFastGlauber.cxx:1315
 AliFastGlauber.cxx:1316
 AliFastGlauber.cxx:1317
 AliFastGlauber.cxx:1318
 AliFastGlauber.cxx:1319
 AliFastGlauber.cxx:1320
 AliFastGlauber.cxx:1321
 AliFastGlauber.cxx:1322
 AliFastGlauber.cxx:1323
 AliFastGlauber.cxx:1324
 AliFastGlauber.cxx:1325
 AliFastGlauber.cxx:1326
 AliFastGlauber.cxx:1327
 AliFastGlauber.cxx:1328
 AliFastGlauber.cxx:1329
 AliFastGlauber.cxx:1330
 AliFastGlauber.cxx:1331
 AliFastGlauber.cxx:1332
 AliFastGlauber.cxx:1333
 AliFastGlauber.cxx:1334
 AliFastGlauber.cxx:1335
 AliFastGlauber.cxx:1336
 AliFastGlauber.cxx:1337
 AliFastGlauber.cxx:1338
 AliFastGlauber.cxx:1339
 AliFastGlauber.cxx:1340
 AliFastGlauber.cxx:1341
 AliFastGlauber.cxx:1342
 AliFastGlauber.cxx:1343
 AliFastGlauber.cxx:1344
 AliFastGlauber.cxx:1345
 AliFastGlauber.cxx:1346
 AliFastGlauber.cxx:1347
 AliFastGlauber.cxx:1348
 AliFastGlauber.cxx:1349
 AliFastGlauber.cxx:1350
 AliFastGlauber.cxx:1351
 AliFastGlauber.cxx:1352
 AliFastGlauber.cxx:1353
 AliFastGlauber.cxx:1354
 AliFastGlauber.cxx:1355
 AliFastGlauber.cxx:1356
 AliFastGlauber.cxx:1357
 AliFastGlauber.cxx:1358
 AliFastGlauber.cxx:1359
 AliFastGlauber.cxx:1360
 AliFastGlauber.cxx:1361
 AliFastGlauber.cxx:1362
 AliFastGlauber.cxx:1363
 AliFastGlauber.cxx:1364
 AliFastGlauber.cxx:1365
 AliFastGlauber.cxx:1366
 AliFastGlauber.cxx:1367
 AliFastGlauber.cxx:1368
 AliFastGlauber.cxx:1369
 AliFastGlauber.cxx:1370
 AliFastGlauber.cxx:1371
 AliFastGlauber.cxx:1372
 AliFastGlauber.cxx:1373
 AliFastGlauber.cxx:1374
 AliFastGlauber.cxx:1375
 AliFastGlauber.cxx:1376
 AliFastGlauber.cxx:1377
 AliFastGlauber.cxx:1378
 AliFastGlauber.cxx:1379
 AliFastGlauber.cxx:1380
 AliFastGlauber.cxx:1381
 AliFastGlauber.cxx:1382
 AliFastGlauber.cxx:1383
 AliFastGlauber.cxx:1384
 AliFastGlauber.cxx:1385
 AliFastGlauber.cxx:1386
 AliFastGlauber.cxx:1387
 AliFastGlauber.cxx:1388
 AliFastGlauber.cxx:1389
 AliFastGlauber.cxx:1390
 AliFastGlauber.cxx:1391
 AliFastGlauber.cxx:1392
 AliFastGlauber.cxx:1393
 AliFastGlauber.cxx:1394
 AliFastGlauber.cxx:1395
 AliFastGlauber.cxx:1396
 AliFastGlauber.cxx:1397
 AliFastGlauber.cxx:1398
 AliFastGlauber.cxx:1399
 AliFastGlauber.cxx:1400
 AliFastGlauber.cxx:1401
 AliFastGlauber.cxx:1402
 AliFastGlauber.cxx:1403
 AliFastGlauber.cxx:1404
 AliFastGlauber.cxx:1405
 AliFastGlauber.cxx:1406
 AliFastGlauber.cxx:1407
 AliFastGlauber.cxx:1408
 AliFastGlauber.cxx:1409
 AliFastGlauber.cxx:1410
 AliFastGlauber.cxx:1411
 AliFastGlauber.cxx:1412
 AliFastGlauber.cxx:1413
 AliFastGlauber.cxx:1414
 AliFastGlauber.cxx:1415
 AliFastGlauber.cxx:1416
 AliFastGlauber.cxx:1417
 AliFastGlauber.cxx:1418
 AliFastGlauber.cxx:1419
 AliFastGlauber.cxx:1420
 AliFastGlauber.cxx:1421
 AliFastGlauber.cxx:1422
 AliFastGlauber.cxx:1423
 AliFastGlauber.cxx:1424
 AliFastGlauber.cxx:1425
 AliFastGlauber.cxx:1426
 AliFastGlauber.cxx:1427
 AliFastGlauber.cxx:1428
 AliFastGlauber.cxx:1429
 AliFastGlauber.cxx:1430
 AliFastGlauber.cxx:1431
 AliFastGlauber.cxx:1432
 AliFastGlauber.cxx:1433
 AliFastGlauber.cxx:1434
 AliFastGlauber.cxx:1435
 AliFastGlauber.cxx:1436
 AliFastGlauber.cxx:1437
 AliFastGlauber.cxx:1438
 AliFastGlauber.cxx:1439
 AliFastGlauber.cxx:1440
 AliFastGlauber.cxx:1441
 AliFastGlauber.cxx:1442
 AliFastGlauber.cxx:1443
 AliFastGlauber.cxx:1444
 AliFastGlauber.cxx:1445
 AliFastGlauber.cxx:1446
 AliFastGlauber.cxx:1447
 AliFastGlauber.cxx:1448
 AliFastGlauber.cxx:1449
 AliFastGlauber.cxx:1450
 AliFastGlauber.cxx:1451
 AliFastGlauber.cxx:1452
 AliFastGlauber.cxx:1453
 AliFastGlauber.cxx:1454
 AliFastGlauber.cxx:1455
 AliFastGlauber.cxx:1456
 AliFastGlauber.cxx:1457
 AliFastGlauber.cxx:1458
 AliFastGlauber.cxx:1459
 AliFastGlauber.cxx:1460
 AliFastGlauber.cxx:1461
 AliFastGlauber.cxx:1462
 AliFastGlauber.cxx:1463
 AliFastGlauber.cxx:1464
 AliFastGlauber.cxx:1465
 AliFastGlauber.cxx:1466
 AliFastGlauber.cxx:1467
 AliFastGlauber.cxx:1468
 AliFastGlauber.cxx:1469
 AliFastGlauber.cxx:1470
 AliFastGlauber.cxx:1471
 AliFastGlauber.cxx:1472
 AliFastGlauber.cxx:1473
 AliFastGlauber.cxx:1474
 AliFastGlauber.cxx:1475
 AliFastGlauber.cxx:1476
 AliFastGlauber.cxx:1477
 AliFastGlauber.cxx:1478
 AliFastGlauber.cxx:1479
 AliFastGlauber.cxx:1480
 AliFastGlauber.cxx:1481
 AliFastGlauber.cxx:1482
 AliFastGlauber.cxx:1483
 AliFastGlauber.cxx:1484
 AliFastGlauber.cxx:1485
 AliFastGlauber.cxx:1486
 AliFastGlauber.cxx:1487
 AliFastGlauber.cxx:1488
 AliFastGlauber.cxx:1489
 AliFastGlauber.cxx:1490
 AliFastGlauber.cxx:1491
 AliFastGlauber.cxx:1492
 AliFastGlauber.cxx:1493
 AliFastGlauber.cxx:1494
 AliFastGlauber.cxx:1495
 AliFastGlauber.cxx:1496
 AliFastGlauber.cxx:1497
 AliFastGlauber.cxx:1498
 AliFastGlauber.cxx:1499
 AliFastGlauber.cxx:1500
 AliFastGlauber.cxx:1501
 AliFastGlauber.cxx:1502
 AliFastGlauber.cxx:1503
 AliFastGlauber.cxx:1504
 AliFastGlauber.cxx:1505
 AliFastGlauber.cxx:1506
 AliFastGlauber.cxx:1507
 AliFastGlauber.cxx:1508
 AliFastGlauber.cxx:1509
 AliFastGlauber.cxx:1510
 AliFastGlauber.cxx:1511
 AliFastGlauber.cxx:1512
 AliFastGlauber.cxx:1513
 AliFastGlauber.cxx:1514
 AliFastGlauber.cxx:1515
 AliFastGlauber.cxx:1516
 AliFastGlauber.cxx:1517
 AliFastGlauber.cxx:1518
 AliFastGlauber.cxx:1519
 AliFastGlauber.cxx:1520
 AliFastGlauber.cxx:1521
 AliFastGlauber.cxx:1522
 AliFastGlauber.cxx:1523
 AliFastGlauber.cxx:1524
 AliFastGlauber.cxx:1525
 AliFastGlauber.cxx:1526
 AliFastGlauber.cxx:1527
 AliFastGlauber.cxx:1528
 AliFastGlauber.cxx:1529
 AliFastGlauber.cxx:1530
 AliFastGlauber.cxx:1531
 AliFastGlauber.cxx:1532
 AliFastGlauber.cxx:1533
 AliFastGlauber.cxx:1534
 AliFastGlauber.cxx:1535
 AliFastGlauber.cxx:1536
 AliFastGlauber.cxx:1537
 AliFastGlauber.cxx:1538
 AliFastGlauber.cxx:1539
 AliFastGlauber.cxx:1540
 AliFastGlauber.cxx:1541
 AliFastGlauber.cxx:1542
 AliFastGlauber.cxx:1543
 AliFastGlauber.cxx:1544
 AliFastGlauber.cxx:1545
 AliFastGlauber.cxx:1546
 AliFastGlauber.cxx:1547
 AliFastGlauber.cxx:1548
 AliFastGlauber.cxx:1549
 AliFastGlauber.cxx:1550
 AliFastGlauber.cxx:1551
 AliFastGlauber.cxx:1552
 AliFastGlauber.cxx:1553
 AliFastGlauber.cxx:1554
 AliFastGlauber.cxx:1555
 AliFastGlauber.cxx:1556
 AliFastGlauber.cxx:1557
 AliFastGlauber.cxx:1558
 AliFastGlauber.cxx:1559
 AliFastGlauber.cxx:1560
 AliFastGlauber.cxx:1561
 AliFastGlauber.cxx:1562
 AliFastGlauber.cxx:1563
 AliFastGlauber.cxx:1564
 AliFastGlauber.cxx:1565
 AliFastGlauber.cxx:1566
 AliFastGlauber.cxx:1567
 AliFastGlauber.cxx:1568
 AliFastGlauber.cxx:1569
 AliFastGlauber.cxx:1570
 AliFastGlauber.cxx:1571
 AliFastGlauber.cxx:1572
 AliFastGlauber.cxx:1573
 AliFastGlauber.cxx:1574
 AliFastGlauber.cxx:1575
 AliFastGlauber.cxx:1576
 AliFastGlauber.cxx:1577
 AliFastGlauber.cxx:1578
 AliFastGlauber.cxx:1579
 AliFastGlauber.cxx:1580
 AliFastGlauber.cxx:1581
 AliFastGlauber.cxx:1582
 AliFastGlauber.cxx:1583
 AliFastGlauber.cxx:1584
 AliFastGlauber.cxx:1585
 AliFastGlauber.cxx:1586
 AliFastGlauber.cxx:1587
 AliFastGlauber.cxx:1588
 AliFastGlauber.cxx:1589
 AliFastGlauber.cxx:1590
 AliFastGlauber.cxx:1591
 AliFastGlauber.cxx:1592
 AliFastGlauber.cxx:1593
 AliFastGlauber.cxx:1594
 AliFastGlauber.cxx:1595
 AliFastGlauber.cxx:1596
 AliFastGlauber.cxx:1597
 AliFastGlauber.cxx:1598
 AliFastGlauber.cxx:1599
 AliFastGlauber.cxx:1600
 AliFastGlauber.cxx:1601
 AliFastGlauber.cxx:1602
 AliFastGlauber.cxx:1603
 AliFastGlauber.cxx:1604
 AliFastGlauber.cxx:1605
 AliFastGlauber.cxx:1606
 AliFastGlauber.cxx:1607
 AliFastGlauber.cxx:1608
 AliFastGlauber.cxx:1609
 AliFastGlauber.cxx:1610
 AliFastGlauber.cxx:1611
 AliFastGlauber.cxx:1612
 AliFastGlauber.cxx:1613
 AliFastGlauber.cxx:1614
 AliFastGlauber.cxx:1615
 AliFastGlauber.cxx:1616
 AliFastGlauber.cxx:1617
 AliFastGlauber.cxx:1618
 AliFastGlauber.cxx:1619
 AliFastGlauber.cxx:1620
 AliFastGlauber.cxx:1621
 AliFastGlauber.cxx:1622
 AliFastGlauber.cxx:1623
 AliFastGlauber.cxx:1624
 AliFastGlauber.cxx:1625
 AliFastGlauber.cxx:1626
 AliFastGlauber.cxx:1627
 AliFastGlauber.cxx:1628
 AliFastGlauber.cxx:1629
 AliFastGlauber.cxx:1630
 AliFastGlauber.cxx:1631
 AliFastGlauber.cxx:1632
 AliFastGlauber.cxx:1633
 AliFastGlauber.cxx:1634
 AliFastGlauber.cxx:1635
 AliFastGlauber.cxx:1636
 AliFastGlauber.cxx:1637
 AliFastGlauber.cxx:1638
 AliFastGlauber.cxx:1639
 AliFastGlauber.cxx:1640
 AliFastGlauber.cxx:1641
 AliFastGlauber.cxx:1642
 AliFastGlauber.cxx:1643
 AliFastGlauber.cxx:1644
 AliFastGlauber.cxx:1645
 AliFastGlauber.cxx:1646
 AliFastGlauber.cxx:1647
 AliFastGlauber.cxx:1648
 AliFastGlauber.cxx:1649
 AliFastGlauber.cxx:1650
 AliFastGlauber.cxx:1651
 AliFastGlauber.cxx:1652
 AliFastGlauber.cxx:1653
 AliFastGlauber.cxx:1654
 AliFastGlauber.cxx:1655
 AliFastGlauber.cxx:1656
 AliFastGlauber.cxx:1657
 AliFastGlauber.cxx:1658
 AliFastGlauber.cxx:1659
 AliFastGlauber.cxx:1660
 AliFastGlauber.cxx:1661
 AliFastGlauber.cxx:1662
 AliFastGlauber.cxx:1663
 AliFastGlauber.cxx:1664
 AliFastGlauber.cxx:1665
 AliFastGlauber.cxx:1666
 AliFastGlauber.cxx:1667
 AliFastGlauber.cxx:1668
 AliFastGlauber.cxx:1669
 AliFastGlauber.cxx:1670
 AliFastGlauber.cxx:1671
 AliFastGlauber.cxx:1672
 AliFastGlauber.cxx:1673
 AliFastGlauber.cxx:1674
 AliFastGlauber.cxx:1675
 AliFastGlauber.cxx:1676
 AliFastGlauber.cxx:1677
 AliFastGlauber.cxx:1678
 AliFastGlauber.cxx:1679
 AliFastGlauber.cxx:1680
 AliFastGlauber.cxx:1681
 AliFastGlauber.cxx:1682
 AliFastGlauber.cxx:1683
 AliFastGlauber.cxx:1684
 AliFastGlauber.cxx:1685
 AliFastGlauber.cxx:1686
 AliFastGlauber.cxx:1687
 AliFastGlauber.cxx:1688
 AliFastGlauber.cxx:1689
 AliFastGlauber.cxx:1690
 AliFastGlauber.cxx:1691
 AliFastGlauber.cxx:1692
 AliFastGlauber.cxx:1693
 AliFastGlauber.cxx:1694
 AliFastGlauber.cxx:1695
 AliFastGlauber.cxx:1696
 AliFastGlauber.cxx:1697
 AliFastGlauber.cxx:1698
 AliFastGlauber.cxx:1699
 AliFastGlauber.cxx:1700
 AliFastGlauber.cxx:1701
 AliFastGlauber.cxx:1702
 AliFastGlauber.cxx:1703
 AliFastGlauber.cxx:1704
 AliFastGlauber.cxx:1705
 AliFastGlauber.cxx:1706
 AliFastGlauber.cxx:1707
 AliFastGlauber.cxx:1708
 AliFastGlauber.cxx:1709
 AliFastGlauber.cxx:1710
 AliFastGlauber.cxx:1711
 AliFastGlauber.cxx:1712
 AliFastGlauber.cxx:1713
 AliFastGlauber.cxx:1714
 AliFastGlauber.cxx:1715
 AliFastGlauber.cxx:1716
 AliFastGlauber.cxx:1717
 AliFastGlauber.cxx:1718
 AliFastGlauber.cxx:1719
 AliFastGlauber.cxx:1720
 AliFastGlauber.cxx:1721
 AliFastGlauber.cxx:1722
 AliFastGlauber.cxx:1723
 AliFastGlauber.cxx:1724
 AliFastGlauber.cxx:1725
 AliFastGlauber.cxx:1726
 AliFastGlauber.cxx:1727
 AliFastGlauber.cxx:1728
 AliFastGlauber.cxx:1729
 AliFastGlauber.cxx:1730
 AliFastGlauber.cxx:1731
 AliFastGlauber.cxx:1732
 AliFastGlauber.cxx:1733
 AliFastGlauber.cxx:1734
 AliFastGlauber.cxx:1735
 AliFastGlauber.cxx:1736
 AliFastGlauber.cxx:1737
 AliFastGlauber.cxx:1738
 AliFastGlauber.cxx:1739
 AliFastGlauber.cxx:1740
 AliFastGlauber.cxx:1741
 AliFastGlauber.cxx:1742
 AliFastGlauber.cxx:1743
 AliFastGlauber.cxx:1744
 AliFastGlauber.cxx:1745
 AliFastGlauber.cxx:1746
 AliFastGlauber.cxx:1747
 AliFastGlauber.cxx:1748
 AliFastGlauber.cxx:1749
 AliFastGlauber.cxx:1750
 AliFastGlauber.cxx:1751
 AliFastGlauber.cxx:1752
 AliFastGlauber.cxx:1753
 AliFastGlauber.cxx:1754
 AliFastGlauber.cxx:1755
 AliFastGlauber.cxx:1756
 AliFastGlauber.cxx:1757
 AliFastGlauber.cxx:1758
 AliFastGlauber.cxx:1759
 AliFastGlauber.cxx:1760
 AliFastGlauber.cxx:1761
 AliFastGlauber.cxx:1762
 AliFastGlauber.cxx:1763
 AliFastGlauber.cxx:1764
 AliFastGlauber.cxx:1765
 AliFastGlauber.cxx:1766
 AliFastGlauber.cxx:1767
 AliFastGlauber.cxx:1768
 AliFastGlauber.cxx:1769
 AliFastGlauber.cxx:1770
 AliFastGlauber.cxx:1771
 AliFastGlauber.cxx:1772
 AliFastGlauber.cxx:1773
 AliFastGlauber.cxx:1774
 AliFastGlauber.cxx:1775
 AliFastGlauber.cxx:1776
 AliFastGlauber.cxx:1777
 AliFastGlauber.cxx:1778
 AliFastGlauber.cxx:1779
 AliFastGlauber.cxx:1780
 AliFastGlauber.cxx:1781
 AliFastGlauber.cxx:1782
 AliFastGlauber.cxx:1783
 AliFastGlauber.cxx:1784
 AliFastGlauber.cxx:1785
 AliFastGlauber.cxx:1786
 AliFastGlauber.cxx:1787
 AliFastGlauber.cxx:1788
 AliFastGlauber.cxx:1789
 AliFastGlauber.cxx:1790
 AliFastGlauber.cxx:1791
 AliFastGlauber.cxx:1792
 AliFastGlauber.cxx:1793
 AliFastGlauber.cxx:1794
 AliFastGlauber.cxx:1795
 AliFastGlauber.cxx:1796
 AliFastGlauber.cxx:1797
 AliFastGlauber.cxx:1798
 AliFastGlauber.cxx:1799
 AliFastGlauber.cxx:1800
 AliFastGlauber.cxx:1801
 AliFastGlauber.cxx:1802
 AliFastGlauber.cxx:1803
 AliFastGlauber.cxx:1804
 AliFastGlauber.cxx:1805
 AliFastGlauber.cxx:1806
 AliFastGlauber.cxx:1807
 AliFastGlauber.cxx:1808
 AliFastGlauber.cxx:1809
 AliFastGlauber.cxx:1810
 AliFastGlauber.cxx:1811
 AliFastGlauber.cxx:1812
 AliFastGlauber.cxx:1813
 AliFastGlauber.cxx:1814
 AliFastGlauber.cxx:1815
 AliFastGlauber.cxx:1816
 AliFastGlauber.cxx:1817
 AliFastGlauber.cxx:1818
 AliFastGlauber.cxx:1819
 AliFastGlauber.cxx:1820
 AliFastGlauber.cxx:1821
 AliFastGlauber.cxx:1822
 AliFastGlauber.cxx:1823
 AliFastGlauber.cxx:1824
 AliFastGlauber.cxx:1825
 AliFastGlauber.cxx:1826
 AliFastGlauber.cxx:1827
 AliFastGlauber.cxx:1828
 AliFastGlauber.cxx:1829
 AliFastGlauber.cxx:1830
 AliFastGlauber.cxx:1831
 AliFastGlauber.cxx:1832
 AliFastGlauber.cxx:1833
 AliFastGlauber.cxx:1834
 AliFastGlauber.cxx:1835
 AliFastGlauber.cxx:1836
 AliFastGlauber.cxx:1837
 AliFastGlauber.cxx:1838
 AliFastGlauber.cxx:1839
 AliFastGlauber.cxx:1840
 AliFastGlauber.cxx:1841
 AliFastGlauber.cxx:1842
 AliFastGlauber.cxx:1843
 AliFastGlauber.cxx:1844
 AliFastGlauber.cxx:1845
 AliFastGlauber.cxx:1846
 AliFastGlauber.cxx:1847
 AliFastGlauber.cxx:1848
 AliFastGlauber.cxx:1849
 AliFastGlauber.cxx:1850
 AliFastGlauber.cxx:1851
 AliFastGlauber.cxx:1852
 AliFastGlauber.cxx:1853
 AliFastGlauber.cxx:1854
 AliFastGlauber.cxx:1855
 AliFastGlauber.cxx:1856
 AliFastGlauber.cxx:1857
 AliFastGlauber.cxx:1858
 AliFastGlauber.cxx:1859
 AliFastGlauber.cxx:1860
 AliFastGlauber.cxx:1861
 AliFastGlauber.cxx:1862
 AliFastGlauber.cxx:1863
 AliFastGlauber.cxx:1864
 AliFastGlauber.cxx:1865
 AliFastGlauber.cxx:1866
 AliFastGlauber.cxx:1867
 AliFastGlauber.cxx:1868
 AliFastGlauber.cxx:1869
 AliFastGlauber.cxx:1870
 AliFastGlauber.cxx:1871
 AliFastGlauber.cxx:1872
 AliFastGlauber.cxx:1873
 AliFastGlauber.cxx:1874
 AliFastGlauber.cxx:1875
 AliFastGlauber.cxx:1876
 AliFastGlauber.cxx:1877
 AliFastGlauber.cxx:1878
 AliFastGlauber.cxx:1879
 AliFastGlauber.cxx:1880
 AliFastGlauber.cxx:1881
 AliFastGlauber.cxx:1882
 AliFastGlauber.cxx:1883
 AliFastGlauber.cxx:1884
 AliFastGlauber.cxx:1885