物理の駅 by onsanai

Physics-station 研究で日々感じたことを忘れないための備忘録

直交座標系のずれ量を飛跡の角度に相当するRadial Lateral空間に変換するコード

放射線の飛跡のように、座標(cx, cy) 角度 (ax, ay)で表すことの出来る量の比較をする時、しばしば飛跡の進行方向(Radial方向)と垂直方向(Lateral方向)に分離して考えることがある。

飛跡1 (cx1, cy1, ax1, ay1)、飛跡2 (cx2, cy2, ax2, ay2) があるとき、角度差は次のように求めることができる。

import math

# 方法1
def xy2radial(cx1,cy1,ax1,ay1,cx2,cy2,ax2,ay2):
    ax = (ax1 + ax2) * 0.5
    ay = (ay1 + ay2) * 0.5
    r = (ax ** 2 + ay ** 2) ** 0.5
    ax /= r
    ay /= r
    dax = ax2 - ax1
    day = ay2 - ay1
    dar = ax * dax + ay * day
    dal = -ay * dax + ax * day
    dcx = cx2 - cx1
    dcy = cy2 - cy1
    dcr = ax * dcx + ay * dcy
    dcl = -ay * dcx + ax * dcy

    return dcr,dcl,dar,dal

# 方法2
def xy2radial2(cx1,cy1,ax1,ay1,cx2,cy2,ax2,ay2):
    ax = (ax1 + ax2) * 0.5
    ay = (ay1 + ay2) * 0.5
    rotation_cos = math.cos(-math.atan2(ay,ax))
    rotation_sin = math.sin(-math.atan2(ay,ax))
    dax = ax2 - ax1
    day = ay2 - ay1
    dar = dax * rotation_cos - day * rotation_sin
    dal = dax * rotation_sin + day * rotation_cos
    dx = cx2 - cx1
    dy = cy2 - cy1
    dcr = dx * rotation_cos - dy * rotation_sin
    dcl = dx * rotation_sin + dy * rotation_cos

    return dcr,dcl,dar,dal

処理時間は若干方法1の方が早かった。三角関数の処理は、CPUで実装されている? らしいのでそんなに遅くもない。

f:id:onsanai:20190808123417p:plain

参考文献 https://dx.doi.org/10.1093/ptep/ptx131

高い飛跡密度における飛跡再構成アルゴリズム (日本語訳)

arxiv.org

の付録Bの日本語訳です。一部に筆者の意訳を含みます。

原子核乾板は、その高分解能により各フィルムに記録された膨大な量のイベント又は10^ 5 tracks/cm^ 2オーダーの飛跡の陽子反応を再構成することができる。 これは、位置分解能0.4 μmと角度分解能2mradのベクトル情報を持つベーストラックの性質にも依っている。

飛跡再構成の基本的概念は、位置及び角度空間における異なるフィルム上のベーストラックの対応付に基づいている。 広く用いられてきたアルゴリズムでは、2つの連続したベーストラックの対応テストが用いられている。 しかし、DsTau実験のように高い飛跡密度環境では再構成が失敗することがある。 特に、数mrad以内の同じ方向に進む2つ以上の飛跡が数ミクロン以内で接近すると、このアルゴリズムは正しい経路を決定できないことがある。

新しいトラッキングアルゴリズムは、高い飛跡密度と狭い角度広がり環境で飛跡を再構成するために開発した。 複数の経路候補が生じた時、まずは全ての可能な経路を生成する。 図21-(a)の場合、Z _ 1からZ _ 2間で2^ 4=16の可能な経路が考えられる。 それぞれの経路に対して、経路に含まれるベーストラックによって作られた面積の平均に基づいた次のテスト変数を評価する。

 \displaystyle
a^{average}=\left(\sum_i^{n-2}a^{pos}+\sum_i^{n-1} a^{angle}\right)/(n-0.5).

ここで、a^ {average}は図21-(b)で示すようにベーストラックの位置a^ {pos}と角度a^ {angle}で作られる平均面積である。nは経路に含まれるベーストラック数で、-0.5は長い経路に重みを付けるための経験値である。最小のa^ {average}を持つ経路が最良の経路として選ばれる。 この手順を、選択した経路に含まれるベーストラックを取り除きながら繰り返す。

 \displaystyle
a^{pos}=\frac{1}{2}\cdot \left|\vec{AB} \times \vec{AC}\right|

 \displaystyle
a^{angle}=\frac{1}{2}\cdot \frac{\left|\vec{AB}\right|}{2} \cdot \left|(\vec{V_A}+\vec{V_B}) \times \frac{\vec{AB}}{2}\right|

f:id:onsanai:20190731104845p:plain f:id:onsanai:20190731104855p:plain

pythonのmatplotlibで2次元ヒストグラム(hist2d)のビンの値を直接操作する

pythonのmatplotlibで2次元ヒストグラム(plt.hist2d)のビンの値 (bin contents)を直接操作する方法はない。なので、 colormesh を使って描画しよう。

import matplotlib.pyplot as plt
import numpy as np

def f(x, y):
    return np.exp(-(x * x) / (2 * 3 ** 2)) * np.exp(-(y * y) / (2 * 3 ** 2))

X, Y = np.mgrid[0:21:1.0, 0:21:1.0]

X-=10.0
Y-=10.0
Y*=2.0

Z = f(X,Y)

fig, ax = plt.subplots()
im = ax.pcolormesh(X, Y, Z, cmap='inferno',vmin=0,vmax=1.1)
ax.set_xlabel("X")
ax.set_ylabel("Y")
cbar = fig.colorbar(im)
cbar.set_label("Z")
plt.show()

このようにしてカラーマップを作成することが出来る。

f:id:onsanai:20190703085418p:plain
二次元ヒストグラムっぽい絵

普通のヒストグラム(一次元ヒストグラム)も、棒グラフを使ってヒストグラムっぽく表示することが出来る。

import matplotlib.pyplot as plt
import numpy as np
 
X = np.array([i for i in range(-10,10)])
Y = np.exp(-(X * X) / (2 * 3 ** 2))
plt.bar(X, Y, width=1.0)
plt.show()

f:id:onsanai:20190430233014p:plain
一次元ヒストグラムっぽい絵

Windows10とVisual Studio 2017でGeant4を動かした

Overview | geant4.web.cern.ch

WindowsでGeant4を入れたお話。

皆さん誤解しているかもしれませんが、WindowsでGeant4は動きます。GUIや、他のLinuxにしか対応していないツールや、マルチスレッド関連を除けば、ちゃんと動きます。 深いことをやろうとすると動かなくなることがありますが、とりあえず構造を作って粒子を打ち込むのであればできる。

と思っていたわけだ。こんなに苦労するとは思わなかった(8時間ほど)。詰まったのは G4LIB_BUILD_DLLプリプロセッサに追加するところだった。一通り導入方法を書きます。質問があればコメントでどうぞ。

Geant4のWindows + Visual Studioでの導入方法

Geant4の入手

まず、ダウンロードページからデータファイルとコンパイル済みのライブラリを入手しましょう。

http://geant4.web.cern.ch/support/download

「compiled using Visual Studio 2017 (version 15.9.3) on Windows 7, 32 bits, executable installer」と書いてありますが、Windows10でも動きます。

コンパイル済みのライブラリは C:\Program Files (x86)\Geant4 10.5 にインストールしたとします。初期設定でそうなっていると思います。 変更する点は次の絵のPATHを通すところ。とりあえず下記のようにしておこう。

f:id:onsanai:20190412153006p:plain:w300
Geant4 Windows10 Install Options

ソースファイルの入手

インストールした先のexamplesにサンプルコードはありますが、分かりにくいです。全然ベーシックじゃねぇです。KEKのGean4初心者講習の実習用ソースコードの方が分かりやすいです。KEK万歳です。 Geant4Tutorial20171129\TutorialMaterials\P01_FirstStep\source にあります。 直下にある *.cc と、srcにある *.cc includeにある *.hh を使います。

KEK wiki https://wiki.kek.jp/display/geant4/Geant4+Japanese+Tutorial+for+Detector+Simulation+2017

ソリューション、プロジェクトを作る

Visual Studio 2017で、適当な名前を付けて作ってください。

先ほどの *.cc *.hh を全部コピーして、ソースファイルとヘッダーファイルに突っ込みます。正しい流儀ではありませんが、気にしたら負けです。

プロジェクトの設定

Visual Studio内で、インクルードパスを設定し、ライブラリファイルを設定します。

  • インクルードパス: プロパティ→構成プロパティ→VC++ディレクトリ→インクルードディレクトC:\Program Files (x86)\Geant4 10.5\include\Geant4
  • ライブラリファイル: プロパティ→構成プロパティ→リンカー→入力→追加の依存ファイル C:\Program Files (x86)\Geant4 10.5\lib\*.lib

適宜インストール先に置きかえてください。

いくつかのプリプロセッサを設定する必要があります。プリプロセッサは プロパティ→構成プロパティ→C/C++プリプロセッサプリプロセッサの定義に、

G4LIB_BUILD_DLL
WIN32
_CRT_SECURE_NO_WARNINGS

を設定してください。

上から、DLLを呼び出す(dllimport を使う)設定、Windowsのヘッダにするための設定、古いCの関数でエラーを出さないための設定です。

環境変数

環境変数のpathに C:\Program Files (x86)\Geant4 10.5\bin を追加しましょう。

データファイルの環境変数は適宜下記のものを設定してください。使う物理リストによっては他の変数を要求されることもあります。

G4ENSDFSTATEDATA=C:\G4ENSDFSTATE2.2
G4LEDATA=C:\G4EMLOW7.7
G4LEVELGAMMADATA=C:\PhotonEvaporation5.3
G4SAIDXSDATA=C:\G4SAIDDATA2.0
G4PARTICLEXSDATA=C:\G4PARTICLEXS1.1

以前は G4NEUTRONXSDATA=C:\G4NEUTRONXS1.4 が必要だったが、G4PARTICLEXSDATAに統合された?

コードの修正

下記のようなエラーが出ると思います。

Geant4\Geometry.cc(43): error C2039: 'Invisible': 'G4VisAttributes' のメンバーではありません。
Geant4\G4VisAttributes.hh(68): note: 'G4VisAttributes' の宣言を確認してください
Geant4\Geometry.cc(43): error C2065: 'Invisible': 定義されていない識別子です。

Geometry.cc の43行目を下記のように修正してください。

   logVol_World->SetVisAttributes (G4VisAttributes::GetInvisible());

動かしてみる

P1というサンプルコードを動かすと、次のようにコンソールに表示されて止まります。ここにコマンドを入れていくと、粒子を買えたりビームを打ったりすることができます。warningが出ていますが、とりあえず気にしません。なぜ出るのかエラい人教えてください。

**************************************************************
 Geant4 version Name: geant4-10-05    (7-December-2018)
                       Copyright : Geant4 Collaboration
                      References : NIM A 506 (2003), 250-303
                                 : IEEE-TNS 53 (2006), 270-278
                                 : NIM A 835 (2016), 186-225
                             WWW : http://geant4.org/
**************************************************************

<<< Geant4 Physics List simulation engine: FTFP_BERT 2.0


 FTFP_BERT : new threshold between BERT and FTFP is over the interval
 for pions :   3 to 12 GeV
 for kaons :   3 to 12 GeV
 for proton :  3 to 12 GeV
 for neutron : 3 to 12 GeV

### Adding tracking cuts for neutron  TimeCut(ns)= 10000  KinEnergyCut(MeV)= 0
Available UI session types: [ Win32, GAG, csh ]

-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : UI0002
      issued by : G4UIExecutive::G4UIExecutive()
Specified session type is not build in your system,
or no session type is specified.
A fallback session type is used.
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Idle>

講義資料の「ユーザ・アプリケーション作成の基礎」に使い方が載っています。コマンドを実際に打ちながら、色々動かしてみましょう。

https://wiki.kek.jp/display/geant4/Geant4+Japanese+Tutorial+for+Detector+Simulation+2017

ソースコードを書く

好きにソースコードを書きましょう。例えば、Geant4でHello World!を書くには次のようにするそうです。

#include <G4coutDestination.hh>

int main() {
    G4cout << "Hello world!" << G4endl;
}

IO部分はとりあえずの実装で良いなら、 std::coutstd::ofstream を使っても動きます。

では、良いシュミレーションライフを。

講演会 川村静児氏「重力波:アインシュタインの奏でる宇宙からのメロディー」質疑応答 @2018年 ぎふサイエンスフェスティバル

会場内でメモしたため間違っているところがあるかもしれません。

Q. インフレーション理論で物理現象(特殊相対性理論)は成立しているのか?
A. 場が広がってるだけなので成立できる。
Q. インフレーションの音はシミュレーションはされているのか?
A. されているが、まだ言わない。
Q. 重さがなくなったのは重力子の質量になったと考えられないのか?
A. 重力波の媒介粒子は質量はないと考えられているので、ほぼ全てがエネルギーになったと考えるのが妥当。
Q. ブラックホール連星では他の情報は出てこないのか?
A. ブラックホール連星の合体では出てこないとされている。中性子星連星の合体では他の情報も出てくる。
Q. 光と重力波は何が違うのか?
A. 媒介する力が違う。また、電磁波はベクトル波、重力波はテンソル波。
Q. 重力子と重力波は同じなのか?
A. 重力波も粒であり波であるから、重力子と考えても良いと思う。

注: インフレーションで発生する重力波はsub-Hz帯なので、実は人間の耳には聞こえないんです。中性子連星やブラックホール連星での合体で発生する重力波は たまたま 人の可聴領域内だっただけなのです。

http://www.city.gifu.lg.jp/17404.htm

川村さん、準備と講演、お疲れさまでした。

SRIM 2013を日本語版 Windows 10 64bitで動かすためのメモ

とても個人的なメモなので、使えなくても泣かないお約束で。

SRIM本体は下記のURLからインストールする。

http://www.srim.org/SRIM/SRIMLEGL.htm

Msvbvm50.dll がないと怒られるので、下記URLからMsvbvm50.exeダウンロードし実行する

https://support.microsoft.com/ja-jp/help/180071/file-msvbvm50-exe-installs-visual-basic-5-0-run-time-files

RICHTX32.OCX がなんとかってエラーが出るので、下記URLから Visual Basic 6.0 Service Pack 6 をダウンロード

Download Visual Basic 6.0 Service Pack 6 from Official Microsoft Download Center

適当なフォルダに展開し、さらにフォルダ内の RichTx32.CAB を展開する。

RICHTX32.OCX というファイルがあるので、 C:\Windows\SysWOW64 にコピーし、管理者権限のコマンドプロンプト

cd C:\Windows\SysWOW64
regsvr32 RICHTX32.OCX

を実行する。同様の処理をさらに4つのファイルで行う

  • COMDLG32.CAB
  • MSFLXGRD.CAB
  • TABCTL32.CAB
  • COMCTL32.CAB
regsvr32 COMDLG32.OCX
regsvr32 MSFLXGRD.OCX
regsvr32 TABCTL32.OCX
regsvr32 COMCTL32.OCX

これで、 SR.exe (静的にdE/dx やStragglingを計算するソフト)や、 TRIM.exe (動的に粒子を打ち込むソフト)が動くはずである。

TRIM.exe を実行するための TRIM.in ファイルを作るための TIN.exe は動かない。 SRIM.exe に特段の機能は無い。

Geant4で原子核乾板 (10の12乗チャンネル) を実装した

ミクロンの空間分解能を持ちながら、ミリメートルからメートルのサイズを持つ放射線検出器は原子核乾板以外に存在しない。通常、チャンネル数はある特定の数になるが、原子核乾板の場合はそういう概念はないに近い。1平方センチあたり10の14乗を程度チャンネル数が存在するとみなしてもよい。

高エネルギー粒子をシミュレーションするときは、乾板の表面に入ってから出ていくまでを直線と仮定して差し支えないが、低エネルギー粒子の飛跡をシミュレーションで再現しようとするときは、PhysicsListが低エネルギー側でどの程度現実の物理を反映しているかどうかは別にして、それなりに細かいStepで飛跡の位置情報を出力したいと誰しも思うだろう。

Geant4でこのような検出器を実装するには、 G4PVReplica という機能を使う。

G4PVReplicaG4VPhysicalVolume から継承されているクラスで、先に作った箱(例えば G4Box)をある方向に分割する機能である。X軸、Y軸、Z軸の3回分割を繰り返せば、細かい要素に箱を区切ることができる。シミュレーションしたい内容に応じて、細かさは調整すればよい。 コメントで指摘があったとおり、最小要素(Elmt)を四角柱(Stick)に詰め、それを平面に伸びる四角柱(Srfc)に詰め、さらに箱(Box)に詰めている、という考え方の方が分かりやすい。

次の例では、簡単のために一辺100ミリメートルの立方体を、一辺10ミクロンの立方体で区切る。計10の12乗のチャンネル数がGeant4にできていることになる。

G4double leng_EmulsionBox = 100.0*mm;     // X-full-length of emulsion box 
G4double leng_EmulsionElmt = 0.01*mm;     // X-full-length of emulsion element
G4int nDiv = std::round(leng_EmulsionBox / leng_EmulsionElmt);

//Emulsionを格納するBOX
auto* solid_EmulsionBox =
    new G4Box("Solid_EmulsionBox", leng_EmulsionBox / 2.0, leng_EmulsionBox / 2.0, leng_EmulsionBox / 2.0);
auto* materi_EmulsionBox = materi_Man->FindOrBuildMaterial("G4_AIR");
auto* logVol_EmulsionBox =
    new G4LogicalVolume(solid_EmulsionBox, materi_EmulsionBox, "LogVol_EmulsionBox", 0, 0, 0);

//XYに伸びる四角柱(平面)
auto* solid_EmulsionSrfc =
    new G4Box("Solid_EmulsionSrfc", leng_EmulsionBox / 2.0, leng_EmulsionBox / 2.0, leng_EmulsionElmt / 2.0);
auto* materi_EmulsionSrfc = materi_Man->FindOrBuildMaterial("G4_AIR");
auto* logVol_EmulsionSrfc =
    new G4LogicalVolume(solid_EmulsionSrfc, materi_EmulsionSrfc, "LogVol_EmulsionSrfc", 0, 0, 0);
new G4PVReplica("PhysVol_EmulsionSrfc", logVol_EmulsionSrfc, logVol_EmulsionBox, kZAxis,
    nDiv, leng_EmulsionElmt);

//X方向に伸びる四角柱(棒)
auto* solid_EmulsionStick =
    new G4Box("Solid_EmulsionStick", leng_EmulsionBox / 2.0, leng_EmulsionElmt / 2.0, leng_EmulsionElmt / 2.0);
auto* materi_EmulsionStick = materi_Man->FindOrBuildMaterial("G4_AIR");
auto* logVol_EmulsionStick =
    new G4LogicalVolume(solid_EmulsionStick, materi_EmulsionStick, "LogVol_EmulsionStick", 0, 0, 0);
new G4PVReplica("PhysVol_EmulsionStick", logVol_EmulsionStick, logVol_EmulsionSrfc, kYAxis,
    nDiv, leng_EmulsionElmt);

//最小単位 (Element)
auto* solid_EmulsionElmt =
    new G4Box("Solid_EmulsionElmt", leng_EmulsionElmt / 2.0, leng_EmulsionElmt / 2.0, leng_EmulsionElmt / 2.0);
auto* materi_EmulsionElmt = materi_Man->FindOrBuildMaterial("EmulsionE373");
auto* logVol_EmulsionElmt =
    new G4LogicalVolume(solid_EmulsionElmt, materi_EmulsionElmt, "LogVol_EmulsionElmt", 0, 0, 0);
new G4PVReplica("PhysVol_EmulsionElmt", logVol_EmulsionElmt, logVol_EmulsionStick, kXAxis,
    nDiv, leng_EmulsionElmt);

auto threeVect_LogV = G4ThreeVector(0, 0, 0);
auto trans3D_LogV = G4Transform3D(G4RotationMatrix(), threeVect_LogV);
new G4PVPlacement(trans3D_LogV, "Emulsion", logVol_EmulsionBox, physVol_World, false, 0);

こういう使い方は、Geant4の開発者が想定する使用方法ではないと思うが、大丈夫ですか? (エラい人へ)

参考文献

https://wiki.kek.jp/display/geant4/Geant4+Japanese+Tutorial+for+Detector+Simulation+2017

原子核乾板の元素組成と密度

更新履歴を含めた情報は、Gitlabで管理しています

https://gitlab.com/yoshimoto/emulsiondensity/blob/master/ja.md

原子核乾板の密度と組成について

原子核乾板 (Nuclear Emulsion) の組成や密度を実測することは難しい。 軽元素から重元素までの測定レンジを持つ組成の測定方法は存在しないため、同じ条件で全ての元素組成を測定することはできない。 また、環境の湿度が変化すると乾板に含まれる水分量も変化し、酸素や水素の組成が変化する。 製造工程から推定する方法については、原料の組成は測定しやすいものの、原子核乳剤を製造する工程に水洗工程があり、この工程で何がどれだけ流出したかの推定は難しい。

よくある測定方法

軽元素である水素、炭素、窒素、硫黄、酸素の質量比(通常の装置では水素、炭素、窒素)は、元素分析装置による測定が可能である1。硫黄についてはその存在比が少ないため、元素分析装置での測定は一般には不可能である。 重元素については、銀、臭素ヨウ素の比はSEM-EDXによる測定が可能である[^8]。 定着(FIX)前後の厚みと質量測定により、AgBrIとバインダーの比及び密度を測定することが可能である[^8]。

アルファ線を用いた密度測定方法

アルファ線の飛程からエネルギーとレンジ(飛程)の関係式を導き出す方法が最も精度の良い密度測定方法である2

アルファ線の飛程は約1μmの標準偏差を持つ分布となり、100イベント使えば飛程につく標準誤差は0.1μmとなる。 密度3.5g/ccのPo-212からのアルファ線の飛程は約50μmなので、飛程につく誤差は0.2%となる。この誤差は、密度の誤差に換算すると0.3%相当になる。

この手法で求められる密度は、アルファ線のエネルギーとレンジ(飛程)の関係式を媒介する、換算密度とも呼ばれるものであり、絶対値としての信頼性ない。

原子核乾板の組成および密度

以下、論文等で紹介されている各原子核乾板の組成を紹介する。

Emulsion Ilford G-5

  • Density 3.907g/cc 3
H C N O S Ag Br I
Mass % ^21 1.4 6.9 1.7 6.9 0.2 47.4 34.8 0.6

50% R.H. 。

[^13]の第3表のIlford G-5 の Nが0.67g/ccとなっているのは0.067g/ccの誤記と思われる。4

Photographic Emulsion (ICRU-215) = Standard emulsion

  • Density 3.815g/cc 5
H C N O S Ag Br I
Mass % [^3] 1.4 7.2 1.9 6.6 0.2 47.4 34.9 0.3

Fuji ET-7B ET-7C/D

  • 密度: 3.73 g/cc (Fuji ET-7B)
  • 密度: 3.60 g/cc (Fuji ET-7C/D)
  • 密度: 3.63g/cc (ET-7D) [^11]
  • 相対湿度: 68%

結晶サイズ

  • ET-7B 240±78 nm 6
  • ET-7C 260 nm ^9
  • ET-7D 176 nm ^9
H C N O S Ag Br I
ET-7D Mass g/cc^9 0.05 0.31 0.1 0.226 0.0067 1.675 1.22 0.034
ET-7D Mass %^9 1.38 8.56 2.76 6.24 0.18 46.25 33.69 0.94
ET-7C/7D Mass %1 1.5 9.3 3.1 6.8 0.2 45.4 33.4 0.3
Binder AgBrI
ET-7D Volume %[^11] 55.0 45.0

Fuji OPERA film

  • 密度: 2.71g/cc 7
  • 密度: 2.84g/cc 乳剤層のみ [^11]
  • 密度: 2.77g/cc 乳剤層と保護層 [^11]

OPERAの論文(シミュレーションも多分)は全ては密度2.71g/cc と表記されている。桑原氏の値となぜズレたのかは不明。

結晶サイズ

  • 200±16 nm [^7]

組成

以下の組成は乳剤層と保護層である。

H C N O S Si Na Sr Ba Ag Br I
Mass % ^6 2.4 13.0 4.81 12.43 0.09 0.08 0.08 0.02 0.01 38.34 27.86 0.81
Binder AgBrI
Mass % 33.0 67.0
Volume %^7 69.0 31.0

組成は原料からの計算値である。水洗工程でバインダーが水溶液側に抜ける割合に大胆な仮定が入っている[^11]。

OPERAタイプ GIF処方[^11]

  • 結晶サイズ 200 nm
  • 密度 3.53g/cc
H C N O S Ag Br I
GIF g/cc[^11] 0.05 0.326 0.11 0.23 - 1.6 1.166 0.033
GIF Mass %[^11] 1.42 9.27 3.13 6.54 - 45.52 33.17 0.94

Nagoya NIT [^8]

  • 乳剤層の密度: 3.44 g/cc
  • AgBrの密度: 6.473g/cc
  • ゼラチンの密度: 1.32g/cc
  • PVAの密度: 1.19g/cc

結晶サイズ

  • NIT 44.2±6.8 nm
  • UNIT 24.8±4.3 nm

組成

H C N O Ag Br I
Mass % 1.6 10.1 2.7 7.4 44.5 31.8 1.9
Atomic % 41.1 21.4 4.9 11.7 10.5 10.1 0.4
Binder AgBrI
Mass % 21.9 78.2
Density g/cc 1.29 6.473

AgBrIの分量は原料からの計算値、HとCとNは実測値、Oはそこからの推定値である。

文責: 吉本

間違いがあればメール or コメントをください

参考文献


  1. 浅田 吉本ら PTEP論文 https://doi.org/10.1093/ptep/ptx076 (2017)

  2. Geant4で伝統的に使われてきたFujiFilm ET-7C/7D emulsion for KEK-PS E373 の組成 元ソースは歳藤氏による

  3. 古関 靖夫 原子核乾板 http://hdl.handle.net/2261/30908 (1959)

  4. 吉田ら NIM論文 https://doi.org/10.1016/j.nima.2016.11.044 (2016)

  5. SRIM ソフトウェア http://www.srim.org/

  6. 桑原謙一氏・吉本 private communication

  7. Photographic Emulsion (ICRU-215) https://physics.nist.gov/cgi-bin/Star/compos.pl?refer=ap&matno=215

  8. 中村琢ら NIM論文 https://dx.doi.org/10.1016/j.nima.2005.08.109 (2006)

  9. L. Patrizii ら JINST論文 https://dx.doi.org/10.1088/1748-0221/3/07/P07002 (2008)

直交座標から球面座標系(極座標)に変換するときの誤差伝搬

Error propagation when converting from rectangular coordinate system to spherical coordinate system

飛跡のベクトルが次のように得られたとする

Assume that a vector of trajectory is obtained as follows

 \displaystyle

V_x = x \pm \delta x \\
V_y = y \pm \delta y \\
V_z = z \pm \delta z \\

Range、Theta、Phiは次のようになる

Range, Theta, Phi are as follows

 \displaystyle
Range = r = \sqrt{x^2+y^2+z^2}\\
Theta = \theta =  acos(z/r)\\
Phi = \phi = atan(y/x)\\

誤差伝搬を計算すると誤差は次のようになる

By calculating the error propagation, the error becomes as follows

 \displaystyle
\delta r = \frac{ \sqrt{  (x\,\delta x)^2+ (y\,\delta y)^2+ (z\,\delta z)^2 }}{r} \\

\delta \theta = \frac { \sqrt{ (x^2 \, \delta x^2+y^2 \, \delta y^2)z^2      + (x^2+y^2)^2\delta z^2    }    }                {r^2\sqrt{x^2+y^2}}\\

\delta \phi =        \frac{\sqrt{(y \, \delta x)^2+(x \, \delta y)^2}}       {x^2+y^2}

検索用

参照

Derivative Calculator with Steps

 \displaystyle
x = r \sin\theta \cos\phi \\
y = r \sin\theta \sin\phi\\
z = r \cos\theta\\

Rangeを1に規格化するとき、

 \displaystyle
\delta x_{norm} = \sqrt{ (\cos\theta \cos\phi \cdot \delta \theta)^2+(\sin\theta\sin\phi\cdot\delta \phi)^2 }\\
\delta y_{norm} = \sqrt{ (\cos\theta \sin\phi \cdot\delta \theta)^2+(\sin \theta\cos \phi \cdot\delta \phi)^2}\\
\delta z_{norm} = \sqrt{(\sin\theta \cdot\delta \phi)^2}\\

tan空間角度と、ラジアン空間の角度の話

3次元角度がtan空間で次のように与えられたとき

ベクトル1

 \displaystyle
  \left( tan\theta _{x1}, tan\theta _{y1}, 1\right)=\overrightarrow {V_1}

ベクトル2

 \displaystyle
  \left( tan\theta _{x2}, tan\theta _{y2}, 1\right)=\overrightarrow {V_2}

角度差は

 \displaystyle
  \Delta tan\theta = \sqrt {\left( tan\theta _{x1} - tan\theta _{x2}\right)^{2}+\left( tan\theta _{y1} - tan\theta _{y2}\right)^{2} }

cosθはドット積を用いて

 \displaystyle
  \overrightarrow {V_1} \cdot  \overrightarrow {V_2} = 
\left| \overrightarrow {V_1}\right| \left| \overrightarrow {V_2}\right| cos\left( \Delta \theta \right)

のように与えられるので、

 \displaystyle
cos\left( \Delta \theta \right) = \dfrac {tan\theta_{x1}tan\theta_{x2}+tan\theta_{y1}tan\theta_{y2}+1} {\sqrt{tan^2\theta_{x1}+tan^{2}\theta_{y1}+1}\cdot \sqrt{tan^2\theta_{x2}+tan^{2}\theta_{y2}+1}}

よってΔθは

 \displaystyle
\Delta \theta = arccos\left( \dfrac {tan\theta_{x1}tan\theta_{x2}+tan\theta_{y1}tan\theta_{y2}+1} {\sqrt{tan^2\theta_{x1}+tan^{2}\theta_{y1}+1}\cdot \sqrt{tan^2\theta_{x2}+tan^{2}\theta_{y2}+1}} \right)

となる

arccosC++だとacosで計算可能。

Windows 10のアップデート後の不本意な再起動を止める方法

f:id:onsanai:20161001184332p:plain:w600

Windows Updateを生かしたまま、アップデート後の不本意な再起動を止める方法は存在しない。Pro版でも、Enterprise版でも。

なので、Windows Updateのサービスを止めるしかない。という結論にたどり着いた。

諦めよう、Microsoftが改心するまで。

blog.bagooon.com

δtanθ - tanθグラフをδθ - tanθグラフに変換する

より正しい記述へ

phst.hateblo.jp

tanθを0.1刻みで、δtanθの値が一律で0.1のときにδθを計算するコードの例。

#include <vector>

void main()
{
    std::vector<double> x, y; //tanθ δtanθ
    for (int i = 0; i < 20; i++) {
        x.emplace_back(i * 0.1 + 0.05);
        y.emplace_back(0.1);
    }

    for (int i = 0; i < 20; i++) {
        double dtheta = atan(x.at(i) + y.at(i) / 2) - atan(x.at(i) - y.at(i) / 2); //δθを計算
        printf("%8.5f %8.5f %8.5f\n", x.at(i), y.at(i), dtheta);
    }
}

f:id:onsanai:20160918044320p:plain