モンテカルロシミュレーションによる雷雲ガンマ線放射モデルの構築
古田 禄大（東大理）
中澤知洋（東大理），和田有希（東大理／理研），楳本大悟（理研），榎戸輝揚（京大宇物），湯浅孝行（理研），土屋晴文（JAEA／理研），牧島一夫（理研），ほか GROWTH コラボレーション
2017/9/15

日本物理学会 2017年秋季大会

1

1. 雷雲起源のガンマ線放射
ロングバースト：数分続く雷雲ガンマ線•

雷雲が近づくとガンマ線量が増大•

必ずしも落雷を伴わない•

スペクトルは 10 MeV に達するべき型
カットオフの徴候も•電場加速された電子の制動放射で解釈可••ガンマ線は地上に届くまでに大気中で多数の散乱・吸収形や強度から加速器の情報が分かるはず
e.g. 放射源高度，広がり，e.g. 加速電子のエネルギー，電場強度
モンテカルロシミュレーションを使い雷雲ガンマ線をモデル化→ 観測データと比較しこれらを推定
2017/9/15＋＋＋＋＋＋ ＋

加速- e-移動

e

制動放射
-1
counts /ssec
103 Counts

地表

3000

検出器
2015/1/25 のロングバースト（柏崎原発）

2000
21分

100000Counts s-1 MeV-1•­ ­ ­­ ­e-4000
4スペクトルから探る電場加速器

­ ­

104
10
4

1:04

01:04:00

1:06

01:06:00

01:08:00
Timetime (JST)

検出器応答を含む

1000
1002
10
10

1 雷雲ガンマ線
1
環境 BGD
0.1
0.1
10
0.1
1110
日本物理学会 2017年秋季大会
Energy (MeV) 2

2. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•入力パラメータ

放射源
高度h

オフセット角度θ

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h

電子ビーム
エネルギー Ee（100‒1000 m）

ガンマ線
検出地表

0.111020.0110−410­2
2017/9/15

15 MeV, 200 m, 0 ‒10

110−3

Normalized
Counts MeV­1 m­2103. オフセット角度 θ （0 ‒45 ）

112351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)3

2. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•

入力パラメータ電子
光子

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h（100‒1000 m）

地表

0.1

1102

0.01

10−410­2
2017/9/15

15 MeV, 200 m, 0 ‒10

110−3

Normalized
Counts MeV­1 m­210

3. オフセット角度 θ （0 ‒45 ）

112

351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)4

2. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•

入力パラメータ電子光子

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h（100‒1000 m）

地表

0.1

1102

0.01

10−410­2
2017/9/15

15 MeV, 200 m, 0 ‒1015 MeV, 200 m, 20 ‒30

110−3Normalized
Counts MeV­1 m­210

3. オフセット角度 θ （0 ‒45 ）11

2351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)52. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•

入力パラメータ電子
光子

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h（100‒1000 m）

地表

0.1

11020.0110−410­2
2017/9/15

15 MeV, 200 m, 0 ‒1015 MeV, 200 m, 20 ‒3015 MeV, 500 m, 0 ‒10110−3Normalized
Counts MeV­1 m­2103. オフセット角度 θ （0 ‒45 ）11

2351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)6

2. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•

入力パラメータ電子光子

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h（100‒1000 m）

地表

0.1

11020.01

10−410­2
2017/9/15

15 MeV, 200 m, 0 ‒1015 MeV, 200 m, 20 ‒3015 MeV, 500 m, 0 ‒1030 MeV, 200 m, 0 ‒10110−3Normalized
Counts MeV­1 m­2103. オフセット角度 θ （0 ‒45 ）

112351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)7

2. Geant4 による大気散乱シミュレーション•

単一エネルギーの加速電子を点源から真下に照射•

地表でガンマ線のオフセット角度とエネルギーを記録•

入力パラメータ電子光子

1. 電子エネルギー Ee （10‒30 MeV）
2. 放射源高度 h（100‒1000 m）

地表

0.1

1102

0.01

10−410­2
2017/9/15

15 MeV, 200 m, 0 ‒1015 MeV, 200 m, 20 ‒30
15 MeV, 500 m, 0 ‒1030 MeV, 200 m, 0 ‒10110−3Normalized
Counts MeV­1 m­210

3. オフセット角度 θ （0 ‒45 ）スペクトルから放射源の情報を引き出せる
112351010

日本物理学会 2017年秋季大会

2030 Energy (MeV)
8•A，Γ を実データと比較•

情報量が少なく一意に決まらず

Counts / m / MeV / Primary

従来の手法（e.g. 古田+17 春）• べき関数 AE で低エネルギー側を近似

CountsMeV­123. 解析的な関数による近似

10−510−2610

べき関数によるフィット

10−7

110−810−910­2100.3−1011

新たな近似関数3 4 5 6 7 8 10310
Energy (MeV)
Energy (MeV)
2↵•

カットオフも近似できる新たな関数を提案： AE [1 (E/Ee ) ]α はカットオフの始まりの位置，β はカットオフの鋭さを反映
10210−310­2 ­2100.1

2017/9/15100.010.11101α 固定
1 べき関数のみ
β = 5.0
β = 3.0
β = 1.5
β = 0.50.1

β 固定
1 べき関数のみ
α = 0.5
α = 1.0
α = 2.0
α = 4.00.01

Counts MeV­110210­1110110­2 ­21010−3•0.1

日本物理学会 2017年秋季大会

10­1
1110

E / Ee9•A，Γ を実データと比較•

情報量が少なく一意に決まらず

Counts / m / MeV / Primary

従来の手法（e.g. 古田+17 春）• べき関数 AE で低エネルギー側を近似

CountsMeV­123. 解析的な関数による近似

10−510−2610べき関数によるフィット

10−7110−810−910­2100.3−1011

新たな近似関数••3 4 5 6 7 8 10310Energy (MeV)
Energy (MeV)
2↵

カットオフも近似できる新たな関数を提案： AE [1 (E/Ee ) ]α はカットオフの始まりの位置，β はカットオフの鋭さを反映
1010­910−620 MeV, 1000 m, 0 ‒101010­710­1110−10

normalized counts s−1 MeV−110−9Γ = ­0.3
α = 1.4
β = 1.610−810­910−9410−1144220−2−4

2017/9/15

Γ = ­1.3
α = 1.0
β = 3.510−1040­420 MeV, 200 m, 30 ‒45−7(data−model)/error

normalized counts s−1 MeV−1−8(data−model)/error

χCountsMeV­1 m­210−71123

Energy (MeV)

510102000−2­4−411

日本物理学会 2017年秋季大会

2310
Energy (MeV)10

Energy (MeV)

51020

•

A，Γ を実データと比較•

情報量が少なく一意に決まらず

Counts / m / MeV / Primary

従来の手法（e.g. 古田+17 春）• べき関数 AE で低エネルギー側を近似

CountsMeV­12

3. 解析的な関数による近似

10−510−2610

べき関数によるフィット

10−7

1

10−810−910­210
0.3−1011

新たな近似関数••3 4 5 6 7 8 10310
Energy (MeV)
Energy (MeV)
2↵

カットオフも近似できる新たな関数を提案： AE [1 (E/Ee ) ]α はカットオフの始まりの位置，β はカットオフの鋭さを反映
1010­910−620 MeV, 1000 m, 0 ‒101010­710­1110−10

normalized counts s−1 MeV−110−9Γ = ­0.3
α = 1.4
β = 1.610−810­910−9410−1144220
解析的な関数でよく近似できた

0−2−4

2017/9/15

Γ = ­1.3
α = 1.0
β = 3.510−1040­420 MeV, 200 m, 30 ‒45−7(data−model)/error

normalized counts s−1 MeV−1−8(data−model)/error

χCountsMeV­1 m­210−71123Energy (MeV)

51010200−2­4−411

日本物理学会 2017年秋季大会

2310
Energy (MeV)11

Energy (MeV)

510204. 関数のパラメータの振る舞い AE [1↵(E/Ee ) ]A は高度 h の 〜 6 乗に逆比例• Ee が大α が小；相対的により低エネルギー側までカットオフ••

θ が大
Γ が小；ソフトなスペクトル（制動放射のビーミング効果）
θ が大
β が大；なだらかなカットオフ
0 deg < θ < 10 deg (Energies in different colors.)
10 deg < θ < 20 deg (Energies in different colors.)−310Normalization

Normalization

A
101010­5−410−4Ee (MeV)
1010 15 20 30 10−5−61010−7−5−610−71010­8−810−8−910−1010−910−10101010­1190100
100−11−310200

300

400

300

10−1190100
1000 高度 (m)1000
Height (m)20 deg < θ < 30 deg (Energies in different colors.)

Normalization

Normalization−31010−4−510−6−310300

400

1000
Height (m)30 deg < θ < 45 deg (Energies in different colors.)

10−4−510−6101010−710−7−810−8−910−101010−910−101010−1120090100

2017/9/15

200

300

400

1000
Height (m)10−1190100

200

日本物理学会 2017年秋季大会

300

400

1000
Height (m)12↵(E/Ee ) ]A は高度 h の 〜 6 乗に逆比例• Ee が大α が小；相対的により低エネルギー側までカットオフ••

θ が大
Γ が小；ソフトなスペクトル（制動放射のビーミング効果）
θ が大
β が大；なだらかなカットオフ
0 deg < θ < 10 deg (Energies in different colors.)
10 deg < θ < 20 deg (Energies in different colors.)−310Normalization

Normalization

A101010­5−4−31010−4Ee (MeV)
1010 15 20 30 10−5−61010−7−5−610−71010­8−810−8−910−1010−910−1010−11−11−3−3Normalization

Normalization
Ratio Index

101010­1190100
200300 400
200300 4001000
1001000
高度 90(m)Height(m)Height (m)100
1000
10 deg < θ < 20 deg 300(Heights in different colors.)
2.4 20 deg < θ < 30 deg (Energies in different colors.)
30 deg < θ < 45 deg (Energies in different colors.)1010
α2.2h (m) 100 20010
102
2.010300 500 100010
1.8−4−5−610
1.6
10−7
1.4
−810
1.2
−910
1−1010
0.8
10−11

0.8

−4−5−610

10−7

1.4

−810

−910−10109 1090100

E (MeV)2030
10−11
Electron Energy 1000(MeV)
40090100
200300eHeight (m)30deg<θ<45deg(Heightsindifferentcolors.)
2017/9/15
日本物理学会 2017年秋季大会
4dexV)

4. 関数のパラメータの振る舞い AE [110200300

20304001000
Height (m)13↵(E/Ee ) ]A は高度 h の 〜 6 乗に逆比例• Ee が大α が小；相対的により低エネルギー側までカットオフ••

Normalization

Normalization

A101010­5−4−31010−4Ee (MeV)
1010 15 20 30 10−5−61010−7−5−610−71010­8−9−1.2−10−1.4­1.4−910−101010

10.50−0.8­1.0−110300 500 1000

−0.6
­0.6−8−810−0.2

Γ−0.4

Photon Index−310

Photon Index

θ が大
Γ が小；ソフトなスペクトル（制動放射のビーミング効果）
θ が大
β が大；なだらかなカットオフ
0 deg < θ < 10 deg (Energies in different colors.)
10 deg < θ < 20 deg (Energiesin different
10 MeV(Heights colors.)in differentcolors.)h (m)100 200−0.5−1−1.5−3−4−5−610
1.6
10−7
1.4
−810
1.2−9101−10100.8
10−11−5−60
3

10−7−2.5
0.5
2.5

10

1.4

0.8

−4CutoffIndex
PhotonIndex−3h (m) 100 200
300 500 1000

−810−910−1010−211.5
1.5−1.5

12030
51015202530
3540
45−11−20100
5 300 10400 15
20251000
30 Angle
35 (degree)
4045
Electron Energy 1000(MeV)
200
300 400
90100
200Offsete
Offset Angle (degree)
Height (m)Height (m)30 deg < θ < 45 deg (Heights in different colors.) 日本物理学会 2017年秋季大会 20 MeV (Heights in different colors.)
2017/9/15
64

1020

30

E (MeV)0dex

9 1090100

03
−0.5

1530

2−11−1.5Oﬀset θ
14dex−11CutoffIndex
PhotonIndex−11

NormalizationNormalization
Ratio Index−1.6−2101010­1105 300 10400 152025100030
354045090100
200300 40090100
200
1000
高度(m)
Offset Angle(degree)
Height(m)Height30(m)15Oﬀset
θ
300
0100
1000
10 deg < θ < 20 deg (Heights in different colors.)
10 MeV (Heights in different colors.)
2.4 20 deg < θ < 30 deg (Energies in different colors.)
620 MeV(Heights colors.)in different colors.)
30 deg < θ <4.545 deg (Energiesin different
111010βα2.245h (m) 100 20010
100.50.522.0
3.5
3.510300 500 100010
41.8dexV)

4. 関数のパラメータの振る舞い AE [10−2006counts / sec

2015/1/25 のロングバースト

4000
4103 Counts s−1

5. 地上観測データへの適用

3000

•

過去最大の光子統計•

原発のモニタリングポストでも検出•

モニタリングポストでの最大光度から放射源の通過位置を推定
00→ 検出器の直上を通過（楳本+ in prep.）

2000
2

•

X 線天文のスペクトル解析ツール
XSPEC（NASA）を使用•

光度曲線の FWHM 内 45 秒間からデータ抽出し近似関数でフィット•θ は十分小 → β = 2 に ﬁx

1000time (JST)100

1010144200−2­4−411

2017/9/15

01:08:00
Timeスペクトル（検出器応答を含む）

102
normalized counts s−1 MeV−1↵(E/Ee ) ]2分
01:06:00
1:06

01:04:00
1:04

103(data−model)/error

AE [1χ

近似関数によるフィット

Counts s­1 MeV­11000

日本物理学会 2017年秋季大会

2310Energy (MeV)
510Energy (MeV)15

counts / sec

2015/1/25 のロングバースト

4000
4103 Counts s−1

5. 地上観測データへの適用

3000

•

過去最大の光子統計•

原発のモニタリングポストでも検出•

モニタリングポストでの最大光度から放射源の通過位置を推定
00→ 検出器の直上を通過（楳本+ in prep.）

光度曲線の FWHM 内 45 秒間からデータ抽出し近似関数でフィットχθ は十分小 → β = 2 に ﬁx

1000

time (JST)100

1010Red. χ2 = 1.21 (430)144200−2­4−4実データも近似関数でよく説明できた
2017/9/15

01:08:00
Timeスペクトル（検出器応答を含む）

102
normalized counts s−1 MeV−1

X 線天文のスペクトル解析ツール
XSPEC（NASA）を使用

2分
01:06:00
1:06

01:04:00
1:04

103(data−model)/error•↵(E/Ee ) ]Counts s­1 MeV­11000

AE [1•（FWHM）
45 秒

2000
2

近似関数によるフィット•

データ抽出時間帯

日本物理学会 2017年秋季大会

1
1

2

310
Energy (MeV)
510

Energy (MeV)16

6. 放射源情報の推定
放射源高度と電子数••

Γ 〜 ­0.8 → 高度 200‒500 m
冬季雷の典型的な雲底高度と無矛盾（200‒800 m；Goto+1992）

AE [1↵(E/Ee ) ]A = 1.50 (0.97‒2.40)
Ee = 21.6 (20.6‒22.7) MeV
Γ = ­0.80 (0.90‒0.56)
α = 0.45 (0.26‒0.61)
0Γ θ 0 ‒10 10 ‒20­0.4 20 ‒30 30 ‒45

シミュレーションとフラックスを比較
A 〜 1.5，Ee 〜 20 MeV，h 〜 300 m­0.8→ 加速電子数 〜 5 108 s­1（cf. 土屋+11 での推定 107‒10 s­1）
­1.2

− 0.2
− 0.4
− 0.6
− 0.8−1− 1.2− 1.4

実データから物理的に妥当な解が得られた ­1.6

− 1.6− 1.8

今後の進展
•

100102300

低統計のデータに適用するには観測との両輪が重要
1. 数十 MeV での有効面積拡大・ダイナミクスレンジの確保
2. 多地点観測でオフセット角度推定・三角測量で高度測定

•

31000
高度 (m)
10

中澤+次公演

電子エネルギーに広がりがある場合 → 今後検証
（原理的にはモデルの畳み込みで再現可能）2017/9/15

日本物理学会 2017年秋季大会

17

7. まとめ•

過去の新潟県柏崎刈羽原発での地上観測により，雷雲ガンマ線の放射源の水平位置は概ねわかってきたが，高さや加速電子のエネルギーは未知。

•

カットオフを含むスペクトル全体を使って放射源の情報を推定するため，モンテカルロシミュレーションを使い，雷雲ガンマ線をモデル化した。

•

カットオフの形状も表現できる解析的な関数で近似し，放射源の高度・電子エネルギー・オフセット角度と関数パラメータの対応関係を得た。

•

過去最大の光子統計があったロングバーストのスペクトルをこの関数でフィットし，そのパラメータから，電子エネルギー 〜 20 MeV，放射源
高度 200‒500 m，加速電子数 〜 5 108 s­1 を推定した。

•

低統計のデータでは推定精度が落ちるが，十分な感度とダイナミックレンジをもった検出器を設置したり，多地点観測網で情報量を拡大することで補える。
2017/9/15

日本物理学会 2017年秋季大会

18Appendix

2017/9/15

日本物理学会 2017年秋季大会

19雷雲電場による電子加速〜 100 keV の種電子と数百 kV/m の電場があれば，電子は数十 MeVまで加速されうる
100

電子の受ける制動力（MeV/m）•〜 100 keV
種電子1010.1

2017/9/15

数十 MeV
逃走電子
過去に雷雲で測定された
典型的な電場強度

電場eE

10­30.1103
1 10電子の運動エネルギー（MeV） 電場から受ける力 ＞ 制動力
日本物理学会 2017年秋季大会

20

Geant4 による大気散乱シミュレーション•

加速電子を模擬し，上空の一点源から真下方向に電子ビームを照射•

異なる電子エネルギー・高度でシミュレーション•

地表でのガンマ線の到達位置とエネルギーを記録し，オフセット角度で分類＜入力パラメータ＞：組み合わせ 80 通り
電子エネルギー（MeV） 10，15，20，30
放射源高度（m）

100，200，300，500，1000

オフセット角度（度）0‒10，10‒20，20‒30，30‒45•

電磁相互作用の Physics List は Livermore を使用•

Cut Range はデフォルト（0.7 mm）
空気中で 1 keV までのトラッキングに相当•

空気密度 1.205

10­3 g cm­3電子ビーム
エネルギー Ee
放射源
高度
h

オフセット
角度θ
ガンマ線
検出地表2017/9/15

日本物理学会 2017年秋季大会

21