Commit 9e12ce0c authored by Vigeesh Gangadharan's avatar Vigeesh Gangadharan
Browse files

cleanup

parent b11ed04a
Pipeline #2973 passed with stage
in 3 minutes and 6 seconds
......@@ -306,7 +306,6 @@ class FitsCardGroup:
self.cardlist = cardlist
self.header = None
def to_header(self):
"""
Constructs the header
......@@ -332,10 +331,12 @@ def get_imap_list(vfisv_params):
"""Get the imap number as a list corresponding to the sorted filenames"""
return (vfisv_params.imap - 1).astype("int")
def get_istep_list(vfisv_params):
"""Get the istep number as a list corresponding to the sorted filenames"""
return (vfisv_params.istep - 1).astype("int")
def get_result_shape_slit(vfisv_params):
"""Get the shape of the data"""
......@@ -345,6 +346,7 @@ def get_result_shape_slit(vfisv_params):
return [nmaps, nlength, nsteps]
def get_result_shape_ifu(vfisv_params):
"""Get the shape of the data"""
......@@ -354,6 +356,7 @@ def get_result_shape_ifu(vfisv_params):
return [nmaps, nlength, nsteps]
def get_stokes_shape(vfisv_params):
"""Get the shape of the data"""
......@@ -364,6 +367,7 @@ def get_stokes_shape(vfisv_params):
return [nmaps, nlength, nsteps, nwave]
def get_stokes_shape_ifu(vfisv_params):
"""Get the shape of the data"""
......@@ -374,6 +378,7 @@ def get_stokes_shape_ifu(vfisv_params):
return [nmaps, nlength, nsteps, nwave]
def get_quantity(data, vfisv_params, index: int, units=Unit.dimensionless_unscaled):
"""
Constructs a numpy array for the requested quantity
......@@ -406,7 +411,7 @@ def get_quantity(data, vfisv_params, index: int, units=Unit.dimensionless_unscal
# reorder according to imap and istep.
# TODO: Check for the flip that is caused while reading the fits file
if vfisv_params.header['INSTRUME'] == 'GRIS':
if vfisv_params.header["INSTRUME"] == "GRIS":
# the lx_range is the range the filenames are ordered.
lx_range = np.arange(vfisv_params.LX)
# Get the shape for the output
......@@ -415,14 +420,20 @@ def get_quantity(data, vfisv_params, index: int, units=Unit.dimensionless_unscal
quantity = np.zeros(data_shape)
# distribute data into maps from the requested index
quantity[imaps, :, isteps] = data[lx_range, ::-1, index]
elif vfisv_params.header['INSTRUME'] == 'GRIS-IFU':
elif vfisv_params.header["INSTRUME"] == "GRIS-IFU":
# Get the shape for the output
data_shape = get_result_shape_ifu(vfisv_params)
# initialize the empty numpy array
quantity = np.zeros(data_shape)
# distribute data into maps from the requested index
for i_map in np.arange(vfisv_params.header['NMAPS']):
quantity[i_map, :, :] = data[i_map*vfisv_params.header['NAXIS1']:vfisv_params.header['NAXIS1']*(i_map+1), :, index].T
for i_map in np.arange(vfisv_params.header["NMAPS"]):
quantity[i_map, :, :] = data[
i_map
* vfisv_params.header["NAXIS1"] : vfisv_params.header["NAXIS1"]
* (i_map + 1),
:,
index,
].T
return quantity * units
......@@ -455,26 +466,33 @@ def get_stokes_fits(data, vfisv_params):
# reorder according to imap and istep.
# TODO: Check for the flip that is caused while reading the fits file
# quantity[imaps, :, isteps, :] = data[lx_range, ::-1, :]
if vfisv_params.header['INSTRUME'] == 'GRIS':
# quantity[imaps, :, isteps, :] = data[lx_range, ::-1, :]
if vfisv_params.header["INSTRUME"] == "GRIS":
lx_range = np.arange(vfisv_params.LX)
# Get the shape for the output
data_shape = get_stokes_shape(vfisv_params)
# initialize the empty numpy array
quantity = np.zeros(data_shape)
quantity[imaps, :, isteps, :] = data[lx_range, ::-1, :]
elif vfisv_params.header['INSTRUME'] == 'GRIS-IFU':
elif vfisv_params.header["INSTRUME"] == "GRIS-IFU":
# print('GRIS-IFU')
data_shape = get_stokes_shape_ifu(vfisv_params)
# initialize the empty numpy array
quantity = np.zeros(data_shape)
for i_map in np.arange(vfisv_params.header['NMAPS']):
quantity[i_map, ...] = np.einsum("ijk->jik",data[
i_map * vfisv_params.header['NAXIS1']:vfisv_params.header['NAXIS1'] * (i_map + 1),
...])
for i_map in np.arange(vfisv_params.header["NMAPS"]):
quantity[i_map, ...] = np.einsum(
"ijk->jik",
data[
i_map
* vfisv_params.header["NAXIS1"] : vfisv_params.header["NAXIS1"]
* (i_map + 1),
...,
],
)
return quantity
def process_vfisv_output(data, sifit, sqfit, sufit, svfit, vfisv_params):
"""
......@@ -1010,7 +1028,7 @@ def create_map(
dpi=100,
cticks=None,
cbar_title=None,
#bbox=None,
# bbox=None,
annotate_starttime=False,
annotate_endtime=False,
):
......@@ -1083,10 +1101,10 @@ def create_map(
fig.suptitle(title)
# set lon_lat axis
if header['INSTRUME']=='GRIS':
set_axis_lon_lat(ax, grid=grid,spacing = 5 * Unit.arcsec)
elif header['INSTRUME']=='GRIS-IFU':
set_axis_lon_lat(ax, grid=grid)#, spacing = 50 * Unit.arcsec)
if header["INSTRUME"] == "GRIS":
set_axis_lon_lat(ax, grid=grid, spacing=5 * Unit.arcsec)
elif header["INSTRUME"] == "GRIS-IFU":
set_axis_lon_lat(ax, grid=grid) # , spacing = 50 * Unit.arcsec)
# additonal axis
set_additional_axis(ax, header)
......@@ -1137,11 +1155,11 @@ def create_map(
# else:
plt.savefig(filename, bbox_inches="tight")
#plt.show()
# plt.show()
# close plot
plt.close("all")
return None#bbox
return None # bbox
def plot_image(
......@@ -1245,7 +1263,7 @@ def plot_image(
annotate_endtime = False
if imap == 0:
annotate_starttime = True
bbox=None
bbox = None
if imap == header["RMAPS"] - 1:
annotate_endtime = True
......@@ -1263,7 +1281,7 @@ def plot_image(
cbar_title=cbar_title,
annotate_starttime=annotate_starttime,
annotate_endtime=annotate_endtime,
#bbox=bbox
# bbox=bbox
)
filelist = glob(os.path.join(td.name, "*"))
......@@ -1614,26 +1632,24 @@ class VFISVpackage:
# Get the pixels that define the boundary
self.get_pixel_bounds()
#get data depending on slit or IFU
if self.header['INSTRUME']=='GRIS':
#print('GRIS')
# get data depending on slit or IFU
if self.header["INSTRUME"] == "GRIS":
# print('GRIS')
self.get_data_slit()
elif self.header['INSTRUME']=='GRIS-IFU':
#print('GRIS-IFU')
elif self.header["INSTRUME"] == "GRIS-IFU":
# print('GRIS-IFU')
self.get_data_ifu()
self.get_map_props()
self.get_line_props()
self.get_rel_wave()
self.nmaps = self.header["NMAPS"] # Number of maps
self.nsteps = self.header["NSTEPS"] # Number of steps per map
self.LX = self.header["NAXIS1"]*len(self.filenames) # width
self.LX = self.header["NAXIS1"] * len(self.filenames) # width
self.LY = self.header["NAXIS2"] # Length of the slit
self.LZ = 16 # VFISV data structure
@staticmethod
def get_line_list():
"""Get the list of available lines"""
......@@ -1684,12 +1700,12 @@ class VFISVpackage:
def get_map_props(self):
"""Create the integer parameters buffer"""
if self.header['INSTRUME']=='GRIS':
#print('GRIS')
if self.header["INSTRUME"] == "GRIS":
# print('GRIS')
LX = len(self.filenames)
elif self.header['INSTRUME']=='GRIS-IFU':
#print('GRIS-IFU')
LX = self.header['NAXIS1']*len(self.filenames)
elif self.header["INSTRUME"] == "GRIS-IFU":
# print('GRIS-IFU')
LX = self.header["NAXIS1"] * len(self.filenames)
self.map_props = np.array(
(
......@@ -1781,13 +1797,13 @@ class VFISVpackage:
self.date_end = datetime.fromisoformat(header_end["DATE-BEG"])
self.headerlist = []
if self.header['INSTRUME'] == 'GRIS-IFU':
for ifiles in range(len(self.filenames)):
_, ifileheader = fits.getdata(self.filenames[ifiles], header=True, memap=True)
if self.header["INSTRUME"] == "GRIS-IFU":
for ifiles in range(len(self.filenames)):
_, ifileheader = fits.getdata(
self.filenames[ifiles], header=True, memap=True
)
self.headerlist.append(ifileheader)
def get_pixel_bounds(self):
"""Get the +/- width Angstrom bounds centered at the line"""
......@@ -1902,17 +1918,17 @@ class VFISVpackage:
# Initialize the Stokes vector and continuum intensity arrays for VFISV
self.SI = np.zeros(
(len(self.filenames)*nx, ny, self.pix_end - self.pix_ini + 1),
(len(self.filenames) * nx, ny, self.pix_end - self.pix_ini + 1),
dtype="float32",
)
self.SQ = np.zeros_like(self.SI, dtype="float32")
self.SU = np.zeros_like(self.SI, dtype="float32")
self.SV = np.zeros_like(self.SI, dtype="float32")
self.ICONT = np.zeros([len(self.filenames)*nx, ny], dtype="float32")
self.ICONT = np.zeros([len(self.filenames) * nx, ny], dtype="float32")
self.imap = np.zeros(len(self.filenames)*nx, dtype="int")
self.istep = np.zeros(len(self.filenames)*nx, dtype="int")
self.time = np.zeros(len(self.filenames)*nx, dtype="datetime64[s]")
self.imap = np.zeros(len(self.filenames) * nx, dtype="int")
self.istep = np.zeros(len(self.filenames) * nx, dtype="int")
self.time = np.zeros(len(self.filenames) * nx, dtype="datetime64[s]")
# # Iterate through the L1 split files/maps
# # for iff in tqdm(
......@@ -1943,19 +1959,16 @@ class VFISVpackage:
):
(
self.SI[nx*iff:nx*(iff+1), ...],
self.SQ[nx*iff:nx*(iff+1), ...],
self.SU[nx*iff:nx*(iff+1), ...],
self.SV[nx*iff:nx*(iff+1), ...],
self.ICONT[nx*iff:nx*(iff+1), :],
self.imap[nx*iff:nx*(iff+1)],
#self.istep[iff],
self.time[nx*iff:nx*(iff+1)],
self.SI[nx * iff : nx * (iff + 1), ...],
self.SQ[nx * iff : nx * (iff + 1), ...],
self.SU[nx * iff : nx * (iff + 1), ...],
self.SV[nx * iff : nx * (iff + 1), ...],
self.ICONT[nx * iff : nx * (iff + 1), :],
self.imap[nx * iff : nx * (iff + 1)],
# self.istep[iff],
self.time[nx * iff : nx * (iff + 1)],
) = self.read_fits_ifu(self.filenames[iff])
# Switch to a FORTRAN array
self.SI = np.asfortranarray(self.SI)
self.SQ = np.asfortranarray(self.SQ)
......@@ -1967,8 +1980,6 @@ class VFISVpackage:
# Truncate unfullfilled maps
self.truncate_maps_ifu()
def truncate_maps(self):
"""Truncate unfullfilled maps"""
......@@ -1986,7 +1997,9 @@ class VFISVpackage:
).astype("int"):
self.header["RMAPS"] = self.header["NMAPS"] - 1
else:
self.header['NSTEPS'] = max(self.istep[self.imap == self.header['NMAPS']]).astype("int")
self.header["NSTEPS"] = max(
self.istep[self.imap == self.header["NMAPS"]]
).astype("int")
def truncate_maps_ifu(self):
"""Truncate unfullfilled maps"""
......@@ -1997,7 +2010,6 @@ class VFISVpackage:
# introduce the RMAPS keyword for realized maps
self.header["RMAPS"] = self.header["NMAPS"]
def read_fits_slit(self, filename):
"""
Read the fits files into numpy array
......@@ -2036,7 +2048,6 @@ class VFISVpackage:
np.datetime64(hdul[0].header["DATE-OBS"]),
)
def read_fits_ifu(self, filename):
"""
Read the fits files into numpy array
......@@ -2057,10 +2068,18 @@ class VFISVpackage:
with fits.open(filename, memmap=True) as hdul:
data = hdul[0].data
si = data[0, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
sq = data[1, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
su = data[2, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
sv = data[3, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
si = (
data[0, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
)
sq = (
data[1, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
)
su = (
data[2, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
)
sv = (
data[3, self.pix_ini : self.pix_end + 1, :, :].copy().astype(np.float32)
)
ic = data[0, self.pix_cont, :, :].copy().astype(np.float32)
return (
si.T,
......@@ -2069,12 +2088,10 @@ class VFISVpackage:
sv.T,
ic.T,
hdul[0].header["IMAP"],
#hdul[0].header["ISTEP"],
# hdul[0].header["ISTEP"],
np.datetime64(hdul[0].header["DATE-OBS"]),
)
def check_fits(self, header):
"""
Check's if the fits file is properly formatted.
......@@ -2090,12 +2107,12 @@ class VFISVpackage:
# print_fail("Aborting.")
# sys.exit()
if not ([s for s in header['HISTORY'] if 'gris_translate_header' in s]):
if not ([s for s in header["HISTORY"] if "gris_translate_header" in s]):
print("Fits format problem: gris_translate_header not applied.")
print_fail("Aborting.")
sys.exit()
if header["OBS_TRGT"] == 'USAF':
if header["OBS_TRGT"] == "USAF":
print("Looks like a target scan")
print_fail(f"Aborting.")
sys.exit()
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment