/* $Id$ */
-//----------------------------------------------------------------------------
-// Implementation of the class to calculate the parton energy loss
-// Based on the "BDMPS" quenching weights by C.A.Salgado and U.A.Wiedemann
-//
-// References:
-// C.A.Salgado and U.A.Wiedemann, Phys.Rev.D68 (2003) 014008 [hep-ph/0302184]
-// A.Dainese, Eur.Phys.J.C, in press, [nucl-ex/0312005]
+
+// Implementation of the class to calculate the parton energy loss
+// Based on the "BDMPS" quenching weights by C.A.Salgado and U.A.Wiedemann
+// References:
+// C.A.Salgado and U.A.Wiedemann, Phys.Rev.D68 (2003) 014008 [hep-ph/0302184]
+// A.Dainese, Eur.Phys.J.C, in press, [nucl-ex/0312005]
//
//
// Origin: C. Loizides constantinos.loizides@cern.ch
//
//=================== Added by C. Loizides 27/03/04 ===========================
//
-// Added support for k-Quenching, where wc=I1*k and R=2I1^2/I0*k
-// (see the AliFastGlauber class for definition of I0/I1)
+// Added support for k-Quenching, where wc=I1*k and R=2I1^2/I0*k
+// (see the AliFastGlauber class for definition of I0/I1)
//-----------------------------------------------------------------------------
#include <Riostream.h>
#include <TGraph.h>
#include <TROOT.h>
#include <TSystem.h>
+#include <TString.h>
#include <TLegend.h>
#include "AliQuenchingWeights.h"
AliQuenchingWeights::AliQuenchingWeights()
- : TObject()
+ : TObject(),
+ fInstanceNumber(fgCounter++),
+ fMultSoft(kTRUE),
+ fECMethod(kDefault),
+ fQTransport(1.),
+ fMu(1.),
+ fK(4.e5),
+ fLengthMax(20),
+ fLengthMaxOld(0),
+ fHistos(0),
+ fHisto(0),
+ fTablesLoaded(kFALSE)
{
//default constructor
- fTablesLoaded=kFALSE;
- fMultSoft=kTRUE;
- fHistos=0;
- SetMu();
- SetQTransport();
- SetK();
- fECMethod=kDefault;
- //SetECMethod();
- SetLengthMax();
- fLengthMaxOld=0;
- fInstanceNumber=fgCounter++;
- Char_t name[100];
- sprintf(name,"hhistoqw_%d",fInstanceNumber);
- fHisto = new TH1F(name,"",fgkBins,0.,fgkMaxBin);
- for(Int_t bin=1;bin<=fgkBins;bin++)
- fHisto->SetBinContent(bin,0.);
}
AliQuenchingWeights::AliQuenchingWeights(const AliQuenchingWeights& a)
- : TObject()
+ : TObject(),
+ fInstanceNumber(fgCounter++),
+ fMultSoft(kTRUE),
+ fECMethod(kDefault),
+ fQTransport(1.),
+ fMu(1.),
+ fK(4.e5),
+ fLengthMax(20),
+ fLengthMaxOld(0),
+ fHistos(0),
+ fHisto(0),
+ fTablesLoaded(kFALSE)
{
// copy constructor
fLengthMax=a.GetLengthMax();
fInstanceNumber=fgCounter++;
Char_t name[100];
- sprintf(name,"hhistoqw_%d",fInstanceNumber);
+ snprintf(name,100, "hhistoqw_%d",fInstanceNumber);
fHisto = new TH1F(name,"",fgkBins,0.,fgkMaxBin);
for(Int_t bin=1;bin<=fgkBins;bin++)
fHisto->SetBinContent(bin,0.);
delete fHisto;
}
+void AliQuenchingWeights::Init()
+{
+// Initialization
+ if (fHisto) return;
+ fHisto = new TH1F(Form("hhistoqw_%d",fInstanceNumber),"",fgkBins,0.,fgkMaxBin);
+ for(Int_t bin=1;bin<=fgkBins;bin++)
+ fHisto->SetBinContent(bin,0.);
+}
+
void AliQuenchingWeights::Reset()
{
//reset tables if there were used
fECMethod=type;
if(fECMethod==kDefault)
Info("SetECMethod","Energy Constraint Method set to DEFAULT:\nIf (sampled energy loss > parton energy) then sampled energy loss = parton energy.");
- else
+ else if(fECMethod==kReweight)
Info("SetECMethod","Energy Constraint Method set to REWEIGHT:\nRequire sampled energy loss <= parton energy.");
+ else Info("SetECMethod","Energy Constraint Method set to REWEIGHTCONT:\nRequire sampled energy loss <= parton energy (only implemented for FAST method.");
}
Int_t AliQuenchingWeights::InitMult(const Char_t *contall,const Char_t *discall)
fMultSoft=kTRUE;
Char_t fname[1024];
- sprintf(fname,"%s",gSystem->ExpandPathName(contall));
+ snprintf(fname,1024, "%s",gSystem->ExpandPathName(contall));
//PH ifstream fincont(fname);
fstream fincont(fname,ios::in);
#if defined(__HP_aCC) || defined(__DECCXX)
}
fincont.close();
- sprintf(fname,"%s",gSystem->ExpandPathName(discall));
+ snprintf(fname,1024, "%s",gSystem->ExpandPathName(discall));
//PH ifstream findisc(fname);
fstream findisc(fname,ios::in);
#if defined(__HP_aCC) || defined(__DECCXX)
fMultSoft=kFALSE;
Char_t fname[1024];
- sprintf(fname,"%s",gSystem->ExpandPathName(contall));
+ snprintf(fname, 1024, "%s",gSystem->ExpandPathName(contall));
//PH ifstream fincont(fname);
fstream fincont(fname,ios::in);
#if defined(__HP_aCC) || defined(__DECCXX)
}
fincont.close();
- sprintf(fname,"%s",gSystem->ExpandPathName(discall));
+ snprintf(fname, 1024, "%s",gSystem->ExpandPathName(discall));
//PH ifstream findisc(fname);
fstream findisc(fname,ios::in);
#if defined(__HP_aCC) || defined(__DECCXX)
Error("CalcSingleHard","Tables are not loaded.");
return -1;
}
- if(!fMultSoft){
+ if(fMultSoft){
Error("CalcSingleHard","Tables are not loaded for Single Hard Scattering.");
return -1;
}
Double_t AliQuenchingWeights::CalcR(Double_t wc, Double_t l) const
{
- //calculate R value and
+ //calculate r value and
//check if it is less then maximum
- Double_t R = wc*l*fgkConvFmToInvGeV;
- if(R>=fgkRMax) {
- Warning("CalcR","Value of R = %.2f; should be less than %.2f",R,fgkRMax);
+ Double_t r = wc*l*fgkConvFmToInvGeV;
+ if(r >= fgkRMax) {
+ Warning("CalcR","Value of r = %.2f; should be less than %.2f", r, fgkRMax);
return fgkRMax-1;
}
- return R;
+ return r;
}
Double_t AliQuenchingWeights::CalcRk(Double_t k, Double_t I0, Double_t I1) const
//calculate R value and
//check if it is less then maximum
- Double_t R = fgkRMax-1;
+ Double_t r = fgkRMax-1;
if(I0>0)
- R = 2*I1*I1/I0*k;
- if(R>=fgkRMax) {
- Warning("CalcRk","Value of R = %.2f; should be less than %.2f",R,fgkRMax);
+ r = 2*I1*I1/I0*k;
+ if(r>=fgkRMax) {
+ Warning("CalcRk","Value of r = %.2f; should be less than %.2f",r,fgkRMax);
return fgkRMax-1;
}
- return R;
+ return r;
}
Double_t AliQuenchingWeights::GetELossRandom(Int_t ipart, Double_t length, Double_t e) const
ret=fHistos[ipart-1][l-1]->GetRandom();
if(++ws==1e6){
Warning("GetELossRandom",
- "Aborted reweighting; maximum loss assigned after 1e6 trials.");
+ "Stopping reweighting; maximum loss assigned after 1e6 trials.");
return e;
}
}
return -1000.;
}
- Double_t R=CalcRk(I0,I1);
- if(R<0.){
+ Double_t r=CalcRk(I0,I1);
+ if(r<0.){
Fatal("GetELossRandomK","R should not be negative");
return -1000.;
}
Fatal("GetELossRandomK","wc should be greater than zero");
return -1000.;
}
- if(SampleEnergyLoss(ipart,R)!=0){
+ if(SampleEnergyLoss(ipart,r)!=0){
Fatal("GetELossRandomK","Could not sample energy loss");
return -1000.;
}
ret=fHisto->GetRandom();
if(++ws==1e6){
Warning("GetELossRandomK",
- "Aborted reweighting; maximum loss assigned after 1e6 trials.");
+ "Stopping reweighting; maximum loss assigned after 1e6 trials.");
return e;
}
}
// all parameters are well within the bounds.
// read-in data tables before first call
- Double_t R=CalcRk(I0,I1);
- if(R<=0.){
+ Double_t r=CalcRk(I0,I1);
+ if(r<=0.){
return 0.;
}
return 0.;
}
- return GetELossRandomKFastR(ipart,R,wc,e);
+ return GetELossRandomKFastR(ipart,r,wc,e);
}
-Double_t AliQuenchingWeights::GetELossRandomKFastR(Int_t ipart, Double_t R, Double_t wc, Double_t e)
+Double_t AliQuenchingWeights::GetELossRandomKFastR(Int_t ipart, Double_t r, Double_t wc, Double_t e)
{
// return DeltaE for new dynamic version
// for given parton type, R and wc value and energy
// all parameters are well within the bounds.
// read-in data tables before first call
- if(R>=fgkRMax) {
- R=fgkRMax-1;
+ if(r>=fgkRMax) {
+ r=fgkRMax-1;
}
-
+
Double_t discrete=0.;
Double_t continuous=0.;
Int_t bin=1;
Double_t xxxx = fHisto->GetBinCenter(bin);
if(fMultSoft)
- CalcMult(ipart,R,xxxx,continuous,discrete);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
else
- CalcSingleHard(ipart,R,xxxx,continuous,discrete);
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
- if(discrete>0.999) {
+ if(discrete>=1.0) {
return 0.; //no energy loss
}
-
+ if (!fHisto) Init();
+
fHisto->SetBinContent(bin,continuous);
- const Int_t kbinmax=fHisto->FindBin(e/wc);
+ Int_t kbinmax=fHisto->FindBin(e/wc);
+ if(kbinmax>=fgkBins) kbinmax=fgkBins-1;
if(kbinmax==1) return e; //maximum energy loss
if(fMultSoft) {
- for(Int_t bin=2; bin<=kbinmax; bin++) {
+ for(bin=2; bin<=kbinmax; bin++) {
xxxx = fHisto->GetBinCenter(bin);
- CalcMult(ipart,R,xxxx,continuous,discrete);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
fHisto->SetBinContent(bin,continuous);
}
} else {
- for(Int_t bin=2; bin<=kbinmax; bin++) {
+ for(bin=2; bin<=kbinmax; bin++) {
xxxx = fHisto->GetBinCenter(bin);
- CalcSingleHard(ipart,R,xxxx,continuous,discrete);
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
fHisto->SetBinContent(bin,continuous);
}
}
- const Double_t kdelta=fHisto->Integral(1,kbinmax);
- Double_t val=discrete*fgkBins/fgkMaxBin;
- fHisto->Fill(0.,val);
-
- if(fECMethod==kReweight)
+ if(fECMethod==kReweight){
fHisto->SetBinContent(kbinmax+1,0);
- else
+ fHisto->Fill(0.,discrete*fgkBins/fgkMaxBin);
+ } else if (fECMethod==kReweightCont) {
+ fHisto->SetBinContent(kbinmax+1,0);
+ const Double_t kdelta=fHisto->Integral(1,kbinmax);
+ fHisto->Scale(1./kdelta*(1-discrete));
+ fHisto->Fill(0.,discrete);
+ } else {
+ const Double_t kdelta=fHisto->Integral(1,kbinmax);
+ Double_t val=discrete*fgkBins/fgkMaxBin;
+ fHisto->Fill(0.,val);
fHisto->SetBinContent(kbinmax+1,(1-discrete)*fgkBins/fgkMaxBin-kdelta);
-
- for(Int_t bin=kbinmax+2; bin<=fgkBins; bin++) {
+ }
+ for(bin=kbinmax+2; bin<=fgkBins; bin++) {
fHisto->SetBinContent(bin,0);
}
-
+ //cout << kbinmax << " " << discrete << " " << fHisto->Integral() << endl;
Double_t ret=fHisto->GetRandom()*wc;
if(ret>e) return e;
return ret;
return e-loss;
}
+Double_t AliQuenchingWeights::GetDiscreteWeight(Int_t ipart, Double_t I0, Double_t I1)
+{
+ // return discrete weight
+
+ Double_t r=CalcRk(I0,I1);
+ if(r<=0.){
+ return 1.;
+ }
+ return GetDiscreteWeightR(ipart,r);
+}
+
+Double_t AliQuenchingWeights::GetDiscreteWeightR(Int_t ipart, Double_t r)
+{
+ // return discrete weight
+
+ if(r>=fgkRMax) {
+ r=fgkRMax-1;
+ }
+
+ Double_t discrete=0.;
+ Double_t continuous=0.;
+ Int_t bin=1;
+ Double_t xxxx = fHisto->GetBinCenter(bin);
+ if(fMultSoft)
+ CalcMult(ipart,r,xxxx,continuous,discrete);
+ else
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
+ return discrete;
+}
+
+void AliQuenchingWeights::GetZeroLossProb(Double_t &p,Double_t &prw,Double_t &prwcont,
+ Int_t ipart,Double_t I0,Double_t I1,Double_t e)
+{
+ //calculate the probabilty that there is no energy
+ //loss for different ways of energy constraint
+ p=1.;prw=1.;prwcont=1.;
+ Double_t r=CalcRk(I0,I1);
+ if(r<=0.){
+ return;
+ }
+ Double_t wc=CalcWCk(I1);
+ if(wc<=0.){
+ return;
+ }
+ GetZeroLossProbR(p,prw,prwcont,ipart,r,wc,e);
+}
+
+void AliQuenchingWeights::GetZeroLossProbR(Double_t &p,Double_t &prw,Double_t &prwcont,
+ Int_t ipart, Double_t r,Double_t wc,Double_t e)
+{
+ //calculate the probabilty that there is no energy
+ //loss for different ways of energy constraint
+ if(r>=fgkRMax) {
+ r=fgkRMax-1;
+ }
+
+ Double_t discrete=0.;
+ Double_t continuous=0.;
+ if (!fHisto) Init();
+ Int_t kbinmax=fHisto->FindBin(e/wc);
+ if(kbinmax>=fgkBins) kbinmax=fgkBins-1;
+ if(fMultSoft) {
+ for(Int_t bin=1; bin<=kbinmax; bin++) {
+ Double_t xxxx = fHisto->GetBinCenter(bin);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
+ fHisto->SetBinContent(bin,continuous);
+ }
+ } else {
+ for(Int_t bin=1; bin<=kbinmax; bin++) {
+ Double_t xxxx = fHisto->GetBinCenter(bin);
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
+ fHisto->SetBinContent(bin,continuous);
+ }
+ }
+
+ //non-reweighted P(Delta E = 0)
+ const Double_t kdelta=fHisto->Integral(1,kbinmax);
+ Double_t val=discrete*fgkBins/fgkMaxBin;
+ fHisto->Fill(0.,val);
+ fHisto->SetBinContent(kbinmax+1,(1-discrete)*fgkBins/fgkMaxBin-kdelta);
+ Double_t hint=fHisto->Integral(1,kbinmax+1);
+ p=fHisto->GetBinContent(1)/hint;
+
+ // reweighted
+ hint=fHisto->Integral(1,kbinmax);
+ prw=fHisto->GetBinContent(1)/hint;
+
+ Double_t xxxx = fHisto->GetBinCenter(1);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
+ fHisto->SetBinContent(1,continuous);
+ hint=fHisto->Integral(1,kbinmax);
+ fHisto->Scale(1./hint*(1-discrete));
+ fHisto->Fill(0.,discrete);
+ prwcont=fHisto->GetBinContent(1);
+}
+
Int_t AliQuenchingWeights::SampleEnergyLoss()
{
// Has to be called to fill the histograms
Char_t meddesc[100];
if(fMultSoft) {
medvalue=(Int_t)(fQTransport*1000.);
- sprintf(meddesc,"MS");
+ snprintf(meddesc, 100, "MS");
} else {
medvalue=(Int_t)(fMu*1000.);
- sprintf(meddesc,"SH");
+ snprintf(meddesc, 100, "SH");
}
for(Int_t ipart=1;ipart<=2;ipart++){
for(Int_t l=1;l<=4*fLengthMax;l++){
Char_t hname[100];
- sprintf(hname,"hDisc-ContQW_%s_%d_%d_%d_%d",meddesc,fInstanceNumber,ipart,medvalue,l);
+ snprintf(hname, 100, "hDisc-ContQW_%s_%d_%d_%d_%d",meddesc,fInstanceNumber,ipart,medvalue,l);
Double_t len=l/4.;
Double_t wc = CalcWC(len);
fHistos[ipart-1][l-1] = new TH1F(hname,hname,fgkBins,0.,fgkMaxBin*wc);
return 0;
}
-Int_t AliQuenchingWeights::SampleEnergyLoss(Int_t ipart, Double_t R)
+Int_t AliQuenchingWeights::SampleEnergyLoss(Int_t ipart, Double_t r)
{
// Sample energy loss directly for one particle type
// choses R (safe it and keep it until next call of function)
Double_t discrete=0.;
Double_t continuous=0;;
Int_t bin=1;
+ if (!fHisto) Init();
Double_t xxxx = fHisto->GetBinCenter(bin);
if(fMultSoft)
- CalcMult(ipart,R,xxxx,continuous,discrete);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
else
- CalcSingleHard(ipart,R,xxxx,continuous,discrete);
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
if(discrete>=1.) {
fHisto->SetBinContent(1,1.);
- for(Int_t bin=2;bin<=fgkBins;bin++)
+ for(bin=2;bin<=fgkBins;bin++)
fHisto->SetBinContent(bin,0.);
return 0;
}
fHisto->SetBinContent(bin,continuous);
- for(Int_t bin=2; bin<=fgkBins; bin++) {
+ for(bin=2; bin<=fgkBins; bin++) {
xxxx = fHisto->GetBinCenter(bin);
if(fMultSoft)
- CalcMult(ipart,R,xxxx,continuous,discrete);
+ CalcMult(ipart,r,xxxx,continuous,discrete);
else
- CalcSingleHard(ipart,R,xxxx,continuous,discrete);
+ CalcSingleHard(ipart,r,xxxx,continuous,discrete);
fHisto->SetBinContent(bin,continuous);
}
Char_t meddesc[100];
if(fMultSoft) {
wc=CalcWC(medval,length);
- sprintf(meddesc,"MS");
+ snprintf(meddesc, 100, "MS");
} else {
wc=CalcWCbar(medval,length);
- sprintf(meddesc,"SH");
+ snprintf(meddesc, 100, "SH");
}
Char_t hname[100];
- sprintf(hname,"hContQWHisto_%s_%d_%d_%d",meddesc,ipart,
+ snprintf(hname, 100, "hContQWHisto_%s_%d_%d_%d",meddesc,ipart,
(Int_t)(medval*1000.),(Int_t)length);
TH1F *hist = new TH1F("hist",hname,fgkBins,0.,fgkMaxBin*wc);
Char_t meddesc[100];
if(fMultSoft) {
wc=CalcWC(medval,length);
- sprintf(meddesc,"MS");
+ snprintf(meddesc, 100, "MS");
} else {
wc=CalcWCbar(medval,length);
- sprintf(meddesc,"SH");
+ snprintf(meddesc, 100, "SH");
}
Char_t hname[100];
- sprintf(hname,"hContQWHistox_%s_%d_%d_%d",meddesc,ipart,
+ snprintf(hname, 100, "hContQWHistox_%s_%d_%d_%d",meddesc,ipart,
(Int_t)(medval*1000.),(Int_t)length);
TH1F *histx = new TH1F("histx",hname,fgkBins,0.,fgkMaxBin);
return histx;
}
-TH1F* AliQuenchingWeights::ComputeQWHistoX(Int_t ipart,Double_t R) const
+TH1F* AliQuenchingWeights::ComputeQWHistoX(Int_t ipart,Double_t r) const
{
// compute P(E) distribution for
// given ipart = 1 for quark, 2 for gluon
Char_t meddesc[100];
if(fMultSoft) {
- sprintf(meddesc,"MS");
+ snprintf(meddesc, 100, "MS");
} else {
- sprintf(meddesc,"SH");
+ snprintf(meddesc, 100, "SH");
}
Char_t hname[100];
- sprintf(hname,"hQWHistox_%s_%d_%.2f",meddesc,ipart,R);
+ snprintf(hname, 100, "hQWHistox_%s_%d_%.2f",meddesc,ipart,r);
TH1F *histx = new TH1F("histx",hname,fgkBins,0.,fgkMaxBin);
histx->SetXTitle("x = #Delta E/#omega_{c}");
if(fMultSoft)
histx->SetYTitle("p(#Delta E/#bar#omega_{c})");
histx->SetLineColor(4);
- Double_t rrrr = R;
+ Double_t rrrr = r;
Double_t continuous=0.,discrete=0.;
//loop on histogram channels
for(Int_t bin=1; bin<=fgkBins; bin++) {
Char_t name[100];
Char_t hname[100];
if(ipart==1){
- sprintf(name,"Energy Loss Distribution - Quarks;E_{loss} (GeV);#");
- sprintf(hname,"hLossQuarks");
+ snprintf(name, 100, "Energy Loss Distribution - Quarks;E_{loss} (GeV);#");
+ snprintf(hname,100, "hLossQuarks");
} else {
- sprintf(name,"Energy Loss Distribution - Gluons;E_{loss} (GeV);#");
- sprintf(hname,"hLossGluons");
+ snprintf(name, 100, "Energy Loss Distribution - Gluons;E_{loss} (GeV);#");
+ snprintf(hname, 100, "hLossGluons");
}
TH1F *h = new TH1F(hname,name,250,0,250);
Char_t name[100];
Char_t hname[100];
if(ipart==1){
- sprintf(name,"Energy Loss Distribution - Quarks;E_{loss} (GeV);#");
- sprintf(hname,"hLossQuarks");
+ snprintf(name, 100, "Energy Loss Distribution - Quarks;E_{loss} (GeV);#");
+ snprintf(hname,100, "hLossQuarks");
} else {
- sprintf(name,"Energy Loss Distribution - Gluons;E_{loss} (GeV);#");
- sprintf(hname,"hLossGluons");
+ snprintf(name,100, "Energy Loss Distribution - Gluons;E_{loss} (GeV);#");
+ snprintf(hname, 100, "hLossGluons");
}
TH1F *h = new TH1F(hname,name,250,0,250);
return h;
}
-TH1F* AliQuenchingWeights::ComputeELossHisto(Int_t ipart,Double_t R) const
+TH1F* AliQuenchingWeights::ComputeELossHisto(Int_t ipart,Double_t r) const
{
// compute energy loss histogram for
// parton type and given R
- TH1F *dummy = ComputeQWHistoX(ipart,R);
+ TH1F *dummy = ComputeQWHistoX(ipart,r);
if(!dummy) return 0;
Char_t hname[100];
- sprintf(hname,"hELossHistox_%d_%.2f",ipart,R);
+ snprintf(hname, 100, "hELossHistox_%d_%.2f",ipart,r);
TH1F *histx = new TH1F("histxr",hname,fgkBins,0.,fgkMaxBin);
for(Int_t i=0;i<100000;i++){
//if(i % 1000 == 0) cout << "." << flush;
return ret;
}
-Double_t AliQuenchingWeights::GetMeanELoss(Int_t ipart,Double_t R) const
+Double_t AliQuenchingWeights::GetMeanELoss(Int_t ipart,Double_t r) const
{
// compute average energy loss over wc
// for parton type and given R
- TH1F *dummy = ComputeELossHisto(ipart,R);
+ TH1F *dummy = ComputeELossHisto(ipart,r);
if(!dummy) return 0;
Double_t ret=dummy->GetMean();
delete dummy;
return ret;
}
-void AliQuenchingWeights::PlotDiscreteWeights(Double_t len) const
+void AliQuenchingWeights::PlotDiscreteWeights(Double_t len,Double_t qm) const
{
// plot discrete weights for given length
c = new TCanvas("cdiscsh","Discrete Weight for Single Hard Scattering",0,0,500,400);
c->cd();
- TH2F *hframe = new TH2F("hdisc","",2,0,5.1,2,0,1);
+ TH2F *hframe = new TH2F("hdisc","",2,0,qm+.1,2,0,1.25);
hframe->SetStats(0);
if(fMultSoft)
hframe->SetXTitle("#hat{q} [GeV^{2}/fm]");
else
hframe->SetXTitle("#mu [GeV]");
- hframe->SetYTitle("Probability #Delta E = 0 , p_{0}");
+ //hframe->SetYTitle("Probability #Delta E = 0 , p_{0}");
+ hframe->SetYTitle("p_{0} (discrete weight)");
hframe->Draw();
- TGraph *gq=new TGraph(20);
+ Int_t points=(Int_t)qm*4;
+ TGraph *gq=new TGraph(points);
Int_t i=0;
if(fMultSoft) {
- for(Double_t q=0.05;q<=5.05;q+=0.25){
+ for(Double_t q=0.05;q<=qm+.05;q+=0.25){
Double_t disc,cont;
CalcMult(1,1.0,q,len,cont,disc);
gq->SetPoint(i,q,disc);i++;
}
} else {
- for(Double_t m=0.05;m<=5.05;m+=0.25){
+ for(Double_t m=0.05;m<=qm+.05;m+=0.25){
Double_t disc,cont;
CalcSingleHard(1,1.0,m,len,cont, disc);
gq->SetPoint(i,m,disc);i++;
}
}
gq->SetMarkerStyle(20);
- gq->Draw("pl");
+ gq->SetMarkerColor(1);
+ gq->SetLineStyle(1);
+ gq->SetLineColor(1);
+ gq->Draw("l");
- TGraph *gg=new TGraph(20);
+ TGraph *gg=new TGraph(points);
i=0;
if(fMultSoft){
- for(Double_t q=0.05;q<=5.05;q+=0.25){
+ for(Double_t q=0.05;q<=qm+.05;q+=0.25){
Double_t disc,cont;
CalcMult(2,1.0,q,len,cont,disc);
gg->SetPoint(i,q,disc);i++;
}
} else {
- for(Double_t m=0.05;m<=5.05;m+=0.25){
+ for(Double_t m=0.05;m<=qm+.05;m+=0.25){
Double_t disc,cont;
CalcSingleHard(2,1.0,m,len,cont,disc);
gg->SetPoint(i,m,disc);i++;
}
}
gg->SetMarkerStyle(24);
- gg->Draw("pl");
+ gg->SetMarkerColor(2);
+ gg->SetLineStyle(2);
+ gg->SetLineColor(2);
+ gg->Draw("l");
TLegend *l1a = new TLegend(0.5,0.6,.95,0.8);
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"L = %.1f fm",len);
+ snprintf(label, 100, "L = %.1f fm",len);
l1a->AddEntry(gq,label,"");
- l1a->AddEntry(gq,"quark","pl");
- l1a->AddEntry(gg,"gluon","pl");
+ l1a->AddEntry(gq,"quark projectile","l");
+ l1a->AddEntry(gg,"gluon projectile","l");
l1a->Draw();
c->Update();
Char_t name[1024];
if(fMultSoft) {
if(itype==1)
- sprintf(title,"Cont. Weight for Multiple Scattering - Quarks");
+ snprintf(title, 1024, "Cont. Weight for Multiple Scattering - Quarks");
else if(itype==2)
- sprintf(title,"Cont. Weight for Multiple Scattering - Gluons");
+ snprintf(title, 1024, "Cont. Weight for Multiple Scattering - Gluons");
else return;
medvals[0]=4;medvals[1]=1;medvals[2]=0.5;
- sprintf(name,"ccont-ms-%d",itype);
+ snprintf(name, 1024, "ccont-ms-%d",itype);
} else {
if(itype==1)
- sprintf(title,"Cont. Weight for Single Hard Scattering - Quarks");
+ snprintf(title, 1024, "Cont. Weight for Single Hard Scattering - Quarks");
else if(itype==2)
- sprintf(title,"Cont. Weight for Single Hard Scattering - Gluons");
+ snprintf(title, 1024, "Cont. Weight for Single Hard Scattering - Gluons");
else return;
medvals[0]=2;medvals[1]=1;medvals[2]=0.5;
- sprintf(name,"ccont-ms-%d",itype);
+ snprintf(name, 1024, "ccont-ms-%d",itype);
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"L = %.1f fm",ell);
+ snprintf(label, 100, "L = %.1f fm",ell);
l1a->AddEntry(h1,label,"");
if(fMultSoft) {
- sprintf(label,"#hat{q} = %.1f GeV^{2}/fm",medvals[0]);
+ snprintf(label, 100, "#hat{q} = %.1f GeV^{2}/fm",medvals[0]);
l1a->AddEntry(h1,label,"pl");
- sprintf(label,"#hat{q} = %.1f GeV^{2}/fm",medvals[1]);
+ snprintf(label, 100, "#hat{q} = %.1f GeV^{2}/fm",medvals[1]);
l1a->AddEntry(h2,label,"pl");
- sprintf(label,"#hat{q} = %.1f GeV^{2}/fm",medvals[2]);
- l1a->AddEntry(h3,label,"pl");
+ snprintf(label, 100, "#hat{q} = %.1f GeV^{2}/fm",medvals[2]);
+ l1a->AddEntry(h3, label,"pl");
} else {
- sprintf(label,"#mu = %.1f GeV",medvals[0]);
+ snprintf(label, 100, "#mu = %.1f GeV",medvals[0]);
l1a->AddEntry(h1,label,"pl");
- sprintf(label,"#mu = %.1f GeV",medvals[1]);
+ snprintf(label, 100, "#mu = %.1f GeV",medvals[1]);
l1a->AddEntry(h2,label,"pl");
- sprintf(label,"#mu = %.1f GeV",medvals[2]);
+ snprintf(label, 100, "#mu = %.1f GeV",medvals[2]);
l1a->AddEntry(h3,label,"pl");
}
l1a->Draw();
Char_t name[1024];
if(fMultSoft) {
if(itype==1)
- sprintf(title,"Cont. Weight for Multiple Scattering - Quarks");
+ snprintf(title,1024, "Cont. Weight for Multiple Scattering - Quarks");
else if(itype==2)
- sprintf(title,"Cont. Weight for Multiple Scattering - Gluons");
+ snprintf(title,1024, "Cont. Weight for Multiple Scattering - Gluons");
else return;
- sprintf(name,"ccont2-ms-%d",itype);
+ snprintf(name,1024, "ccont2-ms-%d",itype);
} else {
if(itype==1)
- sprintf(title,"Cont. Weight for Single Hard Scattering - Quarks");
+ snprintf(title, 1024, "Cont. Weight for Single Hard Scattering - Quarks");
else if(itype==2)
- sprintf(title,"Cont. Weight for Single Hard Scattering - Gluons");
+ snprintf(title, 1024, "Cont. Weight for Single Hard Scattering - Gluons");
else return;
- sprintf(name,"ccont2-sh-%d",itype);
+ snprintf(name, 1024, "ccont2-sh-%d",itype);
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
c->cd();
l1a->SetBorderSize(0);
Char_t label[100];
if(fMultSoft)
- sprintf(label,"#hat{q} = %.1f GeV^{2}/fm",medval);
+ snprintf(label, 100, "#hat{q} = %.1f GeV^{2}/fm",medval);
else
- sprintf(label,"#mu = %.1f GeV",medval);
+ snprintf(label, 100, "#mu = %.1f GeV",medval);
l1a->AddEntry(h1,label,"");
l1a->AddEntry(h1,"L = 8 fm","pl");
c->Update();
}
-void AliQuenchingWeights::PlotAvgELoss(Double_t len,Double_t e) const
+void AliQuenchingWeights::PlotAvgELoss(Double_t len,Double_t qm,Double_t e) const
{
// plot average energy loss for given length
// and parton energy
Char_t title[1024];
Char_t name[1024];
if(fMultSoft){
- sprintf(title,"Average Energy Loss for Multiple Scattering");
- sprintf(name,"cavgelossms");
+ snprintf(title, 1024, "Average Energy Loss for Multiple Scattering");
+ snprintf(name, 1024, "cavgelossms");
} else {
- sprintf(title,"Average Energy Loss for Single Hard Scattering");
- sprintf(name,"cavgelosssh");
+ snprintf(title, 1024, "Average Energy Loss for Single Hard Scattering");
+ snprintf(name, 1024, "cavgelosssh");
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
c->cd();
- TH2F *hframe = new TH2F("avgloss",title,2,0,5.1,2,0,100);
+ TH2F *hframe = new TH2F("avgloss","",2,0,qm+.1,2,0,100);
hframe->SetStats(0);
if(fMultSoft)
hframe->SetXTitle("#hat{q} [GeV^{2}/fm]");
TGraph *gq=new TGraph(20);
Int_t i=0;
- for(Double_t v=0.05;v<=5.05;v+=0.25){
+ for(Double_t v=0.05;v<=qm+.05;v+=0.25){
TH1F *dummy=ComputeELossHisto(1,v,len,e);
Double_t avgloss=dummy->GetMean();
gq->SetPoint(i,v,avgloss);i++;
delete dummy;
}
- gq->SetMarkerStyle(20);
+ gq->SetMarkerStyle(21);
gq->Draw("pl");
- TGraph *gg=new TGraph(20);
+ Int_t points=(Int_t)qm*4;
+ TGraph *gg=new TGraph(points);
i=0;
- for(Double_t v=0.05;v<=5.05;v+=0.25){
+ for(Double_t v=0.05;v<=qm+.05;v+=0.25){
TH1F *dummy=ComputeELossHisto(2,v,len,e);
Double_t avgloss=dummy->GetMean();
gg->SetPoint(i,v,avgloss);i++;
delete dummy;
}
- gg->SetMarkerStyle(24);
+ gg->SetMarkerStyle(20);
+ gg->SetMarkerColor(2);
gg->Draw("pl");
- TGraph *gratio=new TGraph(20);
- for(Int_t i=0;i<20;i++){
+ TGraph *gratio=new TGraph(points);
+ for(i=0;i<points;i++){
Double_t x,y,x2,y2;
gg->GetPoint(i,x,y);
gq->GetPoint(i,x2,y2);
}
gratio->SetLineStyle(4);
gratio->Draw();
- TLegend *l1a = new TLegend(0.5,0.6,.95,0.8);
+ TLegend *l1a = new TLegend(0.15,0.60,0.50,0.90);
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"L = %.1f fm",len);
+ snprintf(label, 100, "L = %.1f fm",len);
l1a->AddEntry(gq,label,"");
- l1a->AddEntry(gq,"quark","pl");
- l1a->AddEntry(gg,"gluon","pl");
+ l1a->AddEntry(gq,"quark projectile","pl");
+ l1a->AddEntry(gg,"gluon projectile","pl");
l1a->AddEntry(gratio,"gluon/quark/2.25*10","pl");
l1a->Draw();
Char_t title[1024];
Char_t name[1024];
if(fMultSoft){
- sprintf(title,"Average Energy Loss for Multiple Scattering");
- sprintf(name,"cavgelossms2");
+ snprintf(title, 1024, "Average Energy Loss for Multiple Scattering");
+ snprintf(name, 1024, "cavgelossms2");
} else {
- sprintf(title,"Average Energy Loss for Single Hard Scattering");
- sprintf(name,"cavgelosssh2");
+ snprintf(title, 1024, "Average Energy Loss for Single Hard Scattering");
+ snprintf(name, 1024, "cavgelosssh2");
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
gg->Draw("pl");
TGraph *gratio=new TGraph(20);
- for(Int_t i=0;i<20;i++){
+ for(i=0;i<20;i++){
Double_t x,y,x2,y2;
gg->GetPoint(i,x,y);
gq->GetPoint(i,x2,y2);
else gratio->SetPoint(i,x,0);
}
gratio->SetLineStyle(4);
- gratio->Draw();
+ //gratio->Draw();
TLegend *l1a = new TLegend(0.5,0.6,.95,0.8);
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"<L> = %.2f fm",hEll->GetMean());
+ snprintf(label, 100, "<L> = %.2f fm",hEll->GetMean());
l1a->AddEntry(gq,label,"");
l1a->AddEntry(gq,"quark","pl");
l1a->AddEntry(gg,"gluon","pl");
- l1a->AddEntry(gratio,"gluon/quark/2.25*10","pl");
+ //l1a->AddEntry(gratio,"gluon/quark/2.25*10","pl");
l1a->Draw();
c->Update();
Char_t name[1024];
Float_t medval;
if(fMultSoft){
- sprintf(title,"Average Energy Loss for Multiple Scattering");
- sprintf(name,"cavgelosslms");
+ snprintf(title, 1024, "Average Energy Loss for Multiple Scattering");
+ snprintf(name, 1024, "cavgelosslms");
medval=fQTransport;
} else {
- sprintf(title,"Average Energy Loss for Single Hard Scattering");
- sprintf(name,"cavgelosslsh");
+ snprintf(title, 1024, "Average Energy Loss for Single Hard Scattering");
+ snprintf(name, 1024, "cavgelosslsh");
medval=fMu;
}
gg->Draw("pl");
TGraph *gratio=new TGraph((Int_t)fLengthMax*4);
- for(Int_t i=0;i<=(Int_t)fLengthMax*4;i++){
+ for(i=0;i<=(Int_t)fLengthMax*4;i++){
Double_t x,y,x2,y2;
gg->GetPoint(i,x,y);
gq->GetPoint(i,x2,y2);
l1a->SetBorderSize(0);
Char_t label[100];
if(fMultSoft)
- sprintf(label,"#hat{q} = %.2f GeV^{2}/fm",medval);
+ snprintf(label, 100, "#hat{q} = %.2f GeV^{2}/fm",medval);
else
- sprintf(label,"#mu = %.2f GeV",medval);
+ snprintf(label, 100, "#mu = %.2f GeV",medval);
l1a->AddEntry(gq,label,"");
l1a->AddEntry(gq,"quark","pl");
l1a->AddEntry(gg,"gluon","pl");
Char_t title[1024];
Char_t name[1024];
if(fMultSoft){
- sprintf(title,"Relative Energy Loss for Multiple Scattering");
- sprintf(name,"cavgelossvsptms");
+ snprintf(title, 1024, "Relative Energy Loss for Multiple Scattering");
+ snprintf(name, 1024, "cavgelossvsptms");
} else {
- sprintf(title,"Relative Energy Loss for Single Hard Scattering");
- sprintf(name,"cavgelossvsptsh");
+ snprintf(title, 1024, "Relative Energy Loss for Single Hard Scattering");
+ snprintf(name, 1024, "cavgelossvsptsh");
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"L = %.1f fm",len);
+ snprintf(label, 100, "L = %.1f fm",len);
l1a->AddEntry(gq,label,"");
l1a->AddEntry(gq,"quark","pl");
l1a->AddEntry(gg,"gluon","pl");
Char_t title[1024];
Char_t name[1024];
if(fMultSoft){
- sprintf(title,"Relative Energy Loss for Multiple Scattering");
- sprintf(name,"cavgelossvsptms2");
+ snprintf(title, 1024, "Relative Energy Loss for Multiple Scattering");
+ snprintf(name, 1024, "cavgelossvsptms2");
} else {
- sprintf(title,"Relative Energy Loss for Single Hard Scattering");
- sprintf(name,"cavgelossvsptsh2");
+ snprintf(title, 1024, "Relative Energy Loss for Single Hard Scattering");
+ snprintf(name, 1024, "cavgelossvsptsh2");
}
TCanvas *c = new TCanvas(name,title,0,0,500,400);
c->cd();
l1a->SetFillStyle(0);
l1a->SetBorderSize(0);
Char_t label[100];
- sprintf(label,"<L> = %.2f fm",hEll->GetMean());
+ snprintf(label, 100, "<L> = %.2f fm",hEll->GetMean());
l1a->AddEntry(gq,label,"");
l1a->AddEntry(gq,"quark","pl");
l1a->AddEntry(gg,"gluon","pl");
Int_t AliQuenchingWeights::GetIndex(Double_t len) const
{
+ //get the index according to length
if(len>fLengthMax) len=fLengthMax;
Int_t l=Int_t(len/0.25);
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff