課題４：高エネルギーイベントのシミュレーション

[1] 概要

高エネルギーの衝突現象は複雑なので、理論に基づいたイベントのシミュレーションプログラムが発達している。ここではその一例を示す。ここで試すPythiaは高エネルギー分野で最もよく使われているプログラムである。

Pythiaのマニュアルはここにある。

[2] Pythia（ピシア）を走らせる

quark/lepton				gauge,H		meson						baryon
1	d	11	e⁻	21	g	111	π⁰	211	π⁺	411	D⁺	2212	p	2112	n
2	u	12	ν_e	22	γ	113	ρ⁰	213	ρ⁺	413	D^*+	1114	Δ^-	2114	Δ⁰
3	s	13	μ^-	23	Z⁰	221	η	223	ω	421	D⁰	2214	Δ⁺	2224	Δ⁺⁺
4	c	14	ν_μ	24	W⁺	311	K⁰	313	K^*0	423	D^*0	3122	Λ⁰	3112	∑^-
5	b	15	τ⁻	25	h⁰	321	K⁺	323	K^*+	431	D_s⁺	3114	∑^*-	3212	∑⁰
6	t	16	ν_τ			331	η'	333	φ	433	Ds^*+	3214	∑^*0	3222	∑⁺

mkdir subject-4でsubject-4/のdirectoryを作りcd subject-4で入る。
PYTHIA 8.1を開きThe current version of the PYTHIA code:（現在は pythia8157.gz, 5481rm -r py KB）をダウンロードして tar xvfz pythia8157.gz で展開するとdirectory pythia8157 ができる。
cd pythia8157の中に入り make でコンパイルする。
cd examplesの中に入りサンプルプログラム main01.cc,...., main84.cc などがあることを確認する。
make main01 で例の01をコンパイルすると main01.exe ができる。
./main01.exe > output01 で実行すると output01.htm が出来る。
最初のイベントの粒子が中間状態でできたクオークやグルーオンを含んで全てリストしてある。数百個の粒子が生成されていることが見える。各粒子の右上の表のようなKF番号（粒子の種類番号）が付いている。

[3] 粒子のラピディティ（rapidity）分布

高エネルギー衝突実験では観測量の変数の一つとしてラピディティ（rapidity）を使う。その定義と性質はRapidityとはに載っている。

LHCアトラス実験の最初の物理結果としてE_CM=900GeVでの荷電粒子のラピディティ分布が出ている（論文と右図）。

cp main01.cc main00.ccで main00.cc をつくり、その内容を以下に従って変更する：
1. pythia.init( 2212, 2212, 14000.);をpythia.init( 2212, 2212, 900.);に変更する (900 GeVに変更）。
2. if (iEvent < 1) {pythia.info.list(); pythia.event.list();} を消す。
3. for (int iEvent = 0; iEvent < 100; ++iEvent) { を　for (int iEvent = 0; iEvent < 1000; ++iEvent) { に変えて1000 eventつくる。
4. if (pythia.event[i].isFinal() && pythia.event[i].isCharged()) のうち　pythia.event[i].isCharged()の条件をはずす。
5. ++nCharged;の代わりに
```
  {
        ++nCharged;
        std::cout<<iEvent<<" "<<pythia.event[i].id()<<" "<<
        pythia.event[i].px()<<" "<<pythia.event[i].py()<<" "<<
        pythia.event[i].pz()<<" "<<pythia.event[i].e()<<endl;
      }
```
  で置き換える。
make main00.cc でコンパイルするとmain00.exeができる。ただしmain00.exeが既に存在するときはrm main00.exeであらかじめ消しておく。
./main00.exe > output00.txt で実行し、output00.txt の前後を不要な部分を削除して
```
0 321 -1.93 0.37 4.02 4.51
0 321 -0.06 -0.05 -42.18 42.18
0 -321 0.45 -0.19 -12.98 13.00
0 2212 0.60 -0.70 -15.83 15.88
0 -2212 1.13 -0.19 -27.10 27.14
0 211 -0.62 0.42 -0.90 1.18
0 211 -0.30 0.23 28.54 28.54
0 211 -0.40 -0.50 0.83 1.06
0 2212 0.18 0.15 1.86 2.09
0 -211 0.88 -1.16 -0.70 1.62
0 211 0.05 0.15 0.02 0.21
0 -211 -0.23 -0.19 0.47 0.58
.
.
```
の部分のみを残す。ただし途中で PYTHIA Warning in MultipleInteractions::pTnext: weight above unityなどの文字列が入っていることがあるので、あらかじめ PY などでそれらの文字行を探してその行を削除削除しておく。
output00.txtを読み込みrapidityを計算してプロットするROOTマクロplot_rapidity.Cをつくると右図のような荷電粒子のrapidity分布ができる。
提出課題4a：ECM=7000 GeVでのrapidity分布を粒子別（π⁺,π^-,）に重ねてプロットするグラフを作る。