はじめに
fits ファイル を扱うのが初めての人向けです．銀河の画像は親しみやすいので，これを例に google Colab だけでも動くサンプルです．もう少し上級者向けに書いた pythonのastroqueryを用いて検索し、skyviewでfits画像を取得する方法 もあります．
ここで紹介する内容は、google Colab 上での実行例 に一覧があります。もし動かない？という場合は、こちらをご覧ください。記事の中では、コードの一部を紹介してますので、もし動かない場合は、Colab上で確認ください。
基本編
サーバーとクライエントの基本
宇宙天文では，大量のデータをやりとりするので，サーバーとクライエント，について知識が必要になります．クライアントサーバ (client server system) とは，宇宙天文でも基本となるコンセプトなので，知っておこう．簡単には，サーバーとは，ネットワーク上のたくさんの機器に情報を提供できる機能が備わったOSやソフトウェア群のことで，クライエントはサーバーからデータを取得するだめのマシンのことです．
OSによってもかなり違いがあります．宇宙天文では，linuxベースのサーバがほぼ100%だと思いますので，基本的なサーバーとのやりとりや仕組みを事前に知っておくと，宇宙天文のデータを取得する際に少しは気が楽になると思います．サーバーOSとは？ などを参考ください．
urllib の例
サーバーとクライエントの仕組みを使う場合は，今ではとても簡単なので，まずは気楽に使ってみよう．例えば，NASAのheasoftという有名なページのファイルを取得してみよう．
import urllib.request
from bs4 import BeautifulSoup
url = 'https://heasarc.gsfc.nasa.gov/docs/software/heasoft/'
req = urllib.request.Request(url)
html = urllib.request.urlopen(req)
soup = BeautifulSoup(html, "html.parser")
topicsindex = soup.find('div', attrs={'class': 'navigation'})
これで，soup に全htmlの情報が取得されている．ブラウザはこれを表示しているのである．この中から，navigation という部分だけ取得したのが，topicsindex であり，
print(topicsindex) # extract sentences begining with div navigation
<div class="navigation">
<div id="nasaHeader">
<div>
....
のように，一部の情報だけ表示される．これは，クライエントからNASAのサーバーにアクセスし，データをもらうという作業で，数行のプログラムで書けることに注目してほしい．データを取得するというのは，ブラウザからマウスでボタンを色々と押してデータを取得する，というだけではないことを認識して欲しいのと，宇宙天文の場合はこのようなコマンドラインでデータを取得することが多い．
astroquery の設定
pip install astroquery
で astroquery をインストールしよう．これはskyviewやvizierなど，多々のカタログから，同じコーディングでデータを取得するためのツールです．
vizier の使い方
vizierは論文ごとにカタログを抽出したものである．例えば，blackhole で検索すると，大量のカタログを含む論文のリストが現れる．vizier上でも簡単なプロットなどがすぐにできる．論文とも紐付けられているので，もと情報にもすぐに辿れるのも便利．astroquery からは，
# download galaxy data via Vizier
from astroquery.vizier import Vizier
v = Vizier(catalog="J/ApJS/80/531/gxfluxes",columns=['Name',"logLx","Bmag","dist","type"],row_limit=-1)
data = v.query_constraints()
sname = data[0]["Name"]
namelist = []
olist=["HRI","DSS"]
def save(p,name,obs): # save file into fits
for onep,oneo in zip(p,obs):
onep.writeto(name+"_"+oneo+".fits",overwrite=True)
for one in sname: # change "N" to NGC, "I" to "IC"
name=one.strip().split()[0]
name=name.replace("N","NGC").replace("I","IC")
namelist.append(name)
J/ApJS/80/531/gxfluxes は銀河のカタログで，ここから，'Name',"logLx","Bmag","dist","type" のデータを抽出して，sname = data[0]["Name"] として， Name だけ取り出している．
skyview の使い方
skyview は，天体の画像を取得するのに便利なカタログです．今回はここから銀河の画像を取得してみましょう．
# get image from skyview
from astroquery.skyview import SkyView
name=namelist[0]
paths = SkyView.get_images(position=name, survey=olist)
save(paths,name,olist)
ここでは，position="天体名", survey=["DSS",..]などを指定して，
skyview から画像を取得し，save という関数で fits ファイルに保存する．
fitsファイルへの保存方法
def save(p,name,obs): # save file into fits
for onep,oneo in zip(p,obs):
onep.writeto(name+"_"+oneo+".fits",overwrite=True)
obs の中に書かれたカテゴリリストの名前をつけて，fits ファイルに保存する．
matplotlib を用いた天空座標でのプロット例
画像をプロットするのに，必要なモジュールをインポートする．
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
from astropy.wcs import WCS
from astropy.io import fits
import glob
import numpy as np
次に，fits ファイルを astropy で読んで，header から WCS 情報を引き出す．
WCS は，World Coordinate System (WCS) の意味で，天体画像の座標変換に必要な情報で，具体的には，検出器座標から天空座標や，銀経銀緯座標への変換に必要な情報です．座標変換は，matplotlib の projection というキーワードに指定するだけで変換してくれます．
# set IC342 imaga of DSS
dssname = "IC0342_DSS.fits"
dsshdu = fits.open(dssname)[0]
dsswcs = WCS(dsshdu.header)
dssdata = dsshdu.data
matplotlib でプロットする．
# plot image
F = plt.figure(figsize=(8,8))
plt.subplot(111, projection=dsswcs)
plt.imshow(dssdata, origin='lower', norm=LogNorm())
plt.colorbar()
このような絵が表示されれば成功です．
画像解析の例
画像の回転，転置，差分，など基本的な操作をして画像で遊んでみよう．
例えば，差分画像も簡単に表示できます．
fits ファイルとは，テキストで書かれた基礎情報と，バイナリで書かれたデータが合わさったファイル，というだけで，本当にそれだけなので，まずは苦手意識を持たずに扱おう．画像データは，openCV など，python にはたくさんのモジュールがあるので，まずはいろいろと遊んでみましょう．画像も，銀河でも波長が違えば見え方が違うし，超新星残骸，銀河団など，他の画像で遊んでみるのもよいです．画像の検索方法は，vizier でもよいですが，意外と欲しいものが見つからないこともあるので，適当にgoogleで検索した方が早いときもあります．．
モルフォロジー変換を用いた画像解析の例
専門的な天体解析では，望遠鏡の応答関数(像の広がり方，感度など)を用いて，画像中にある"点"が点源なのか，広がった構造なのかを調べます．例えば，Chandra衛星には，点源の解析のためのパターンマッチングを行う wavdetect というツールが用意されています．
ここではもう少し簡単に天体画像を遊べる方法を紹介します．モルフォロジー変換というのもので，画像を２値化(0 or 1)して，適当な小さい領域との and/or 演算を繰り返して画像の形状変換を行う方法です．例えば，銀河の面積が100ピクセルで繋がった構造をしていて，点源は一ピクセルだけだとすると，収縮(Erosion)処理(=適当な小さい領域に1ピクセルでも0があれば0とする演算)を行うと，ノイズは消えます．その次に，膨張(Dilation)処理(==適当な小さい領域に1ピクセルでも1があれば1とする演算)を行うと元に戻りますが，一旦0になったノイズは出てきません．このように，収縮(Erosion)と膨張(Dilation)を繰り返して，ノイズを処理する方法を，オープニング処理と呼びます．cv2.morphologyEx(binarydssdata,cv2.MORPH_OPEN, kernel) を用いて実行できます．モルフォロジー変換が参考になりました．大元はC言語で書かれていて，Morphology_2_8cpp-example の実行例があります．
使い方は簡単で，画像を，
# convert the data into binary image
binarydssdata = 1.0 * (dssdata > np.median(dssdata))
のように，median よりも大小で２値化します．int だと動かない関数があるようなので，float の 1.0 にします．
次に，getStructuringElement(cv2.MORPH_ELLIPSE,(8,8)) で，8x8の画像で，楕円形の 1, 0 の行列を生成して，これを kernel として，cv2.morphologyEx として変換します．
# 1.0 is used because morphologyEx may not work for int32
import cv2
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(8,8))
opening = cv2.morphologyEx(binarydssdata,cv2.MORPH_OPEN, kernel)
コードはこれだけです．実行例はこちらです．
左から順に，「元画像」「２値化画像」「モルフォロジー変換後」で，一番右は，「元画像」 x 「モルフォロジー変換後」をしたものです．実際の天体解析でも，点源除去は同じような操作をします．
google Colab での実行例
google Colab 上で実行するには，astroquery がデフォルトではインストールされてないので，astroquery の install が必要です．google Colab 上での実行例 を参照ください．
また，fits ファイルを google drive 上に保存しないと，データは外に出せないので，必要ならgoogle drive をマウントして，そこに fits ファイルを書き出しましょう．