Geant4でピクセル検出器を作って画像認識したい（ジオメトリ編）

はじめに　　

東京大学理学部物理学科ではありがたいことに、理論系の大学院に進む人でも、卒論の代わりに特別実験として、実験系の研究室で修業しないといけないことになっています。
その一環で相原研究室で「Geant4でピクセル検出器作って画像認識してみたら」と言われたので、自分のメモも兼ねてqiita初投稿です。
先生にはあんまり教えてもらってないし、調べると講習会の難しめの資料しかなかったので、自分でとりあえず動かせるくらいの説明をつくりました。
コード全部やっているとものすごく長くなりそうだったので

• ジオメトリ編
• 有感検出器(Sensitive Volume)編
• データ処理編
• 画像認識編
  に分割して投稿します。

目的

Geant4でピクセル検出器を作り、何か見分けのつきにくい反応同士を画像認識で見分けたい。
とりあえず10keVくらいのγ線とβ線は見分けにくとのことなのでチャレンジ。見分けられなかったら別の反応を探してみよう。

環境構築

大学のパソコンにCentOS 8を入れてしまったので、
https://wiki.kek.jp/display/geant4/Geant4+Virtual+Machine から仮想環境をダウンロードしてきて、VirtualBoxでマウントして動かしてます。

元にするファイル

[Geant4初心者講習会のウェブページ]
(https://wiki.kek.jp/display/geant4/Geant4+Japanese+Tutorial+for+Detector+Simulation+2019)から、実習用ソースコードをダウンロードしてくる。
P02_GeometryとP05_Scoringのソースコードを参考にして、検出器とスコアリングプログラムを書く。

ジオメトリの作成

Geant4で物体を定義して配置する時には3ステップあります。

1. 物体の形や寸法から、実体積（Solid Volume）定義する。
2. 物体の素材を決めて、論理物体（Logical Volume）と実体積を関連付けて定義する。
3. 物体の座標や回転を決めて、物理的実体（Physical Volume）と論理物体を関連付けて配置する。物体のコピーナンバーは後でわかりやすいように付けると便利。

ピクセル検出器を作るとき、ピクセルを一個一個配置してもいいんですが面倒なので
4. 以上で作った物理的実体を複製して一定間隔で並べる。

また、当たった粒子の情報をとれるような有感検出器にする場合は
5. 検出器となる物体にSensitiveな属性を持たせる。

注意点として、これはほぼGeant4上で動かすためG4doubleとかを型指定しないといけないんですが、大抵は関数の方でうまくやってくれるのでauto型でも問題ないかと思います。

実体積（Solid Volume）

// ステップ A-1　物体の形と寸法
   G4double leng_X_PixEnvG = 1.*mm;     // X-full-length 
   G4double leng_Y_PixEnvG = 1.*mm;     // Y-full-length
   G4double leng_Z_PixEnvG = 10.*cm;     // Z-full-length 

ここでピクセル検出器全体の寸法を決めています。
一つのピクセルのはX方向とY方向では、微細な電子の動きを取りたいので、現在の検出器の最小の大きさに近い10ミクロン（100分割）としました。
Z方向はガンマ線を止めてほしいので長めに適当に設定しました。

   auto solid_PixEnvG = 
     new G4Box("Solid_PixEnvG", leng_X_PixEnvG/2.0, leng_Y_PixEnvG/2.0, leng_Z_PixEnvG/2.0);

直方体の検出器なのでG4Boxを使って定義しました。（円柱形だったらG4Tubsなはず）
Geant4では中心から測った長さを引数として要求するので、それぞれ/2.0しています。

論理物体（Logical Volume）

   // ステップ A-２　logical volume の定義
   auto materi_PixEnvG = materi_Man->FindOrBuildMaterial("G4_Galactic");

物体の素材をここで決めます。自分で密度や原子番号など定義して新しい物質を作ることもできますが、ここではひとまずGeant4標準装備の真空（G4_Galactic）を使います。
ピクセルを構成する物質については後でやるレプリケーションで定義するので、この段階ではひとまず真空で大丈夫なはずです。

   auto logVol_PixEnvG = 
     new G4LogicalVolume(solid_PixEnvG, materi_PixEnvG, "LogVol_PixEnvG");

実体積、素材と適当な名前（粒子がぶつかった時に表示する名前）を引数として、論理物体を定義します。

物理的実体（Physical Volume）

   // ステップ A-３　Placement of the global envelop to the world
   G4double pos_X_LogV_PixEnvG = 0.0*cm;       // X-location LogV_PixEnvG 
   G4double pos_Y_LogV_PixEnvG = 0.0*cm;       // Y-location LogV_PixEnvG
   G4double pos_Z_LogV_PixEnvG = 0.0*cm;       // Z-location LogV_PixEnvG
   auto threeVect_LogV_PixEnvG = 
   G4ThreeVector(pos_X_LogV_PixEnvG, pos_Y_LogV_PixEnvG, pos_Z_LogV_PixEnvG);
   auto rotMtrx_LogV_PixEnvG = G4RotationMatrix();
   auto trans3D_LogV_PixEnvG = G4Transform3D(rotMtrx_LogV_PixEnvG, threeVect_LogV_PixEnvG);

   G4int copyNum_LogV_PixEnvG = 1000;                // Set ID number of LogV_PixEnvG
   new G4PVPlacement(trans3D_LogV_PixEnvG, "PhysVol_PixEnvG", logVol_PixEnvG, physVol_World, 
                     false, copyNum_LogV_PixEnvG);

検出器を配置する座標とその回転を決めて、G4Transform3Dで渡せる形にしてもらったのちに、コピーナンバー（後で当たった検出器の判別などに使う番号）を割り振って、あらかじめ定義しておいたワールドに配置してます。
回転をさせたい場合は、例えばZ軸を軸として45度回転させたいときは、
rotMtrx_LogV_PixEnvG.RotateZ(45.0 *deg)を追加します。

物体のレプリカを作成して並べる

// 段取りB Local Envelop (the global envelop をY軸方法へ分割)
   // ステップ B-１　Define the shape of the local envelop
   G4int nDiv_Y = 1000;                                 // Number of divisions in Y-direction 
   auto solid_PixEnvL = 
     new G4Box("Solid_PixEnvL", leng_X_PixEnvG/2.0, (leng_Y_PixEnvG/nDiv_Y)/2.0, leng_Z_PixEnvG/2.0);

   // ステップ B-２　Define logical volume of the local envelop
   auto materi_PixEnvL = materi_Man->FindOrBuildMaterial("G4_Galactic");
   auto logVol_PixEnvL = 
     new G4LogicalVolume(solid_PixEnvL, materi_PixEnvL, "LogVol_PixEnvL");

他の物体と同様に、ここではY軸方向に分割した後の一つ一つのパーツの直方体の論理物体を定義しておく。

   // ステップ B-３　the local envelop を the global envelop の中へ Replicaを使って配置する
   new G4PVReplica("PhysVol_PixEnvL", logVol_PixEnvL, logVol_PixEnvG, kYAxis,  
                   nDiv_Y, leng_Y_PixEnvG/nDiv_Y);

物理的実体とその論理物体、分割元となる論理物体、レプリカ方向、分割数、一分割の長さを引数として、レプリケーションさせてます。

// 段取りC　 Pixel Element ( the local envelop をX軸方向へ分割)
   // ステップ C-１　Define the shape of the pixel element
   G4int nDiv_X = 1000;                                 // Number of divisions in X-direction 
   auto solid_PixElmt = 
     new G4Box("Solid_PixElmt", (leng_X_PixEnvG/nDiv_X)/2.0, (leng_Y_PixEnvG/nDiv_Y)/2.0, leng_Z_PixEnvG/2.0);

   // ステップ C-２　Define logical volume of the pixel element
   auto materi_PixElmt = materi_Man->FindOrBuildMaterial("G4_Si");
   auto logVol_PixElmt = 
     new G4LogicalVolume(solid_PixElmt, materi_PixElmt, "LogVol_PixElmt");

   // ステップ C-３　Placement of pixel elements to the local envelop using Replica
   new G4PVReplica("PhysVol_PixElmt", logVol_PixElmt, logVol_PixEnvL, kXAxis,  
                   nDiv_X, leng_X_PixElmt);

Y軸方向の分割と同じくX軸方向も分割する。ここで定義した論理物体がピクセル検出器のピクセル一つ一つを構成する。

今論理物体の階層は、
logVol_PixEnvG > LogVol_PixEnvL > LogVol_PixElmt
になっていて、Geant4ではより内側の階層のものが外側の階層との共通部分では優先されるので、この段階でLogVol_PixEnvG全体をSiのピクセルで敷き詰めたものが構成できた。

検出器を有感（Sensitive）化

/////// 有感検出器
    auto aSV = new SensitiveVolume("SensitiveVolume");
    logVol_PixElmt->SetSensitiveDetector(aSV);         // Add sensitivity to the logical volume
    auto SDman = G4SDManager::GetSDMpointer();
    SDman->AddNewDetector(aSV);

次の有感検出器編で書く"SenstiveVolume.CC"のクラスからインクルードしているので、現時点では呪文ですが、ピクセル一つ一つの論理物体である"logVol_PixElmt"に有感の属性を付与する感じです。

ジオメトリ編のおわりに

本当はヘッダファイルやクラスの別の関数もあるのですが、ジオメトリの主要な部分はここでは省略したワールドを定義すれば終わりです。
これで検出器の実体はできたので、次は検出器の機能を持たせてみましょう。
Sensitive Volume編につづく。