Logo ROOT   6.18/05
Reference Guide
quasirandom.C File Reference

Detailed Description

View in nbviewer Open in SWAN Example of generating quasi-random numbers

Time for gRandom Real time 0:00:00, CP time 0.000
Time for Sobol Real time 0:00:00, CP time 0.010
Time for Niederreiter Real time 0:00:00, CP time 0.000
number of empty bins for pseudo-random = 31139
number of empty bins for sobol = 30512
number of empty bins for niederreiter-base-2 = 30512
(int) 0
#include "Math/Random.h"
#include "TH2.h"
#include "TCanvas.h"
#include "TStopwatch.h"
#include <iostream>
using namespace ROOT::Math;
int quasirandom(int n = 10000, int skip = 0) {
TH2D * h0 = new TH2D("h0","Pseudo-random Sequence",200,0,1,200,0,1);
TH2D * h1 = new TH2D("h1","Sobol Sequence",200,0,1,200,0,1);
TH2D * h2 = new TH2D("h2","Niederrer Sequence",200,0,1,200,0,1);
// quasi random numbers need to be created giving the dimension of the sequence
// in this case we generate 2-d sequence
// generate n random points
double x[2];
TStopwatch w; w.Start();
for (int i = 0; i < n; ++i) {
r0.RndmArray(2,x);
h0->Fill(x[0],x[1]);
}
std::cout << "Time for gRandom ";
w.Print();
w.Start();
if( skip>0) r1.Skip(skip);
for (int i = 0; i < n; ++i) {
r1.Next(x);
h1->Fill(x[0],x[1]);
}
std::cout << "Time for Sobol ";
w.Print();
w.Start();
if( skip>0) r2.Skip(skip);
for (int i = 0; i < n; ++i) {
r2.Next(x);
h2->Fill(x[0],x[1]);
}
std::cout << "Time for Niederreiter ";
w.Print();
TCanvas * c1 = new TCanvas("c1","Random sequence",600,1200);
c1->Divide(1,3);
c1->cd(1);
h0->Draw("COLZ");
c1->cd(2);
// check uniformity
h1->Draw("COLZ");
c1->cd(3);
h2->Draw("COLZ");
gPad->Update();
// test number of empty bins
int nzerobins0 = 0;
int nzerobins1 = 0;
int nzerobins2 = 0;
for (int i = 1; i <= h1->GetNbinsX(); ++i) {
for (int j = 1; j <= h1->GetNbinsY(); ++j) {
if (h0->GetBinContent(i,j) == 0 ) nzerobins0++;
if (h1->GetBinContent(i,j) == 0 ) nzerobins1++;
if (h2->GetBinContent(i,j) == 0 ) nzerobins2++;
}
}
std::cout << "number of empty bins for pseudo-random = " << nzerobins0 << std::endl;
std::cout << "number of empty bins for " << r1.Name() << "\t= " << nzerobins1 << std::endl;
std::cout << "number of empty bins for " << r2.Name() << "\t= " << nzerobins2 << std::endl;
int iret = 0;
if (nzerobins1 >= nzerobins0 ) iret += 1;
if (nzerobins2 >= nzerobins0 ) iret += 2;
return iret;
}
#define gPad
Definition: TVirtualPad.h:286
The Canvas class.
Definition: TCanvas.h:31
virtual Int_t GetNbinsY() const
Definition: TH1.h:293
virtual Int_t GetNbinsX() const
Definition: TH1.h:292
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3258
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2981
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4882
2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:289
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:292
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH2.h:82
Stopwatch class.
Definition: TStopwatch.h:28
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:58
void Print(Option_t *option="") const
Print the real and cpu time passed between the start and stop events.
Definition: TStopwatch.cxx:219
return c1
Definition: legend1.C:41
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
TH1F * h1
Definition: legend1.C:5
Random< ROOT::Math::GSLRngMT > RandomMT
Definition: GSLRandom.h:11
QuasiRandom< ROOT::Math::GSLQRngNiederreiter2 > QuasiRandomNiederreiter
Definition: QuasiRandom.h:179
QuasiRandom< ROOT::Math::GSLQRngSobol > QuasiRandomSobol
Definition: QuasiRandom.h:178
Author
Lorenzo Moneta

Definition in file quasirandom.C.