令和2年度　修士論文次世代ガンマ線天文台CTA大口径望遠鏡の波形記録チップDRS4におけるサンプリング時間幅の較正茨城大学大学院理工学研究科理学専攻宇宙物理学コース学籍番号19NM156N野上　優人2021年３月31日
概要
次世代ガンマ殺天文台Cherenkov Telescope Array（CTA）計画は、口径の異なる3種類の解像型大気チェレンコフ望遠鏡を南北のサイトに合計で約100台設置することで20 GeVから300 TeVのエネルギーのガンマ線を観測することを目的とする観測計画である。CTAの中で20 GeVから3 TeVの低エネルギー側に感度をもつ口径23 mの大口径望遠鏡（Large Sized Telescope, LST）は、光電子増倍管（Photo Multiplier Tube, PMT）と７つと波形読み出し回路からなる光検出器265本で構成しているカメラを持つ。2018年10月にスペイン領カナリア諸岛ラパルマ島にてLST初号機（LST1）が完成し、現在は試験観測を行っている。また、LST2-4号機（LST2-4）は、LST1と同じ場所に建設予定であり準備が進められている。 LSTカメラに用いている波形読み出し回路には、数nsに1回程度の頻度で混入する夜光由来のバックグラウンドを除去するためにGHzでの高速波形サンプリングを可能とする波形記録チップDomino Ring Sampler version 4（DRS4）を搭載している。DRS4は、1入力チャンネルあたり1024個のキャパシタアレイを8チャンネル有しており、約1 GHzの周期でスイッチ切り替えを行うことで波形情報をキャパシタに電荷として保存し、波形読み出しは約30 MHzほどの低速で行われる。DRS4では、スイッチを切り替える間隔（サンプリング時間幅）が不均一であることが知られており、この特性を補正することでピクセル毎の電荷分解能を良くすることができ、LSTのエネルギー関値や宇宙線・ガンマ線イベントの識別などに効果があると考えられる。LSTの検出器について、電荷分解能が1000光電子以上で5%以下という目標性能を設定しているため、較正によってDRS4のサンプリング時間幅による寄与を5%以下に必要がある。LST1において検出器を取り外すことなくこの特性を補正するために、読み出し回路の前段にあるスローコントロールボード（Slow Control Board, SCB）に搭載されたテストパルス生成回路から読み出し回路に約100万発のテストパルスを入力し、DRS4の各キャパシタでパルスのピークを捉えた回数で重み係数を求めることで補正する方法が考案された。また、これから建設するLST2-4の波形読み出し回路では、LST1のものからいくつか変更点があり、そのうちの一つがサンプリング時間幅の較正用に搭載されたサイン波生成回路である。これにより、前述したテストパルスを用いた較正に比べて精度の良い較正ができるようになった。回路上で生成した周波数30 MHzのサイン波を測定し得られた波形に対しサイン関数でフィッティングを行い、フィット関数と測定点との時間差を測定することでサンプリング時間幅の揺らぎを較正する。そこで本研究では、LST1の読み出し回路とLST2-4の読み出し回路それぞれにおけるサンプリング時間幅の較正をそれぞれ検証した。 １つ目の検証では、テストパルスを用いた較正方法によってLST1カメラの全ピクセルでの電荷分解能が改善するかを行った。始めに2019年12月と2020年6月に取得された約200万イベントを用いて重み係数を求め、 DRS4のサンプリング時間幅が不均一であることを確かめた。前述で用いた2019年12月の測定データに対して、先ほど得た重み係数を使った補正をする前後でのカメラ平均の電荷分解能を求め比較したところ、高ゲイン系統(High gain, HG)で7.3±0.7%から3.7±0.5%、低ゲイン系統（Low gain, LG）では5.4±0.6%から2.9±0.3%まで改善し、サンプリング時間幅の寄与を5%以下まで較正できることを確認した。また、サンプリング時間幅の時間変化を調査するために、2019年5月、同年12月、2020年6月に取得した時間幅補正用のデータから重み係数を計算し、それぞれ元のデータに対して補正したときの各電荷分解能のカメラ平均を比較した結果、年間の変動率はHGで約1.9%、LGで約4%であることが確かめられた。2019年12月と2020年6月のデータから計算した重み係数が異なる時期のデータでも用いて問題ないかの検証として、2019年5月から2020年9月の間に取得されたテストパルスデータに対してサンプリング時間輻補正をしたときの電荷分解能のカメラ平均を比較し、年間の変動率はHGで3.52±0.03%、LGで3.02±0.02%という結果が得られた。これらの結果より、サンプリング時間幅の時間変化によって約10年までは電荷分解能が5%以下であることが確認できた。 2つ目の検証では、LST2号機以降で用いる読み出し回路でサイン被の生成および測定を行い前述の方法でサンプリング時間幅の較正ができているかを確認するため、始めに1万イベント分測定したサイン波の周期の平均を算出したところ、33.32 nsと理想値である33.33 nsから1 ns短い値となったが標準偏差が0.04 nsであり数十psオーダーの精度で較正できていることを確認した。次に、テストパルスを用いた較正とサイン波を用いた較正の精度を、電荷分解能と時間分解能の観点で比較を行った。その結果、電荷分解能ではテストパルスを用いた較正ではHGとLGでそれぞれ2.4%と2.7%、サイン波を用いた較正ではHGとLGでそれぞれ1.9%と2.7%となり、本研究の目標である5%以下に電荷分解能が抑えられ、テストパルスを用いた方法よりも較正に用いるイベント数を削減しつつ同等以上の精度で較正できていることを確認した。また時間分解能について、テストパルスを用いた較正では0.3%、サイン波を用いた較正では0.07%となり、サイン波を用いた方が較正精度が良いことを確認した。
The Cherenkov Telescope Array (CTA), the next generation gamma-ray observatory, will detect the gamma rays with energies from 20 GeV to 300 TeV to set three kinds of telescopes with different sizes of reflectors on the north and south site.  The Large-Sized Telescope (LST), is the largest telescope in the CTA and sensitives to the gamma-ray energies from 20 GeV to 3 TeV, has the camera consisted of 265 photo detector modules that consists 7 photomultiplier tubes (PMTs) and one read-out board.  The prototype LST (LST1) built on La Palma Island in Canary Islands in October 2018 and is in the commissioning phase.  The second to forth LST (LST2-4) is preparing for the next construction on La Palma Islands. The readout boards of LSTs has the Domino Ring Sampler version 4 (DRS4), the analog memory can sample the waveforms in GHz, due to remove the incidents of the night sky background in once every several nanoseconds. DRS4 has 8 input channels with 1024 capacitors per channel, and by switching at a cycle of several GHz the waveforms are saved the electric charges and they are read out with low speed at 30 MHz.  DRS4 has the non- uniform characteristic of the sampling interval, by correcting this characteristic the electric charge resolution of each camera pixel improves and this correction is effective for lower energy threshold of LST and the identification of cosmic-ray and gamma-ray events. LST is required to reduce the systematic error due to the read-out board to less than 2%. Therefore we should correct the non-uniform of sampling interval of DRS4.
v
宇宙ではGeV（= 10° eV）領域を超えるような超高エネルギーのガンマ線が飛び交い、その一部は地球に飛来する。超高エネルギーのガンマ線は、同等あるいはそれ以上のエネルギーを持った宇宙線陽子やガンマ線バーストのような短時間に膨大なエネルギーを電磁波の放射によって解放する爆発的な現象から生じていると考えられている。しかし、宇宙線陽子の加速機構やガンマ線バーストの放射機構など宇宙で生じている現象は現在も謎に包まれたままである。これらを、解明する鍵の一つが地球に飛米する超高エネルギーガンマ線である。ガンマ線は宇宙空間の磁場によって軌道が曲げられることがないため到来方向の情報が失われにくいことから、どのような天体ないし場所でガンマ線が放射されたかを知ることができる。また、検出したガンマ線のエネルギーが分かればどのような過程で放射されたかを推定することもできる。超高エネルギーガンマ線は、地球大気に入射すると大気中の原子核との相互作用により空気シャワーを生じる。シャワーで生じた荷電粒子が大気中の光速を超える速度で走ることでチェレンコフ光と呼ばれる可視光を放射し、ガンマ線の入射軸から約１度の広がりを持って地上に降り注ぐ。このシャワーによって生じたチェレンコフ光を集光することで間接的にガンマ線を検出するのが大気チェレンコフ望遠鏡である。超高エネルギーガンマ線は到来類度が少ないことに加え、地球に飛来する超高エネルギー宇宙線のほとんどを占める陽子もまた大気中の原子核との相互作用により空気シャワーを生じチェレンコフ光を放射する。そのため、ガンマ線観測をするうえでバックグラウンドとなる宇宙線陽子由米のイベントを除去するのが重要となる。宇宙線陽子とガンマ線が形成する空気シャワーの形状は異なることが知られており、この造いを利用することで宇宙線陽子とガンマ線の識別が可能となる。現行の大気チェレンコフ望遠鏡は、焦点面に設置した複数の光検出器から成るカメラヘチェレンコフ光を反射鏡で集光し、カメラ面上での分布および光子数を記録しており、このような望遠鏡は特に解像型大気チェレンコフ望遠鏡（Imaging Atmospheric Cherenkov Telescope, IACT）と呼ばれている。現行のIACTでは、複数の望遠鏡を並べて同時に観測することで検出有効面積を稼ぎつつガンマ線の到米方向を一意に推定できるステレオ観測が主流となっている。 Cherenkov Telescope Array（CTA）計画は、口径の異なる3種類のIACTを南北に合計約100台を設置することにより、20 GeVから300 TeVと広いエネルギー領域のガンマ線を従来の望遠鏡の10倍の感度で検出することにより次世代のガンマ線天文台となることを目指す国際計画である。口径23 mとCTAのなかで最大のIACTである大口径望遠鏡（Large-Sized Telesope, LST）は、 20 GeVから3 TeVとCTAの低エネルギー領域のガンマ線に感度をもち、南北に4台ずつ設置する予定である。 L S Tに用いているカメラは、光電子増倍管（Photo Multiplier Tube, PMT）7本と波形読み出し回路1枚からなる光検出器モジュール265組で構成している。LSTでは、低エネルギー領域のガンマ線に由来するチェレンコフ光を観測するため、数百MHzの頻度でランダムに混入する星の明かりなどに由来する夜光バックグラウンドの影響が無視できない。この夜光バックグラウンドを除去するには、 GHzでの高速波形サンプリングが可能な波形読み出し回路が要求される。そこでLSTでは、1 GHzで波形サンプリングがができるDomino Ring Sampler version 4（DRS4）と呼ばれるアナログメモリチップを搭載した波形読み1出し回路を開発し用いている。2018年10月には、CTAの北サイトであるスペインカナリア諸島ラ・パルマ島にLST初号機（LST1）の建設が完了し、現在は安定観測中である。また、LST2から4号機に関しては検出器モジュールの組み立てと品質管理試験と同時に建設も進められている。 DRS4チップは、入力チャンネル1つあたり1024個のキャパシタをもつチャンネルを9つ持っており、複数のインバータをつないだ回路でスイッチ切り替えを行い入力チャンネル内のキャパシタに電荷を随時替えることで読み出し回路に入力された波形情報保持し、約30 MHzの低速で波形読み出すことができる。DRS4では、複数の入力チャンネルを繋げて１つの入力チャンネルとすることができチャンネルカスケードと呼ばれ、LSTでは、４つのチャンネルをカスケード接続させることで約4μsのメモリ深さを確保している。また、入力チャンネル内で読み出すキャパシタの範囲を指定することもできこの範囲をRegion of Interst（ROI）と呼ぶ。読み出すキャパシタの数はスライスと表現され、LSTではROIを40スライスにしている。DRS4チップにはいくつかの特性をもっていることが知られており、各キャパシタ毎に固有のペデスタル波高値を持っていること（ペデスタル特性）、一度電荷を酱えたキャパシタが放電しきる前に再度電荷を酱えてしまうことにより本来よりも電荷が大きくなってしまうこと(残留電荷特性)、キャパシタに蓄えた電荷がまわりのキャパシタよりも高くなるものが存在すること（スパイク特性）、キャパシタの切り替えを行うスイッチ回路に用いている半導体の個性により波形サンプリングの時間間隔にばらつきが生じること（サンプリング特性）が挙げられる。いずれの特性も、読み出された信号の電荷量を変化させるため望遠鏡で捉えた光子量の推定を誤ってしまう。各特性の補正方法は、DRS4チップを開発したスイスポール・シェラー研究所や現行のIACTで波形読み出し回路に用いているMAGICで研究されてきた。そして、CTAでもベデスタル特性、残留電荷特性、スパイク特性の補正方法を確立しLSTの解析で実際に行っている。サンプリング特性については、波形の電荷量だけでなく時間情報も歪めてしまうため望遠鏡への入射光量およびガンマ線の到来方向の推定に影響を及ぼすと考えられる。信号の到来時間情報に関しては、MAGICでの補正方法を参考にしてLSTの解析にすでに取り入れている。信号波形の電荷量の補正についても、補正方法の確立を目指して波形読み出し回路の開発段階の時から日本チームで検証してきた。具体的には、2011年の佐々木修論（茨城大学）では読み出し回路に搭載しているField Programmable Gate Array（FPGA）から出力できるクロック信号(矩形波)を用いた方法、2016年の谷川修論（京都大学）ではLST1のスローコントロールボード(Slow Control Board, SCB)で生成できるようになったテストパルスを使った方法、2017年の野崎修論（京都大学）ではLST2号機以降用に改良された読み出し回路で生成できるサイン波を使った方法が挙げられる。本研究では、既に建設が完了し稼働しているLST1および建設が進められているLST2以降で、望遠鏡の稼働を止めずにDRS4のサンプリング特性を補正する方法の検証および実装を目指し、LST1に関しては谷川修論（京都大学）で行われたSCBで生成できるテストパルスを使った較正方法を用いてLST1カメラの全1855ピクセルに対して較正を行った。またLST2号機以降については、野崎修論（京都大学）で行った2号機以降用の読み出し回路をで生成できるサイン波を用いた補正およびLST1で検証したテストパルスを用いた補正を行い、電荷分解能と時間分解能の観点で補正精度の比較を行った。これらの検証結果について本論文で述べていく。

2.1 ガンマ線天文学
宇宙からのガンマ線を確実に検出したのは、1967年に打ち上げられたOrbiting Solar Observatory 3衛星が初めてである[1] [2]。この検出のあと、いくつもの衛星が打ち上げられたことにより次々とガンマ線天体が発見されていった。1973年には地上核実験の監視のため打ち上げられたVela衛星の観測データから、短時間に大量のガンマ線が放射されるガンマ線バーストを検出していたことが報告された[3]。ガンマ線天文学が飛躍的に発展していったのは、1989年に打ち上げられ100 keVから1 MeVのガンマ線を観測したGranat衛星や1991年に打ち上げられ数十keVから10 GeVのガンマ線を観測したコンプトンガンマ線衛星の登場以降のことである。そして2008年には、20 MeVから300 GeVのガンマ線に感度をもつフェルミガンマ線宇宙望遠鏡（フェルミ衛星）が打ち上げられたことにより5000以上の天体を検出している[4]。しかし、エネルギーが数十GeV以上のガンマ線になると飛米する数が減ってくるため大きな検出面を有する検出器が要求されるが、検出器に伴って人工衛星自体も大きくする必要があり費用の増加などの理由から従来の検出器よりもさらに大きいものを新たに搭載した人工衛星を打ち上げるのは困難である。そこでエネルギーが数十GeV以上のガンマ線を観測する際は、宇宙より飛来するガンマ線と地球大気の原子核との相互作用によって生成される電子から起こる電磁カスケードシャワーを観測することによって間接的にガンマ線の検出をする地上に配置された大気チェレンコフ望遠鏡が用いられる。宇宙より飛来する陽子を始めとした原子核（宇宙線）もまた大気の原子核と相互作用することでガンマ線と同様に空気シャワーを形成する。そのため、宇宙線とガンマ線の識別が大気チェレンコフ望遠鏡においても重要となる。宇宙線とガンマ線それぞれから形成されるシャワー形状の造いを使った識別方法の有用性は、1989年のWhipple望遠鏡によるかに星雲の観測[5]によって示された。その後は、ステレオ観測技術や高感度の光検出器、超高速電子回路などを用いている現行の大気チェレンコフ望遠鏡であるH.E.S.S. [6]、MAGIC [7]、VERITAS [8]がガンマ線天文学の発展に貢献をしている。
2.2 ガンマ線天体
ここでは、ガンマ線天文学で観測対象となる天体の中でガンマ線バーストと活動銀河核について紹介する。
2.2.1 ガンマ線バースト
ガンマ線バースト（Gamma Ray Burst, GRB）とは、光速の99.9999%もの速さで放出されるジェット流に伴って数秒から数十秒の間の短時間に1053 ergと非常に大きなエネルギーをガンマ線をはじめとした様々な波長の電磁波の放射によって解放する宇宙で最大の爆発現象である[9]。GRB放射のイメージを図2.1に示す。典型的には地球から数十億光年と宇宙論的距離でGRBが発生することが分かっており、Swift衛星の観測では地球から約132億光年の距離で発生するものが発見されている[10]。GRBには放射継続時間が2秒以下の短い種族と2秒以上の長い種族が存在しており、放射を起こす中心駆動天体の有力候補として短い種族では中性子星やプラックホールが、長い種族では大質量量の重力崩壊が挙げられているものの正体は未だ不明である[11]。GRBの放射には種族に関わらず、ミリ秒以下での不規則な時間変動で放射する即時放射および数日にわたって徐々に減光していく残光放射がある。即時放射では主にMeV領域と少し遅れてGeV領域のガンマ線が、残光放射では電波からGeV領域のガンマ線までと広い波長域の電磁波が放射される。フェルミ衛星とSwift衛星で観測されたGRBの光度曲線を示す図2.2を見ると、数十keVから10 MeVのガンマ線に比べて100 MeV以上のものは数秒遅れて放射されていることがわかる[12]。GRBのMeV領域のガンマ線に関しては、人工衛星が主として観測できているがGeV領域のガンマ線は検出有効面積の制限から困難となるため、広い検出有効面積により数十GeV以上のガンマ線への感度を有することが可能な大気チェレンコフ望遠鏡との連携が必要となる。大気チェレンコフ望遠鏡におけるGRBの観測は、フェルミ衛星やSwift衛星などからのアラートに応じて追尾観測を行うことが主となる。Cherenenkov Telescope Array（CTA）の大口径望遠鏡では、人工衛星からのアラートを受けてから20秒で180度という速さで回転させることが可能なため、放射継続時間の長いGRBであれば即時放射を途中から観測でき、20 GeV以上のガンマ線に感度を持つことからGeV領域の統計数を稼ぐことによりGRBの放射機構などの解明に寄与できると考えられる。
2.2.2 活動銀河核
活動銀河核（Active Galactic Nuclei, AGN）とは、銀河中心にあるプラックホールへガスが降着する際に解放される重力エネルギーによって輝く現象であると考えられている。図2.3に現在考えられているAGNの構造を示す。活動銀河核の約１割は電波銀河と呼ばれる非熱的放射が優勢の天体であり、相対論的なジェット流が観測されている。特にこのジェット流が観測者に対して正面を向き、相対論的ピーミングにより非常に明るく見える天休はプレーザーと呼ばれ、電波からガンマ線と広い波長域で観測されている。プレーザーでは、数分から数日の間にガンマ線の放射強度が大きく変動することが知られており、図2.4はH.E.S.Sが2006年7月に観測したブレーザー PKS 2155–304のガンマ線放射強度の時間変動を示している[14]。また、プレーザーの１つであるMarkarian 421(Mrk 421)の多波長観測によって得られたエネルギースペクトルを図2.5に示す。これまでに観測されているプレーザーの多くは図2.5のように2つのピークを持っており、低エネルギー側のピークは電子のシンクロトロン放射、高エネルギー側のピークは高エネルギーの電子による逆コンプトン散乱に対応すると考えられている[11]。 2
この章では、はじめに超高エネルギーガンマ線が大気に突入した際に生じる空気シャワー現象および空気シャワーの観測によって間接的にガンマ線の観測を行う解像型大気チェレンコフ望遠鏡（Imaging Atmospheric Cherenkov Telescope、IACT）について述べる。そして、複数台のIACTを用いる次世代ガンマ線天文台であるCherenkov Telescope Array（CTA）およびCTAの中で最大である大口径望遠鏡(Large-sized Telescope、LST)、LSTの焦点面カメラに用いている波形読み出し回路について述べる。
3.1 空気シャワー現象
地球大気には宇宙より飛来する高エネルギーの粒子（1次粒子）が突入し、大気中の原子核との相互作用から多数の2次粒子が生じる。 1次粒子が超高エネルギーのガンマ線である場合は、対生成によって電子と陽電子が生成される[17]。 γ +原子核→ e¯ + e† +原子核′対生成で生じた電子や陽電子は、大気中の原子核の近傍を通過する際に制動放射によってガンマ線を放射し再び対生成が生じる。このように対生成と割動放射の繰り返しによって荷電粒子の数が増加していく現象を空気シャワー、特に１次粒子がガンマ線の時は電磁シャワーという。空気シャワーによって生成された荷電粒子のエネルギーが制動放射の臨界エネルギーである約84 MeVまで低下すると、制動放射によるガンマ線の放射がなくなることで対生成が生じず、代わりに電離損失が卓越し荷電粒子がエネルギーを失うことで大気に吸収され荷電粒子の数は減少に転じていく[2]。 ー次粒子が宇宙線陽子の場合は、大気中の原子核との相互作用によって陽子や中性子、π 0 、π±、k±が生成される。 p +原子核→ p' +原子核' + π', π ± π⁰の寿命は0.84 × 10−16秒であり、崩壊するとガンマ線を放射し電磁カスケードシャワーを生じる。 π → γ + γ
π±の場合は、寿命である2.6 × 10− 8秒を迎えるとμ粒子とニュートリノへ崩壊し、生じたμ粒子は寿命2.2 × 10− 6秒で電子とニュートリノに崩壊したのち、電子が制動放射によって放射されるガンマ線より電磁カスケードシャワーが発生する。
このように宇宙線陽子由来の空気シャワーの発生過程は、超高エネルギーガンマ線のときと比べて少し異なり、電磁シャワーと区別してハドロンシャワーとも呼ばれる。ハドロンシャワーは、図3.1のように電磁シャワーに比べて形状が横方向に広がる。空気シャワーで生じた荷電粒子が大気中の光速を超える速度で走ると荷電粒子の後方で原子・分子で分極が生じ、粒子が通過した後に平衡状態へ戻ろうとする際に電磁波が放射される。この電磁波が作る放射面はよく音の衝撃波面に例えられ、放射面で原子・分子の分極で生じた電磁波が強め合うことでチェレンコフ光と呼ばれる可視領域の光子を荷電粒子の進行方向に放射する。このチェレンコフ光の放射角θ。は、大気中での荷電粒子の速度vと真空中での光速c及び大気の屈折率nからcos θc = 品で与えられる。1 TeVのガンマ線から生じる電磁カスケードシャワーで荷電粒子の数が最大となる高度（最大発達高度）はおよそ8 kmであり[2]、対応する大気の屈折率は地上より小さくなり１と近似すると、放射角はθε〜１度となり、地上では半径約140 mの円の範囲にチェレンコフ光が降り注ぐことになる。
3.2解像型大気チェレンコフ望遠鏡
解像型大気チェレンコフ望遠鏡（Imaging Atmospheric Cherenkov Telescope, IACT）は、前節で述べた空気シャワーによって生成された荷電粒子より放射されるチェレンコフ光を反射鏡で焦点面にある検出部へ集光し、検出部でのイメージを記録することにより高エネルギーガンマ線を間接的に観測する手法である。実際の観測では、１次粒子のほとんどが宇宙線であるためガンマ線観測をする上ではバックグラウンドとなり宇宙線とガンマ線を識別することが必要となる。そこでIACTでは、前節で述べたガンマ線と宇宙線が作る空気シャワーの形状の造いを利用しており、検出部で捉えた空気シャワーのイメージは楕円形となるため、シャワー形状を楕円近似することで得られるwidthやlengthといった楕円の形状を特徴づけるパラメータを用いて識別をすることが可能であり、この手法はA. M. Hillasによって考案されたためHillasパラメーターと呼ばれており[20]、図3.3に概念図を示す。IACTでの観測では、検出部で捉えたシャワーイメージでの光量からガンマ線のエネルギーを、シャワーイメージの形状からガンマ線の到来方向をそれぞれ推定している。また、図3.3のように望遠鏡を複数台用いたステレオ観測をすることで、複数のシャワーイメージが取得できそれらの軸の交点からガンマ線の到来方向を推定ができるため、1台の望遠鏡の時に比べて到来方向の推定を精度良く行うことができる。
Cherenkov Telascope Array（CTA）計画
3.3
Cherenkov Telascope Array (CTA)計画は、20 GeVから300 TeVと幅広いエネルギー領域のガンマ線の観測を実現するために口径23 mの大口径望遠鏡(Large-Sized Telescope, LST)と口径12 mの中口径望遠鏡(Middle-Sized Telescope, MST）、口径4.3 mの小口径望遠鏡（Small-Sized Telescope, SST）の3種類の異なる大きさの望遠鏡を北半球のスペイン領カナリア諸岛ラパルマ島と南半球のチリパラナルに合計で約100基を配置することにより現行の望遠鏡の10倍の感度で観測することを目指す国際計画である。CTAの北サイトの完成予想図を図3.4に、各望遠鏡の仕様は表3.1に示す[9]。

3.4
大口径望遠鏡（LST）
LSTは、口径23 mとCTAでは最大の大きさとなる望遠鏡であることから、低エネルギーのガンマ線より発生する空気シャワーによって放射される僅かなチェレンコフ光を集光することが可能であり、低エネルギーである20 GeV以上のガンマ線を検出することが期待される。LSTの反射鏡や焦点面検出器の支持構造の主要な部分には軽量で頑丈なカーボンファイバー強化チュープを用いており、望遠鏡全体の質量を100トンまで軽量化することにより、20秒で180度の高速回転を行いガンマ線バーストのような突発天体の観測に対応する[9]。口径が23 mあるLSTの主鏡は、チェレンコフ光を焦点面カメラへ集光するために198枚の球面分割鏡で構成される放物面鏡を実現している。焦点面カメラには、光検出器として1855本の光電子増倍管（Photo Multiplier Tube, PMT）を用いており、PMT7本で1組とした光検出器モジュール265台でカメラを構成している。2018年10月にはLST初号機（LST1）が北サイトであるラ・パルマ岛に完成し、完成したLST1の様子を図3.5に示す。また、LST1と同じ場所に建設するLTS2号機から4号機（LST2-4）については、現在までに光検出器の組み立てを800モジュールの組み立てが完了し、うち657モジュールまでは品質管理試験を行い問題なく望遠鏡に用いられると判定したのは589モジュールであり、LST2台分の準備は完了している。2021年度にはLST2-4の建設が始まる予定である。
3.5  LST波形読み出し回路
IACTが空気シャワー由来のチェレンコフ光を捉える際、光検出器より出力される信号の幅は数nsと非常に短く、加えて実際の観測で得られるデータには星などに由来する夜光が数百MHzでランダムに混入するため、チェレンコフ光の信号幅に合わせた時間幅で積分した電荷量の測定が必要となる。そのため、高速で信号を読み出すことができ数GHzでの高速波形サンプリングによって信号の波形情報を記録できるような波形読み出し回路が要求される。また、LSTでは複数台の望遠鏡を使ったアレイ観測を行い、波形読み出しを行うためのトリガー信号は各望遠鏡間でコインシデンスを取ることにより生成されるため、チェレンコフ光を捉えてから数μs後に読み出し回路へ入力される。そのため、数μ sのあいだ波形情報を保持できるメモリ深さが要求される。
3.5.1    Dragon Board
LSTの波形読み出し回路上での信号の流れを概説し、図3.6にプロック図を示す。空気シャワーより生じたチェレンコフ光が焦点面カメラに集光されると、検出器モジュール前面にあるPMTによって光信号から電気信号へと変換され、プリアンプ（PreAmplifier for the CTA camera, PACTA）で信号を増幅するとともに、0.2光電子から2000光電子のダイナミックレンジを確保するため高ゲイン系統（High Gain, HG）と低ゲイン系統（Low Gain, LG）の2系統に分けたのち差動信号として出力する。検出器モジュール1台はPMT7本から成るので、HGとLG合わせて14チャンネル分の差動信号がスローコントロールボード(Slow Control Board, SCB)を介して1枚の波形読み出し回路へ入力される。読み出し回路に入るとメインアンプで各系統の信号が増幅され、特にHGでは信号をさらに二つに分けトリガー生成用の信号も生成される。増幅した信号はアナログメモリ茶子Domino Ring Sampler version 4（DRS4）に入力され、LSTでは約1 GHzでサンプリングされる。そしてトリガーがかかると、サンプリングした信号はAnalog Digital Converter (ADC)でA/D変換され、Field Pragrammable Gate Array (FPGA)を介してギガビットイーサネット通信によって外部の観測データストレージへ転送される。LSTでは、夜光各カメラピクセルからの入力信号を隣接する複数のピクセル間で足し合わせたのちに関値にかけることでトリガーを生成するアナログサム方式を採用している[9]。この方式によって、短時間に集中して観測されるチェレンコフ光信号を大きくし、各ピクセルにランダムに入射する夜光でトリガーが生成されるのを抑えることができる。
3.5.2プリアンプ
PMTより出力される電流信号は、PACTAと呼ばれるCTAの光検出器用に開発されたプリアンプによって電圧信号への変換および増幅されたのち差動信号として出力される。PACTAでは、増幅率の異なるゲインを用意しており、入力された信号はHGでは24倍、LGでは1.6倍に増幅するように設計されている。PACTAで増幅された信号が波形読み出し回路へ入力されるとメインアンプでさらに増幅される。これによりLSTの読み出し回路では2000光電子までのダイナミックレンジを確保することが可能となっている[24]。
3.5.3    Domino Ring Sampler version 4（DRS4）
概要
まずアナログメモリについて述べる。アナログメモリとは複数のキャパシタを並べて、入出力部分にスイッチの役割をするインバータを多数接続したものであり、図3.7にアナログメモリの回路の例を示す。入力部分のスイッチを高速で切り替え、各キャパシタに入力信号の電荷を蓄えることで波形情報を保持することができる。メモリ深さは、アナログメモリに用いるキャパシタの数とサンプリング周波数によって決まるため、例えば1 GHzで波形サンプリングをする場合は数千個のキャパシタを用いることで数μsのメモリ深さを確保することが可能である。そして、シフトレジスターの切り替えによりキャパシタに替えた電荷を読み出すタイミングを制御し、ADCでデジタル信号に変換することで入力信号の波形搭報を得ることができる。シフトレジスターは、入力部のスイッチとは独立であるので数十MHzほどの低速で動作させることできる。そのため、アナログメモリでは低速での波形読み出しによる低消費電力を実現できる。 LSTの波形読み出し回路では、GHzでの高速波形サンプリング、数μsのあいだの波形情報の保持および低消費電力を実現するために、スイスのポール・シェラー研究所がミューオンの電子・光子崩壊の探索を目的として2008年より行っているMEG実験[26]のために開発されたアナログメモリであるDRS4というASICを用いている。 DRS4は現行のIACTの一つであるMAGIC望遠鏡の波形読み出し回路でも用いられている[27]。 DRS4は、700 MHzから5 GHzの範囲でサンプリング周波数を変更することができるアナログメモリであり、 9つの差動入力チャンネルを持っている。DRS4のプロック図は図3.8に示す。1チャンネルあたり1024個のキャパシタから成り、８つのチャンネルはPMTからの入力信号の波形サンプリングに用いられ、残る１チャンネルはサンプリング時間の揺らぎによる波形の歪みを補正する測定を行うために確保されている。DRS4での入力部分のスイッチ切り替えは、Domino wave circuitと呼ばれる多数のインバータを連結した回路上をパルス信号が伝搬することによって行われる。DRS4のサンプリング速度は、Domino wave circuit上をパルス信号が一周する速度で決まり、phase locked loop（PLL）回路で制御している。DRS4外部からリファレンスとして周波数freferenceのクロックを入力すると2048 × freference Hzで波形のサンプリングが行われる[28]。DRS4に入力された信号はこのようにして1024個のキャパシタに電荷を保存し、1024番目のキャパシタまで保存されるとチャンネルの始めのキャパシタの電荷を放電し１番目のキャパシタから再び電荷を保存していく。DRS4のwrite shift registerを書き換えることで、あるチャンネルの1024番目のキャパシタまで電荷を保存したあと、同チャンネル内の1番目のキャパシタでなく次のチャンネルへ信号を入力・記録させることが可能であり、チャンネルカスケードと呼ばれる。これにより、1チャンネルあたりのメモリ深さを増やすことができ、LSTでは各ゲイン系統のPMT1本あたり4チャンネルを用いることで1 GHzサンプリング時で約4μsのメモリ深さを確保している。そして、DRS4に波形読み出しトリガーがかかると、33 MHzの低速クロックによって波形読み出しが行われる。波形読み出しを行う際、ファイルへ転送するデータ量を抑えるためにチャンネル内で読み出すキャパシタの範囲を設定することができ、この範囲をRegion of Interest（ROI）と呼ばれる。LSTではROI = 40で設定しており、これは1 GHzであれば40 ns分の波形槽報に相当する。DRS4が出力する信号には、いくつかの特性があるため観測データを解析する際には補正する必要がある。これらの特性について紹介していく。
ペデスタル特性
DRS4には図3.9で示すように、DRS4の各チャンネルの511番目のキャパシタを境にペデスタルの電圧値が異なることが分かっている。この特性は各ゲイン系統の各PMTチャンネル（ピクセル）毎に、波形読み出し回路にのみ電源をいれた状態にして電気ノイズを測定し（ペデスタル測定）、キャパシタ毎のノイズ電圧の平均値を解析に用いる観測データに対して差し引くことにより補正が可能である。
dt特性
DRS4では、入力信号は入力チャンネル内の複数のキャパシタによって電荷として一時的に蓄えることで波形情報を記録している。そしてトリガーがかかると各コンデンサに蓄えていた電荷のうち、ROIの数のキャパシタだけ電荷の読み出しが行われる。そして読み出しが終わると再度各キャパシタに電荷を保存し始める。電荷の読み出しをされたキャパシタが再び読み出されるまでの時間をdtとしたとき、このdtに応じて読み出される電荷量が変化することが知られている[30]。これは、波形読み出し後にキャパシタで蓄えていた電荷を放電したあとも残留電荷が僅かに残っていることが原因と考えられている。図3.10はdt（横軸）と読み出された電荷量（縦軸）の関係を表しており、dtが大きくなるにつれて読み出される電荷量が一定になっている。LST波形読み出し回路では、dtが100 msを超えると残留電荷が放電され読み出される電荷量が一定になることが分かっている[23]。この特性は、dtに応じた電荷量を差し引くことにより補正することができる。
スパイク特性
DRS4には、キャパシタから読み出した電荷が局所的に高くなることで偽のパルスが１つのキャパシタもしくは隣接する二つのキャパシタに現れるスパイク特性があり、図3.11にLST1でのペデスタルの測定データにスパイクが混入した波形を示す。このスパイクには、発生する条件は１イベント前に読み出されたキャパシタのidに依存し、さらに発生パターンに対応して一定の電荷量を持つことが知られている[23]。この特性は、スパイクの発生条件を満たすキャパシタidを見つけ出し、スパイクを記録しているキャパシタの両隣のキャパシタの電荷量で線形補間するもしくはペデスタルの測定ごとに各パターンのスパイクの電荷量の平均を導出し差し引くことによって補正できる。
サンプリング時間幅の不均一性
DRS4では、多数のインバータを繋げたDomino wave circuitでスイッチの切り替えを行っているが、このインバータに用いている半導体の特性によりスイッチの切り替わる時間間隔が各キャパシタ毎に少し異なることが知られている。サンプリング時間間隔を一定としてしまうと、DRS4によってサンプリングされた信号波形からチェレンコフ光信号の到来時間や入力電荷量を推定する際、元の入力波形から歪んだ波形を再構成してしまうため、入力信号に対する到来時間や電荷量が正しく求められない。図3.12にサンプリング時間幅の揺らぎにより入力波形が歪むイメージを示す。LSTでは、夜光のバックグラウンドが0.125 p.e./nsの時の入力電荷に対する電荷分解能の要求仕様および最終目標性能を図3.13に示すように設定している。この電荷分解能は、ポアソン統計、夜光の揺らぎ、読み出し回路の電気ノイズの揺らぎ、PMTゲインの揺らぎによる成分からなる。数百光電子以上になるとポアソン統計による揺らぎの寄与が小さくなりおよそ5%で頭打ちになることから、LSTの目標性能を達成するには統計以外の寄与を5%以下にする必要がある。そこで、本研究の目標は、サンプリング時間幅による寄与を5%以下に抑
えることである。
この特性について、DRS4を搭載した波形読み出し回路を用いているMAGIC望遠鏡ではDRS4の各キャパシタに対応した信号の到来時間に対して補正を行っており[30]、D.Striker-Shaver et al. [31]ではサイン波を用いてサンプリング時間幅のそのものの較正が行われている。またCTAにおいては、2011年の佐々木修論（茨城大学） [32]では矩形波を、2016年の谷川修論（京都大学） [33]ではテストパルス、2017年の野崎修論（京都大学）ではサイン波をそれぞれ用いたサンプリング時間幅の較正を行ってきた。本研究で行ったテストパルスとサイン波を用いた補正方法については、第４章と第５章でそれぞれ紹介する。
3.5.4サイン波生成回路
サイン波生成回路は、DRS4のサンプリング時間幅の不均一性をテストパルスを用いた方法よりも精度良く較正するために、LST2号機から4号機の波形読み出し回路上に新たに搭載された[23]。生成回路上にあるオシレーターから30 MHzで生成される矩形波を同じく回路上に実装したアナログ回路のローパスフィルターを用いて高周波成分を除去し、30 MHzのサイン波を取り出している。このサイン波は、読み出し回路上にある8つのDRS4チップのサンプリング時間幅較正用の9チャンネル目に入力し、測定で得られた波形情報から較正を行う。
3.5.5周辺回路
スローコントロールボード（Slow Control Board, SCB） SCBは、波形読み出し回路とPMTをつなぐ治具としての役割の他にも、PMTに印加する高電圧の設定、PMTのアノードカレントや基板上の温湿度のモニター、テストパルスの生成などの機能を持っており、図3.14にSCBの写真を示す。テストパルスは、読み出し回路上のFPGAから矩形波をSCBへと送り、SCB上にあるテストパルス生成回路でPMT信号を横した擬似信号を生成する。この生成された信号は、PMT信号と同じようにPACTAとメインアンプで増幅されたのちDRS4へ入力される。
バックプレーン
バックプレーン（Back Plane, BP）は、波形読み出し回路の後段に接続されており、主な役割としては読み出し回路へ24 Vの電源を供給すること、カメラサーバーとのイーサネット通信、トリガー信号、クロック信号、Busy信号の伝撒及び分配が挙げられる。また、BPでは10 MHzのクロック信号を生成することができ、本研究ではテストパルス生成回路のトリガー周波数のリファレンスとして用いている。
LSTの波形読み出し回路に搭載している波形記録チップDRS4にはサンプリング時間幅に揺らぎがあり、この特性について較正を行わないと電荷分解能が悪くなることを第３章で述べた。この章では、現在稼働しているLSTIでの波形サンプリング時間幅の不均一性に較正について、まず4.1で目標について述べる。そして4.2節ではLST1でテストパルスを測定したデータを用いた補正係数の導出及び較正方法、4.3節で解析環境について述べる。4.4節では、LST1で測定したテストパルスの電荷を算出する際に4.2節で述べた方法で導出した補正係数を適用し、適用の前後で電荷分解能を比較した結果について述べる。
4.1
現在稼働しているLST1の読み出し回路では、SCBより生成できるテストパルスを用いてサンプリング時間幅の軟正を行う手法が考案されているが、CTAの解析ツールでは確立していなかったためLST1が完成してからこれまでサンプリング時間幅の較正は行われていなかった。LSTの検出器では、電荷分解能について最終目標性能を設定しており、本研究全体としてはサンプリング時間幅による寄与を5%以下に抑えることを目標とすることを第３章で述べた。よってこの章での目標は、サンプリング時間幅による寄与を5%以下に抑えるとともにテストパルスを用いた較正をCTAの解析ツールに確立することである。
4.2サンプリング時間幅の較正方法
LSTIでは、波形読み出し回路とPMTの間にあるSCB上でテストパルス生成でき、これを用いて較正する方法が宇宙線研究所の齋藤隆之氏によって考案され、2016年の谷川修論（京都大学） [33]では実際にこの方法を用いた検証を行っている。そして、本研究ではこの較正方法を現在稼働しているLST1に対して適用した。そこで、以降では較正方法について述べる。較正方法のイメージを図4.2に示す。まず、読み出し回路の前段にあるSCBでテストパルスを一定の周期で生成し読み出し回路へ100万回ほど入力する。テストパルスの生成周期を一定にすることで、波形読み出しが行われるキャパシタが一定問隔で遷移していき各キャパシタにテストパルスが一様に記録させるためである。そして、入力されたテストパルスの最大波高値を記録した回数をDRS4内のキャパシタ毎に数えていく。入力したテストパルスの総数をN total 、DRS4チップのキャパシタid（＝ i）でテストパルスの最大波高値を記録した回数をN i、DRS4の入力チャンネルあたりのキャパシタ数をNcap（= 1024）とすると、以下の式のように重みづけを行うことで各キャパシタのサンプリング時間幅に対応する補正保数wiが与えられる。
この操作をゲイン系統毎、ピクセル毎に行い2 × 1855 × 1024のテープルにまとめる。第３章で述べたように、LSTではDRS4の入力チャンネルを1ピクセルあたり4つ、つまりキャパシタ4096個を用いているがキャパシタの切り替えは1024個のスイッチで行っているため、求めるべき補正係数は1ピクセルあたり1024個である。今回の軟正ではテストパルスを100万イベントほど用いると前述したのは、テストパルスが各キャパシタに記録される確率がポアソン統計に従うとすれば、DRS4のキャパシタ１個あたりおよそ1000イベント分の統計があれば約3%の精度で較正できるためである。また、高ゲイン系統（High Gain、HG）と低ゲイン系統（Low Gain、LG）では信号の増幅率が異なるため較正に用いるテストパルスの波高値を各系統で変えている。ちなみに、ここまでで述べたテストパルスの生成に関して、波高値や生成周期は読み出し回路のファームウェアを書き換えることにより設定することができる。最後に、今回の方法を少し簡潔に述べると、サンプリング時間幅の大きさに比例してテストパルスの最大波高値を捉える確率が上がることを利用して較正している。
4.3解析環境およびツール
ここでは、この章で行った解析の際の環境およびツールについて述べる。
4.3.1 解析環境
今回の解析を行った場所は、La Palma島に配置しているLST1の傍にあるコンテナ内に設置された計算機サーバーであり、2000 CPU core、3.4 PBのディスク容量を有している。望遠鏡の観測データはこのサーバー内に保存されており、サーバー上で観測データの解析を行うことができる。
4.3.2解析ツール
この章での解析には、プログラミング言語の一種であるpythonをベースにしたCTAの解析ソフトウェアであるctapipe[35]、同じくpythonをベースにしたLST用解析ソフトウェアのテストベンチであるcta-lstchain (lstchain) [36]を用いており、当時用いていたのはctapipeでバージョンv.0.7.0、lstchainでバージョンv.0.4.5であった。また、本章で行ったサンプリング時間幅の補正係数を算出するソースコードは付録に載せた。
4.4  LST初号機テストデータを用いた時間幅較正の結果
ここでは、4.2節で述べた方法を用いてLSTIで取得したテストパルスデータからサンプリング時間幅の補正係数を導出し、LST1で測定したテストパルスの電荷量推定にサンプリング時間幅の較正を適用させた時の結果を紹介していく。
4.4.1　サンプリング時間幅の選出
ここでは、4.4.2節で用いるサンプリング時間幅の補正係数テープルを異なる2つの時期の測定データから作成した経緯を述べる。本章で用いたテストパルスの測定データおよびペデスタルの測定データの全てを表4.1に示す。このうち本節および次節で用いたデータセットはRun01732からRun02066の計6つである。4.2節で述べた方法でサンプリング時間幅を算出するにあたり、LST1で測定したテストパルスが全1855ピクセルのROI内に入力されているかを確認するために2019年12月に測定したHG用の鞍正データであるRun01732およびLG用の鞍正データであるRun01734について全イベントの積算電荷を確認した。電荷の算出にはctapipeに実装された電荷算出方法の一つで、ROI内の入力信号波形の最大波高値となるキャパシタ位置（ピーク位置）を探し出し、その位置を含めた任意の範囲を波高値g （iはキャパシタid）の時間積分によって電荷量Qを推定するLocalPeakWindowSumを用いた。 DRS4のキャパシタid iのサンプリング時間幅をΔtとすると、電荷量Qは以下式で表せる。
この時の設定としては、ピーク位置の３つ前のキャパシタ（セル）から７つ先のキャパシタ（セル）までを積分範囲とプログラム上で設定した。図4.3にRun01732およびRun01734の全イベントに対する積算電荷を示す。図4.3のカラーバーはLST1の各カメラピクセルの積算電荷量を示しており、図4.3(a)で暗くなっている部分は該当の検出器モジュール（7ピクセル分）のBP上での遅延時間が他のモジュールとずれていたためである。この問題は、ずれてしまっているモジュールのBPの遅延時間の設定をファームウェア上で変更することによって2020年2月以降の測定では解消している•!。この結果から、HGのデータセットでは全ピクセル分のサンプリング時間幅の較
正ができないことになる。そこで、新たに2020年6月にテストパルスを測定したデータセットであるRun02064 (HG用）およびRun02065（LG用）についても先ほどと同様にして全イベントに対する積算電荷を算出し、カメラマップに描画した結果を図4.4に示す。2020年6月のデータセットでは、HGとLGともに同じ検出器モジュールでテストパルスがROI内に入っていないが、これはテストパルス生成のリファレンスクロックとして用いている10 MHzクロックのカウント数をリセットする際に用いているPulse Per Second（PPS）クロック信号の位相がカメラ中心とカメラ端のモジュールで100 nsずれていることに起因している。この問題はバックプレーンのファームウェアを書き換えることで解消できる。前述のピクセル以外に、テストパルスがDRS4の各セルにまんべんなく入力できなかったことにより4.2節の軟正がうまくできないピクセルが複数存在していた。5万3千イベント毎（subrun+²）にDRS4の各セルでテストパルスの最大波高値を記録した回数のヒストグラムをRun1734およびRun2065でのピクセルid 1487について200万イベント分プロットしたものを図4.5に示す。図4.5の縦軸はDRS4のセルid、横軸はsubrun id、カラーバーは各セルでテストパルスのピークを記録した回数である。縦軸にとったセルidに対するテストパルスを記録した回数の分布がsubrun idに依らず一様であることが期待されるが、Run1734（LG）のpixel id 1487の結果である図4.5(a)ではsubrun id 2から6の間をはじめとしてプロット上では大きく3箇所でテストパルスの記録回数に偏りが見えている。一方、別の時期に同じように時間輻較正用にテストパルスを測定したRun2065（LG）のpixel id 1487の結果である図4.5(b)ではsubrunに対するテストパルスの記録回数の偏りが見えておらず較正に用いることができる。この現象はRun1734だけでなくRun2065、HG用のデータであるRun1732、Run2064でも確認されており、原因は調査中である。ここまでの結果から、各時期のデータセットのみではサンプリング時間幅の較正ができない場合も生じうることがあることが分かった。そこで、今回は2019年12月（Run1732とRun1734）と2020年6月（Run2064と2065）の異なる時期のデータセットを用いて較正を行うことにした。手順としては、まず各データセットからサンプリング時間幅の補正係数のテープルを4.2節で述べたように作成する。次に、各データセットから作成した補正係数のテープルをそれぞれの作成元のデータセットに対して適用してテストパルスの電荷量を算出する。サンプリング時刻もおよびね+1とサンプリング時間幅ΔもにはΔも=ね+1 −ねの関係があることを用いて、式??を以下の(b) (a) (b)ように変形できる。
式4.5の最後の項の△年年年が本章で算出している補正係数にあたる。このようにして補正係数を適用した場合の電荷量を算出する。そして、電荷分解能がより良くなっているものをピクセル毎に選び出すことで最終的なテーブルを作成する。以後の結果では、Run1732とRun2064、Run1734とRun2065に対してはこの節で述べた方法で作成したテープルを用い、そのほかのデータセットについては各Runのみを用いてテープルを作成した。
4.4.2テストパルスに対する時間幅較正の検証
前節で説明した方法を実証するために、サンプリング時間輻較正用のデータを取得しLST1カメラの全ピクセルの補正係数の導出を行った。サンプリング時間幅が較正できているかの検証として、今回は補正係数の導出に用いた較正用データから電荷を補正係数の適用前後で算出し電荷分解能を比較した。今回の較正用データは、2019年の12月および2020年の6月にLST1の各検出器モジュールのSCBからテストパルスを平均1 kHzでそれぞれ約200万発を読み出し回路に入力したときのものを用い、高ゲイン系統と低ゲイン系統それぞれで読み出し回路のペデスタルと区別できるように、PMTから50光電子、500光電子それぞれの信号が入力されたときと同等の波高値に調整した。通常は各系統に対する測定データが１セットずつあれば十分であるが、バックプレーンでの信号の遅延時間の調整ができていなったことやカメラ端のピクセルでテストパルス生成に用いているPulse Per Second (a) (a) (PPS)クロック信号の位相がずれてしまったことによって読み出し範囲内にテストパルスが入らなかったピクセルが各データセットで存在しており片方のデータでは全ピクセル分の補正係数を取得できないため今回は4.4.1節で述べたように２つのデータセットを用いている。補正係数の算出および電荷の計算をする際にはペデスタル補正をした上で行った。また、波形の読み出し範囲である40スライス内で最大の波高値をもつキャパシタを中心とした7スライス分の波高の積分値を電荷量とした。まず、補正係数は各データで200万イベントを用いて前節の式4.1より算出をし、各データに算出した補正係数を用いて電荷を計算して電荷分解能が良い補正係数をピクセル毎に選ぶことにより全ピクセル分の補正係数を用意した。ここでは、テストパルスのピークの３スライス前から７スライス分の波高値を足しわせた値を電荷量とし、ヒストグラムに詰めてガウス分布でフィッティングして得られたパラメータμとσの比を電荷分解能とする。
これにより得られた1ピクセルの補正係数の分布が図4.6である。次にサンプリング時間幅の較正ができているかの検証として、2019年12月の較正用データから補正係数を使った補正前後の電荷分解能の比較を行い、その結果を図4.7で示している。電荷及び電荷分解能の計算は前述と同様である。ただし、高ゲイン系統のデータでは先ほど述べたテストパルスが読み出し範囲内に入らなかったピクセルがいたためプロットからは除いた。高ゲイン系統では、ほとんどのピクセルでサンプリング時間幅の較正によって電荷分解能が改善し、平均と標準偏差は7.3± 0.7%から3.7± 0.5%まで改善することが確認できた。低ゲイン系統では較正によって全ピクセルで改善し、平均と標準偏差は5.4± 0.6%から2.9± 0.3%まで改善し、本研究の目標であるサンプリング時間幅による寄与を5%以下に抑えることができた。
4.4.3サンプリング時間幅の時間変化の調査
前節では4.2節で述べた方法でサンプリング時間幅の補正係数を導出し、鞍正前後の電荷分解能を比較することで較正ができていることを確認した。この節では、解析によって導出した補正係数が時間変化するのかを調査した。 LST1で異なる時期に取得されたテストパルスが入力された較正用データから補正係数を導出し、時期的に独立な測定データに適用した時の電荷分解能の年間の変動率を指標として評価を行った。この調査では、補正係数の導出には前節で用いた2019年12月と2020年6月の較正用データに加えて2019年5月に取得した同様の較正用デー (a) (b)
タを、また電荷分解能の算出には2019年7月に取得したテストパルスの測定データを用い、使用したデータをまとめたものを表4.2に示す。ただし、2019年10月に行った4つの検出器モジュールの交換によって読み出し回路のDRS4チップが変わりサンプリング時間幅も変わってしまっているため交換したモジュール及び前節で述べた理由で読み出し範囲内にテストパルスが入っていないモジュールを除いたうえで評価を行った。この調査で得られた結果のプロットは図4.8に示す。補正係数の導出および電荷の計算には、各データセットの全イベントを用いた。プロットした2点（2019年12月および2020年6月の値）に対して直線のフィッティングを行うことで得られた傾きおよび切片は、高ゲイン系統では−6.35 × 10−5 ± 3.12 × 10−5と4.09 × 10−2 ± 2.67 × 10−4、低ゲイン系統では−1.07 × 10− 4 ± 2.49 × 10− 5と3.25 × 10− 2 ± 2.25 × 10− 4となった。ここでは最小二乗法によるフィッティングの誤差を各パラメータの誤差としている。年間の変動率fyearは、フィッティングを行うことで得られる傾きaと切片なから式4.7のように見積もった。図4.8での傾きは月間での値なので年間の値に直すために12倍する。
LST1のカメラの平均の電荷分解能の年間の変動率とフィッティングで得た傾きと切片の各誤差を考慮した誤差は、高ゲイン系統では約1.9±0.1%、低ゲイン系統では約4±0.1%であることが確かめられた。今回の検証では、サンプリング時間幅の較正に用いるデータを変え、電荷分解能の算出に用いたデータを共通にしていることから、前述の変動率はサンプリング時間幅の変動率を示していると考えられる。ここで、サンプリング時間幅の較正によって得られた電荷分解能をσρ、LSTの目標性能である電荷分解能をσ goal（=5%）とすると、σ goalに至るまでの時間は、簡易的には次の式で与えられる。
4.4.2項での結果では、σQは高ゲイン系統で3.7%であり低ゲイン系統で2.9%であった。誤差を大きく見積(a) (b)もったとき、上式より高ゲイン系統では変動率が1.91％となりt,earは約18年、低ゲイン系統では変動率が4.01％となりtyearは約18年とそれぞれ求まり、サンプリング時間幅は10年間は問題なく較正に用いることができることを確認した。
4.4.4サンプリング時間幅の汎用性の調査
DRS4のサンプリング時間幅はベデスタルのように測定ごとの変動がほとんどないことが知られているが、補正係数を導出してからどのくらいの時間までは問題なく較正に用いることができるかは知っておく必要がある。この節では、4.4.2項において2019年12月及び2020年6月のテストパルスの測定データから作成した補正係数が、県なる時期の測定データに対してもサンプリング時間幅の較正ができるかを調査した。今回の調査では、表4.3に示すように2019年5月、7月、12月、2020年6月、7月、9月の6つの異なる時期のテストパルスを測定したデータに対して補正係数を用いて較正をしたうえで電荷分解能を算出し、補正係数の導出に用いたデータの取得月からの変動率を求めた。電荷分解能の算出には、2020年7月のデータのみ5万イベント、その他のデータでは100万イベント用い、各データにおける全ピクセルの電荷分解能の平均をプロットしたものが図4.9である。このプロットに対しても直線でフィッティングを行い得られた傾きおよび切片は、高ゲイン系統で1.09 × 10− 4 ± 1.05 × 10− 5と3.72 × 10−² ± 8.01 × 10− 5、低ゲイン系統で7.40 × 10− 5 ± 6.03 × 10− 6と2.94 × 10− 2 ± 4.31 × 10− 5のようになった。ここでの誤差も、最小二乗法によるフィッティングパラメータの誤差を用いた。年間の変動率は4.4.3項と同様に式4.7で算出した。この調査から、高ゲイン系統では年間での変動率が3.52±0.03%、低ゲイン系統では3.02±0.02%という結果が得られた。また、4.4.2項で得たLST1の平均の電荷分解能、誤差を大きく見積もったときの変動率および式4.8を用いて、本研究の目標である電荷分解能が5%に至るまでの時間を求めると、高ゲイン系統では変動率を3.55%としたときでtyearは約10年、低ゲイン系統では変動率を3.04%としたときでtyearは約24年となり、各系統ともに10年間は問題なく較正に用いることができることを確認した。
4.5 議論
LST1の読み出し回路では、DRS4チップのサンプリング時間幅の揺らぎを約100万イベントのテストパルスを用いて較正を行うことができ、4.4.2項ではLST1のカメラ全体の電荷分解能が高ゲイン系統で3.7%、低ゲイン系
(a) (b)統で2.9%となり、LSTの目標性能である5%以下を満たしていることが確認できた。しかし、LSTの電荷分解能の目標性能は0.125 p.e./nsの夜光がバックグラウンドを加味した場合であるので、今後はLST1のカメラに実際に光を入射させた観測データからLST1のカメラの電荷分解能について評価する必要がある。4.4.3項および4.4.4項で得た年間の変動率から、カメラモジュールの交換が行われない限りサンプリング時間幅のテープルは約10年間は問題ないと考える。
LST2号機以降用の波形読み出し回路は、LST1で用いているものから改良されている。主な変更点としては、読み出し回路上の各要素に対する温度依存性や入力電圧依存性をモニターするための温湿度計・電圧モニター回路およびDRS4チップのサンプリング時間間隔較正用サイン波生成回路の導入である[23]。サイン波生成回路を導入したことにより、LST1で行っているテストパルスを用いた方法に比べて精度良くサンプリング時間幅の較正ができ、さらに較正に用いるイベント数の削減によって計算時間が短縮され、測定したあとすぐに較正結果を確認できるようになった。5.1節ではこの章での目標、5.2節でサイン波生成回路を用いたDRS4のサンプリング特性の較正方法について、5.3節では実験室でLST2号機以降用の読み出し回路を用いた較正用データの測定方法などについて、 5.4節ではサイン波生成回路を使った較正の検証結果および第４章で紹介したSCBで生成できるテストパルスを用いた方法との較正精度の比較した結果について述べる。
5.1 目標
第３章で述べたように、LSTでは検出器の電荷分解能について最終目標性能を設定しており、本研究の共通の目標としてはLSTの電荷分解能についてサンプリング時間幅に依る寄与を5%以下に抑えることである。この章で扱うLST2号機以降の読み出し回路では、テストパルスを用いた方法よりも精度良く較正を行うためにサイン波生成回路を搭載している。よってこの章での目標は、サンプリング時間幅の寄与を5%以下に抑えるとともにサイン波を用いた較正によってテストパルスのときよりも精度良く較正し、較正にかかる時間を抑えることである。
5.2サンプリング時間幅の較正方法
LST2号機以降の波形読み出し回路では、回路上に新たに搭載したサイン波生成回路で生成したサイン波をDRS4チップのサンプリング時間間隔較正用である9つ目のチャンネルに入力している。生成回路では周波数30 MHzのサイン波を生成し、ROIは100個のキャパシタ（100スライス）に設定して1万イベント分の測定を行った。今回はサンプリング周波数を1.034 GHzに設定しているためサンプリング時間間隔は等間隔の場合であれば0.967 nsである。また、サイン波１周期の時間は33.33 nsである。サイン波を使った較正には、D. Striker-Shaver et. alの較正手法[31]を元に開発した野崎修論の較正手法[23]を用いており、隣り合う2つのキャパシタ間で補正するローカルキャリプレーションを行い、そのあとに入力サイン波１周期ごとのぶれを補正するグローバルキャリプレーションを行う。
5.2.1 ローカルキャリブレーション
ローカルキャリプレーションは、図5.1に示すように隣り合う2つの各データ点とフィット関数との時間方向の差からサンプリング時間間隔を較正をする方法である。まず、測定したサイン波について以下のようなサイン関数でフィッティングをROI全体に対して行う。
このフィット関数のパラメータのうち入力周波数に対応する|p[1]|は、入力周波数が既知であるため30 MHzで図定している。較正には、フィッティングで得た波高値のパラメータ|p[0]|(1 ± 0.75)のサイン波の立ち上がり・立ち下がりの部分に含まれるデータ点を用いる。キャパシタi、i + 1がサイン波をサンプリングした時刻をt[i]、t[i + 1]とすると、較正前のサンプリング時間幅は△t[i] = t[i + 1] − t[i]となる。サンプリングされた波形に対するフィット関数が各キャパシタの波高値と同じであるときの時刻をそれぞれt′[i]、t′[i + 1]とする。この際ペデスタル特性や残留電荷の補正によって各キャパシタの波高値は正しく揺らぎがないとしている。このとき、キャパシタi、i+1でのサンプリング時刻とフィット関数上でのサンプリング時刻の差はdt[i] = t ′ [i] − t [i]、dt[i + 1] = t ′ [i + 1] − t [i + 1]とすれば、予想されるキャパシタi、i + 1のサンプリング時間幅△ねc[i]は以下の式で表せる。
式5.3を各イベントのキャパシタ毎に行い、全イベントNeventで平均をとることで各キャパシタに対応するサンプリング時間幅としていく。
今回はサンプリング周波数が約1 GHz、ROIが100スライス、入力サイン波の周波数が30 MHzであるのでサイン波1周期は33.33 nsとなり、ROI内には3周期分のサイン波が入ることになる。
5.2.2 グローバルキャリブレーション
隣り合う2点間のデータ点から局所的な補正をするのがローカルキャリプレーションであったが、グローバルキャリプレーションでは、ローカルキャリプレーションで得たサンプリング時間幅を用いて図5.2に示すようにサイン波1周期に含まれる複数のキャパシタに対して補正を行っていく。流れとしては、初めにROI内のサイン波1周期分を見つけ出すために測定波形が波高値0の直線との交点（crossing point）を探す。このcrossing pointは、波高値が0になる前後のキャパシタの波高値を用いて線形補問することにより決定する。2点のcrossing pointからサイン波の周期を推定できこの時間をTest、入力したサイン波の理想的な周期をTerne(＝ 33.33 ns)とする。着目したサイン波１周期に含まれたキャパシタのサンプリング時間幅△tʰɛ[i]を測定から得た周期と理想値との比を用いて以下のように補正する。
式5.5を各イベント・各サイン波周期内の各キャパシタ毎に行い、全イベントNeventで平均したものを最終的なサンプリング時間幅とする。

5.3
測定・解析環境
ここでは、実験室でLST2号機以降用の波形読み出し回路を用いた測定をする際の環境や、実際の測定や解析に用いたツールについて触れる。
5.3.1
測定環境
図5.3に実際の測定をする際のセットアップを示す。測定には、実際にLST2号機以降の焦点面カメラに用いるものと同じ波形読み出し回路、テストパルスの生成をするためのスローコントロールボード（Slow Control Board, SCB）および信号増幅をするプリアンプ（Preamplifer for CTA, PACTA）、読み出し回路への電源やテストパルス生成に用いるクロック信号の供給および測定した波形情報やコマンドの送受信を担っているバックプレーン（Back Plane, BP)、読み出し回路への電源の供給には24 Vの電圧と1 Aの電流が流せる直流安定化電源◆I、読み出し回路の制御や測定データの解析に用いるデスクトップPC+2を用いている。波形読み出し回路、SCBおよびPACTA、 BPの概要については第３章の通りである。波形読み出し回路にはの直流電圧が必要であるため、電源装置の正出力端子とGND端子にBPからの銅線を繋ぎ、出力電圧を+24 Vに設定している*3。

5.3.2解析ツール
この辛での解析では大きく分けて2種類の解析ツールを用いた。1つは、LSTで用いる波形読み出し回路の品質管理用のC、C++言語ベースのソフトウェアdragonqc**で、読み出し回路上の操作や測定したサイン波の解析に用いた。もう一つは、CTAのpython言語ベースの解析ソフトウェアであるctapipe[35]およびLST用の解析ソフトウェアのテストベンチであるcta-lstchain（lstchain） [36]を第４章と同様にテストパルスの解析に用いている。また、実験室での測定データはバイナリ形式、LST1の観測データファイルはfits.fz形式とファイル形式が異なり、 ctapipeおよびlstchainでのファイル読み出しは実験客データで取得したバイナリ形式を想定していないため、実験室データファイルをLST1のファイル形式に変換するソフトウェアConverToZfitsも用いている。表5.1に解析に用いたツールを示す。
5.3.3解析に用いた波形データ
今回の検証では、各測定データの較正に用いる波形読み出し回路のペデスタル、DRS4サンプリング時間幅の較正に用いる読み出し回路上で生成したサイン波およびSCBで生成したテストパルスの3種類の信号を測定し解析を行った。テストパルスおよびテストパルスの補正用ペデスタルの測定時は、LST1と同様にROIは40スライスとした。また、テストパルスの測定では、トリガー周波数を2 kHzで固定することにより、DRS4の読み出し始めの位置を約７スライスの間隔でずらし、各キャパシタにテストパルスのピークが一様に記録されるようにしている。実験室で測定したテストパルスの波形を図5.4に示す。また、テストパルスの測波高値はデスクトップPCからのコマンドを介して読み出し回路上のファームウェアを書き換えることによりdB単位で設定することができ、LST1での設定と同じく高ゲイン系統では40 dB、低ゲイン系統では60 dBと設定することによりそれぞれ50光電子と500光電子に対応する波高値を実現している。光電子数と電圧値の対応は、高ゲイン系統を基準にするとベデスタル、残留電荷、スパイクの較正後でおよそ40 ADC counts/光電子である♂ 。テストパルスに対するペデスタルの較正用のデータは、LST1でのペデスタル測定を参考にROIは40スライス、トリガー周波数はランダムに平均で2 kHzとなるように設定して20 万イペント分を測定した。ランダムなトリガーの生成は読み出し回路のファームウェア上で実装した線形帰還シフトレジスターによって行っている。次に、読み出し回路で生成できるサイン波について述べる。図5.5にペデスタル補正後のサイン波の測定波形を示す。サイン波の測定では、通常の観測やテストパルスの測定のときに比べてDRS4にかけるオフセット電圧をデスクトップPCからコマンドを送り読み出し回路のファームウェアを書き換えることで高く設定している。これは、生成サイン波の波高値が± 450 mVであるのに対しDRS4に入力できる電圧の範囲である−50 mVから+950 mVを超えてしまうためである。また、サイン波の測定時は、LST2号機以降の読み出し回路に対して行った品質管理測定を参考にROIを100スライス、トリガー周波数を300 Hzに設定している。サイン波用のペデスタルの測定でもオフセット電圧およびトリガー周波数はサイン波の測定時と同じであるが、ROIは1024スライスに設定している。また、DRS4のサンプリング周波数はDRS4外部からのリファレンスクロック信号の周波数f…およびクロック数Nclockを用いて以下の式5.7を用いて表される。
今回の測定では、LST1でのサンプリング周波数1 GHzを実現するために、FPGA内の66.666 MHzクロック信号をリファレンスとして用い、クロック数をNclock = 66とすることでサンプリング周波数を1.034 GHzに設定した。
5.4実験室でのサンプリング時間幅較正の結果
前節で説明したサイン波の測定データを用いたDRS4のサンプリング時間幅の較正を行った結果を5.4.1項で示す。また、第４章で行ってきたテストパルスの測定データを用いたサンプリング時間幅の較正を２号機以降の波形読み出し回路でも行い、サイン波データを用いた方法との較正精度を電荷分解能と時間分解能の2点での比較をした結果を5.4.2項で示す。
5.4.1サイン波を用いたサンプリング時間幅較正の検証
ここでは、LST2号機以降の波形読み出し回路で生成したサイン波を測定したデータを用いてサンプリング時間幅の較正ができているかの検証を行った。較正に用いた測定データについては表5.2に示す。ペデスタルやサイン波の測定後の解析の流れとしては、まず測定したベデスタルデータからオフセットの平均の波高値をDRS4チップ毎、キャパシタ毎に算出し7[pixels] × 4096[capacitors]テープルにまとめる。次に測定サイン波に対して先ほど算出したベデスタルテープルの値を差し引き、5.2節で述べた方法でサンプリング時間幅の較正を読み出し回路上のDRS4チップ毎に行う。表5.2に示したdata2のサイン波を用いてローカルキャリプレーションおよびグローバルキャリプレーションによって得られた読み出し回路上の１チップの全キャパシタのサンプリング時間幅の分布を
図5.6に示す。今回の測定では、サンプリング周波数fsamplingは1.034 GHzと設定しているため理想的なサンプリング時間幅は0.967 nsとなる。図5.6の平均値を改めて見ると、ローカルキャリプレーションのみだと平均のサンプリング時間幅が0.974 nsと理想値に比べて約0.007 ns大きくなっているがグローバルキャリプレーションまで行うと平均は理想値と一致しサンプリング時間幅の較正ができていることが確認できた。また、グローバルキャリプレーションまで行った時の標準偏差が0.086 nsであることからサンプリング時間幅にはおよそ9%の幅で各キャパシタ固有のサンプリング時間幅を持っていることが分かる。次に、先ほど述べたサンプリング時間幅の較正が正しくできているかを確かめるためにD. Stricker-Shaver et al.や2017年の野崎修論（京都大学）で行っているperiod time testを行った[31] [23]。 period time testは、入力したサイン波に対してcrossing point（波高値が0になる点）を探すことでサイン波１周期に含まれるキャパシタを求め（図5.2の左側)、これらのキャパシタに対応したサンプリング時間幅の足し上げによって周期Ｔを、ROI内のサイン波周期毎、イベント毎に算出し、ヒストグラムから平均の周期を求める。今回は、先ほどdata2で鞍正したサンプリング時間幅を用いて、data2とは別の時期に読み出し回路上で生成したサイン波を測定したdata4に対してperiod time testを行い、サンプリング時間幅の較正前、ローカルキャリプレーション後およびグローバルキャリプレーション後の結果を図5.7に示す。入力したサイン波の周波数は30 MHzであるので理想的な周期は33.33 nsであり、まず較正前（図5.11(a)）に関しては平均値は理想値と一致するが標準偏差が0.28 nsと揺らぎが大きいことが分かる。ローカルキャリプレーション後（図5.7(b)）に関しては、算出した周期の平均が理想値に比べて0.01 nsほど小さくなっているが標準偏差は較正前に比べてまの0.05 nsまで抑えられている。グローバルキャリプレーション後（図5.7(c)）でも周期の平均は理想値より0.01 ns小さくなっているが、標準偏差は0.04 nsとさらに改善し、数十psオーダーとなっていることは2017年の野崎修論（京都大学）の結果とも整合していることが確認できた。
5.4.2テストパルスでの較正との比較
第４章では、LST1カメラの全ピクセルに対してSCBで生成できるテストパルスを用いたサンプリング時間幅の較正によって電荷分解能が改善することを述べた。そして、前項まではLST2号機以降の読み出し回路で生成できるサイン波を用いてサンプリング時間幅を算出し、period time testの結果から較正ができていることを確認した。そこで、今度はテストパルスを用いた較正も行い、サイン波を用いた方法との較正精度を比較し、 5.1項で述べた
ようにサイン波を用いた較正によって電荷分解能が5%以下に抑えることおよびテストパルスを用いた方法よりも較正精度が良いかを検証した。この検証で用いた測定データについては表5.3にまとめた。テストパルスを用いた較正については、第４章で述べたようにDRS4の各キャパシタに一様にテストパルスが入るようにおよそ100万イベント分を読み出し回路に入力し、各キャパシタでテストパルスを記録した回数で重みづけすることで各キャパシタのサンプリング時間幅に対応した補正係数をから算出する。今回は、表5.3のdata3およびdata4から補正係数の算出を行い、各系統の１ピクセルでの補正係数の分布を図5.8に示す。また、サイン波を用いた較正については5.2節で述べたようにローカルキャリプレーションの後にグローバルキャリプレーションを施すことで得られるサンプリング時間幅を算出する。表5.3のdata2を用いて算出した。初めに、各方法で求めたサンプリング時間幅および補正係数を用いてテストパルスの電荷の算出を行い電荷分解能の比較を行った。電荷の算出には、補正係数の算出にも用いた表5.3のdata3とdata4を用いており、電荷は第4章と同じようにキャパシタidに対する波高値をadc i 、補正係数をviとしたとき∑ adc i ×  w;のように計算した。また、サイン波から求めたサンプリング時間幅Δねは、テストパルスから求めた補正係数viと以下の関係があるため補正係数に換算したうえで電荷計算に用いた。
高ゲイン系統について、サンプリング時間幅の較正前および各方法で較正した後の電荷を1000イベント分のヒストグラムを図5.9に示す。各ヒストグラムに対してガウス分布（式（ 4.6)）でフィッティングを行い、得られた平均値と標準偏差のパラメータの比を電荷分解能とすると、較正前は3.3%、テストパルスを用いた補正後で2.4%、サイン波を用いた較正後で2%となった。また、低ゲイン系統でも同様に1000イベント分の電荷を詰めたヒストグラムを図5.10に示す。低ゲイン系統での電荷分解能は、較正前で3.1%、テストパルスおよびサイン波を用いた各較正ではどちらも2.7%であった。この結果から、サイン波を用いたサンプリング時間幅の較正によって今回の検証の目標である電荷分解能を5%以下に抑え、さらにテストパルスを用いた方法と同等以上の精度で較正でき、使用するイベント数が100分の1まで削減できることにより測定した後すぐにサンプリング時間幅の較正ができることが確かめられた。次に、各較正方法によって求めたサンプリング時間幅および補正係数を用いて5.4節で行ったperiod time testによって時間の観点での比較も行った。この検証で用いたのは読み出し回路上で生成し測定したサイン波データである表5.3のdata1であり、period timeを1万イベント分算出しヒストグラムにしたものを図5.11に示す。ここで、 period timeのヒストグラムの平均値と標準偏差の比を時間分解能とすると、較正前は0.8%、テストパルスを用いた較正では0.3%、サイン波を用いた較正では0.07%となりサイン波を用いた較正が最も稍度良く較正ができることが確かめられた。
5.5 議論
LST2号機以降の波形読み出し回路に搭載されたDRS4のサンプリング時間幅較正用サイン波生成回路で生成したサイン波を用いてサンプリング時間幅の較正を行った結果、5.4.1項で述べたようにDRS4の各キャパシタのサンプリング時間幅は、設定したサンプリング周波数に対する理想的な時間幅を中心に9%程度のばらつきをもつこ(a)
とを確認できた。
テストパルスを用いた較正との比較では、サイン波を用いた方法によって各系統の電荷分解能を5%以下に抑え、さらにテストパルスを用いた方法よりも用いるイベント数を削減しつつ同等以上の精度で較正するという本研究の目標を達成できた。しかし、今回行ったサイン波を用いたサンプリング時間幅の較正方法はまだCTAの解析ツールに実装していないこと、ROIをLSTIと同じ40セル分のサイン波データで較正ができるかの検証ができていないこと、読み出し回路の電荷分解能がLSTの性能に及ぼす影響についての調査ができていないことから、今後の課恩と考える。
本修士論文では、LSTの波形読み出し回路に搭載したアナログメモリDRS4チップの特性の一つであるサンプリング時間幅の揺らぎの較正について述べた。まず、LSTIではSCBで生成できるテストパルスを用いて試験観測中のLST1カメラの全1855ピクセルについてサンプリング時間幅の較正を行い、同様に測定したテストパルスの電荷を較正前と較正後とで電荷分解能の比較を行ったところ、高ゲイン系統ではカメラ平均で7.3%から3.7%、低ゲイン系統では5.4%から2.9%まで改善し、サンプリング時間幅の寄与を5%以下に抑えることができた。また、サンプリング時間幅の時間変化の調査から、電荷分解能が5%を超えるまでには約10年かかることが分かった。 LST2号機以降の読み出し回路では、新たに搭載したサイン波生成回路を用いて、測定したサイン波に対してフィッティングを行い、フィット関数と測定点との時間差から較正するローカルキャリプレーションおよびサイン波の周期ごとに波形サンプリングによる時間の伸縮を補正するグローバルキャリプレーションによってサンプリング時間幅の較正を行った。その後、period time testを行い入力サイン波の周期を数十psのオーダーで再現する結果が得られ較正ができていることが確認できた。また、同じ2号機以降の読み出し回路で、LST1で検証したテストパルスを用いた較正と前述のサイン波を用いた較正の較正精度を電荷分解能と時間分解能の2つの観点で比較した。その結果、電荷分解能ではテストパルスを用いた場合は高ゲイン系統で2.4%、低ゲイン系統では2.7%、サイン波を用いた場合は高ゲイン系統で1.9%、低ゲイン系統で2.7%となり目標性能である5%以下を満たし、テストパルスを用いた方法と同等以上の精度で較正できることさらには、較正に用いるイベント数を100分の1まで抑えられたことで測定後すぐに較正ができることが確認できた。また、時間分解能に関しては、テストパルスを用いた場合は0.3%でサイン波を用いた場合では0.07%となった。これらの結果から、サイン波を用いた方がサンプリング時間幅の較正の精度が良いことが確認できた。今後の課題として、読み出し回路の電荷分解能がLSTの性能にどのように影響するかの調査をしていく必要がある。また、LST2-4用の読み出し回路に新たに搭載したサイン波生成回路を用いたサンプリング時間幅の較正方法については、LST1と同じようにROIを40セルにした場合でも同等の較正精度が得られるかの検証およびCTAの解析ツールへの実装も課題として挙げられる。

A　テストパルスを使ったサンプリング時間幅較正のスクリプト
テストパルスを使ったサンプリング時間幅較正に用いたpythonスクリプトを以下に記す。なお、ctapipeおよびcta-lstchainに依存したスクリプトであり、現在も各解析ツールの開発が進められているため仕様変更が行われた際には修正が必要となるので注意が必要である。 CTAの解析ツールであるctapipeおよびcta-lstchainはGitHub上で開発が行われており、以下のスクリプトは、 2021年3月時点でcta-lstchainへpull requestをしているところである*¹。現状としては、サンプリング時間幅のテープル作成に１日以上かかってしまうため、計算時間が改善することが必要である。すると、4.4.1項で述べたような問題があるのかを測定後すぐに確認でき、測定当日に再測定を行うことができる。ソースコード6.1  lstchain/script/sampling_interval_creation.py This script creates the sampling intervals for each DR$4 cell of all pixels only high gain or low gain, and the result outputs as fits file. Note: 1) The interval table creation is needed ~10^6 events (or more) of test pulse data. 2) This script subtracts pedestal and corrects the spikes & time–lapse before timing–calibration. 3) The sampling interval table should be re–created if the LST camera modules are changed. This script was validated with ctapipe(v0.8.0), cta–lstchain(v0.6.3) and ctapipe_io_lst(v0.5.3) in 2020/12/30. Authors: Yuto Nogami (Ibaraki Univ.), Hideaki Katagri (Ibaraki Univ.) Last update: January 18, 2021 import os import time import h5py import argparse import numpy as np import scipy optimize as so from ctapipe.io import event_source from ctapipe.image import LocalPeakWindowSum from ctapipe.version import get_version as ctapipe_ver from lstchain import get_version as lstchain_ver from lstchain.calib.camera.r0 import LSTR0Corrections from traitlets.config.loader import Config parser = argparse.ArgumentParser( usage='$> _python _sampling_interval_creation.py _--input_file _/.../LST.1.1.Run02064.0000. fits.fz _--output_file_sampling_interval_20201222.hdf5_--pedestal_file_ pedestal_Run2065.fits_--coefficient_file_sampling_interval_20201222.hdf5_-- gain_selection_0_--num_events_100000,_--num_charges_1000', description=’This script creates DRS4 sampling interval for each capacitor of all pixels.הwercedה10^6הtestהpulseこeventsהtocorrectהthe sampling interval. ', add_help=True, ) # required arguments parser.add_argument('-i','--input_file',help='Input_event_file',action='store', type=str, dest='path_to_data' , required=True) parser.add_argument( '-o' , '--output_file' , help='Output_file, file_format_is__fits' , action='store' , type=str, dest='output_file' , required=True) parser.add_argument('-ped','--pedestal_file',help='Pedestal_file_for_the_subtraction', action='store' , type=str, dest='pedestal_file' , required=True) parser.add_argument('-gain','--gain_selection',help='Select_gain_to_create. High-gain_is_0,0, low-gainىisى1)’, action='store' , type=int, choices=[0,1], dest='gain' , required=True) # optinal arguments parser.add_argument('-coef','--coefficient_file',help='please_add_if_the_coefficient_file_was created_and_error_pixels_were_existed.', action='store', type=str, dest='coefficient_file') parser.add_argument(' -ne', ' --num_events', help='Set_the_number_of_events_you_want_to_use. The_ default cis "None", ci.e. callevents are used. ' action='store' , type=str, default=None, dest='n_events') parser.add_argument('-nc','--num_charges',help='Set_the_number_of_events_for_charge_ caluculation. Default: 1000', action='store' , type=int, default=1000, dest='n_charges') args = parser.parse_args() # set the number of events to calculate the sampling interval. if args.n_events != None: args.n_events = int(args.n_events) print(' {}こeventsسwillとbeسused create create create create cinterval.'.format(args.n_events)) else: print('All events will be used to create new sampling interval.') # serach the event file for the sampling interval calculation. # Example filename was as LST–1.1.Run00001.0000.fits.fz. path, name = os.path.split(os.path.abspath(args.path_to_data)) if 'Run' in name: # ex) divide to "LST-1.1." and "0001.0000.fits.fz" stream, run = name.split( ' Run' , 1) # ex) divide to "0001" and "0000.fits.fz" run_num, run_sub = run.split( ' . ' ,1) run_name = ' Run ' + run_num run_num = stream + run_name # for log command = 'python_sampling_interval_creation.py'\ + ' _-i_' + args.path_to_data\
+ ' _-o_' + args.output_file\ + '_-ped_' + args.pedestal_file\ + '_-gain_{}'.format(args.gain)\ + ' _-ne _' + str(args.n_events)\ + 'こ-nc _{} ' .format(args.n_charges) if args.coefficient_file != None: command = command + 'こ-coefこ' + args.coefficient_file method_name = 'counting_method' date = time.strftime( ' %Y-%m-%d ' ) ls = os.listdir(path) file_list = [] for file_name in ls: if run_num in file_name: full_name = os.path.join(path,file_name) file_list.append(full_name) file_list = sorted(file_list) # setting for reading and analyzing events n_gain = 2 # number of gain channels n_cell = 1024 # number of DRS4 capacitors per 1 domino ring. n_pixels = 7 # number of pixels per 1 module. n_module = 265 # number of modules. roi = 40 # region of interesta num_pixels = n_pixels*n_module # number of all pixels for LST. offset_adc = 400 # offset adc counts of LST1. # fill the counts of test pulse peaks for each drs4 cell, each pixel. peak_counts = np.zeros((num_pixels, n_cell), dtype=np.float32) sampling_interval = np.zeros((n_gain, num_pixels, n_cell), dtype=np.float32) # high gain if args.gain == 0: # for drs chip id chan = 0 # low gain elif args.gain == 1: # for drs chip id chan = 4 config = Config( { "LSTR0Corrections": { "pedestal_path":args.pedestal_file, "offset": 0, "tel_id": 1, "r1_sample_start": None, "r1_sample_end": None } }) lst_r0 = LSTR0Corrections(config=config) # read events and counting test pulse peak for sampling interval creation. print("----sampling_interval_calculation_is_started-----") # this flag controls the number of events to readout. flag = False # control the number of events to create the interval. num_used_events = 0 start_time = time.time() for filename in file_list: print(filename) with event_source(filename) as source: for event in source: tel_id = event.lst.tels_with_data[0] # convert first capacitor ids from 0–4095 to 0–1023 fc_all = event.lst.tel[tel_id].evt.first_capacitor_id % 1024 pid = event.lst.tel[tel_id].svc.pixel_ids # DR$4 correction lst_r0.calibrate(event) # The counts of first 3 cells & last 1 cell become strange. event.r1.tel[tel_id].waveform[args.gain, :, :3] = 0 event.r1.tel[tel_id].waveform[args.gain, :, 39] = 0 for n_mod in range(n_module): for n_pix in range(n_pixels): drs_chip_id = 8*n_mod+n_pix//2 + chan fc = fc_all[drs_chip_id] pix_id = pid[n_mod*7+n_pix] adc = event.r1.tel[tel_id].waveform[args.gain, pix_id] max_cell = adc.argmax() + fc true_max = max_cell # for the cell ids of pulse peak > 1023. if true_max >= n_cell: true_max = max_cell — n_cell peak_counts[pix_id, true_max] += 1 num_used_events += 1 # when event_id reach to number of events that user set, the counting process finishes. if args.n_events != None: if num_used_events == args.n_events: flag = True break if flag == True: break # convert from the counts to the sampling intervals sampling_interval[args.gain] = (peak_counts / num_used_events) * n_cell end_time_1 = time.time() print('{}こeventscareused ctocobtaincthecsamplingcinterval.'.format(num_used_events)) print('-----sampling_interval_calculation_was_finished----') time1 = (end_time_1 — start_time) / 3600 print('this_process_took_{}_hours_to_calculate.' .format(time1)) # validate the sampling interval coefficients print('----Chargecalculationciscstarted----') num_loop = 0 # number of events we calculate the charge. num_events = args.n_charges print('{}Jevents_were_used_to_obtain_the_charge_resolution'.format(num_events)) # histogram setting num_bins = 1000 pix_arr = np.arange(num_pixels) x_min = 0 x_max = num_pixels y_min = 0 y_max = 1000 # charge histogram without the coefficients hist_without_coef = np.zeros((num_pixels, num_bins), dtype=np.int16) # with previous coefficients hist_old_coef = np.zeros((num_pixels, num_bins), dtype=np.int16) # with new coefficients hist_new_coef = np.zeros((num_pixels, num_bins), dtype=np.int16) # read old sampling interval constant if args.coefficient_file == None: old_coef = np.ones((num_pixels, n_cell), dtype=np.float32) else: with h5py.File(args.coefficient_file, 'r') as f: old_coef = f[' sampling_interval_coefficient'][args.gain] print(old_coef) # charge calculation process flag = False for filename in file_list: print(filename) with event_source(filename) as source: integrator = LocalPeakWindowSum(source.subarray) for event in source: tel_id = event.lst.tels_with_data[0] # convert ids 0–4095 to 0–1023 fc_all = event.lst.tel[tel_id].evt.first_capacitor_id%1024 pid = event.lst.tel[tel_id].svc.pixel_ids # DR$4’s properties correction process lst_r0.calibrate(event) # The counts of first 3 cells & last 1 cell become strange. event.r1.tel[tel_id].waveform[args.gain, :, :3] = 0 event.r1.tel[tel_id].waveform[args.gain, :, 39] = 0 raw_r1 = event.r1.tel[tel_id].waveform[args.gain] # fill the charge to histogram raw_charge = integrator(raw_r1, tel_id, args.gain)[0] count, xedge, yedge = np.histogram2d(pix_arr, raw_charge, bins=(num_pixels, num_bins), range=(( x_min, num_pixels), (y_min, y_max))) hist_without_coef += np.int16(count) old_coef_r1 = event.r1.tel[tel_id].waveform[args.gain] # generate the independent array new_coef_r1 = event.r1.tel[tel_id].waveform[args.gain].copy() for n_mod in range(n_module): for n_pix in range(n_pixels): drs_chip_id = 8*n_mod+n_pix//2 + chan fc = fc_all[drs_chip_id] pix_id = pid[n_mod*7+n_pix] residual = n_cell — fc # apply DR$4 sampling interval coefficients # treatments for some capacitor ids of ROI over capacitor id=1023. if residual < roi: old_coef_r1[pix_id, :residual] = old_coef_r1[pix_id, :residual]*old_coef[pix_id, fc:] old_coef_r1[pix_id, residual:] = old_coef_r1[pix_id, residual:]*old_coef[pix_id, :roi – residual] new_coef_r1[pix_id, :residual] = new_coef_r1[pix_id, :residual]*sampling_interval[args. gain, pix_id, fc:] new_coef_r1[pix_id, residual:] = new_coef_r1[pix_id, residual:]*sampling_interval[args. gain, pix_id, :roi – residual] else: old_coef_r1[pix_id, :] = old_coef_r1[pix_id, :]*old_coef[pix_id, fc:fc+roi] new_coef_r1[pix_id, :] = new_coef_r1[pix_id, :]*sampling_interval[args.gain, pix_id, fc:fc +roi] # fill the charge to histogram old_coef_charge = integrator(old_coef_r1, tel_id, args.gain)[0] count, xedge, yedge = np.histogram2d(pix_arr, old_coef_charge, bins=(num_pixels, num_bins), range=(( x_min, num_pixels), (y_min, y_max))) hist_old_coef += np.int16(count) new_coef_charge = integrator(new_coef_r1, tel_id, args.gain)[0] count, xedge, yedge = np.histogram2d(pix_arr, new_coef_charge, bins=(num_pixels, num_bins), range=(( x_min, num_pixels), (y_min, y_max))) hist_new_coef += np.int16(count) num_loop += 1 # when event_id reach to number of events that user set, the counting process finishes. if num_loop == num_events: flag = True break if flag == True: break end_time_2 = time.time() time2 = (end_time_2 — end_time_1) / 60 print(yedge) print(' {}_events_were_used' .format(num_loop)) print('----Charge_calculation_was_finished----') print('this process _ print('-----sampling_interval_selection_is_started-----') def gaus_fit(x, amp, mu, sigma): return amp * np.exp( − (x − mu)**2 / (2 * sigma**2)) # number of pixels considered to have the best coefficients of the three cases. pixel_ratio = { 'without_coef': 0, 'old_coef': 0, 'new_coef': 0} # charge resolution reso_without_coef = np.zeros(num_pixels) reso_old_coef = np.zeros(num_pixels) reso_new_coef = np.zeros(num_pixels) y_axis = np.array([sum(yedge[i:i+2])/2 for i in range(num_bins)]) for pix in range(num_pixels): # estimate the max amplitude of histogram peak_index = np.argmax(hist_new_coef[pix]) estimate_amp = hist_new_coef[pix, peak_index] # convert the amplitude of histogram to charge values array charge_arr = [] for i, n_events in enumerate(hist_new_coef[pix]): for j in range(n_events): charge_arr.append(y_axis[i]) estimate_mu = np.mean(charge_arr) estimate_sigma = np.std(charge_arr) # initial parameters for fitting param_init = np.array([estimate_amp, estimate_mu, estimate_sigma]) # the case of pedestal event. if estimate_mu < offset_adc: print('pixel offset_adc)) reso_without_coef[pix] = 1 reso_old_coef[pix] = 1 reso_new_coef[pix] = 1 pixel_ratio['old_coef'] += 1 sampling_interval[args.gain, pix, :] = old_coef[pix] # the normal case. else: try: popt, pcov = so.curve_fit(gaus_fit, y_axis, hist_without_coef[pix], p0=param_init) reso_without_coef[pix] = popt[2] / popt[1] popt, pcov = so.curve_fit(gaus_fit, y_axis, hist_old_coef[pix], p0=param_init) reso_old_coef[pix] = popt[2] / popt[1] popt, pcov = so.curve_fit(gaus_fit, y_axis, hist_new_coef[pix], p0=param_init) reso_new_coef[pix] = popt[2] / popt[1] # for comparison of charge resolution tmp_reso = { 'without_coef' : reso_without_coef[pix]} tmp_reso[' old_coef'] = reso_old_coef[pix] tmp_reso['new_coef'] = reso_new_coef[pix] except (so.OptimizeWarning, RuntimeError, ZeroDivisionError) as err: print('pixel_{}:'.format(pix), err) reso_without_coef[pix] = 1 reso_old_coef[pix] = 1 reso_new_coef[pix] = 1 # for comparison of charge resolution tmp_reso = { 'without_coef' : reso_without_coef[pix]} tmp_reso['old_coef'] = reso_old_coef[pix] tmp_reso['new_coef'] = reso_new_coef[pix] # the coefficient selection using charge resolutions min_key = min(tmp_reso, key=tmp_reso.get) if min_key == 'new_coef': pixel_ratio[min_key] += 1 elif min_key == ' old_coef ' : sampling_interval[args.gain, pix] = old_coef[pix] pixel_ratio[min_key] += 1 elif min_key == 'without_coef': sampling_interval[args.gain, pix, :] = 1 pixel_ratio[min_key] += 1 print('pixel_id_{}_was_bad!'.format(pix)) end_time_3 = time.time() time3 = (end_time_3 — end_time_2) / 60 print('this _process print('----sampling_interval_selection_was_finished----') # set the key, high gain or low gain if args.gain == 0: gain_key = ' hg ' elif args.gain == 1: gain_key = 'lg' # output as hdf5 format. with h5py.File(args.output_file, 'a') as f: # the case if the constants file was already created. if 'sampling_interval_coefficient' in f: dataset = f['sampling_interval_coefficient'] dataset[args.gain] = sampling_interval[args.gain] # meta data dataset.attrs[ ' creation_date '] = np.append(dataset.attrs[ ' creation_date ' ], date) # add to existing attribute if gain_key + '_dataset' in dataset.attrs: dataset.attrs[gain_key + '_dataset'] = np.append(dataset.attrs[gain_key + '_dataset'] , run_name) # create the new attribute for dataset else: dataset.attrs[gain_key + '_dataset'] = np.array([run_name], dtype=object) dataset.attrs[ ' shell ' ] = np.append(dataset.attrs[ ' shell ' ], command) dataset.attrs[ 'num_events' ] = np.append(dataset.attrs[ 'num_events' ], num_used_events) dataset.attrs[gain_key + '_num_pixels_wo_coef'] = pixel_ratio['without_coef'] dataset.attrs[gain_key + '_num_pixels_old_coef'] = pixel_ratio[' old_coef' ] dataset.attrs[gain_key + '_num_pixels_new_coef' ] = pixel_ratio[ 'new_coef' ] dataset.attrs[gain_key + '_charge_resolution_wo_coef'] = reso_without_coef dataset.attrs[gain_key + '_charge_resolution_old_coef'] = reso_old_coef dataset.attrs[gain_key + '_charge_resolution_new_coef'] = reso_new_coef # for the first time else: dataset = f.create_dataset( ' sampling_interval_coefficient ' , data=sampling_interval) # meta data dataset.attrs['creation_date'] = np.array([date], dtype=object) dataset.attrs[gain_key + '_dataset'] = np.array([run_name], dtype=object) dataset.attrs[ ' shell ' ] = np.array([command], dtype=object) dataset.attrs['method'] = method_name # number of events to create sampling interval coefficients. dataset.attrs['num_events'] = np.array([num_used_events], dtype=np.int) # number of pixels selected as best coefficients of the three. dataset.attrs[gain_key + '_num_pixels_wo_coef'] = pixel_ratio[' without_coef'] dataset.attrs[gain_key + '_num_pixels_old_coef'] = pixel_ratio['old_coef'] dataset.attrs[gain_key + '_num_pixels_new_coef'] = pixel_ratio['new_coef'] dataset.attrs[gain_key + '_charge_resolution_wo_coef'] = reso_without_coef dataset.attrs[gain_key + '_charge_resolution_old_coef'] = reso_old_coef dataset.attrs[gain_key + '_charge_resolution_new_coef'] = reso_new_coef dataset.attrs['ctapipe_version'] = ctapipe_ver() dataset.attrs[ ' lstchain_version' ] = lstchain_ver() total_time = (end_time_3 – start_time) / 3600 print('Total _time_is_{}_hours' .format(total_time)) print('Finished!!') %\end{ verbatim }
本修士論文を作成するにあたり、指導教員である片桐秀明准教授には学部4年次の研究室配蔵からの3年間、研究をする姿勢や研究内容の相談・議論、ミーティングでの報告の作法など多方面で手厚いご指導を賜りました。研究究の相談や議論時には、私の理解力・説明力の至らなさで長時間付き合わせてしまうこともありましたが、納得するまでお付き合い頂いたこと感謝しております。また、学会発表やテネリフェへの出張の打診があったときに積極的になれなず中々決断できない私に、良い経験になるのではと後押しをしてくださったこと感謝しております。また、自分の至らなさで数えきれないほどご迷惑をおかけしたこと、大変申し訳なく思っております。それでも私のことを指導してくださったことには感謝してもしきれません、本当にありがとうございました。最後に、昼食に毎回誘っていただきありがとうございました。吉田龍生教授には、学内でのミーティングで進捗報告やゼミで多くの助言を頂いたこと、本当に感謝しております。また修論執筆時には心配していただいていたことを小原君から伝え聞きました、ご心配をおかけしました。柳田昭平名誉教授には、ゼミで多くの助言やご指導をいただきましたこと感謝申し上げます。 CTA-japanグループの皆様、特にAS/ODAグループおよびFPI/CAL/ELECグループの皆様には、各定期ミーティングでの進捗報告に対するアドバイスや研究を進める上でのサポートなど多くのことをしていただき心から感謝しております。FPI/CAL、ELECのメンバーの皆様、特に京都大学の笛秀利教授には研究の進捗報告をする場を与えてくださったこと、研究を進めるうえで数多くの助言を頂きましたこと、感謝申し上げます。AS/ODAグループの皆様、特に宇宙線研究所の吉越貴紀准教授、大石理子助教には、まだ正式にCTAメンバーでなかった学部4年次のから研究の進捗報告の場とさせて頂いたこと、報告内容に対して質問やご指摘を数多くしていただいたこと、CTAの解析等に関する情報共有および学外のCTA-japanメンバーとの交流する場としてF2F meetingを開催してくださったことなど多岐にわたってお世話になりました、感謝しております。宇宙線研究所の齋藤隆之助教には、2018年度の学部4年次のころからお世話になりました。本研究をする上では欠かせないLSTの読み出し回路周辺について多くのご指導を賜り、定期ミーティングでは進捗報告のたびに質問や助言、議論などを欠かさずして頂いたことには、心から感謝しております。修士論文に向けた研究内容や学会発表、テネリフェ出張の打診をしてくださったことにも感謝しております。また。ヨーロッパの方々とのミーティングでは、相手から質問やコメント等をいただいた際に英会話能力が乏しい私に代わって伝えていただくなどサポートいただきましたこと、感謝しております。京都大学の野崎誠也さんには、LSTの読み出し回路周辺のことやLSTの実データ解析の仕方、解析結果などをCTA-japanの各定期ミーティングやメール等を通じて数々の助言や議論、サポートなど多岐にわたってお世話になりました、本当に感謝しております。また、対面する機会が2019年度のF2F meetingでの1度切りであったことがとても残念であり、コミュニケーション能力が低いですがそれでも対面でお話しする機会があったらとミーティングやメールの文面から野崎さんの優しそうな雰囲気を垣間見るたびに思っておりました。宇宙線研究所の小林志鳳さんには、LSTの実データ解析の方法に関しての情報共有や定期ミーティングでは多くの助言を頂いたきました、本当にありがとうございました。宇宙線研究所およびテネリフェで行ったLST2号機以降の光検出器モジュールの組立て・品質管理では、宇宙線研究所の高橋光成特任研究員、学術支援専門職員の猪目祐介さん、櫻井駿介さん、埼玉大学の砂田裕志さん、京都大学の岡知彦さん、甲南大学の川村孔明君、塚本友祐君にはお世話に