Logo ROOT   6.10/09
Reference Guide
pirndm.C
Go to the documentation of this file.
1 //Test program for random number generators (spped and quality)
2 //The program get n random number pairs x and y i [0,1]
3 //It counts the ratio of pairs in the circle of diameter 1
4 //compared to the total number of pairs.
5 //This ratio must be Pi/4
6 //The test shows the graph of the difference of this ratio compared to PI.
7 //To test only the speed, call pirndm(1) (default)
8 //To test quality and speed call pirndm(50)
9 // root pirndm.C+
10 //or
11 // root "pirndm.C+(10)"
12 //
13 //Author: Rene Brun
14 
15 #include "TROOT.h"
16 #include "TStopwatch.h"
17 #include "TMath.h"
18 #include "TRandom1.h"
19 #include "TRandom2.h"
20 #include "TRandom3.h"
21 #include "TCanvas.h"
22 #include "TH2.h"
23 #include "TGraph.h"
24 #include "TSystem.h"
25 #include "TLegend.h"
26 #include "TPaveLabel.h"
27 #include "TSystem.h"
28 #ifdef USE_MATHMORE
29 #include "Math/Random.h"
30 #include "Math/GSLRndmEngines.h"
31 #endif
32 #include <vector>
33 #include <iostream>
34 
36 TCanvas *c1 = 0;
39 
40 std::vector<TH2D *> vh2;
41 
42 
43 //____________________________________________________________________
44 template<class Random>
45 void piRandom(const char *name, Random *r, Long64_t n, Int_t color) {
46 
47  TH2D * h2 = new TH2D("h2",name,300,0,0.000003,300,0.,1.);
48 
49  timer.Start();
50  TGraph *gr = new TGraph();
51  gr->SetMarkerStyle(20);
52  gr->SetMarkerSize(0.9);
53  gr->SetMarkerColor(color);
54  gr->SetLineColor(color);
55  gr->SetLineWidth(2);
56 
57  Int_t k = 0;
58  Double_t diffpi;
59  Long64_t npi = 0;
60  Double_t pi = TMath::Pi();
61  const Int_t NR = 20000;
62  const Int_t NR2 = NR/2;
63  Double_t rn[NR];
64  Long64_t i = 0;
65  //double r1,r2;
66  while (i<=n) {
67  i += NR2;
68  r->RndmArray(NR,rn);
69  for (Int_t j=0;j<NR;j+=2) {
70  if (rn[j]*rn[j]+rn[j+1]*rn[j+1] <= 1) npi++;
71  if (rn[j] < 0.001) h2->Fill(rn[j],rn[j+1]);
72  }
73 // r1 = r->Rndm();
74 // r2 = r->Rndm();
75  // if (r1*r1+r2*r2 <= 1) npi++;
76  if (i && i % (n/10) == 0) {
78  Double_t norm = 4./Double_t(i);
79  diffpi = norm*npi - pi;
80  gr->SetPoint(k,i,diffpi);
81  if (k ==0) gr->Draw("lp");
82  else {
83  c1->Modified();
84  c1->Update();
85  }
86  k++;
87  }
88  }
89  timer.Stop();
90  Double_t cpu = timer.CpuTime();
91  cputot += cpu;
92  Double_t nanos = 1.e9*cpu/Double_t(2*n);
93  legend->AddEntry(gr,Form("%-14s: %6.1f ns/call",name,nanos),"lp");
94  c1->Modified();
95  c1->Update();
96  printf("RANDOM = %s : RT=%7.3f s, Cpu=%7.3f s\n",name,timer.RealTime(),cpu);
97 
98 // TCanvas * c2 = new TCanvas();
99 // h2->Draw();
100 // c2->Update();
101 
102 // c1->SetSelected(c1);
103  vh2.push_back(h2);
104 
105 }
106 
107 //________________________________________________________________________
109  Int_t np = 40;
110  TGraph *g = new TGraph(2*np+2);
111  Double_t xmax = Double_t(n)/Double_t(np);
112  for (Int_t i=1;i<=np;i++) {
113  Double_t x = i*xmax;
114  //Double_t e = 1./TMath::Sqrt(x);
115  Double_t e = TMath::Sqrt( 2 * TMath::Pi() * (4 - TMath::Pi() )/x );
116  g->SetPoint(i,x,e);
117  g->SetPoint(2*np-i+1,x,-e);
118  }
119  Double_t x0 = 0.1*xmax;
120  Double_t e0 = 1./TMath::Sqrt(x0);
121  g->SetPoint(0,x0,e0);
122  g->SetPoint(2*np+1,0,-e0);
123  g->SetPoint(2*np+2,0,e0);
124  g->SetFillColor(1);
125  g->SetFillStyle(3002);
126  g->Draw("f");
127 }
128 
129 
130 //________________________________________________________________________
131 void pirndm(Long64_t n1=1, unsigned int seed = 0) {
132  Long64_t n = n1*20000000;
133  c1 = new TCanvas("c1");
134  c1->SetLeftMargin(0.12);
135  c1->SetFrameFillColor(41);
136  c1->SetFrameBorderSize(6);
137  c1->SetGrid();
138  Double_t dy = 10/TMath::Sqrt(n);
139  // Double_t dy = 1.5e-3;
140  //if (n1 < 4) dy *=2;
141  TH2F *frame = new TH2F("h","",100,0,1.1*n,100,-dy,dy);
142  frame->GetXaxis()->SetTitle("Number of TRandom calls");
143  frame->GetYaxis()->SetTitle("Difference with #pi");
144  frame->GetYaxis()->SetTitleOffset(1.3);
145  frame->GetYaxis()->SetDecimals();
146  frame->SetStats(0);
147  frame->Draw();
148  legend = new TLegend(0.6,0.7,0.88,0.88);
149  legend->Draw();
150 
151 
152  ErrorBand(n);
153  std::cout << "seed is " << seed << std::endl;
154 
155  piRandom("TRandom",new TRandom(seed),n,kYellow);
156  piRandom("TRandom2",new TRandom2(seed),n,kBlue);
157  piRandom("TRandom3",new TRandom3(seed),n,kRed);
158  piRandom("TRandom1",new TRandom1(seed),n,kGreen);
159 
160 #ifdef USE_MATHMORE
161 #define GSL2
162 #ifdef GSL1
163  piRandom("TRandom()",new TRandom(),n,kYellow);
164  piRandom("TRandom3(0)",new TRandom3(0),n,kBlack);
169 #endif
170 #ifdef GSL2
171  piRandom("TRandom(1)",new TRandom(1),n,kYellow);
172  piRandom("TRandom(2^30)",new TRandom(1073741824),n,kRed);
173  piRandom("TRandom(256)",new TRandom(256),n,kMagenta);
174  piRandom("TRandom(2)",new TRandom(2),n,kBlue);
176  // piRandom("RANMAR",new ROOT::Math::Random<ROOT::Math::GSLRngRanMar>(),n,kBlue);
177  //piRandom("MINSTD",new ROOT::Math::Random<ROOT::Math::GSLRngMinStd>(),n,kBlack);
178  piRandom("new TRandom2",new TRandom2(),n,kCyan);
179  //piRandom("TRandom3",new TRandom3(),n,kCyan);
180 // piRandom("CMRG",new ROOT::Math::Random<ROOT::Math::GSLRngCMRG>(),n,kRed);
181 // piRandom("MRG",new ROOT::Math::Random<ROOT::Math::GSLRngMRG>(),n,kMagenta);
182 #else
183  piRandom("new TRandom2",new TRandom2(std::rand()),n,kYellow);
184  piRandom("new TRandom2",new TRandom2(std::rand()),n,kRed);
185  piRandom("new TRandom2",new TRandom2(std::rand()),n,kMagenta);
186  piRandom("new TRandom2",new TRandom2(std::rand()),n,kBlack);
187  piRandom("new TRandom2",new TRandom2(std::rand()),n,kBlue);
188  piRandom("new TRandom2",new TRandom2(std::rand()),n,kRed);
189  piRandom("new TRandom2",new TRandom2(std::rand()),n,kGreen);
190 #endif
191 #endif
192 
193  // reftime calculated on MACOS Intel dualcore 2GHz
194  // time for TRandom + TRandom2 + TRandom3 + TRandom1 for n = 10**7 (n1=5000)
195  Double_t reftime = (4629.530 + 5358.100 + 5785.240 + 26012.17)/5000.;
196  const Double_t rootmarks = 900*Double_t(n1)*reftime/cputot;
197  TPaveLabel *pl = new TPaveLabel(0.2,0.92,0.8,0.98,Form("cpu time = %6.1fs - rootmarks = %6.1f",cputot,rootmarks),"brNDC");
198  pl->Draw();
199  printf("******************************************************************\n");
200  printf("* ROOTMARKS =%6.1f * Root%-8s %d/%d\n",rootmarks,gROOT->GetVersion(),
201  gROOT->GetVersionDate(),gROOT->GetVersionTime());
202  printf("******************************************************************\n");
203 
204  printf("Time at the end of job = %f seconds\n",cputot);
205  c1->Print("pirndm.root");
206  c1->Print("pirndm.gif");
207 
208 
209  // draw 2D histos
210  TCanvas * c2 = new TCanvas();
211  int nx = 0;
212  int ny = vh2.size();
213  for ( nx = 1; nx < ny; nx++) {
214  double r = double(vh2.size())/nx;
215  ny = int(r - 0.01) + 1;
216  }
217  nx--;
218  std::cout << nx << " " << ny << " " << vh2.size() << std::endl;
219 
220  c2->Divide(ny,nx);
221  for (unsigned int i = 0; i < vh2.size(); ++i) {
222  c2->cd (i+1);
223  vh2[i]->Draw();
224  }
225  c2->Update();
226 
227 }
const int nx
Definition: kalman.C:16
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title Offset is a correction factor with respect to the "s...
Definition: TAttAxis.cxx:262
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
virtual Bool_t ProcessEvents()
Process pending events (GUI, timers, sockets).
Definition: TSystem.cxx:423
Random number generator class based on M.
Definition: TRandom3.h:27
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
Definition: TStopwatch.cxx:110
long long Long64_t
Definition: RtypesCore.h:69
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:23
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:58
const double pi
Definition: Rtypes.h:56
Random number generator class based on the maximally quidistributed combined Tausworthe generator by ...
Definition: TRandom2.h:27
Double_t cputot
Definition: pirndm.C:38
virtual void Draw(Option_t *option="")
Draw this legend with its current attributes.
Definition: TLegend.cxx:452
Definition: Rtypes.h:55
TVirtualPad * cd(Int_t subpadnumber=0)
Set current canvas & pad.
Definition: TCanvas.cxx:679
#define gROOT
Definition: TROOT.h:375
Double_t CpuTime()
Stop the stopwatch (if it is running) and return the cputime (in seconds) passed between the start an...
Definition: TStopwatch.cxx:125
int Int_t
Definition: RtypesCore.h:41
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
Definition: Rtypes.h:56
Definition: Rtypes.h:56
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:745
TLegend * legend
Definition: pirndm.C:35
TStopwatch timer
Definition: pirndm.C:37
void Stop()
Stop the stopwatch.
Definition: TStopwatch.cxx:77
Double_t x[n]
Definition: legend1.C:17
The Ranlux Random number generator class.
Definition: TRandom1.h:27
This is the base class for the ROOT Random number generators.
Definition: TRandom.h:27
const int ny
Definition: kalman.C:17
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)
Definition: TPad.h:323
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
constexpr Double_t Pi()
Definition: TMath.h:40
A Pave (see TPave) with a text centered in the Pave.
Definition: TPaveLabel.h:20
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
void ErrorBand(Long64_t n)
Definition: pirndm.C:108
TRandom2 r(17)
R__EXTERN TSystem * gSystem
Definition: TSystem.h:539
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2851
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
tomato 2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:249
char * Form(const char *fmt,...)
std::vector< TH2D * > vh2
Definition: pirndm.C:40
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TAxis * GetYaxis()
Definition: TH1.h:301
float xmax
Definition: THbookFile.cxx:93
virtual void Print(const char *filename="") const
Save Pad contents in a file in one of various formats.
Definition: TPad.cxx:4581
const long long NR
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
TGraphErrors * gr
Definition: legend1.C:25
void piRandom(const char *name, Random *r, Long64_t n, Int_t color)
Definition: pirndm.C:45
The Canvas class.
Definition: TCanvas.h:31
return c2
Definition: legend2.C:14
virtual void Draw(Option_t *option="")
Draw this pavelabel with its current attributes.
Definition: TPaveLabel.cxx:77
double Double_t
Definition: RtypesCore.h:55
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition: TLegend.cxx:359
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
THist< 2, double, THistStatContent, THistStatUncertainty > TH2D
Definition: THist.hxx:316
Definition: Rtypes.h:56
void SetDecimals(Bool_t dot=kTRUE)
Sets the decimals flag By default, blank characters are stripped, and then the label is correctly ali...
Definition: TAxis.h:203
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2156
virtual void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0)
Automatic pad generation by division.
Definition: TPad.cxx:1135
void SetFrameBorderSize(Width_t size=1)
Definition: TAttPad.h:78
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
TCanvas * c1
Definition: pirndm.C:36
Definition: Rtypes.h:56
void SetFrameFillColor(Color_t color=1)
Definition: TAttPad.h:73
Double_t Sqrt(Double_t x)
Definition: TMath.h:591
virtual void Update()
Update canvas pad buffers.
Definition: TCanvas.cxx:2208
void pirndm(Long64_t n1=1, unsigned int seed=0)
Definition: pirndm.C:131
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:292
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:155
THist< 2, float, THistStatContent, THistStatUncertainty > TH2F
Definition: THist.hxx:317
double norm(double *x, double *p)
Definition: unuranDistr.cxx:40
const Int_t n
Definition: legend1.C:16
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8103
void Modified(Bool_t flag=1)
Definition: TPad.h:410
TAxis * GetXaxis()
Definition: TH1.h:300
Documentation for the Random class.
Definition: Random.h:39
virtual void SetLeftMargin(Float_t leftmargin)
Set Pad left margin in fraction of the pad width.
Definition: TAttPad.cxx:110
tomato 2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:290
Stopwatch class.
Definition: TStopwatch.h:28