5-2a 測定器応答の実装例(トラッカーとカロリメータ)


  1. トラッカータイプの特徴
  2. トラッカーのSDクラス
  3. トラッカーのHitクラス
  4. カロリメータタイプの特徴
  5. カロリメータのSDクラス
  6. 参考資料

1. トラッカータイプの特徴

トラッカータイプの検出器では、一般に次のような情報が要求される。

これらの情報が、

Hit情報として残されるものが、優秀なトラッカーであると言える。勿論、実際の実験ではトラックが時間的・空間的に分離出来ない場合も多々あるわけだが、まず完全に理想的な環境で実験を行った時にどのようなデータが得られるべきかをシミュレートするのが、シミュレータの仕事の大本命である。ここでは、このような理想的な測定器で測定した場合に得られるhit情報(Monte-Carlo Truthと呼ぶ)を書き出す方法を考える(※1)。

(※1)この段階が終わると、次は実際の検出器の性能をシミュレータに組み込みたくなるのがソフトウェア屋と言われる人々の習性であろう。Geant4は自由度の高いプログラムであるから、勿論このようなシミュレーションをGeant4本体に組み込むことも可能であろうが、Geant4本体は決してそれほど軽いプログラムではない。検出器の性能評価の為だけに何度もモンテカルロ計算を行う事が非常な時間ロスになる場合には、上のようにしてMonte-Carlo Truthを一度書き出し、検出器の性能評価の部分はこのMonte-Carlo Truthを用いて別のプログラムに計算させる方が得策であると思われる。

2. トラッカーのSDクラス(例:T00CounterSD)

まず、トラッカータイプのSensitiveDetectorを書くことから始める。
SensitiveDetectorは、必ずG4VSensitiveDetectorを継承して作成しなければならない。

■コンストラクタ

まずコンストラクタで行うことは、ヒットの集合であるHitsCollectionの名前を定義することが要求される。

いきなりcollectionNameなどと出てくるが、これはG4VSensitiveDetectorのprotectedフィールドである。
Geant4では、いかにもローカル変数のような名前のデータメンバが割合頻繁に出てくる。手元にご先祖クラスのヘッダーファイルも置いて、気合いを入れてデータメンバ探しをしよう。(※2)

(※2)残念ながら、Geant4には変数名の命名規則がないため、コードを読む際に大変混乱する。特に、データメンバとローカル変数の区別がない点が非常に苦しい。命名規則さえちゃんとしていれば、ヘッダファイルを目を皿のようにして何度も確認する必要などないのだが。
Geant4は確かにパワフルなプログラムだが、Geant4をきっかけにC++を勉強する方は、この点だけは見習わない方がいいのではないか、というのが私の個人的な意見である。自分で命名規則を考えるのが面倒な人は、ROOTのコーディング規則の中にあるNaming conventionが非常に美しく体系化されていてお勧めである。こちらはコードの書き方からコメント文の付け方までヒステリックなまでにうるさいが、これを完璧に守っておくと、あとでプログラムをシェルスクリプトなどで一括変換する羽目に陥った時に大変助けになる。

■ Initialize関数

コンストラクタの次は、Initialize関数である。ここでは1イベント分のHitを格納するための入れ物の準備、変数のイニシャライズを行う。この関数は、1イベント毎にイベントの最初に呼ばれる

G4THitsCollectionは1イベント分のHitを格納しておく入れ物(配列)を提供するが、この配列につめられるHitの型は、そのHitを吐き出す検出器によって異なるはずである。そこで、どんな型のHitが入るのかを指定してインスタンスを生成する。ここでは、T00CounterHit型のオブジェクトをつめることを想定しているため、ユーザーは必ずT00CounterHitクラスも同時に定義しなければならない(Hitクラスについては後述)。

HitCollectionをnewで作ったら、次はこれを1イベントにまとめる配列(G4HCofThisEvent)に格納する。G4HCofThisEventは、全ての検出器の1イベント分のHitCollectionを集めた配列である。

ほとんどおまじないだと思ってコピー&ペーストしても問題ないような部分だが、テンプレートの引数など正しいクラス名に直すのを忘れないように。また、根性のある人はG4THitsCollectionを用いず、その先祖クラスであるG4VHitsCollectionから直接継承して独自のHitsCollectionクラスを作ることも出来る。

■ ProcessHits関数

ProcessHits関数は、粒子がSensitiveDetectorに差し掛かった瞬間から飛び出す瞬間(あるいはその中で止まる瞬間)まで、1ステップごとに呼ばれる関数であり、測定器反応をシミュレートした場合にTrueを返す。トラッカータイプの場合、ここで行う作業は以下の通りである。

まず、反応を書き出す粒子の選択を行う。粒子の情報は、ProcessHits関数の引数であるG4Step型のポインタから得られる。

この例では、最初のif文で粒子の電荷を見て、荷電粒子でなければそのままfalseを返して関数を抜ける。
その次のif文は、検出器の表面に粒子が達した時のみヒットを書き出すための条件である。ここで用いられている preStepPointは、G4Stepのポインタ(ここではaStep)から得られる。

const G4StepPoint* preStepPoint = aStep->GetPreStepPoint();

このように、特にジオメトリに関する情報(今バウンダリに居るか、今居るジオメトリのコピーナンバーはいくつか、等)は、必ずpre step pointから得なければならない。

次に、書き出したい情報を取得する。

シンチレーション前面での位置が知りたいので、pre step pointから位置を得る。この他G4Stepがどのような情報を与えてくれるかは、G4Stepのヘッダファイルを参照のこと。

トラッカーの場合には、ここですぐに情報を保存したい(イベント毎でなくトラック毎にデータを書き出すため)。そこで、次にデータをHitに格納し、HitをHitsCollectionに格納する。

かならずHitはnewで作ることを忘れないように。次節で、このデータを保存するHitクラスの作り方に簡単にふれる。

■ 注意点

[top]

3.トラッカーのヒットクラス(例:T00CounterHit)

Hit情報を格納するユーザーHitクラスは、必ずG4VHitを継承して作成する。ヒットクラスの役目は、書き出したい情報をプライベートデータメンバに保存しておくことである(Getメソッドを作る手間を面倒がってデータメンバを公開するのはお勧めしない)。また、データを標準出力に書き出す関数、出力の為の関数などを定義しておくと、何かと便利である。

T00ではその他にプログラムの実行速度を上げるため、new演算子とdelete演算子の上書きをしている。以下はT00CounterHit.hhの抜粋である。

先頭のexternは「どこかで定義されている」という意味。従って、定義の本体が別の場所にあるはずで、この場合T00CounterHit.ccに書かれている。

一般にnew演算子は、newで作る型が要求するメモリ領域を確保する(Memory Allocationという)際に時間がかかる。HitオブジェクトはStepごとに非常に頻繁にnewされる可能性があるため、このMemory Allocationを逐一行わないように変更したものが上のコードである。こうしておくと、deleteが呼ばれても確保したメモリを解放せず、次のnewで再び同じメモリ領域を使う。おまじないのようなコードが嫌なら削っても良いが、そのかわり実行速度は遅くなる。

[top]

4.カロリメータタイプの場合

カロリメータタイプの検出器は、トラッカーと異なり、1トラック毎の出力はしない。最終的にHit情報として書き出されるのは、1イベントが終わった時の全トラック分の重ねあわせ情報である。代表的なHit情報としては、

がある。これらを正しく記録するためには、以下の手順が必要である。

[top]

5. カロリメータのSDクラス(例:T00CalSD)

■ コンストラクタ

トラッカータイプのコンストラクタと同様。ただしcollectionNameに与える名前は別のものにすること。

■ Initialize関数

HitCollectionを生成し、G4HCofThisEventにつめるところまではトラッカータイプと同様。その後、1イベント分のエネルギー損失を保存しておくバッファの0クリアを行う。

edepbufはデータメンバで、T00CalSD.hhに定義されている。

カロリメータタイプの検出器では、最初の粒子の通過時間を知りたい場合もある。このような場合には、G4doubleでもう一つTOFの為のバッファをデータメンバに用意する。

# in T00CalSD.hh

private:
G4double fTof;

# in T00CalSD::Initialize()

// initialize TOF-buffer by huge value
fTof = DBL_MAX;

DBL_MAXはG4doubleで表現できる最大値である。最も早く来た粒子の到達時間を知りたいので、この値よりも早い到達時間のhitが来たら、fTofを上書きすれば良い。

■ ProcessHits関数

カロリメータタイプの検出器では、ステップごとにhitを作る必要はないので、バッファのデータを更新していくだけで良い。ここでは、エネルギー損失を足し上げている。エネルギー損失をとってくるメソッドには、GetTotalEnergyDeposit関数を使う。似たような関数にGetDeltaEnergy関数(そのステップで落ちる運動エネルギーを返す)があるため混乱するが、GetTotalEnergyDepositの方には、cut値にかかった為生成されなかったセカンダリ粒子が本来持ち去るべきエネルギー損失も含まれており、より現実のエネルギー損失に近い値を返すのはGetTotalEnergyDeposit関数の方である。

最初の粒子の到達時間を更新するには、return trueの前に次の二行を書けばよい。

G4double tof = aStep->GetPreStepPoint()->GetGlobalTime();
if (tof < fTof) fTof = tof;

■ EndOfEvent関数

ProcessHits関数の中でHitオブジェクトを作らなかったので、この場合にはイベントの最後にHitを作ってやる必要がある。

作り方、つめ方はトラッカータイプと同じ。tofの値も保存するなら、T00CalHitクラスにデータメンバをその分増やす。

Hitクラスの作り方は、トラッカータイプと同様である。

[top]

6. 参考資料

[top]


[ top | 2-2a 物質定義 | 2-2b 測定器記述 | 5-2a 測定器応答実装例 | 5-2b アプリケーション作成 | 5-2c ROOTによる結果処理例]

このページに関するご質問、ご叱責はこちらへ。