【C++/ROOT】C++11でMersenne Twisterを使う / ROOT6でhistogramを描画する

C言語ではメルセンヌツイスターになじみがあったのですが、C++ではなかったのでちょっとやってみます。

はじめに

今回のコードはこちらに置いてあります。
ダウンロードして使うほどのコードではないですが、GitHubを覚えるためなので大目に見てください。
github.com

メルセンヌツイスター

まず以下のライブラリをincludeする

#include <random>

乱数の生成器は以下の２行で設定します。

random_device seed_gen;
mt19937_64 engine(seed_gen()); //メルセンヌツイスター法の乱数生成器

最後に乱数の分布について、以下のクラスが利用できます。

std::uniform_int_distribution  //一様整数乱数
std::uniform_real_distribution  //一様実数乱数
std::normal_distribution  //正規乱数
std::discrete_distribution  //確率を指定した乱数

一様乱数

/*
C++11のメルセンヌツイスターで一様乱数を生成する
*/
#include <random>
#include <iostream>
#include <iomanip>
using namespace std;

int main(int, char**){

  random_device seed_gen;
  mt19937_64 engine(seed_gen());
  uniform_real_distribution<double> dist_r(0., 1.); //dist_r(最小値, 最大値)

  for (int i = 0; i < 10000; ++i){
    double uni_real_num = dist_r(engine);
    cout << fixed << setprecision(5) << uni_real_num << endl;
  }

}

ROOT6でhistogramを描画する

値だけではそのデータの傾向なんてわかりません。
なので、このC++のコードをROOTでhistogramするために以下のように書き換える。

/*
C++11のメルセンヌツイスターで一様乱数を生成->ROOTで描画
*/
#include <random>
#include <iostream>
#include <iomanip>
using namespace std;

void mt_uni_real_dist(){

  random_device seed_gen;
  mt19937_64 engine(seed_gen());
  uniform_real_distribution<double> dist_r(-0.5, 0.5); //dist_r(最小値, 最大値)

  TCanvas *c1 = new TCanvas("c1", "c1");
  TH1D* h1 = new TH1D("mt_uniform_real_dist","", 200, -1, 1); // [-1, 1]に200個binを用意

  for (int i = 0; i < 10000; ++i){
    double uni_real_num = dist_r(engine);
//    cout << fixed << setprecision(5) << uni_real_num << endl;

    h1->Fill(uni_real_num);  //bin詰め
  }

  h1 ->Draw("");  //output histgram
  c1 ->Print(".pdf");  //output pic
  c1 ->Print(".png");  //output pic
}

ROOTで出力したhistogramがこちら

メルセンヌツイスターから[-0.5, 0.5]の一様乱数を生成しROOTで描画したhistogram

正規乱数

#include <random>
#include <iostream>
#include <iomanip>
using namespace std;

int main(int, char**){
  
  random_device seed_gen;
  mt19937_64 engine(seed_gen());
  normal_distribution<double> dist_r(0, 1);

  for (int i = 0; i < 1000; ++i){
    double real_num = dist_r(engine);

    cout << fixed << setprecision(5) << real_num << endl;

  }
}

実行する。

$ g++ -std=c++11 mt_nomal.cpp
$ ./a.out 
0.69976
0.85291
-1.19238
-1.69977
-0.40224
0.33892
-0.23855
0.04096
0.60685
0.70412
....

ランダムなseedを設定しているので、値は合わないはずです。

ROOT6でhistogramを描画する

値だけではそのデータの傾向なんてわかりません。
なので、このC++のコードをROOTでhistogramするために以下のように書き換える。

/*
C++11のメルセンヌツイスターで正規乱数を生成->ROOTで描画
*/
#include <random>
#include <iostream>
#include <iomanip>
using namespace std;

void mt_nomal_dist(){
  
  random_device seed_gen;
  mt19937_64 engine(seed_gen());
  normal_distribution<double> dist_r(0., 1.); //dist_r(平均値, 標準偏差)

  TCanvas *c1 = new TCanvas("c1", "c1");
  TH1D* h1 = new TH1D("mt_nomal_dist","",100, -5, 5); // [-5, 5]に100個binを用意

  for (int i = 0; i < 10000; ++i){
    double real_num = dist_r(engine);
//    cout << fixed << setprecision(5) << real_num << endl;

    h1->Fill(real_num);  //bin詰め
  }

  h1 ->Draw("");  //output histgram
  c1 ->Print(".pdf");  //output pic
}

ROOTで出力したhistogramがこちら

メルセンヌツイスターからmean=0.0, sigma=1.0で生成した正規乱数のhistogram

おわりに

最初はgistでコードを張り付けながら書いていたのですが、あまり見やすくなかった気がしたのではてな記法を使いました。

Gistを使ったほうが効率的ですが、notebookファイルの時にのみ利用するのがbetterみたいです。

jupyterにC++のカーネルを入れるのもありかなぁ。

相変わらずPDFをpngに直すのが面倒ですね。

googleドライブ使えばいいっぽいけど、ちょっとあれはまだよくわからんです。