From f05f6b789c0b2d82ce3f732756f47a552cc8fa01 Mon Sep 17 00:00:00 2001 From: Thibault Barnouin Date: Thu, 4 Sep 2025 11:17:14 +0200 Subject: [PATCH] update scripts and plots to new fits data shape --- package/lib/plots.py | 44 ++++++------- package/src/Combine.py | 109 +++++++++++++++++--------------- package/src/emission_center.py | 8 +-- package/src/overplot_IC5063.py | 10 +-- package/src/overplot_MRK463E.py | 46 +++++++------- 5 files changed, 114 insertions(+), 103 deletions(-) diff --git a/package/lib/plots.py b/package/lib/plots.py index 2cfcdee..62a15d4 100755 --- a/package/lib/plots.py +++ b/package/lib/plots.py @@ -181,9 +181,9 @@ def plot_Stokes(Stokes, savename=None, plots_folder=""): Defaults to current folder. """ # Get data - stkI = Stokes["I_stokes"].data.copy() - stkQ = Stokes["Q_stokes"].data.copy() - stkU = Stokes["U_stokes"].data.copy() + stkI = Stokes["stokes"].data[0].copy() + stkQ = Stokes["stokes"].data[1].copy() + stkU = Stokes["stokes"].data[2].copy() data_mask = Stokes["Data_mask"].data.astype(bool) for dataset in [stkI, stkQ, stkU]: @@ -288,10 +288,10 @@ def polarization_map( The figure and ax created for interactive contour maps. """ # Get data - stkI = Stokes["I_stokes"].data.copy() - stkQ = Stokes["Q_stokes"].data.copy() - stkU = Stokes["U_stokes"].data.copy() - stk_cov = Stokes["IQU_COV_MATRIX"].data.copy() + stkI = Stokes["STOKES"].data[0].copy() + stkQ = Stokes["STOKES"].data[1].copy() + stkU = Stokes["STOKES"].data[2].copy() + stk_cov = Stokes["STOKES_COV"].data.copy() pol = Stokes["Pol_deg_debiased"].data.copy() pol_err = Stokes["Pol_deg_err"].data.copy() pang = Stokes["Pol_ang"].data.copy() @@ -314,7 +314,7 @@ def polarization_map( for j in range(3): stk_cov[i][j][np.logical_not(data_mask)] = np.nan - wcs = WCS(Stokes[0]).deepcopy() + wcs = WCS(Stokes[0]).celestial.deepcopy() pivot_wav = Stokes[0].header["photplam"] convert_flux = Stokes[0].header["photflam"] @@ -968,10 +968,10 @@ class overplot_radio(align_maps): self.wcs_UV = self.map_wcs # Get Data obj = self.Stokes_UV[0].header["targname"] - stkI = self.Stokes_UV["I_STOKES"].data - stkQ = self.Stokes_UV["Q_STOKES"].data - stkU = self.Stokes_UV["U_STOKES"].data - stk_cov = self.Stokes_UV["IQU_COV_MATRIX"].data + stkI = self.Stokes_UV["STOKES"].data[0] + stkQ = self.Stokes_UV["STOKES"].data[1] + stkU = self.Stokes_UV["STOKES"].data[2] + stk_cov = self.Stokes_UV["STOKES_COV"].data pol = deepcopy(self.Stokes_UV["POL_DEG_DEBIASED"].data) pol_err = self.Stokes_UV["POL_DEG_ERR"].data pang = self.Stokes_UV["POL_ANG"].data @@ -1187,10 +1187,10 @@ class overplot_chandra(align_maps): self.wcs_UV = self.map_wcs # Get Data obj = self.Stokes_UV[0].header["targname"] - stkI = self.Stokes_UV["I_STOKES"].data - stkQ = self.Stokes_UV["Q_STOKES"].data - stkU = self.Stokes_UV["U_STOKES"].data - stk_cov = self.Stokes_UV["IQU_COV_MATRIX"].data + stkI = self.Stokes_UV["STOKES"].data[0] + stkQ = self.Stokes_UV["STOKES"].data[1] + stkU = self.Stokes_UV["STOKES"].data[2] + stk_cov = self.Stokes_UV["STOKES_COV"].data pol = deepcopy(self.Stokes_UV["POL_DEG_DEBIASED"].data) pol_err = self.Stokes_UV["POL_DEG_ERR"].data pang = self.Stokes_UV["POL_ANG"].data @@ -1398,10 +1398,10 @@ class overplot_pol(align_maps): self.wcs_UV = self.map_wcs # Get Data obj = self.Stokes_UV[0].header["targname"] - stkI = self.Stokes_UV["I_STOKES"].data - stkQ = self.Stokes_UV["Q_STOKES"].data - stkU = self.Stokes_UV["U_STOKES"].data - stk_cov = self.Stokes_UV["IQU_COV_MATRIX"].data + stkI = self.Stokes_UV["STOKES"].data[0] + stkQ = self.Stokes_UV["STOKES"].data[1] + stkU = self.Stokes_UV["STOKES"].data[2] + stk_cov = self.Stokes_UV["STOKES_COV"].data pol = deepcopy(self.Stokes_UV["POL_DEG_DEBIASED"].data) pol_err = self.Stokes_UV["POL_DEG_ERR"].data pang = self.Stokes_UV["POL_ANG"].data @@ -1702,8 +1702,8 @@ class align_pol(object): def single_plot(self, curr_map, wcs, v_lim=None, ax_lim=None, P_cut=3.0, SNRi_cut=3.0, savename=None, **kwargs): # Get data - stkI = curr_map["I_STOKES"].data - stk_cov = curr_map["IQU_COV_MATRIX"].data + stkI = curr_map["STOKES"].data[0] + stk_cov = curr_map["STOKES_COV"].data pol = deepcopy(curr_map["POL_DEG_DEBIASED"].data) pol_err = curr_map["POL_DEG_ERR"].data pang = curr_map["POL_ANG"].data diff --git a/package/src/Combine.py b/package/src/Combine.py index d5b19bb..77f2e98 100755 --- a/package/src/Combine.py +++ b/package/src/Combine.py @@ -1,5 +1,6 @@ #!/usr/bin/env python3 # -*- coding:utf-8 -*- +from copy import deepcopy from pathlib import Path from sys import path as syspath @@ -26,7 +27,7 @@ def same_reduction(infiles): except KeyError: pass test_IQU = True - for look in ["I_stokes", "Q_stokes", "U_stokes", "IQU_cov_matrix"]: + for look in ["STOKES", "STOKES_COV"]: test_IQU *= look in datatype params["IQU"].append(test_IQU) # test for orientation and pixel size @@ -88,73 +89,78 @@ def same_obs(infiles, data_folder): def combine_Stokes(infiles): """ - Combine I, Q, U from different observations of a same object. + Combine Stokes matrices from different observations of a same object. """ from astropy.io.fits import open as fits_open from lib.reduction import align_data, zeropad + from lib.utils import remove_stokes_axis_from_header from scipy.ndimage import shift as sc_shift - I_array, Q_array, U_array, IQU_cov_array, data_mask, headers = [], [], [], [], [], [] + Stokes_array, Stokes_cov_array, Stokes_cov_stat_array, data_mask, headers = [], [], [], [], [] shape = np.array([0, 0]) for file in infiles: with fits_open(file) as f: headers.append(f[0].header) - I_array.append(f["I_stokes"].data) - Q_array.append(f["Q_stokes"].data) - U_array.append(f["U_stokes"].data) - IQU_cov_array.append(f["IQU_cov_matrix"].data) + Stokes_array.append(f["stokes"].data) + Stokes_cov_array.append(f["stokes_cov"].data) + Stokes_cov_stat_array.append(f["stokes_cov_stat"].data) data_mask.append(f["data_mask"].data.astype(bool)) - shape[0] = np.max([shape[0], f["I_stokes"].data.shape[0]]) - shape[1] = np.max([shape[1], f["I_stokes"].data.shape[1]]) + shape[0] = np.max([shape[0], f["stokes"].data[0].shape[0]]) + shape[1] = np.max([shape[1], f["stokes"].data[0].shape[1]]) exposure_array = np.array([float(head["EXPTIME"]) for head in headers]) shape += np.array([5, 5]) data_mask = np.sum([zeropad(mask, shape) for mask in data_mask], axis=0).astype(bool) - I_array = np.array([zeropad(I, shape) for I in I_array]) - Q_array = np.array([zeropad(Q, shape) for Q in Q_array]) - U_array = np.array([zeropad(U, shape) for U in U_array]) - IQU_cov_array = np.array([[[zeropad(cov[i, j], shape) for j in range(3)] for i in range(3)] for cov in IQU_cov_array]) + Stokes_array = np.array([[zeropad(stk[i], shape) for i in range(4)] for stk in Stokes_array]) + Stokes_cov_array = np.array([[[zeropad(cov[i, j], shape) for j in range(4)] for i in range(4)] for cov in Stokes_cov_array]) + Stokes_cov_stat_array = np.array([[[zeropad(cov_stat[i, j], shape) for j in range(4)] for i in range(4)] for cov_stat in Stokes_cov_stat_array]) - sI_array = np.sqrt(IQU_cov_array[:, 0, 0]) - sQ_array = np.sqrt(IQU_cov_array[:, 1, 1]) - sU_array = np.sqrt(IQU_cov_array[:, 2, 2]) + I_array = deepcopy(Stokes_array[:, 0]) + sI_array = deepcopy(np.sqrt(Stokes_cov_array[:, 0, 0])) - _, _, _, _, shifts, errors = align_data(I_array, headers, error_array=sI_array, data_mask=data_mask, ref_center="center", return_shifts=True) + heads = [remove_stokes_axis_from_header(head) for head in headers] + _, _, _, _, shifts, errors = align_data( + I_array, heads, error_array=sI_array, background=sI_array[:, 0, 0], data_mask=data_mask, ref_center="center", return_shifts=True + ) data_mask_aligned = np.sum([sc_shift(data_mask, s, order=1, cval=0.0) for s in shifts], axis=0).astype(bool) - I_aligned, sI_aligned = ( - np.array([sc_shift(I, s, order=1, cval=0.0) for I, s in zip(I_array, shifts)]), - np.array([sc_shift(sI, s, order=1, cval=0.0) for sI, s in zip(sI_array, shifts)]), + Stokes_aligned = np.array([[sc_shift(stk[i], s, order=1, cval=0.0) for i in range(4)] for stk, s in zip(Stokes_array, shifts)]) + Stokes_cov_aligned = np.array( + [[[sc_shift(cov[i, j], s, order=1, cval=0.0) for j in range(4)] for i in range(4)] for cov, s in zip(Stokes_cov_array, shifts)] ) - Q_aligned, sQ_aligned = ( - np.array([sc_shift(Q, s, order=1, cval=0.0) for Q, s in zip(Q_array, shifts)]), - np.array([sc_shift(sQ, s, order=1, cval=0.0) for sQ, s in zip(sQ_array, shifts)]), + Stokes_cov_stat_aligned = np.array( + [[[sc_shift(cov_stat[i, j], s, order=1, cval=0.0) for j in range(4)] for i in range(4)] for cov_stat, s in zip(Stokes_cov_stat_array, shifts)] ) - U_aligned, sU_aligned = ( - np.array([sc_shift(U, s, order=1, cval=0.0) for U, s in zip(U_array, shifts)]), - np.array([sc_shift(sU, s, order=1, cval=0.0) for sU, s in zip(sU_array, shifts)]), - ) - IQU_cov_aligned = np.array([[[sc_shift(cov[i, j], s, order=1, cval=0.0) for j in range(3)] for i in range(3)] for cov, s in zip(IQU_cov_array, shifts)]) - I_combined = np.sum([exp * I for exp, I in zip(exposure_array, I_aligned)], axis=0) / exposure_array.sum() - Q_combined = np.sum([exp * Q for exp, Q in zip(exposure_array, Q_aligned)], axis=0) / exposure_array.sum() - U_combined = np.sum([exp * U for exp, U in zip(exposure_array, U_aligned)], axis=0) / exposure_array.sum() + Stokes_combined = np.zeros((4, shape[0], shape[1])) + for i in range(4): + Stokes_combined[i] = np.sum([exp * stk for exp, stk in zip(exposure_array, Stokes_aligned[:, i])], axis=0) / exposure_array.sum() - IQU_cov_combined = np.zeros((3, 3, shape[0], shape[1])) - for i in range(3): - IQU_cov_combined[i, i] = np.sum([exp**2 * cov for exp, cov in zip(exposure_array, IQU_cov_aligned[:, i, i])], axis=0) / exposure_array.sum() ** 2 - for j in [x for x in range(3) if x != i]: - IQU_cov_combined[i, j] = np.sqrt( - np.sum([exp**2 * cov**2 for exp, cov in zip(exposure_array, IQU_cov_aligned[:, i, j])], axis=0) / exposure_array.sum() ** 2 + Stokes_cov_combined = np.zeros((4, 4, shape[0], shape[1])) + Stokes_cov_stat_combined = np.zeros((4, 4, shape[0], shape[1])) + for i in range(4): + Stokes_cov_combined[i, i] = np.sum([exp**2 * cov for exp, cov in zip(exposure_array, Stokes_cov_aligned[:, i, i])], axis=0) / exposure_array.sum() ** 2 + Stokes_cov_stat_combined[i, i] = ( + np.sum([exp**2 * cov_stat for exp, cov_stat in zip(exposure_array, Stokes_cov_stat_aligned[:, i, i])], axis=0) / exposure_array.sum() ** 2 + ) + for j in [x for x in range(4) if x != i]: + Stokes_cov_combined[i, j] = np.sqrt( + np.sum([exp**2 * cov**2 for exp, cov in zip(exposure_array, Stokes_cov_aligned[:, i, j])], axis=0) / exposure_array.sum() ** 2 ) - IQU_cov_combined[j, i] = np.sqrt( - np.sum([exp**2 * cov**2 for exp, cov in zip(exposure_array, IQU_cov_aligned[:, j, i])], axis=0) / exposure_array.sum() ** 2 + Stokes_cov_combined[j, i] = np.sqrt( + np.sum([exp**2 * cov**2 for exp, cov in zip(exposure_array, Stokes_cov_aligned[:, j, i])], axis=0) / exposure_array.sum() ** 2 + ) + Stokes_cov_stat_combined[i, j] = np.sqrt( + np.sum([exp**2 * cov_stat**2 for exp, cov_stat in zip(exposure_array, Stokes_cov_stat_aligned[:, i, j])], axis=0) / exposure_array.sum() ** 2 + ) + Stokes_cov_stat_combined[j, i] = np.sqrt( + np.sum([exp**2 * cov_stat**2 for exp, cov_stat in zip(exposure_array, Stokes_cov_stat_aligned[:, j, i])], axis=0) / exposure_array.sum() ** 2 ) header_combined = headers[0] header_combined["EXPTIME"] = exposure_array.sum() - return I_combined, Q_combined, U_combined, IQU_cov_combined, data_mask_aligned, header_combined + return Stokes_combined, Stokes_cov_combined, Stokes_cov_stat_combined, data_mask_aligned, header_combined def main(infiles, target=None, output_dir="./data/"): @@ -190,21 +196,24 @@ def main(infiles, target=None, output_dir="./data/"): infiles = new_infiles - I_combined, Q_combined, U_combined, IQU_cov_combined, data_mask_combined, header_combined = combine_Stokes(infiles=infiles) - I_combined, Q_combined, U_combined, IQU_cov_combined, data_mask_combined, header_combined = rotate_Stokes( - I_stokes=I_combined, Q_stokes=Q_combined, U_stokes=U_combined, Stokes_cov=IQU_cov_combined, data_mask=data_mask_combined, header_stokes=header_combined + Stokes_combined, Stokes_cov_combined, Stokes_cov_stat_combined, data_mask_combined, header_combined = combine_Stokes(infiles=infiles) + Stokes_combined, Stokes_cov_combined, data_mask_combined, header_combined, Stokes_cov_stat_combined = rotate_Stokes( + Stokes=Stokes_combined, + Stokes_cov=Stokes_cov_combined, + Stokes_cov_stat=Stokes_cov_stat_combined, + data_mask=data_mask_combined, + header_stokes=header_combined, ) P, debiased_P, s_P, s_P_P, PA, s_PA, s_PA_P = compute_pol( - I_stokes=I_combined, Q_stokes=Q_combined, U_stokes=U_combined, Stokes_cov=IQU_cov_combined, header_stokes=header_combined + Stokes=Stokes_combined, Stokes_cov=Stokes_cov_combined, Stokes_cov_stat=Stokes_cov_stat_combined, header_stokes=header_combined ) filename = header_combined["FILENAME"] figname = "_".join([target, filename[filename.find("FOC_") :], "combined"]) - Stokes_combined = save_Stokes( - I_stokes=I_combined, - Q_stokes=Q_combined, - U_stokes=U_combined, - Stokes_cov=IQU_cov_combined, + Stokes_c = save_Stokes( + Stokes=Stokes_combined, + Stokes_cov=Stokes_cov_combined, + Stokes_cov_stat=Stokes_cov_stat_combined, P=P, debiased_P=debiased_P, s_P=s_P, @@ -219,7 +228,7 @@ def main(infiles, target=None, output_dir="./data/"): return_hdul=True, ) - pol_map(Stokes_combined, **kwargs) + pol_map(Stokes_c, **kwargs) return "/".join([data_folder, figname + ".fits"]) diff --git a/package/src/emission_center.py b/package/src/emission_center.py index 4c421c5..6443199 100755 --- a/package/src/emission_center.py +++ b/package/src/emission_center.py @@ -20,13 +20,13 @@ def main(infile, P_cut=0.99, target=None, display="pf", output_dir=None): output = [] Stokes = fits_open(infile) - stkI = Stokes["I_STOKES"].data - s_I = np.sqrt(Stokes["IQU_COV_MATRIX"].data[0, 0]) + stkI = Stokes["STOKES"].data[0] + s_I = np.sqrt(Stokes["STOKES_COV"].data[0, 0]) SNRi = np.zeros(stkI.shape) SNRi[s_I > 0.0] = stkI[s_I > 0.0] / s_I[s_I > 0.0] QN, UN, QN_ERR, UN_ERR = np.full((4, stkI.shape[0], stkI.shape[1]), np.nan) for sflux, nflux in zip( - [Stokes["Q_STOKES"].data, Stokes["U_STOKES"].data, np.sqrt(Stokes["IQU_COV_MATRIX"].data[1, 1]), np.sqrt(Stokes["IQU_COV_MATRIX"].data[2, 2])], + [Stokes["STOKES"].data[1], Stokes["STOKES"].data[2], np.sqrt(Stokes["STOKES_COV"].data[1, 1]), np.sqrt(Stokes["STOKES_COV"].data[2, 2])], [QN, UN, QN_ERR, UN_ERR], ): nflux[stkI > 0.0] = sflux[stkI > 0.0] / stkI[stkI > 0.0] @@ -41,7 +41,7 @@ def main(infile, P_cut=0.99, target=None, display="pf", output_dir=None): Stokescentconf, Stokescenter = CenterConf(Stokesconf > P_cut, Stokes["POL_ANG"].data, Stokes["POL_ANG_ERR"].data) else: Stokescentconf, Stokescenter = CenterConf(Stokessnr > P_cut, Stokes["POL_ANG"].data, Stokes["POL_ANG_ERR"].data) - Stokespos = WCS(Stokes[0].header).pixel_to_world(*Stokescenter) + Stokespos = WCS(Stokes[0].header).celestial.pixel_to_world(*Stokescenter) if target is None: target = Stokes[0].header["TARGNAME"] diff --git a/package/src/overplot_IC5063.py b/package/src/overplot_IC5063.py index bea2fbe..2862d3d 100755 --- a/package/src/overplot_IC5063.py +++ b/package/src/overplot_IC5063.py @@ -10,7 +10,7 @@ from astropy.io import fits from lib.plots import overplot_pol, overplot_radio from matplotlib.colors import LogNorm -Stokes_UV = fits.open("./data/IC5063/5918/IC5063_FOC_b0.10arcsec_c0.20arcsec.fits") +Stokes_UV = fits.open("./data/IC5063/5918/IC5063_5918_F502M_FOC_b0.10arcsec_c0.15arcsec.fits") # Stokes_18GHz = fits.open("./data/IC5063/radio/IC5063_18GHz.fits") # Stokes_24GHz = fits.open("./data/IC5063/radio/IC5063_24GHz.fits") # Stokes_103GHz = fits.open("./data/IC5063/radio/IC5063_103GHz.fits") @@ -47,10 +47,12 @@ Stokes_IR = fits.open("./data/IC5063/IR/u2e65g01t_c0f_rot.fits") G = overplot_pol(Stokes_UV, Stokes_IR, cmap="inferno") G.plot( - P_cut=0.99, - SNRi_cut=1.0, + P_cut=3, + SNRi_cut=10, savename="./plots/IC5063/IR_overplot.pdf", - scale_vec=None, + scale_vec=5, norm=LogNorm(Stokes_IR[0].data.max() * Stokes_IR[0].header["photflam"] / 1e3, Stokes_IR[0].data.max() * Stokes_IR[0].header["photflam"]), cmap="inferno", + disptype="pf", + levels="Default", ) diff --git a/package/src/overplot_MRK463E.py b/package/src/overplot_MRK463E.py index 505e39d..6f25f4c 100755 --- a/package/src/overplot_MRK463E.py +++ b/package/src/overplot_MRK463E.py @@ -11,9 +11,9 @@ from lib.plots import overplot_chandra, overplot_pol from matplotlib.colors import LogNorm Stokes_UV = fits.open("./data/MRK463E/5960/MRK463E_FOC_b0.05arcsec_c0.10arcsec.fits") -# Stokes_IR = fits.open("./data/MRK463E/WFPC2/IR_rot_crop.fits") +Stokes_IR = fits.open("./data/MRK463E/WFPC2/IR_rot_crop.fits") # Stokes_Xr = fits.open("./data/MRK463E/Chandra/X_ray_crop.fits") -Radio = fits.open("./data/MRK463E/EMERLIN/Voorwerpjes_1356+1822_1356+1822_uniform-image.fits") +# Radio = fits.open("./data/MRK463E/EMERLIN/Voorwerpjes_1356+1822_1356+1822_uniform-image.fits") # levels = np.geomspace(1.0, 99.0, 7) @@ -21,26 +21,26 @@ Radio = fits.open("./data/MRK463E/EMERLIN/Voorwerpjes_1356+1822_1356+1822_unifor # A.plot(levels=levels, P_cut=0.99, SNRi_cut=1.0, scale_vec=5, zoom=1, savename="./plots/MRK463E/Chandra_overplot.pdf") # A.write_to(path1="./data/MRK463E/FOC_data_Chandra.fits", path2="./data/MRK463E/Chandra_data.fits", suffix="aligned") -# levels = np.array([0.8, 2, 5, 10, 20, 50]) / 100.0 * Stokes_UV[0].header["photflam"] -# B = overplot_pol(Stokes_UV, Stokes_IR, norm=LogNorm()) -# B.plot(levels=levels, P_cut=0.99, SNRi_cut=1.0, scale_vec=5, norm=LogNorm(8.5e-18, 2.5e-15), savename="./plots/MRK463E/IR_overplot.pdf") -# B.write_to(path1="./data/MRK463E/FOC_data_WFPC.fits", path2="./data/MRK463E/WFPC_data.fits", suffix="aligned") +levels = np.array([0.8, 2, 5, 10, 20, 50]) / 100.0 * Stokes_UV[0].header["photflam"] +B = overplot_pol(Stokes_UV, Stokes_IR, norm=LogNorm()) +B.plot(levels=levels, P_cut=0.99, SNRi_cut=1.0, scale_vec=5, norm=LogNorm(8.5e-18, 2.5e-15), savename="./plots/MRK463E/IR_overplot.pdf") +B.write_to(path1="./data/MRK463E/FOC_data_WFPC.fits", path2="./data/MRK463E/WFPC_data.fits", suffix="aligned") # levels = np.array([0.8, 2, 5, 10, 20, 50]) / 100.0 * Stokes_UV[0].header["photflam"] -levels = np.array([5, 10, 20, 50]) -C = overplot_pol(Stokes_UV, Radio, norm=LogNorm()) -C.other_im.set(norm=LogNorm(1e-4, 2e-2)) -C.plot( - levels=levels, - P_cut=0.99, - SNRi_cut=1.0, - step_vec=0, - scale_vec=3, - norm=LogNorm(1e-4, 2e-2), - cmap="inferno_r", - width=0.5, - linewidth=0.5, - disptype="snri", - savename="./plots/MRK463E/EMERLIN_snri_overplot.pdf", -) -C.write_to(path1="./data/MRK463E/FOC_data_EMERLIN.fits", path2="./data/MRK463E/EMERLIN_data.fits", suffix="aligned") +# levels = np.array([5, 10, 20, 50]) +# C = overplot_pol(Stokes_UV, Radio, norm=LogNorm()) +# C.other_im.set(norm=LogNorm(1e-4, 2e-2)) +# C.plot( +# levels=levels, +# P_cut=0.99, +# SNRi_cut=1.0, +# step_vec=0, +# scale_vec=3, +# norm=LogNorm(1e-4, 2e-2), +# cmap="inferno_r", +# width=0.5, +# linewidth=0.5, +# disptype="snri", +# savename="./plots/MRK463E/EMERLIN_snri_overplot.pdf", +# ) +# C.write_to(path1="./data/MRK463E/FOC_data_EMERLIN.fits", path2="./data/MRK463E/EMERLIN_data.fits", suffix="aligned")