- 入門編
- インストール
- 起動~タスクを使ってみる
- イメージングのチュートリアル (CASA Guides の日本語訳)
- ALMAデータ・アーカイヴ
- FAQ
イメージング
このチュートリアルの目的
このチュートリアルでは、CASA を初めて使う人向けに、ALMA データのイメージングの初歩を解説しています。この解説は、CASA guides "First Look at Imaging CASA 6.4" (本家のサイト:英語) の日本語訳(と多少の追加説明)となっています。ALMA では、観測でデータが取得された後、パイプラインや ARC スタッフの作業によってキャリブレーション (大気や装置の影響を取り除く作業) が行われます。 そのキャリブレーションで使われたスクリプトは、ALMA アーカイヴ上で生データと一緒に配布されますので、ユーザーはキャリブレーション済みデータをすぐに入手できます。 ALMA の品質保証 (Quality Assurance 2 (QA2)) では、ビームサイズや到達感度の確認のためにイメージングも行われますが、サイエンスの目的に応じて、自分でイメージングをやり直したいと思うことも多いでしょう。 また、イメージングで用いられる CLEAN は、ALMA データに限らず、干渉計データの解析で頻繁に用いられる手法です。 このページと、他の一連のページ (セルフキャリブレーション、ラインのイメージング、イメージの解析) で、イメージングの基本的な流れを示していきます。 なお、解析における基本的な概念、開口合成やデコンボリューションの解説、また、
tclean
実行時に出てくるマスク、clean component、残差 (residual) などの意味については、例えば干渉計サマースクール2005 教科書や Imaging Interferometric Data (英語)、その他、リンク集にある資料を参照してください。
CASA のインストール
まだ CASA をインストールしていなければ、こちらのページを参考にダウンロード、インストールしてください。このチュートリアルはバージョン 6.4 を対象として書かれていますので、そのバージョンをインストールして下さい。バージョンが多少違って動作する可能性は高いですが、多少見た目やオプションが違う可能性があります。チュートリアルで用いるデータについて
例として用いるデータは、原始惑星系円盤を持つ若い星、うみへび座 TW 星 (TW Hydra) です。プロジェクトコード 2011.0.00340.S の "Searching for H2D+ in the disk of TW Hya v1.5" で得られたものです。2012 年 11 月にサブミリ波の連続光放射、H2D+ と N2H+ のラインをターゲットとして Band 7 で観測されました。 ここでイメージングするデータのバンド幅は 234.375 MHz で、610 kHz ごとに並んだ 384 個のチャンネルを含んでおり、12-m アンテナ 12 台を用いて Cycle 0 観測で得られたデータです。サンプルデータのダウンロード
チュートリアル用のパッケージ (4.1 GB) は以下の場所から取得できます。https://bulk.cv.nrao.edu/almadata/public/working/sis14_twhya.tgz
パッケージは tar コマンドで解凍できます。
# In bash
tar -xvzf sis14_twhya.tgz
このイメージングのチュートリアルでは、キャリブレーション後のデータを使います。上記を実行すると、"sis14_twhya_calibrated_flagged.ms.tar" というファイルができるので、これを更に展開し、データ "sis14_twhya_calibrated_flagged.ms" (measurement set, MS) を作業ディレクトリにコピーします。 そのディレクトリから CASA を起動するようにします (CASA ログは起動したディレクトリに保存されます)。
# In bash
cd working
tar -xvf sis14_twhya_calibrated_flagged.ms.tar
cd ..
mkdir MyTutorial
cd MyTutorial
cp -r ../working/sis14_twhya_calibrated_flagged.ms .
本ページに書かれているイメージングチュートリアル部分だけのデータ (600 MB) が欲しい場合は以下のように部分的に取得することもできます。
# In bash
mkdir MyTutorial
cd MyTutorial
wget -r -np -nH --cut-dirs=4 --reject "index.html*" https://bulk.cv.nrao.edu/almadata/public/working/sis14_twhya_calibrated_flagged.ms.tar
CASA を起動する
CASA を起動する前に、ls
コマンドで sis14_twhya_calibrated_flagged.ms
がカレントディレクトリにあるかどうかを確認しましょう。# In bash
ls
CASA を起動します。
# In bash
casa
CASA を起動すると、CASA のロガー・ウィンドウが新たに立ち上がります。また、
casa
コマンドを打ったターミナルには CASA プロンプトが表示されます。なお、CASA は Python で動いていますので、Python の文法が使えます。バージョン 6 系は Python 3 になります。CASA での作業を終えて CASA を終了するときは、
exit()
又は quit()
と入力します。
CASA の基本操作
CASA の操作に関するドキュメントとしては、CASA Guides の他の様々なページや cookbook があります (リンク集も参考にしてください)。ここでは、CASA を使い始めるにあたって知っておいた方が良いことを、少しだけ紹介します。 一部の Unix コマンドは、CASA でも利用できます。例えば、ls
や pwd
などです。
# In CASA
ls
pwd
次のように
os.system
を用いれば、すべての Unix/Linux コマンドを使うことができます。
# In CASA
os.system('ls')
また、CASA での解析操作を Python スクリプトとして記述しておき、それを一度に実行することもできます。スクリプトは、次のようにして実行します。
# In CASA
execfile('my_script.py',globals())
このチュートリアルは CASA ver. 6.4 を前提として書かれていますので、以下を実行して、CASA のバージョンを確認しておきます。
# In CASA
import casalith
version = casalith.version_string()
print ("You are using {}".format(version))
if (version < '6.4'):
print ("YOUR VERSION OF CASA IS TOO OLD FOR THIS GUIDE.")
print ("PLEASE UPDATE IT BEFORE PROCEEDING.")
else:
print ("Your version of CASA is appropriate for this guide.")
タスクの使用とデータの確認
データ解析の手始めはデータのヘッダー情報やデータ構造を調べることです。listobs
タスクを使用します。CASA でタスクを実行する方法は 2 つあります。 パラメータを
inp
で予め設定しておき go
で実行させる方法と、直接パラメータ込みのコマンドを実行する方法です。# In CASA
inp listobs
vis='sis14_twhya_calibrated_flagged.ms'
go
# In CASA
listobs(vis='sis14_twhya_calibrated_flagged.ms')
タスクの使い方やパラメータがわからないときは、
help
で調べます。# In CASA
help listobs
このヘルプは "q" で抜けることができます。
listobs
の出力は以下のようになります。
================================================================================
MeasurementSet Name: /lustre/naasc/sciops/qa2/dkim/casa_testing/MyTutorial/sis14_twhya_calibrated_flagged.ms MS Version 2
================================================================================
Observer: cqi Project: uid://A002/X327408/X6f
Observation: ALMA
Computing scan and subscan properties...
Data records: 80563 Total elapsed time = 5647.68 seconds
Observed from 19-Nov-2012/07:36:57.0 to 19-Nov-2012/09:11:04.7 (UTC)
ObservationID = 0 ArrayID = 0
Date Timerange (UTC) Scan FldId FieldName nRows SpwIds Average Interval(s) ScanIntent
19-Nov-2012/07:36:57.0 - 07:39:13.1 4 0 J0522-364 4200 [0] [6.05] [CALIBRATE_BANDPASS#ON_SOURCE,CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
07:44:45.2 - 07:47:01.2 7 2 Ceres 3800 [0] [6.05] [CALIBRATE_AMPLI#ON_SOURCE,CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
07:52:42.0 - 07:53:47.6 10 3 J1037-295 1900 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
07:56:23.5 - 08:02:11.3 12 5 TW Hya 8514 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
08:04:36.3 - 08:05:41.9 14 3 J1037-295 1900 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
08:08:09.6 - 08:13:57.3 16 5 TW Hya 10360 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
08:16:20.6 - 08:17:26.2 18 3 J1037-295 2100 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
08:19:53.9 - 08:25:41.7 20 5 TW Hya 10321 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
08:28:17.1 - 08:29:22.6 22 3 J1037-295 2100 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
08:32:00.5 - 08:37:48.2 24 5 TW Hya 10324 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
08:40:11.9 - 08:41:17.4 26 3 J1037-295 2100 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
08:43:45.6 - 08:49:33.4 28 5 TW Hya 9462 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
08:51:57.1 - 08:53:02.6 30 3 J1037-295 1900 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
08:58:12.0 - 09:00:28.1 33 6 3c279 3402 [0] [6.05] [CALIBRATE_BANDPASS#ON_SOURCE,CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
09:01:35.7 - 09:02:41.2 34 3 J1037-295 1900 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
09:05:15.6 - 09:07:31.6 36 5 TW Hya 4180 [0] [6.05] [OBSERVE_TARGET#ON_SOURCE]
09:09:59.1 - 09:11:04.7 38 3 J1037-295 2100 [0] [6.05] [CALIBRATE_PHASE#ON_SOURCE,CALIBRATE_WVR#ON_SOURCE]
(nRows = Total number of rows per scan)
Fields: 5
ID Code Name RA Decl Epoch SrcId nRows
0 none J0522-364 05:22:57.984648 -36.27.30.85128 J2000 0 4200
2 none Ceres 06:10:15.950590 +23.22.06.90668 J2000 2 3800
3 none J1037-295 10:37:16.079736 -29.34.02.81316 J2000 3 16000
5 none TW Hya 11:01:51.796000 -34.42.17.36600 J2000 4 53161
6 none 3c279 12:56:11.166576 -05.47.21.52464 J2000 5 3402
Spectral Windows: (1 unique spectral windows and 1 unique polarization setups)
SpwID Name #Chans Frame Ch0(MHz) ChanWid(kHz) TotBW(kHz) CtrFreq(MHz) BBC Num Corrs
0 ALMA_RB_07#BB_2#SW-01#FULL_RES 384 TOPO 372533.086 610.352 234375.0 372649.9688 2 XX YY
Sources: 5
ID Name SpwId RestFreq(MHz) SysVel(km/s)
0 J0522-364 0 - -
1 Ceres 0 - -
2 J1037-295 0 - -
3 TW Hya 0 - -
4 3c279 0 - -
Antennas: 21:
ID Name Station Diam. Long. Lat. Offset from array center (m) ITRF Geocentric coordinates (m)
East North Elevation x y z
1 DA42 A050 12.0 m -067.45.16.2 -22.53.29.3 43.0352 -744.9713 21.6702 2225079.880016 -5440041.377534 -2481724.598031
2 DA44 A068 12.0 m -067.45.20.6 -22.53.25.7 -82.4232 -631.7828 23.5810 2224981.097784 -5440131.250387 -2481621.066374
3 DA45 A070 12.0 m -067.45.11.9 -22.53.29.3 166.1833 -743.4934 19.8811 2225193.450167 -5439993.764157 -2481722.540534
4 DA46 A067 12.0 m -067.45.12.7 -22.53.27.2 142.4097 -678.7318 20.1280 2225181.070532 -5440026.290790 -2481662.975103
5 DA48 A046 12.0 m -067.45.17.0 -22.53.29.3 21.4267 -742.7987 21.6757 2225060.202580 -5440050.344436 -2481722.598651
6 DA49 A029 12.0 m -067.45.18.2 -22.53.25.8 -12.9134 -636.4552 22.1350 2225044.239583 -5440102.022535 -2481624.808405
7 DA50 A045 12.0 m -067.45.17.9 -22.53.30.1 -5.4183 -767.4398 22.6034 2225032.051652 -5440052.426015 -2481745.660003
9 DV02 A077 12.0 m -067.45.10.1 -22.53.25.9 217.6299 -637.5333 15.8376 2225255.259272 -5440008.987869 -2481623.352052
11 DV05 A082 12.0 m -067.45.08.3 -22.53.29.2 269.0433 -740.9521 15.7832 2225287.593766 -5439952.243679 -2481718.605314
12 DV06 A037 12.0 m -067.45.17.5 -22.53.28.8 6.7403 -727.3003 21.2086 2225048.729287 -5440061.085777 -2481708.139136
14 DV08 A021 12.0 m -067.45.17.2 -22.53.27.0 14.3196 -672.8108 21.3420 2225063.814715 -5440077.948261 -2481657.992572
15 DV10 A071 12.0 m -067.45.19.9 -22.53.23.5 -60.7887 -563.2541 23.3799 2225011.141945 -5440147.560932 -2481557.855663
16 DV13 A072 12.0 m -067.45.12.6 -22.53.24.0 147.1742 -580.5887 18.1825 2225199.254375 -5440058.161494 -2481571.803699
17 DV15 A074 12.0 m -067.45.12.1 -22.53.32.0 161.8159 -828.6196 18.7688 2225176.483514 -5439963.820451 -2481800.529842
18 DV16 A069 12.0 m -067.45.21.3 -22.53.30.2 -101.4797 -770.1047 23.2972 2224942.993176 -5440088.421459 -2481748.384855
19 DV17 A138 12.0 m -067.45.17.1 -22.53.34.4 19.1461 -901.2603 26.0137 2225036.269025 -5439997.853009 -2481870.267607
20 DV18 A053 12.0 m -067.45.17.3 -22.53.31.2 12.5939 -802.9941 21.5281 2225043.111690 -5440031.889497 -2481777.995870
21 DV19 A008 12.0 m -067.45.15.4 -22.53.26.8 67.5592 -667.6872 20.9574 2225113.709955 -5440059.310545 -2481653.122797
22 DV20 A020 12.0 m -067.45.17.8 -22.53.28.0 -2.9649 -703.4389 21.6629 2225043.419055 -5440073.737929 -2481686.333574
24 DV22 A011 12.0 m -067.45.14.4 -22.53.28.4 95.9131 -716.5005 21.0898 2225132.810230 -5440031.115405 -2481698.143589
25 DV23 A007 12.0 m -067.45.15.1 -22.53.27.3 74.0152 -681.2926 21.3231 2225117.809276 -5440052.280005 -2481665.79904
UV データをグラフィカルに確認するタスクは、
plotms
です。まず、データの UV カバレッジをプロットしてみます。# In CASA
plotms(vis='sis14_twhya_calibrated_flagged.ms', xaxis='u', yaxis='v', avgchannel='10000', avgspw=False, avgtime='1e9', avgscan=False, coloraxis="field", showgui=True)
図 1 のように GUI が表示されます。UV カバレッジのプロットの色は、異なる観測ターゲットを示しています。
他にも様々なプロットを見ることができます。例えば、Amplitude (振幅) vs. UV distance を見るためには、"Axes" タブを選択し、"X Axis" を "UVDist" に、"Y Axis" を "Amp" に設定します。 "Page" タブで、"Iteration Axis" を "Field" に設定します。そして、Plot ボタンを押すと、プロットが表示されます。
plotms GUI の下部にある緑色の矢印ボタンを使って、様々なキャリブレーター (較正天体) とサイエンスターゲットに進むことができます。 ケレス (Celes) を含むいくつかの天体では、図 2 のように UV 距離の関数として振幅の変化を見ることができます。これらの天体は分解されており、下でおこなうイメージングではっきりと見ることができます。
tclean
を用いたイメージング
CASA 4.7 から、clean
のイメージング機能は、tclean
というタスクに改良されて置き換えられました。このガイドでは tclean
のみを使用します。
TCLEAN と ALMA のガイド (英語) では、clean
と tclean
の構文の違いが詳細に書かれています。
参考までに、イメージングに tclean
ではなく clean
を使用する古いバージョンの CASA Guide は、
Archived ALMA CASAguides に保管されています (英語のみ)。このデータセットでは、位相較正天体はフィールド 3 で識別されます。この較正天体を "phase_cal" という画像ファイルにイメージ化することにします。 まず、以前このタスクを実行したことがある場合に備えて、この名前を使った古いバージョンの画像を削除しておきます。(画像化すると、拡張子が違うだけの同じルート名のファイルがいくつかできるので、".*" が必要です)。
# In CASA
os.system('rm -rf phase_cal.*')
それでは、
tclean
タスクを使ってイメージングをおこないます。tclean
オプションの完全なリストはドキュメントを参照して下さい。
- 最初に位相較正天体をイメージングし、イメージ名を "phase_cal" に設定します。MS 内の位相較正天体を特定するには、
listobs
の出力の "ScanIntent" 列を見ます。 位相較正天体は、"CALIBRATE_PHASE#ON_SOURCE" で示されています。この MS では、位相較正天体は J1037-295 であり、field = "3" で識別されます。 - 多周波合成 (specmode = mfs) を使用して、単一の連続波イメージを作成します。多周波合成は、選択されたすべてのスペクトルチャンネルからのデータを 1 つの連続波イメージに合成します。
観測が広い周波数範囲をカバーする場合、周波数によって天体の振幅や構造が大きく変化する可能性があります。これは、バンド幅の割合 (Δν/ν_center) が 10 % より大きい場合にのみ考慮する必要があります。
このデータセットでは問題にはなりません。そのため、deconvolver = hogbom に設定し、nterms = 1 を使用し、
tclean
にデコンボリューションされた各成分がすべての周波数で単一の振幅を持つように設定します。 - ポインティングは 1 つだけなので、gridder='standard' とします。複数の連続したポインティングがある場合、あるいは 7 [m] と 12 [m] のアレイデータを同時に撮像する場合 (それぞれ 1 ポインティングだけであっても)、gridder='mosaic' に設定する必要があります。
- セルサイズは 0.1 [arcsec] とし、ビーム全体で 5 ピクセルとした。経験則として、楕円ビームの最小方向で 5 ピクセル程度が望ましいということが知られているためです。
小規模なデータセットの場合、クイックなイメージングで適切なセルサイズを決定し、
tclean
で計算されたビームを記録することができるかもしれません。 しかし、多くの ALMA プロジェクトでは、これにはかなりの時間がかかります。uv-coverage をもう一度見ることで、セルサイズを推定することができます。 x 軸を "uvwave" に変えると、セルサイズはおおよそ 206265/(波長で最も長いベースライン)/(ビーム全体のセル数) となります。 この MS の場合、これは 0.09 [arcsec] になり、0.1 [arcsec] に切り上げることにします。tclean
は 0.58 [arcsec] × 0.51 [arcsec] の合成ビームサイズを出力しており、このセルサイズは適当であると言えるでしょう。 - 画像サイズは 128 × 128 に設定しました (ただし CASA では 2 の累乗は特殊な数ではないので、別の大きさでも構いません)。
これは観測周波数における主ビームの大部分をカバーするのに十分ですが、非点源に対してはもっと広い視野を撮像したいかもしれません。画像サイズは任意ですが、対称であるべきです。
もし設定した値を
tclean
があまりよくないと判断したならば、ログウィンドウでもっと良い値を提示するでしょう。
tclean
タスクは対話的モードで開始され、閾値、メジャーサイクル、マスクを手動でコントロールできます。
以前の clean
に慣れ親しんだユーザーにとっては、clean
タスクとは対照的に、niter (最大反復回数) を設定しないままにしておくと、
niter = 0 がデフォルトとなり、tclean
はクリーニングをおこなわないことに注意してください。
モデル内にクリーンなコンポーネントがない場合、例えば、与えられた MS に対して、このフィールドとスペクトルウィンドウに対して tclean
を最初に起動した場合、作成されるのはダーティキューブまたはイメージです。ここでは、まず niter を設定しない状態で tclean
を実行してみます。# In CASA
tclean(vis='sis14_twhya_calibrated_flagged.ms',
imagename='phase_cal',
field='3',
spw='',
specmode='mfs',
deconvolver='hogbom',
gridder='standard',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='briggs',
threshold='0.0mJy',
interactive=True)
CASA は次のような、クリーニングをしていないことを示す警告を出すはずです。
WARN task_tclean Restoring with an empty model image. Only residuals will be processed to form the output restored image.
もう一度コマンドを実行し、今度は niter を大きな数字に設定します。
# In CASA
os.system('rm -rf phase_cal.*')
tclean(vis='sis14_twhya_calibrated_flagged.ms',
imagename='phase_cal',
field='3',
spw='',
specmode='mfs',
deconvolver='hogbom',
gridder='standard',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='briggs',
threshold='0mJy',
niter=5000,
interactive=True)
tclean
タスクを実行すると、図 2 のような GUI が表示されます。
tclean
のビューワーで、ボタンが新しい楕円形のマスク領域を追加するように設定されていることを確認します。
楕円の中に "R" が表示されているアイコンをクリックする必要があるかもしれません。
点光源である二次較正天体が選択されているので、この時点でフィールドの中央に点光源が見えるはずです。
シグナルがある部分 (中央の点だけ) の周りに楕円形のマスクを描きます。楕円の内側をダブルクリックして、白くなるのを見ます (図 3 参照)。
マスクを設定するときは、シグナルの部分のみを選択し、それ以外はあまり含まないようにします。GUI の "Next Action" セクションに、幾つかのコントロールボタンがあります。
緑色の丸いボタンを押すと、クリーンが開始されます。これでクリーンのサイクルが実行され、その後戻ります。最初の clean の後、残差が表示されます。
残差を光源と比較し、十分だと判断したら(または閾値 (デフォルトでは 0 [mJy]) を満たしたら、つまり残差が負になったら)、赤い "X" をクリックすると tclean
は終了します。
この例では、2 回のクリーンでうまくいくはずです。
より複雑なターゲットの場合何回もクリーンを行う必要があり、各サイクルの後に、残差の見た目に基づいて新しいマスク領域を追加することが可能です。
tclean
でできたファイルを確認します。# In CASA
ls
"phase_cal" で始まる、拡張子だけが違ったファイルが幾つかできています。
- .image :
tclean
した最終的なイメージ - .mask : (最後のクリーンのサイクルで) 使用したマスクのイメージ
- .model :
tclean
で用いた (抽出した) clean components (単位は [Jy/pixel]) - .pb : 主ビームのレスポンス
- .residual : 最終的な残差マップ (最終イメージの "dirty" 部分に相当)
- .psf : 合成ビームのイメージ (dirty beam)
- .sumwt : 重みの和を含む単一ピクセルのイメージ
viewer()
又は、imview()
でおこないます。
UNIX / Linux のコマンドプロンプト (CASA 外) からは casaviewer
でスタンドアローンのビューワーを起動することもできます。# In CASA
imview("phase_cal.image")
ビューワーを使って対話的に他の画像も読み込んでみましょう。 ビューワーにある "Images" 用のテープデッキ・コントローラ (進む、戻る、止まるなどのボタン) を使って変更できます。
tclean
のパラメーターを変更する
tclean
にはたくさんのオプションがあります。inp tclean
と入力すれば、タスクの入力リストを見ることができます。
よく変更されるオプションの 1 つは、UV データをフーリエ平面画像に展開するのに使われる重み付けのパラメーターです。
この重み付けは、最初の例では robust=0.5
でした (デフォルト)。robust パラメータの値を -2 (uniform) と 2 (natural) の間でいくつか変えてみてください。
ビームサイズがどのように変化するかと、また最終画像のノイズに注目して下さい。
ノイズレベルは、ビューワーを使ってチェックできます。
画像上にボックスを描き、"Regions" ボックスの "Statistics" タブでその領域の画像の統計値を見てください ('Regions' ボックスが表示されない場合は、ビューワー GUI 上部の "View" ドロップダウンメニューからオンにしてください)。
以前に実行したことがある場合は、古いバージョンの画像を削除してください。# In CASA
os.system('rm -rf phase_cal_robust.*')
briggs weighting と robust = -1 を指定して
tclean
を呼び出します。
クリーンマスクを作成し、残差のレベルに満足するまで、クリーンを数サイクル実行します。# In CASA
tclean(vis='sis14_twhya_calibrated_flagged.ms',
imagename='phase_cal_robust',
field='3',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='briggs',
robust=-1.0,
threshold='0mJy',
niter=5000,
interactive=True)
結果を確認します。
# In CASA
imview("phase_cal_robust.image")
他の較正天体 (フィールド 0 と 2 。
listobs
出力で "ScanIntent" をもう一度確認してください) をイメージングして、画像サイズとセルサイズを大きくしたり小さくしたりしてみてください。
例えば、わずかに分解された一次 (フラックス) 較正天体である Ceres (フィールド 2; "ScanIntent" = "CALIBRATE_AMPLI#ON_SOURCE") を見てみましょう:# In CASA
os.system('rm -rf amp_cal_robust.*')
tclean(vis='sis14_twhya_calibrated_flagged.ms',
imagename='amp_cal_robust',
field='2',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='briggs',
threshold='0mJy',
niter=5000,
interactive=True)
imview("amp_cal_robust.image")
上の
plotms
で見たように、ケレスはある程度分解されています。
ピクセルサイズを非常に大きくしてみると、画像が壊れるのがわかるでしょう。
tclean
(tclean
はデコンボリューションをピクセル単位で量子化する) の目的から、ピクセルサイズは合成ビームに比べて小さくすることをお勧めします。
ピクセルサイズが合成ビームに比べて大きい場合、tclean
の方法とは関係なく、一般的に画像は劣化します。# In CASA
os.system('rm -rf amp_cal_bigpix.*')
tclean(vis='sis14_twhya_calibrated_flagged.ms',
imagename='amp_cal_bigpix',
field='2',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[32,32],
cell=['0.5arcsec'],
weighting='briggs',
threshold='0mJy',
niter=5000,
interactive=True)
imview("amp_cal_bigpix.image")
設定したピクセルサイズと最終的なビームサイズを比較して、どのような画像になるかを確認してみましょう。
較正とフラグの効果
ここでは、各種較正とフラグがどのような効果を与えるかを説明します。 イメージングの工程だけを見たい場合にはスキップできます。ここで扱っているデータは、様々な較正過程を経ています。これらの較正はどのような効果があったのでしょうか? 較正をおこなった場合とおこなわなかった場合、フラグを立てた場合と立てなかった場合で、二次較正天体の最終画像がどのように変わってくるかを試してみます。 まず、未較正のデータを作業ディレクトリーからコピーする必要があります。 これはデータを保存している場所によって、どこからコピーしてくるかが変わりますが、最終的に未較正のデータが現在の作業ディレクトリーにあれば問題ありません。 完全なデータパッケージをダウンロードしている場合は、以下のコマンドによって未較正データを作業ディレクトリーにコピーします。
# In CASA
os.system("tar -xvf ../working/sis14_twhya_uncalibrated.ms.tar")
os.system("rm -rf sis14_twhya_uncalibrated.ms")
os.system("cp -r ../working/sis14_twhya_uncalibrated.ms .")
較正済みデータセットのみをダウンロードした場合は、以下のように
wget
を使用して未較正データをダウンロードすることができる。# In bash
% wget -r -np -nH --cut-dirs=4 --reject "index.html*" https://bulk.cv.nrao.edu/almadata/public/working/sis14_twhya_uncalibrated.ms/
又は、このリンクからブラウザ経由でファイルをダウンロードすることもできます。
未較正データを
tclean
してみます。
二次較正天体 (フィールド 3) に対して、前と同じコマンドを実行してみます。# In CASA
os.system('rm -rf phase_cal_uncalibrated.*')
tclean(vis='sis14_twhya_uncalibrated.ms',
imagename='phase_cal_uncalibrated',
field='3',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='natural',
threshold='0mJy',
niter=5000,
interactive=True)
クリーンマスクをかけようとする天体が見つかればよいのですが…。
# In CASA
imview("phase_cal_uncalibrated.image")
生データ (それでも Tsys と WVR は補正されています)では、フィールド全体に較正天体のエコーが見えます。
次はフラグありなしによる違いです。フラグのついていないデータを作業ディレクトリーからコピーしてきます。
# In CASA
os.system("tar -xvf ../working/sis14_twhya_calibrated.ms.tar")
os.system("rm -rf sis14_twhya_calibrated.ms")
os.system("cp -r ../working/sis14_twhya_calibrated.ms .")
wget
でダウンロードすることもできます。# In bash
% wget -r -np -nH --cut-dirs=4 --reject "index.html*" https://bulk.cv.nrao.edu/almadata/public/working/sis14_twhya_calibrated.ms/
ブラウザ経由ならば、このリンクからダウンロードできます。 前と同じパラメーターで二次較正天体の、較正済み、かつ、フラグ処理をしていないデータを
tclean
してみます。# In CASA
os.system('rm -rf phase_cal_unflagged.*')
tclean(vis='sis14_twhya_calibrated.ms',
imagename='phase_cal_unflagged',
field='3',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[128,128],
cell=['0.1arcsec'],
weighting='natural',
threshold='0mJy',
niter=5000,
interactive=True)
imview("phase_cal_unflagged.image")
較正されていないデータとは違って、フラグ処理をしていないデータは幾分まともですが、残差に明らかに人工的な模様が見えます。 フラグ処理をすることによって最終的なイメージの質は向上しますが、今回のような全体に質のよいデータでは、フラグ処理をしていなくてもターゲットの目的天体も一応見えています。
目的天体のイメージング
様々な較正の過程は、目的天体のデータ (サイエンスデータ) を較正することにあります。ここでは、TW Hydra の観測はフィールド 5 がサイエンスデータになります。 基本的なイメージングのチュートリアルの最後のステップは、TW Hydra の連続波イメージングです。 まず、CASA のタスクsplit
を使って、サイエンスデータをデータセットに "分割" します。
これは必ずしも必要ではないですが、これはデータの管理を容易にする一般的なステップです。
同時に、width=8 を使用してデータを平滑化 (平均化) します。この時点では連続波イメージングが対象のため、多くの情報を失うことなくデータ量を減らすことができます。
このデータセットはすでに管理しやすいように設計されていますが、ALMA は非常に大きなデータセットを生成することがあるので、このスムージングのトリックは役に立つでしょう。# In CASA
os.system('rm -rf twhya_smoothed.ms')
split(vis='sis14_twhya_calibrated_flagged.ms', field='5', width='8', outputvis='twhya_smoothed.ms', datacolumn='data')
listobs('twhya_smoothed.ms')
次に、分割したデータの連続波イメージを作成します。 TW Hydra は、新しいデータセットではフィールド 0 として再ラベル化されています。 ここでも多周波合成モード ("mfs") を使い、ピクセルサイズを上よりもやや小さく (=イメージサイズをやや大きく) します。 これは TW Hydra は広がっている天体であるため、そして、"briggs" 重み付けを使うためビームがやや小さくなるためです。 ここでは対話的モードを指定し、閾値はとりあえず設定しないことにします。
# In CASA
os.system('rm -rf twhya_cont.*')
tclean(vis='twhya_smoothed.ms',
imagename='twhya_cont',
field='0',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[250,250],
cell=['0.08arcsec'],
weighting='briggs',
robust=0.5,
threshold='0mJy',
niter=5000,
interactive=True)
マスクツールを使って天体の周りにクリーンマスクを描き、TW Hydra ディスクからの放射がその周りの残差以下か同等になるまでクリーンを実行します。 この対話的なクリーニングモードでは、いつクリーンを停止させるかを決め、その時点で赤い X のボタンを押します。 代わりに、閾値まで自動的にクリーニングしたい場合は、上記の
tclean
の呼び出しの際に、或いはビューワーウィンドウに入力して、RMS ノイズの数倍の閾値を指定できます。
画像を見てください。TW Hydra は非常に明るく、ビームの大きさに対して大きいのがわかります。残差は完璧ではありませんが、この後の過程で改善していきます。# In CASA
imview("twhya_cont.image")
しかし、待ってください。 このスペクトルウィンドウにはスペクトル線 (輝線) があると言いませんでしたか?
はい、今回連続光として撮像したスペクトルウィンドウには、N2H+ の輝線も含まれています。 この場合、N2H+ の輝線は暗いので、イメージング前に輝線のチャンネルにフラグを立てなくても、最終的な連続体画像にほとんど違いはありません。 このチュートリアルと他の連続光ファーストルック・チュートリアルでは、イメージングとセルフキャリブレーションの基本的なステップに重点をおいているため、輝線を無視しています。 しかし、一般的には、連続体をイメージングする前に、自分のデータで輝線を含むチャンネルにフラグを立てるべきです。 スペクトルウィンドウで連続光をイメージングする前に、輝線にフラグを立てる方法については、 IRAS16293_Band9_-_Imaging (英語) の上級チュートリアルを参照してください。
非対話的モードの tclean
tclean
は、対話的なガイダンスなしで実行するように実行することもできます。
指定する 3 つの主なパラメーターは、停止する閾値 (クリーンする領域の最大残差がこの閾値より低いとき、クリーンのサイクルがは停止します)、
マスク (クリーンの工程がシグナルを識別する領域)、および最大反復回数です。
最大反復回数は厳密には必要ではなく、一般的にはフェイルセーフとして使用することが推奨されます。
つまり、もしクリーンのサイクルがその回数に到達したら何かが間違っている、というくらいの数値に設定するのがよいでしょう。まず、クリーンマスクを決定する方法を試してみましょう。 先ほど作った画像を見てください。天体からのものと思われる明らかな放射は全て、(100,100) から (150,150) のような範囲で囲まれたボックスに含まれています。
tclean
の呼び出しで "mask "パラメータを使って、このボックスをマスクに設定します。
また、ファイル (たとえば、以前の対話的バージョンの tclean
で作成されたもの) を指定して設定することもできます。次に、
tclean
の停止閾値を決めましょう。
もう一度、imview
を使って前の画像を開き、ノイズを推定するために、天体からかなり離れたところにボックスを描きます。
先ほどと同じように、ビューワーの "Regions" ボックスの "Statistics" タブを見て、RMS ノイズを求めます。
7 [mJy/beam] のようなものが見えます。閾値をこの約 2 倍、~15 [mJy/beam] に設定します。
通常、デコンボリューションされた画像に偽の天体が追加されるのを避けるために、rms ノイズの数倍の閾値が推奨されます。
つまり、tclean
がランダムなノイズスパイクを天体として扱い、画像からそれをデコンボリューションすることを避けたいからです。
これは、セルフキャリブレーションを行っている場合に特に問題となります。
最後に niter=10000 を設定します。この回数に達する前にクリーンが終了することを期待します。
これはクリーンが永遠に実行されないようにするための単なる大きな数字です。# In CASA
os.system('rm -rf twhya_cont_auto.*')
tclean(vis='twhya_smoothed.ms',
imagename='twhya_cont_auto',
field='0',
spw='',
specmode='mfs',
gridder='standard',
deconvolver='hogbom',
imsize=[250,250],
cell=['0.08arcsec'],
mask='box [ [ 100pix , 100pix] , [150pix, 150pix ] ]',
weighting='briggs',
robust=0.5,
threshold='15mJy',
niter=10000,
interactive=False)
imview('twhya_cont_auto.image')
この非対話的モードは、多くの時間を節約でき、非常に再現性が高いという利点があります。 interactive=True でクリーンプロセスを開始し、右上の青い矢印ボタンをクリックすると、"ハイブリッド" モードも利用できます。 これは
tclean
に最大反復回数または閾値に達するまで処理を進めるように指示する。
この組み合わせモードは、tclean
に使用するマスクを手動で描画できるところが利点です。
ビューワー GUI で、閾値と最大反復回数 (メジャーサイクル数とサイクルあたりの反復回数の積) の両方を手動で設定することもできます。
ただし、較正が不確かな画像、特に光源が明るい画像については、少なくとも初回は対話的にクリーニングするのがよいでしょう。
画像が "ダイナミックレンジに制限がある" (すなわち、品質が較正とデコンボリューションの精度によって設定される) 場合、正しいしきい値を予測することは困難です。
クリーニングのサイクルが進むにつれて、画像の残差は減少するはずです。
ロガーの tclean
に関する部分でこれを確認できます。残差が増加し始めたら、停止しきい値を変更する必要があるかもしれません。
プライマリービーム補正
tclean
の重要な点は、デフォルトでは tclean
によって生成されたイメージが、アレイ内の個々のディッシュのプライマリービーム (視野) に対して補正されないことです。
プライマリービームの応答は、通常、視野の中心で値 1 のガウスです。
天文学的に正しい空の画像を形成するには、tclean
の出力をこのプライマリービーム (モザイクの場合は、モザイクを作成するために使用されるプライマリービームパターンの組み合わせ) で分割する必要があります。tclean
の場合、これを行うには 2 つの方法があります。まず、
tclean
の実行時にパラメータ pbcor = True を設定する方法です。
これにより、.image.pbcor という拡張子を持つ追加画像が生成され、これがプライマリービーム用に補正されたクリーニング画像となります。もう一つの方法は、CASA はこの補正に必要なプライマリービーム情報が保存された拡張子 .pb のファイルを用いる方法です。 CASA のタスク
impbcor
を使用して、.pb 画像と tclean
からの出力画像を結合し、プライマリービーム補正された画像を生成することができます。まず、古いプライマリービーム補正画像があれば、それを削除します。
# In CASA
os.system('rm -rf twhya_cont.pbcor.image')
次にイメージを補正します。
# In CASA
impbcor(imagename='twhya_cont.image',
pbimage='twhya_cont.pb',
outfile='twhya_cont.pbcor.image')
出力されたイメージをチェックします。
# In CASA
imview('twhya_cont.pbcor.image')
プライマリービーム補正前のイメージで作業することは、ノイズがフィールド全体で同じであるため、非常に便利です (例えば、これはシグナルを探索するためのクリーンなデータセットです)。 しかし、サイエンスのためにフラックスや強度を計算する前に、この補正を適用することを忘れないことが非常に重要です。
ページ先頭に戻る
CASA guides "First Look at Imaging CASA 6.4" へ (本家のサイト:英語)
Last Update: 2023.09.30