'''
Created on 20 mars 2018
based on 'Tools.py' python program of Cyril Lachaud
Generic functions to manipulate fits file
@author: Catalano Camille, APC
'''
import os
import numpy as np
import matplotlib.pyplot as plt
from astropy.io import fits
import astropy.time as astrotime
import common.instru.array_convention as ac
from common.mission.time import get_mjd_ref
##############################
[docs]def common_keyword(fits_header, tstart, tstop, creator='ECLAIRs_pipeline', version='0.1', proc_id="01"):
"""add to the fits header the common keywords used in all the pipeline products
.. seealso:: doc : ECL-GP-FITS-KEYWORDS.pdf ECLAIRs GP fits keywords
https://forge.in2p3.fr/dmsf/files/6731/view
:param fits_header: header data unit to add the keyword (ex: fits.PrimaryHDU)
:param creator: name of the program that has created the fits file
:param version: fits file format version
:param proc_id: id of the pipeline process
:type fits_header: astropy.io.fits.*HDU
:type creator: string (default='ECLAIRs_pipeline')
:type version: string (default='0.1')
:type proc_id: string (default="01")
:param tstart: start time of the observation in s from MJDREF
:type tstart: float
:param tstop: end time of the observation in s from MJDREF
:type tstop: float
"""
fits_header.header['TELESCOP'] = ('SVOM', 'Telescop/satellite name')
fits_header.header['INSTRUME'] = ('ECLAIRs', 'Instrument/detector')
fits_header.header['TIMESYS'] = ('TT', 'Time frame system')
fits_header.header['TIMEUNIT'] = ('s', 'Time unit')
fits_header.header['TIMEREF'] = ('LOCAL', 'Time reference frame')
fits_header.header['MJDREF'] = (get_mjd_ref(), 'Modified Julian Date of origin')
fits_header.header['REVOL'] = ('', 'Revolution number')
fits_header.header['OBSID'] = ('', 'Observation Identifier (type_number)')
fits_header.header['OBS_TYPE'] = ('', 'Slew and pointing status')
fits_header.header['OBSBOUND'] = ('', 'Reason for Observation ending')
fits_header.header['OBTSTART'] = ('', 'OBT of the start of the Observation') #obsolete TBC
fits_header.header['OBTEND'] = ('', 'OBT of the end of the Observation') #obsolete TBC
fits_header.header['TSTART'] = (tstart, 'Start time of the Observation')
fits_header.header['TSTOP'] = (tstop, 'End time of the Observation')
fits_header.header['EXTREL'] = (version, 'FSC release number')
fits_header.header['CREATOR'] = (creator, 'software name and version')
fits_header.header['CONFIGUR'] = ('', 'Software system configuration')
fits_header.header['DATE'] = (astrotime.Time.now().tt.isot, 'file creation date YYYY-MM-DDThh:mm:ss UT')
fits_header.header['STAMP'] = ('', '')
fits_header.header['FSCLEVL'] = ('', 'FSC level of data processing')
fits_header.header['ORIGIN'] = ('FSC', 'Origin of FITS file')
fits_header.header['PROCID'] = (proc_id, 'process id')
fits_header.header['EXTREL'] = ('', 'FSC release number')
fits_header.header['EXTVER'] = ('1')
fits_header.header['GRPID1'] = ('-1', 'EXTVER of Group containing this HDU')
fits_header.header['GRPLC1'] = ('', 'URL of file containing Group')
fits_header.header['LONGSTRN'] = ('OGIP 1.0', 'The HEASARC Long String Convention may be used')
fits_header.header['HDUCLASS'] = ('OGIP', 'Format conforms mostly to OGIP standards')
fits_header.header['HDUVERS'] = ('1.1.0', 'Version of format')
fits_header.header['CHANTYPE'] = ('PI', 'Type of detector channels')
##############################
[docs]def grouphdu_keyword(fits_header, grp_table_name):
"""
add to the fits header the specific keywords for a group HDU
:param fits_header: header data unit to add the keyword (ex: fits.PrimaryHDU)
:type fits_header: astropy.io.fits.*HDU
:param grp_table_name: grouping table name
:type grp_table_name: string
"""
fits_header.header['EXTNAME'] = ('GROUPING', 'HDU contains a Grouping Table')
fits_header.header['GRPNAME'] = (grp_table_name, 'Grouping Table name')
### possibly other keywords to add
##############################
[docs]def energy_keywords(fits_header, chanmin, chanmax, e_min, e_max):
"""add to fits header the keywords related to energy
:param fits_header: header data unit to add the keyword (ex: fits.PrimaryHDU)
:type fits_header: astropy.io.fits.*HDU
:param chanmin: lowest channel of the energy range
:type chanmin: float
:param chanmax: highest channel of the energy range
:type chanmax: float
:param e_min: lower energy limit in keV
:type e_min: float
:param e_max: upper energy limit in keV
:type e_max: float
"""
fits_header.header['chanmin'] = (chanmin, 'lowest channel of the energy range')
fits_header.header['chanmax'] = (chanmax, 'highest channel of the energy range')
fits_header.header['e_min'] = (e_min, 'lower energy limit in keV')
fits_header.header['e_max'] = (e_max, 'upper energy limit in keV')
##############################
[docs]def wcs_keyword(fits_header, attitude):
"""add to the fits header the wcs keywords
.. seealso:: doc : ECL-GP-FITS-KEYWORDS.pdf ECLAIRs GP fits keywords
https://forge.in2p3.fr/dmsf/files/6731/view
.. warning:: deprecated. Use of SkyImages.wcs.header() instead.
:param fits_header: header data unit to add the keyword (ex: fits.PrimaryHDU)
:type fits_header: astropy.io.fits.*HDU
:param attitude: [ra, dec, ori]
:type attitude: [float, float, float]
"""
cdelt1 = -0.5633
cdelt2 = 0.5633
fits_header.header['CTYPE1'] = ('RA---TAN', 'first parameter RA, TANgential projection')
fits_header.header['CTYPE2'] = ('DEC--TAN', 'second parameter DEC, TANgential projection')
fits_header.header['CRPIX1'] = (100, 'reference pixel')
fits_header.header['CRPIX2'] = (100, 'reference pixel')
fits_header.header['EQUINOX'] = (2000, 'equinox of coordinates')
fits_header.header['RADECSYS'] = ('FK5', 'stellar reference frame')
fits_header.header['CDELT1'] = (cdelt1, 'angular scale of the ref pixel in degrees/pixel')
fits_header.header['CDELT2'] = (cdelt2, 'angular scale of the ref pixel in degrees/pixel')
fits_header.header['LONGPOLE'] = (180, 'native longitude of celestial pole')
fits_header.header['LATPOLE'] = (0, 'native latitude of celestial pole')
fits_header.header['CRVAL1'] = (attitude[0], 'right ascension in degrees of ref pixel')
fits_header.header['CRVAL2'] = (attitude[1], 'declination in degrees of ref pixel')
fits_header.header['CROTA2'] = (-90-attitude[2], 'rotation in degrees')
fits_header.header['CD1_1'] = (np.cos(np.radians(-90-attitude[2])) * cdelt1, 'coord. transf. matrix in degrees/pixel')
fits_header.header['CD2_1'] = (-np.sin(np.radians(-90-attitude[2])) * cdelt2, 'coord. transf. matrix in degrees/pixel')
fits_header.header['CD1_2'] = (np.sin(np.radians(-90-attitude[2])) * cdelt1, 'coord. transf. matrix in degrees/pixel')
fits_header.header['CD2_2'] = (np.cos(np.radians(-90-attitude[2])) * cdelt2, 'coord. transf. matrix in degrees/pixel')
##############################
[docs]def write_fits(filename,data):
"""Save an array into a fits file
data are saved in the primaryHDU
:param filename: PATH/filename of the fits file
:param data: data to write in the fits file
:type data: array
:type filename: string (PATH)
"""
try:
os.remove(filename)
except OSError:
pass
hdu = fits.PrimaryHDU(data)
hdu.writeto(filename)
##############################
[docs]def read_fits(name, ugts_standard=0):
"""Read an array from a fits file
.. warning:: The structure of the fits file must be simple: data into hdu[0].data
:param name: PATH/name of the fits file
:type name: string
:param ugts_standard: flag to read from ugts Idet,Jdet orientation and make the transformation to Yecl,Zecl
:type ugts_standard: bool
:return: the data from the fits file
:rtype: array
"""
hdu = fits.open(name)
tdata=hdu[0].data
if ugts_standard:
tdata = ac.ijdet_to_yzegp_array(tdata)
return tdata
##############################
[docs]def set_array_from_file(filename, ugts_standard=0):
"""read a fits file (using read_fits)
:param filename: PATH/filename of the fits file
:type filename: string
:param ugts_standard: flag to read from ugts Idet,Jdet orientation and make the transformation to Yecl,Zecl
:type ugts_standard: bool
:return: array of the data, len of the first dim (=size of the square array)
:rtype: array, int
"""
the_array = read_fits(filename, ugts_standard)
if the_array.shape[0]!=the_array.shape[1]:
print("ERROR : arrays have to have same size in both dimensions")
return the_array,the_array.shape[0]
##############################
[docs]def set_array_and_header_from_file(filename, ugts_standard=0):
"""read a fits file
:param filename: PATH/filename of the fits file
:type filename: string
:param ugts_standard: flag to read from ugts Idet,Jdet orientation and make the transformation to Yecl,Zecl
:type ugts_standard: bool
:return: array of the data, len of the first dim (=size of the square array), header object
:rtype: array, int, HDUList
"""
hdu = fits.open(filename)
the_array=hdu[0].data
if the_array.shape[0]!=the_array.shape[1]:
print("ERROR : arrays have to have same size in both dimensions")
if ugts_standard:
the_array = ac.ijdet_to_yzegp_array(the_array)
return the_array,the_array.shape[0], hdu[0].header
##############################
[docs]def image_plot(image,save=None,show=False,title=None,clim=None,scale_legend=None, **kwargs):
"""plot/show/save the image of the array
the function uses imshow
:param image: data array to plot
:param save: PATH/name of the file to save the fig. If None, no saving.
:param show: show the plot
:param title: title of the plot
:param clim: color limit
:param scale_legend: legend for the color scale
:param **kwargs: other arguments directly passed to plt.imshow
:type image: float array
:type save: string (default=None)
:type show: bool (default=False)
:type title: string (default=None)
:type clim: (float, float) (default=None)
:type scale_legend: string (default=None)
.. warning:: clim has not been tested.
"""
plt.imshow(image,interpolation="nearest",origin='lower',**kwargs)
if title!=None:
plt.title(title)
if clim!=None :
plt.clim(clim)
cbar=plt.colorbar()
if scale_legend!=None:
cbar.set_label(scale_legend)
if save!=None:
plt.savefig(save)
if show==True:
plt.show()