fitsをアタマから作りたい

# はじめに
天文データでよく使われるfitsと仲良くなれたと思ったときに記録を残していきます。
これはfitsとはなんぞやというのを勉強するために1から作ってみました。

まず、fitsは天文データによく使われるファイル形式であり、座標、周波数といった情報を格納するヘッダー部と、観測から得られた情報を格納するデータ部の2から構成される。詳しくはwikipediaで...(逃げ)
[https://ja.wikipedia.org/wiki/FITS](https://ja.wikipedia.org/wiki/FITS)

#fitsの作成
なんだかんだで悩んだfitsの作成。本当にこの方法でいいのか。。。いいと信じたい。
とりあえず[100, 100]の配列を持つfitsを作るスクリプトを示します。

```ruby:make_fits_1.py
import numpy as np
from astropy.io import fits

x = 100      #x_axis pixel number
y = 100      #y_axis pixel number
z = 100      #z_axis pixel number

###配列の作成###
#3次元の時
data = np.zeros((z, y, x), dtype=np.float32)
#2次元の時
#data = np.zeros((y, x), dtype=np.float32)

hdu = fits.PrimaryHDU(data = data)

#作ったfitsの保存
hdu.writeto('make_fits.fits',overwrite=True)
```

3次元cubeのfitsを作るときは「**配列の作成**」の部分を変更します。
この状態ではまだヘッダーを入れていないのでds9などで表示しても味気ないのしか出てこない。
![スクリーンショット 2019-11-28 0.05.13.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/316070/40ce94f0-833c-057b-f081-25c698c039c0.png)

ヘッダーを入れていきましょう。
ヘッダーの追加はこれでできる。
### ヘッダーの追加

```ruby:
hdu.header['keyword'] = 'something' #keywordの時
hdu.header['keyword'] = something   #座標といった数値情報を入れる時
```

お試しに、100x100 の配列をもち、
WCS: Galactic
center position 00:00:00.00 00:00:00.00 gal
beamsize 20.0 arcsec
gridsize 10.0 arcsec
のfitsを作ってみる。

```ruby:make_fits_2.py
import numpy as np
from astropy.io import fits

x = 100      #x_axis pixel number
y = 100      #y_axis pixel number
z = 100      #z_axis pixel number

beamsize = 20.0 #arcsec

###配列の作成###
#3次元の時
data = np.zeros((z, y, x), dtype=np.float32)
#2次元の時
data = np.zeros((y, x), dtype=np.float32)

hdu = fits.PrimaryHDU(data = data)

#入れたいヘッダ一覧                                      
hdu.header['BITPIX']  =                  -64                                                 
hdu.header['NAXIS']   =                    2                                                 
hdu.header['NAXIS1']  =                  x                                                  
hdu.header['NAXIS2']  =                  y
hdu.header['BPA']     = 0.000000
hdu.header['BMAJ']    = beamsize/3600                                                  
hdu.header['BMIN']    = beamsize/3600
hdu.header['EPOCH']   =                 2000                                                                                                    
hdu.header['BUNIT']  = 'K       '                                                                                                    
hdu.header['CTYPE1']  = 'GLON-GLS'                                                            
hdu.header['CRVAL1']  =                  0.0                                                  
hdu.header['CDELT1']  =     -beamsize/2/3600                                                  
hdu.header['CRPIX1']  =                  x/2                                                  
hdu.header['CROTA1']  =                  0.0                                                  
hdu.header['CUNIT1']  = 'deg     '                                                            
hdu.header['CTYPE2']  = 'GLAT-GLS'                                                            
hdu.header['CRVAL2']  =                  0.0                                                  
hdu.header['CDELT2']  =      beamsize/2/3600                                                  
hdu.header['CRPIX2']  =                  y/2                                                  
hdu.header['CROTA2']  =                  0.0                                                  
hdu.header['CUNIT2']  = 'deg     '          

#作ったfitsの保存
hdu.writeto('make_fits_2.fits',overwrite=True)
```

これでちゃんと座標情報が入ったfitsが作れました。
![スクリーンショット 2019-11-28 0.48.52.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/316070/cf3a4a64-9170-24c0-4005-5ad83bd20010.png)

#fitsを作る注意点(自分への戒め)
ヘッダーを間違えると間違ったデータになっちゃうので確認はしっかりしましょう。