MATLAB Runtimeもインストールします。インストール先ディレクトリをユーザーのワークスペースにして大丈夫です。
pulse_searcher は、指定した範囲の DM (dispersion measure) で de-dispersion を行いパルスを探査するソフトで、FRB (fast radio burst) や GRP (giant radio pulse) の検出に利用できます。
パルサーの本格的解析に必要な folding や acceleration 補正などはできません。
また DM を振って FRB 探査もできるものの、処理時間がかかりすぎて実用的ではありません。
ホームディレクトリなど適当な場所にワークスペース (e.g., /home/USER/pulsar) を作っておくことをおすすめします。
USER は自身のユーザー名に読み替えてください。
$ mkdir /home/USER/pulsar
適当なディレクトリにソフトをダウンロードし、実行権限を与え、実行してください。
$ wget http://astro.sci.yamaguchi-u.ac.jp/t-aoki/tools/pulse_searcher/pulse_searcher_cent7gnome.install $ chmod 555 pulse_searcher_cent7gnome.install $ ./pulse_searcher_cent7gnome.install
少し時間がかかるかもしれませんが以下のような画面が表示されるはずなので、指示に沿ってインストールしてください。 インストール先ディレクトリは先に作成したワークスペースを指定すると良いと思います。
MATLAB Runtimeもインストールします。インストール先ディレクトリをユーザーのワークスペースにして大丈夫です。
これですべてのファイルはワークスペース /home/USER/pulsar 以下に入ったはずです。
下記の情報を記したテキストファイル crab-6g.conf (6GHz用) と crab-8g.conf (8GHz用) を用意してください。
ファイル名は任意ですが、ターゲットと観測周波数を入れておくとわかりやすいです。
$ vi crab-6g.conf # Set Dispersion Measure. # In the program, there is a loop such as for(dm = min; dm <= max; dm += step). survey_DM_pc_per_cm3_min 56.7450 survey_DM_pc_per_cm3_max 56.7450 survey_DM_pc_per_cm3_step 56.7450 # Set FFT Window Width in milli-second. fft_width_ms_desired 100 # Set Width of Short Time Fourier Transform in mili-second, # relating to a data size transferred to a GPU at once. # e.g., 10 ms stft_width_ms_desired 100 # Set Integration Time in micro-second. # e.g., 2 us. If 0 then no integration. integration_time_us_desired 2 # Set Signal-to-Noise Ratio Threshold for Pulse Detection. snr_threshold 7 # Set Observation Parameters. sampling_freq_MHz 1024 # 1024 mega-sample/s adc_RF_min_MHz 6600 # Minimum RF/MHz of the A/D converted band; e.g., 6600 or 8192. eff_RF_min_MHz 6606 # Minimum RF/MHz of the effective band for analysis; e.g, 6606 or 8198. eff_RF_max_MHz 7106 # Maximum RF/MHz of the effective band for analysis; e.g, 7106 or 8698. band_inverted true # 'true' when inverted by Higher Order Sampling, or 'false' datum_size_bit 2 data_format ADS1K # Linear1B, ADS1K (= ADS3KP_VSISEL11 via OCTAVIA) use_gpu true # true or false save_images true # true or false partition_width_byte 256000000 # Data size/B loaded from a file at once.
各パラメータの意味はおおよそ予想できるかと思います。要点だけ列挙すると下記の通りです。
survey_DM_...: DM を 56.745 pc/cm^3 に固定してデディスパージョンしています。DMを振ることも可能ですが解析に時間がかかります。fft_window_width...: FFT 窓幅を約 100 ms に設定しています。integration_time...: ジャイアントパルス探査には積分時間 2 us 程度が適切です。snr_threshold: S/N 閾値は 7 がマンパワー的に妥当です。 S/N 5 だと偽陽性検出確率がたかだか1000万分の1なので、パルサー自転位相に乗らない偽陽性検出が大量に発生します。adc_RF_...: A/D変換の対象となる記録帯域が 6600 MHz から始まることを示しています。eff_RF_...: データの有効帯域を 6606 - 7106 MHz (幅 500 MHz) としバンドパスフィルタリングしています。 エイリアシング等による帯域端の荒れを除去するためです。band_inverted: 高次サンプリングによって、アナログとディジタルで帯域が反転している場合に true に設定します。data_format: ADS1K としていますが実際には ADS3KP (VSISEL 1 1) + OCTAVIA によるビット並びを想定しており、実際の ADS1K には使えません。pulse_searcher は四つの引数を取り、例えば次のように実行します。
$ ./pulse_searcher 2018001020000_hit32cr_crab.raw crab-6g.conf u18001-band.txt /home/USER/pulsar
各引数は次の通りです。
2018001020000_hit32cr_crab.raw)YYYYDOYhhmmss_(any).(any) という書式とし、そのファイルのデータの記録開始時刻を UTC で {year 4 digits}{day-of-year 3}{hour 2}{minute 2}{second 2} と示す必要があります。
文字列 any は任意で、上記例の hit32cr_crab は Hitachi32m C-band RHCP による Crab pulsar のデータであることを示しています。crab-6g.conf)power_spectrum の引数 sideband に DSB を指定すれば、そのデータを得られます。
データ点数 (FFT点数) は任意ですが、帯域特性を十分表現できるような点数であるべきです (昔のアジレントスペアナの401点データでは不十分です)。
/home/USER/pulsar)ただし上記のように実行できるのはシステムにMATLAB環境変数を設定している場合です。
あまり使わない環境変数をシステムに設定するのは利口ではないので、環境変数を含めたスクリプト pulse_searcher_interface.sh (下記) を用意しました。
$ vi pulse_searcher_interface.sh
#!/bin/bash
set -u
LANG=C
###
### EDIT BELOW: MATLAB Environment
###
# mcr_root = the absolute path of the directory where the MATLAB Runtime is installed.
# version = the version of the MATLAB Runtime
mcr_root=/home/oper/MATLAB/MATLAB_Runtime
version=v94
export LD_LIBRARY_PATH=${mcr_root}/${version}/runtime/glnxa64:${mcr_root}/${version}/bin/glnxa64:${mcr_root}/${version}/sys/os/glnxa64:${mcr_root}/${version}/sys/opengl/lib/glnxa64
execut=/home/oper/MATLAB/pulse_searcher/application/pulse_searcher
outdir=/home/oper/workspace/pulsar
###
### MAIN
###
EXEC=$(basename $0)
if [ $# -ne 3 ]; then
echo "USAGE:"
echo "$ ${EXEC} {directory including *.raw} {config} {band}"
echo "e.g.,"
echo "$ ${EXEC} /mnt/data/p19001a/ crab-6GHz.txt double-side-band.txt"
exit 1
fi
datadir=$1
config=$2
band=$3
# Lists data files
echo "Conf: ${config}"
echo "Band: ${band}"
echo "Data:"
filelist=$(find ${datadir} -name *.raw | sort)
for file in ${filelist}; do
echo " ${file}"
done
# Confirm to continue
read -p "Continue? [y/n]: " yn
if [ ${yn} != y ]; then
echo "Process terminated."
exit 2
fi
# Analysis
for file in ${filelist}; do
${execut} ${file} ${config} ${band} ${outdir}
done
exit 0
実行方法は以下。
$ ./pulse_searcher_interface.sh ← 使い方が表示される $ ./pulse_searcher_interface.sh /mnt/raid/p19001a crab-6g.conf 2019001010203_spectrum.txt
3つの引数はそれぞれ以下の通りです。
power_spectrum で引数 sideband にDSBを指定すれば得られる。出力ファイルは以下の四種類です。
dd を列挙しただけ)。COUNT の値で抽出範囲を指定可能。bs=256000 としているので、1024 MS/s, 2 bit/S (= 256 000 000 B/s) のデータの場合 bs = 1 ms であり、COUNT の単位は ms となります。MEAN, STDEV の値は外れ値 (解析設定ファイルで指定された snr_threshold) を除くノイズの平均値と偏差ですが、自動処理なので信号成分を含んでいます。よってそれから計算される SNR も厳密には S/N を表しません。
正しく S/N を出すには *_dd.sh で取り出したパルスデータを逐一再解析する必要があります。