/* $Id$ */
//
// Utility class to make simple Glauber type calculations
-// for SYMMTRIC collision geometries (AA):
+// for SYMMETRIC collision geometries (AA):
// Impact parameter, production points, reaction plane dependence
//
// The SimulateTrigger method can be used for simple MB and hard-process
// from AliRoot
#include "AliFastGlauber.h"
// from root
-#include <TStyle.h>
-#include <TH1F.h>
-#include <TH2F.h>
+#include <TCanvas.h>
#include <TF1.h>
#include <TF2.h>
+#include <TFile.h>
+#include <TH1F.h>
+#include <TH2F.h>
#include <TLegend.h>
-#include <TCanvas.h>
+#include <TMath.h>
#include <TRandom.h>
-#include <TFile.h>
-#include <TROOT.h>
-#include <Riostream.h>
+#include <TStyle.h>
ClassImp(AliFastGlauber)
TF2* AliFastGlauber::fgWKParticipants = NULL;
TF1* AliFastGlauber::fgWParticipants = NULL;
TF2* AliFastGlauber::fgWAlmondCurrent = NULL;
-TF2 AliFastGlauber::fgWAlmondFixedB[40];
+TF2* AliFastGlauber::fgWAlmondFixedB[40];
const Int_t AliFastGlauber::fgkMCInts = 100000;
-Int_t AliFastGlauber::fgCounter = 0;
-
-AliFastGlauber::AliFastGlauber() : fName()
-{
- // Default Constructor
- //
- fgCounter++;
- if(fgCounter>1)
- Error("AliFastGlauber","More than more instance (%d) is not supported, check your code!",fgCounter);
-
+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()
{
- fgCounter--;
- //if(fgCounter==0) Reset();
+// Destructor
+ for(Int_t k=0; k<40; k++) delete fgWAlmondFixedB[k];
}
void AliFastGlauber::SetAuAuRhic()
fgWSb->SetParameter(2, fWSw);
fgWSb->SetParameter(3, fWSn);
- fgWSbz = new TF2("WSbz", WSbz, 0, fgBMax, 4);
+ fgWSbz = new TF2("WSbz", WSbz, 0, fgBMax, 0, fgBMax, 4);
fgWSbz->SetParameter(0, fWSr0);
fgWSbz->SetParameter(1, fWSd);
fgWSbz->SetParameter(2, fWSw);
for(Int_t k=0; k<40; k++) {
sprintf(almondName,"WAlmondFixedB%d",k);
fgWAlmondCurrent = (TF2*)ff->Get(almondName);
- new(&fgWAlmondFixedB[k]) TF2(*fgWAlmondCurrent);
+ fgWAlmondFixedB[k] = fgWAlmondCurrent;
}
delete ff;
}
fgWPathLength->SetParameter(2, 0); //Pathlength definition
}
-void AliFastGlauber::Reset()
+void AliFastGlauber::Reset() const
{
//
// Reset dynamic allocated formulas
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);
+ 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}]");
fgWIntRadius->Draw("same");
}
-Double_t AliFastGlauber::WSb(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WSb(const Double_t* x, const Double_t* par)
{
//
// Woods-Saxon Parameterisation
return y; //fm^-3
}
-Double_t AliFastGlauber::WSbz(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WSbz(const Double_t* x, const Double_t* par)
{
//
// Wood Saxon Parameterisation
return y; //fm^-3
}
-Double_t AliFastGlauber::WSz(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WSz(const Double_t* x, const Double_t* par)
{
//
// Wood Saxon Parameterisation
return y; //fm^-3
}
-Double_t AliFastGlauber::WSta(Double_t* x, Double_t* /*par*/)
+Double_t AliFastGlauber::WSta(const Double_t* x, const Double_t* /*par*/)
{
//
// Thickness function T_A
return y; //fm^-2
}
-Double_t AliFastGlauber::WStarfi(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WStarfi(const Double_t* x, const Double_t* par)
{
//
// Kernel for overlap function: T_A(s)*T_A(s-b)
return y; //fm^-3
}
-Double_t AliFastGlauber::WStaa(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WStaa(const Double_t* x, const Double_t* par)
{
//
// Overlap function
// 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);
return y;
}
-Double_t AliFastGlauber::WKParticipants(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WKParticipants(const Double_t* x, const Double_t* par)
{
//
// Kernel for number of participants
return y; //fm^-1
}
-Double_t AliFastGlauber::WParticipants(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WParticipants(const Double_t* x, const Double_t* par)
{
//
// Number of Participants as
return y; //no units
}
-Double_t AliFastGlauber::WSgeo(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WSgeo(const Double_t* x, const Double_t* par)
{
//
// Geometrical Cross-Section
return y; //fm
}
-Double_t AliFastGlauber::WSbinary(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WSbinary(const Double_t* x, const Double_t* par)
{
//
// Number of binary hard collisions
return y; //fm
}
-Double_t AliFastGlauber::WSN(Double_t* x, Double_t* /*par*/)
+Double_t AliFastGlauber::WSN(const Double_t* x, const Double_t* /*par*/)
{
//
// Number of hard processes per event
return y; //no units
}
-Double_t AliFastGlauber::WEnergyDensity(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WEnergyDensity(const Double_t* x, const Double_t* par)
{
//
// Initial energy density
return y; //fm^-4
}
-Double_t AliFastGlauber::WAlmond(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WAlmond(const Double_t* x, const Double_t* par)
{
//
// Almond shaped interaction region
return y; //fm^-4
}
-Double_t AliFastGlauber::WIntRadius(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WIntRadius(const Double_t* x, const Double_t* par)
{
//
// Average interaction density over radius
return y; //fm^-3
}
-Double_t AliFastGlauber::WPathLength0(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WPathLength0(const Double_t* x, const Double_t* par)
{
//
// Path Length as a function of phi
} // 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);
+ 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(Double_t* x, Double_t* par)
+Double_t AliFastGlauber::WPathLength(const Double_t* x, const Double_t* par)
{
//
// Path Length as a function of phi
return fgWSbinary->Integral(b1, b2)/fgWSbinary->Integral(0., 100.);
}
-Double_t AliFastGlauber::NHard(Double_t b1, Double_t b2) const
+Double_t AliFastGlauber::NHard(const Double_t b1, const Double_t b2) const
{
//
// Number of binary hard collisions
//
// Return number of binary hard collisions normalized to 1 at b=0
//
- if(b==0) b=1e-4;
+ 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==0) b=1e-4;
+ if(b<1.e-4) b=1e-4;
return fgWSN->Eval(b);
}
//
// Return the number of participants normalized to 1 at b=0
//
- if(b==0) b=1e-4;
+ if(b<1.e-4) b=1e-4;
return (fgWParticipants->Eval(b)/fgWParticipants->Eval(1e-4));
}
//
// Return the number of participants for impact parameter b
//
- if(b==0) b=1e-4;
+ if(b<1.e-4) b=1e-4;
return (fgWParticipants->Eval(b));
}
//
// Return the number of collisions for impact parameter b
//
- if(b==0) b=1e-4;
+ if(b<1.e-4) b=1e-4;
return (fgWStaa->Eval(b)*fSigmaNN);
}
b = fgWSbinary->GetRandom();
Int_t bin = 2*(Int_t)b;
if( (b-(Int_t)b) > 0.5) bin++;
- fgWAlmondCurrent = &fgWAlmondFixedB[bin];
+ fgWAlmondCurrent = fgWAlmondFixedB[bin];
return;
}
return;
}
-void AliFastGlauber::PlotLengthDistr(Int_t n,Bool_t save,Char_t *fname)
+void AliFastGlauber::PlotLengthDistr(Int_t n,Bool_t save,const char *fname)
{
//
// Plot length distribution
return;
}
-void AliFastGlauber::PlotLengthB2BDistr(Int_t n,Bool_t save,Char_t *fname)
+void AliFastGlauber::PlotLengthB2BDistr(Int_t n,Bool_t save,const char *fname)
{
//
// Plot lengths back-to-back distributions
gStyle->SetPalette(1,0);
c->Divide(2,2);
c->cd(1);
- fgWAlmondFixedB[0].Draw("cont1");
+ fgWAlmondFixedB[0]->Draw("cont1");
c->cd(2);
- fgWAlmondFixedB[10].Draw("cont1");
+ fgWAlmondFixedB[10]->Draw("cont1");
c->cd(3);
- fgWAlmondFixedB[20].Draw("cont1");
+ fgWAlmondFixedB[20]->Draw("cont1");
c->cd(4);
- fgWAlmondFixedB[30].Draw("cont1");
+ fgWAlmondFixedB[30]->Draw("cont1");
return;
}
}
void AliFastGlauber::PlotI0I1Distr(Int_t n,Double_t ellCut,
- Bool_t save,Char_t *fname)
+ Bool_t save,const char *fname)
{
//
// Plot I0-I1 distribution
}
void AliFastGlauber::PlotI0I1B2BDistr(Int_t n,Double_t ellCut,
- Bool_t save,Char_t *fname)
+ Bool_t save,const char *fname)
{
//
// Plot I0-I1 back-to-back distributions
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");
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff