Page MenuHomeHEPForge

No OneTemporary

Index: trunk/macros/tools/runArea.C
===================================================================
--- trunk/macros/tools/runArea.C (revision 86)
+++ trunk/macros/tools/runArea.C (revision 87)
@@ -1,122 +1,127 @@
#if !defined(__CINT__) || defined(__MAKECINT__)
#include "../runglauber.C"
#include <TProfile2D.h>
#endif
void runArea(const Int_t n,
const char *sysA,
const char *sysB,
const Double_t signn,
const Double_t mind,
const char *fname)
{
TGlauberMC *mcg = new TGlauberMC(sysA,sysB,signn);
mcg->SetMinDistance(mind);
if (0) {
cerr << "warning: this is for central events" << endl;
mcg->SetBmin(0);
mcg->SetBmax(0.);
}
const Double_t d2 = (Double_t)signn/(TMath::Pi()*10);
+ const Double_t r2 = d2/4.;
- const Int_t nbins=100;
+ const Int_t nbins=200;
const Double_t ell=10;
const Double_t da=2*ell*2*ell/nbins/nbins;
TH2D *areaA = new TH2D("hAreaA",";x (fm);y (fm)",nbins,-ell,ell,nbins,-ell,ell);
areaA->SetStats(0);
TH2D *areaB = new TH2D("hAreaB",";x (fm);y (fm)",nbins,-ell,ell,nbins,-ell,ell);
areaB->SetStats(0);
TFile *out = TFile::Open(fname,"recreate",fname,9);
if (!out)
return;
TNtuple *nt = new TNtuple("nt","nt",
- "Npart:Ncoll:B:Area:S");
+ "Npart:Ncoll:B:Ais:Aov:S");
nt->SetDirectory(out);
const TGlauNucleus *nucA = mcg->GetNucleusA();
const Int_t AN = nucA->GetN();
const TGlauNucleus *nucB = mcg-> GetNucleusB();
const Int_t BN = nucB->GetN();
for (Int_t ev=0; ev<n; ++ev) {
while (!mcg->NextEvent());
if (!(ev%100)) {
cout << "Event # " << ev << " out of " << n << " \r" << flush;
}
areaA->Reset();
areaB->Reset();
Double_t mx = mcg->GetMeanXParts();
Double_t my = mcg->GetMeanYParts();
Double_t npart = mcg->GetNpart();
Double_t bmc = mcg->GetB();
const TObjArray *nucleonsA = nucA->GetNucleons();
for (Int_t i = 0; i<AN; ++i) {
TGlauNucleon *nucleonA=(TGlauNucleon*)(nucleonsA->At(i));
if (!nucleonA->IsWounded())
continue;
Double_t x = nucleonA->GetX()-mx;
Double_t y = nucleonA->GetY()-my;
for (Int_t xi=1; xi<=nbins; ++xi) {
for (Int_t yi=1; yi<=nbins; ++yi) {
Int_t bin = areaA->GetBin(xi,yi);
Double_t val=areaA->GetBinContent(bin);
if (val>0)
continue;
Double_t dx=x-areaA->GetXaxis()->GetBinCenter(xi);
Double_t dy=y-areaA->GetYaxis()->GetBinCenter(yi);
- if (dx*dx+dy*dy<d2)
+ if (dx*dx+dy*dy<r2)
areaA->SetBinContent(bin,1);
}
}
}
const TObjArray *nucleonsB = nucB->GetNucleons();
for (Int_t i = 0; i<BN; ++i) {
TGlauNucleon *nucleonB=(TGlauNucleon*)(nucleonsB->At(i));
if (!nucleonB->IsWounded())
continue;
Double_t x = nucleonB->GetX()-mx;
Double_t y = nucleonB->GetY()-my;
for (Int_t xi=1; xi<=nbins; ++xi) {
for (Int_t yi=1; yi<=nbins; ++yi) {
Int_t bin = areaB->GetBin(xi,yi);
Double_t val=areaB->GetBinContent(bin);
if (val>0)
continue;
Double_t dx=x-areaB->GetXaxis()->GetBinCenter(xi);
Double_t dy=y-areaB->GetYaxis()->GetBinCenter(yi);
- if (dx*dx+dy*dy<d2)
+ if (dx*dx+dy*dy<r2)
areaB->SetBinContent(bin,1);
}
}
}
- Double_t overlap=0;
+ Double_t overlap1=0;
+ Double_t overlap2=0;
for (Int_t xi=1; xi<=nbins; ++xi) {
for (Int_t yi=1; yi<=nbins; ++yi) {
Int_t bin = areaA->GetBin(xi,yi);
Double_t vA=areaA->GetBinContent(bin);
Double_t vB=areaB->GetBinContent(bin);
if (vA>0&&vB>0)
- ++overlap;
+ ++overlap1;
+ if (vA>0||vB>0)
+ ++overlap2;
}
}
Double_t ap = TMath::Sqrt(mcg->GetSx2()*mcg->GetSy2()-mcg->GetSxy()*mcg->GetSxy());
- Float_t v[5]; Int_t f=0;
+ Float_t v[6]; Int_t f=0;
v[f++] = mcg->GetNpart();
v[f++] = mcg->GetNcoll();
v[f++] = mcg->GetB();
- v[f++] = overlap*da;
+ v[f++] = overlap1*da;
+ v[f++] = overlap2*da;
v[f++] = ap;
nt->Fill(v);
}
delete areaA;
out->Write();
out->Close();
delete out;
}

File Metadata

Mime Type
text/x-diff
Expires
Sat, May 3, 6:43 AM (1 d, 8 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4983121
Default Alt Text
(4 KB)

Event Timeline