[PYTHON] I want to make fits from my head
Introduction
I will keep a record when I think that I have become friends with fits, which is often used in astronomical data. I made this from scratch to learn what fits are.
First, fits is a file format often used for astronomical data, and consists of two parts: a header part that stores information such as coordinates and frequencies, and a data part that stores information obtained from observations. For more information on wikipedia ... (escape) https://ja.wikipedia.org/wiki/FITS
Creating fits
Creating fits that I was worried about for some reason. Is this really the way to go? .. .. I want to believe it is good. Here is a script that creates fits with an array of [100, 100] for the time being.
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
###Creating an array###
#3D time
data = np.zeros((z, y, x), dtype=np.float32)
#Two-dimensional time
#data = np.zeros((y, x), dtype=np.float32)
hdu = fits.PrimaryHDU(data = data)
#Save the fits you made
hdu.writeto('make_fits.fits',overwrite=True)
When making fits of 3D cube, change the part of "** Create array **".
In this state, the header has not been inserted yet, so even if it is displayed with ds9 etc., it will only come out as dull.
Let's put in the header. You can now add headers.
Add header
python
hdu.header['keyword'] = 'something' #At the time of keyword
hdu.header['keyword'] = something #When entering numerical information such as coordinates
As a trial, it has a 100x100 array and WCS: Galactic center position 00:00:00.00 00:00:00.00 gal beamsize 20.0 arcsec gridsize 10.0 arcsec Let's make fits.
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
###Creating an array###
#3D time
data = np.zeros((z, y, x), dtype=np.float32)
#Two-dimensional time
data = np.zeros((y, x), dtype=np.float32)
hdu = fits.PrimaryHDU(data = data)
#List of headers you want to put
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 '
#Save the fits you made
hdu.writeto('make_fits_2.fits',overwrite=True)
Now you have a fits that contains the coordinate information.
Precautions for making fits (commandments to yourself)
If you make a mistake in the header, the data will be incorrect, so be sure to check it.
Recommended Posts