From 97479f90fb3b44e76c0704ddef1730b4c7943323 Mon Sep 17 00:00:00 2001 From: Thibault Barnouin Date: Tue, 4 Jun 2024 11:59:33 +0200 Subject: [PATCH] polarisation -> polarization everywhere --- package/FOC_reduction.py | 24 +++++------ package/lib/fits.py | 8 ++-- package/lib/plots.py | 82 ++++++++++++++++++------------------- package/lib/reduction.py | 54 ++++++++++++------------ package/overplot_MRK463E.py | 18 +++----- 5 files changed, 90 insertions(+), 96 deletions(-) diff --git a/package/FOC_reduction.py b/package/FOC_reduction.py index 565ff7a..ecf94c2 100755 --- a/package/FOC_reduction.py +++ b/package/FOC_reduction.py @@ -165,7 +165,7 @@ def main(target=None, proposal_id=None, infiles=None, output_dir="./data", crop= I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, SNRi_cut=None) I_bkg, Q_bkg, U_bkg, S_cov_bkg, _, _ = proj_red.rotate_Stokes(I_bkg, Q_bkg, U_bkg, S_cov_bkg, np.array(True).reshape(1, 1), headers, SNRi_cut=None) - # Compute polarimetric parameters (polarisation degree and angle). + # Compute polarimetric parameters (polarization degree and angle). P, debiased_P, s_P, s_P_P, PA, s_PA, s_PA_P = proj_red.compute_pol(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers) P_bkg, debiased_P_bkg, s_P_bkg, s_P_P_bkg, PA_bkg, s_PA_bkg, s_PA_P_bkg = proj_red.compute_pol(I_bkg, Q_bkg, U_bkg, S_cov_bkg, headers) @@ -194,28 +194,28 @@ def main(target=None, proposal_id=None, infiles=None, output_dir="./data", crop= I_bkg[0, 0]*headers[0]['photflam'], np.sqrt(S_cov_bkg[0, 0][0, 0])*headers[0]['photflam'], 2, out=int))) print("P_bkg = {0:.1f} ± {1:.1f} %".format(debiased_P_bkg[0, 0]*100., np.ceil(s_P_bkg[0, 0]*1000.)/10.)) print("PA_bkg = {0:.1f} ± {1:.1f} °".format(princ_angle(PA_bkg[0, 0]), princ_angle(np.ceil(s_PA_bkg[0, 0]*10.)/10.))) - # Plot polarisation map (Background is either total Flux, Polarization degree or Polarization degree error). + # Plot polarization map (Background is either total Flux, Polarization degree or Polarization degree error). if px_scale.lower() not in ['full', 'integrate'] and not interactive: - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname]), plots_folder=plots_folder) - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "I"]), plots_folder=plots_folder, display='Intensity') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "P_flux"]), plots_folder=plots_folder, display='Pol_Flux') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "P"]), plots_folder=plots_folder, display='Pol_deg') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "PA"]), plots_folder=plots_folder, display='Pol_ang') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "I_err"]), plots_folder=plots_folder, display='I_err') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "P_err"]), plots_folder=plots_folder, display='Pol_deg_err') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "SNRi"]), plots_folder=plots_folder, display='SNRi') - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim, step_vec=step_vec, vec_scale=vec_scale, savename="_".join([figname, "SNRp"]), plots_folder=plots_folder, display='SNRp') elif not interactive: - proj_plots.polarisation_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, + proj_plots.polarization_map(deepcopy(Stokes_test), data_mask, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, savename=figname, plots_folder=plots_folder, display='integrate') elif px_scale.lower() not in ['full', 'integrate']: proj_plots.pol_map(Stokes_test, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, flux_lim=flux_lim) diff --git a/package/lib/fits.py b/package/lib/fits.py index e966f43..811afa5 100755 --- a/package/lib/fits.py +++ b/package/lib/fits.py @@ -144,10 +144,10 @@ def save_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, P, debiased_P, s_P, header['targname'] = (ref_header['targname'], 'Target name') header['orientat'] = (ref_header['orientat'], 'Angle between North and the y-axis of the image') header['filename'] = (filename, 'Original filename') - header['P_int'] = (ref_header['P_int'], 'Integrated polarisation degree') - header['P_int_err'] = (ref_header['P_int_err'], 'Integrated polarisation degree error') - header['PA_int'] = (ref_header['PA_int'], 'Integrated polarisation angle') - header['PA_int_err'] = (ref_header['PA_int_err'], 'Integrated polarisation angle error') + header['P_int'] = (ref_header['P_int'], 'Integrated polarization degree') + header['P_int_err'] = (ref_header['P_int_err'], 'Integrated polarization degree error') + header['PA_int'] = (ref_header['PA_int'], 'Integrated polarization angle') + header['PA_int_err'] = (ref_header['PA_int_err'], 'Integrated polarization angle error') # Crop Data to mask if data_mask.shape != (1, 1): diff --git a/package/lib/plots.py b/package/lib/plots.py index 0ea1354..83b26d4 100755 --- a/package/lib/plots.py +++ b/package/lib/plots.py @@ -9,8 +9,8 @@ prototypes : - plot_Stokes(Stokes, savename, plots_folder) Plot the I/Q/U maps from the Stokes HDUList. - - polarisation_map(Stokes, data_mask, rectangle, SNRp_cut, SNRi_cut, step_vec, savename, plots_folder, display) -> fig, ax - Plots polarisation map of polarimetric parameters saved in an HDUList. + - polarization_map(Stokes, data_mask, rectangle, SNRp_cut, SNRi_cut, step_vec, savename, plots_folder, display) -> fig, ax + Plots polarization map of polarimetric parameters saved in an HDUList. class align_maps(map, other_map, **kwargs) Class to interactively align maps with different WCS. @@ -22,13 +22,13 @@ prototypes : Class inherited from align_maps to overplot chandra data as contours. class overplot_pol(align_maps) - Class inherited from align_maps to overplot UV polarisation vectors on other maps. + Class inherited from align_maps to overplot UV polarization vectors on other maps. class crop_map(hdul, fig, ax) Class to interactively crop a region of interest of a HDUList. class crop_Stokes(crop_map) - Class inherited from crop_map to work on polarisation maps. + Class inherited from crop_map to work on polarization maps. class image_lasso_selector(img, fig, ax) Class to interactively select part of a map to work on. @@ -37,7 +37,7 @@ prototypes : Class to interactively simulate aperture integration. class pol_map(Stokes, SNRp_cut, SNRi_cut, selection) - Class to interactively study polarisation maps making use of the cropping and selecting tools. + Class to interactively study polarization maps making use of the cropping and selecting tools. """ from copy import deepcopy @@ -204,10 +204,10 @@ def plot_Stokes(Stokes, savename=None, plots_folder=""): return 0 -def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_cut=3., +def polarization_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_cut=3., flux_lim=None, step_vec=1, vec_scale=2., savename=None, plots_folder="", display="default"): """ - Plots polarisation map from Stokes HDUList. + Plots polarization map from Stokes HDUList. ---------- Inputs: Stokes : astropy.io.fits.hdu.hdulist.HDUList @@ -230,12 +230,12 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c Defaults to None, limits are computed on the background value and the maximum value in the cut. step_vec : int, optional - Number of steps between each displayed polarisation vector. + Number of steps between each displayed polarization vector. If step_vec = 2, every other vector will be displayed. Defaults to 1 vec_scale : float, optional - Pixel length of displayed 100% polarisation vector. - If vec_scale = 2, a vector of 50% polarisation will be 1 pixel wide. + Pixel length of displayed 100% polarization vector. + If vec_scale = 2, a vector of 50% polarization will be 1 pixel wide. Defaults to 2. savename : str, optional Name of the figure the map should be saved to. If None, the map won't @@ -246,8 +246,8 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c be saved. Not used if savename is None. Defaults to current folder. display : str, optional - Choose the map to display between intensity (default), polarisation - degree ('p', 'pol', 'pol_deg') or polarisation degree error ('s_p', + Choose the map to display between intensity (default), polarization + degree ('p', 'pol', 'pol_deg') or polarization degree error ('s_p', 'pol_err', 'pol_deg_err'). Defaults to None (intensity). ---------- @@ -295,12 +295,12 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c poldata[np.logical_not(mask)] = np.nan pangdata[np.logical_not(mask)] = np.nan - # Look for pixel of max polarisation + # Look for pixel of max polarization if np.isfinite(pol).any(): p_max = np.max(pol[np.isfinite(pol)]) x_max, y_max = np.unravel_index(np.argmax(pol == p_max), pol.shape) else: - print("No pixel with polarisation information above requested SNR.") + print("No pixel with polarization information above requested SNR.") # Plot the map plt.rcParams.update({'font.size': 10}) @@ -325,7 +325,7 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c print("Total flux contour levels : ", levelsI) ax.contour(stkI*convert_flux, levels=levelsI, colors='grey', linewidths=0.5) elif display.lower() in ['pol_flux']: - # Display polarisation flux + # Display polarization flux display = 'pf' if flux_lim is None: if mask.sum() > 0.: @@ -340,19 +340,19 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c print("Polarized flux contour levels : ", levelsPf) ax.contour(stkI*convert_flux*pol, levels=levelsPf, colors='grey', linewidths=0.5) elif display.lower() in ['p', 'pol', 'pol_deg']: - # Display polarisation degree map + # Display polarization degree map display = 'p' vmin, vmax = 0., 100. im = ax.imshow(pol*100., vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno', alpha=1.) fig.colorbar(im, ax=ax, aspect=50, shrink=0.75, pad=0.025, label=r"$P$ [%]") elif display.lower() in ['pa', 'pang', 'pol_ang']: - # Display polarisation degree map + # Display polarization degree map display = 'pa' vmin, vmax = 0., 180. im = ax.imshow(princ_angle(pang), vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno', alpha=1.) fig.colorbar(im, ax=ax, aspect=50, shrink=0.75, pad=0.025, label=r"$\theta_P$ [°]") elif display.lower() in ['s_p', 'pol_err', 'pol_deg_err']: - # Display polarisation degree error map + # Display polarization degree error map display = 's_p' if (SNRp > SNRp_cut).any(): vmin, vmax = 0., np.max([pol_err[SNRp > SNRp_cut].max(), 1.])*100. @@ -384,7 +384,7 @@ def polarisation_map(Stokes, data_mask=None, rectangle=None, SNRp_cut=3., SNRi_c im = ax.imshow(SNRi, aspect='equal', cmap='inferno', alpha=1.) fig.colorbar(im, ax=ax, aspect=50, shrink=0.75, pad=0.025, label=r"$I_{Stokes}/\sigma_{I}$") elif display.lower() in ['snrp']: - # Display polarisation degree signal-to-noise map + # Display polarization degree signal-to-noise map display = 'snrp' vmin, vmax = 0., np.max(SNRp[np.isfinite(SNRp)]) if vmax*0.99 > SNRp_cut: @@ -710,7 +710,7 @@ class overplot_radio(align_maps): self.fig_overplot, self.ax_overplot = plt.subplots(figsize=(10, 10), subplot_kw=dict(projection=self.wcs_UV)) self.fig_overplot.subplots_adjust(hspace=0, wspace=0, bottom=0.1, left=0.1, top=0.8, right=1) - # Display UV intensity map with polarisation vectors + # Display UV intensity map with polarization vectors vmin, vmax = stkI[np.isfinite(stkI)].max()/1e3*self.map_convert, stkI[np.isfinite(stkI)].max()*self.map_convert for key, value in [["cmap", [["cmap", "inferno"]]], ["norm", [["norm", LogNorm(vmin, vmax)]]]]: try: @@ -728,7 +728,7 @@ class overplot_radio(align_maps): self.cbar = self.fig_overplot.colorbar(self.im, ax=self.ax_overplot, aspect=50, shrink=0.75, pad=0.025, label=r"$F_{{\lambda}}$ [{0:s}]".format(self.map_unit)) - # Display full size polarisation vectors + # Display full size polarization vectors if vec_scale is None: self.vec_scale = 2. pol[np.isfinite(pol)] = 1./2. @@ -738,7 +738,7 @@ class overplot_radio(align_maps): self.X, self.Y = np.meshgrid(np.arange(stkI.shape[1]), np.arange(stkI.shape[0])) self.U, self.V = pol*np.cos(np.pi/2.+pang*np.pi/180.), pol*np.sin(np.pi/2.+pang*np.pi/180.) self.Q = self.ax_overplot.quiver(self.X[::step_vec, ::step_vec], self.Y[::step_vec, ::step_vec], self.U[::step_vec, ::step_vec], self.V[::step_vec, ::step_vec], units='xy', angles='uv', scale=1./self.vec_scale, - scale_units='xy', pivot='mid', headwidth=0., headlength=0., headaxislength=0., width=0.5, linewidth=0.8, color='white', edgecolor='black', label="{0:s} polarisation map".format(self.map_observer)) + scale_units='xy', pivot='mid', headwidth=0., headlength=0., headaxislength=0., width=0.5, linewidth=0.8, color='white', edgecolor='black', label="{0:s} polarization map".format(self.map_observer)) self.ax_overplot.autoscale(False) # Display other map as contours @@ -750,7 +750,7 @@ class overplot_radio(align_maps): self.ax_overplot.set_xlabel(label="Right Ascension (J2000)") self.ax_overplot.set_ylabel(label="Declination (J2000)", labelpad=-1) - self.fig_overplot.suptitle("{0:s} polarisation map of {1:s} overplotted with {2:s} {3:.2f}GHz map in {4:s}.".format( + self.fig_overplot.suptitle("{0:s} polarization map of {1:s} overplotted with {2:s} {3:.2f}GHz map in {4:s}.".format( self.map_observer, obj, self.other_observer, other_freq*1e-9, self.other_unit), wrap=True) # Display pixel scale and North direction @@ -770,7 +770,7 @@ class overplot_radio(align_maps): self.cr_other, = self.ax_overplot.plot(*(self.other_wcs.celestial.wcs.crpix-(1., 1.)), 'g+', transform=self.ax_overplot.get_transform(self.other_wcs)) handles, labels = self.ax_overplot.get_legend_handles_labels() - handles[np.argmax([li == "{0:s} polarisation map".format(self.map_observer) for li in labels]) + handles[np.argmax([li == "{0:s} polarization map".format(self.map_observer) for li in labels]) ] = FancyArrowPatch((0, 0), (0, 1), arrowstyle='-', fc='w', ec='k', lw=2) labels.append("{0:s} contour".format(self.other_observer)) handles.append(Rectangle((0, 0), 1, 1, fill=False, lw=2, ec=other_cont.collections[0].get_edgecolor()[0])) @@ -829,7 +829,7 @@ class overplot_chandra(align_maps): self.fig_overplot, self.ax_overplot = plt.subplots(figsize=(11, 10), subplot_kw=dict(projection=self.wcs_UV)) self.fig_overplot.subplots_adjust(hspace=0, wspace=0, bottom=0.1, left=0.1, top=0.8, right=1) - # Display UV intensity map with polarisation vectors + # Display UV intensity map with polarization vectors vmin, vmax = stkI[np.isfinite(stkI)].max()/1e3*self.map_convert, stkI[np.isfinite(stkI)].max()*self.map_convert for key, value in [["cmap", [["cmap", "inferno"]]], ["norm", [["norm", LogNorm(vmin, vmax)]]]]: try: @@ -847,7 +847,7 @@ class overplot_chandra(align_maps): self.cbar = self.fig_overplot.colorbar(self.im, ax=self.ax_overplot, aspect=50, shrink=0.75, pad=0.025, label=r"$F_{{\lambda}}$ [{0:s}]".format(self.map_unit)) - # Display full size polarisation vectors + # Display full size polarization vectors if vec_scale is None: self.vec_scale = 2. pol[np.isfinite(pol)] = 1./2. @@ -857,7 +857,7 @@ class overplot_chandra(align_maps): self.X, self.Y = np.meshgrid(np.arange(stkI.shape[1]), np.arange(stkI.shape[0])) self.U, self.V = pol*np.cos(np.pi/2.+pang*np.pi/180.), pol*np.sin(np.pi/2.+pang*np.pi/180.) self.Q = self.ax_overplot.quiver(self.X[::step_vec, ::step_vec], self.Y[::step_vec, ::step_vec], self.U[::step_vec, ::step_vec], self.V[::step_vec, ::step_vec], units='xy', angles='uv', scale=1./self.vec_scale, - scale_units='xy', pivot='mid', headwidth=0., headlength=0., headaxislength=0., width=0.5, linewidth=0.8, color='white', edgecolor='black', label="{0:s} polarisation map".format(self.map_observer)) + scale_units='xy', pivot='mid', headwidth=0., headlength=0., headaxislength=0., width=0.5, linewidth=0.8, color='white', edgecolor='black', label="{0:s} polarization map".format(self.map_observer)) self.ax_overplot.autoscale(False) # Display other map as contours @@ -870,7 +870,7 @@ class overplot_chandra(align_maps): self.ax_overplot.set_xlabel(label="Right Ascension (J2000)") self.ax_overplot.set_ylabel(label="Declination (J2000)", labelpad=-1) - self.fig_overplot.suptitle("{0:s} polarisation map of {1:s} overplotted\nwith {2:s} contour in counts.".format( + self.fig_overplot.suptitle("{0:s} polarization map of {1:s} overplotted\nwith {2:s} contour in counts.".format( self.map_observer, obj, self.other_observer), wrap=True) # Display pixel scale and North direction @@ -889,7 +889,7 @@ class overplot_chandra(align_maps): self.cr_map, = self.ax_overplot.plot(*(self.map_wcs.celestial.wcs.crpix-(1., 1.)), 'r+') self.cr_other, = self.ax_overplot.plot(*(other_wcs.celestial.wcs.crpix-(1., 1.)), 'g+', transform=self.ax_overplot.get_transform(other_wcs)) handles, labels = self.ax_overplot.get_legend_handles_labels() - handles[np.argmax([li == "{0:s} polarisation map".format(self.map_observer) for li in labels]) + handles[np.argmax([li == "{0:s} polarization map".format(self.map_observer) for li in labels]) ] = FancyArrowPatch((0, 0), (0, 1), arrowstyle='-', fc='w', ec='k', lw=2) labels.append("{0:s} contour in counts".format(self.other_observer)) handles.append(Rectangle((0, 0), 1, 1, fill=False, lw=2, ec=other_cont.collections[0].get_edgecolor()[0])) @@ -944,7 +944,7 @@ class overplot_pol(align_maps): self.ax_overplot.set_xlabel(label="Right Ascension (J2000)") self.ax_overplot.set_ylabel(label="Declination (J2000)", labelpad=-1) - self.fig_overplot.suptitle("{0:s} observation from {1:s} overplotted with polarisation vectors and Stokes I contours from {2:s}".format( + self.fig_overplot.suptitle("{0:s} observation from {1:s} overplotted with polarization vectors and Stokes I contours from {2:s}".format( obj, self.other_observer, self.map_observer), wrap=True) # Display "other" intensity map @@ -965,7 +965,7 @@ class overplot_pol(align_maps): self.cbar = self.fig_overplot.colorbar(self.im, ax=self.ax_overplot, aspect=80, shrink=0.75, pad=0.025, label=r"$F_{{\lambda}}$ [{0:s}]".format(self.other_unit)) - # Display full size polarisation vectors + # Display full size polarization vectors if vec_scale is None: self.vec_scale = 2. pol[np.isfinite(pol)] = 1./2. @@ -976,7 +976,7 @@ class overplot_pol(align_maps): self.X, self.Y = np.meshgrid(np.arange(stkI.shape[1]), np.arange(stkI.shape[0])) self.U, self.V = pol*np.cos(np.pi/2.+pang*np.pi/180.), pol*np.sin(np.pi/2.+pang*np.pi/180.) self.Q = self.ax_overplot.quiver(self.X[::step_vec, ::step_vec], self.Y[::step_vec, ::step_vec], self.U[::step_vec, ::step_vec], self.V[::step_vec, ::step_vec], units='xy', angles='uv', scale=px_scale/self.vec_scale, scale_units='xy', pivot='mid', - headwidth=0., headlength=0., headaxislength=0., width=2.0, linewidth=1.0, color='white', edgecolor='black', transform=self.ax_overplot.get_transform(self.wcs_UV), label="{0:s} polarisation map".format(self.map_observer)) + headwidth=0., headlength=0., headaxislength=0., width=2.0, linewidth=1.0, color='white', edgecolor='black', transform=self.ax_overplot.get_transform(self.wcs_UV), label="{0:s} polarization map".format(self.map_observer)) # Display Stokes I as contours if levels is None: @@ -1009,7 +1009,7 @@ class overplot_pol(align_maps): self.legend_title = r"{0:s} image".format(self.other_observer) handles, labels = self.ax_overplot.get_legend_handles_labels() - handles[np.argmax([li == "{0:s} polarisation map".format(self.map_observer) for li in labels]) + handles[np.argmax([li == "{0:s} polarization map".format(self.map_observer) for li in labels]) ] = FancyArrowPatch((0, 0), (0, 1), arrowstyle='-', fc='w', ec='k', lw=2) labels.append("{0:s} Stokes I contour".format(self.map_observer)) handles.append(Rectangle((0, 0), 1, 1, fill=False, ec=cont_stkI.collections[0].get_edgecolor()[0])) @@ -1358,7 +1358,7 @@ class crop_map(object): class crop_Stokes(crop_map): """ - Class to interactively crop a polarisation map to desired Region of Interest. + Class to interactively crop a polarization map to desired Region of Interest. Inherit from crop_map. """ @@ -1440,10 +1440,10 @@ class crop_Stokes(crop_map): 2 + Q_diluted**2*U_diluted_err**2 - 2.*Q_diluted*U_diluted*QU_diluted_err) for dataset in self.hdul_crop: - dataset.header['P_int'] = (P_diluted, 'Integrated polarisation degree') - dataset.header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarisation degree error') - dataset.header['PA_int'] = (PA_diluted, 'Integrated polarisation angle') - dataset.header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarisation angle error') + dataset.header['P_int'] = (P_diluted, 'Integrated polarization degree') + dataset.header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarization degree error') + dataset.header['PA_int'] = (PA_diluted, 'Integrated polarization angle') + dataset.header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarization angle error') self.fig.canvas.draw_idle() @property @@ -1738,7 +1738,7 @@ class aperture(object): class pol_map(object): """ - Class to interactively study polarisation maps. + Class to interactively study polarization maps. """ def __init__(self, Stokes, SNRp_cut=3., SNRi_cut=3., flux_lim=None, selection=None): @@ -1769,7 +1769,7 @@ class pol_map(object): # Display selected data (Default to total flux) self.display() - # Display polarisation vectors in SNR_cut + # Display polarization vectors in SNR_cut self.pol_vector() # Display integrated values in ROI self.pol_int() @@ -2474,7 +2474,7 @@ if __name__ == "__main__": parser = argparse.ArgumentParser(description='Interactively plot the pipeline products') parser.add_argument('-f', '--file', metavar='path', required=False, help='the full or relative path to the data product', type=str, default=None) - parser.add_argument('-p', '--snrp', metavar='snrp_cut', required=False, help='the cut in signal-to-noise for the polarisation degree', type=float, default=3.) + parser.add_argument('-p', '--snrp', metavar='snrp_cut', required=False, help='the cut in signal-to-noise for the polarization degree', type=float, default=3.) parser.add_argument('-i', '--snri', metavar='snri_cut', required=False, help='the cut in signal-to-noise for the intensity', type=float, default=3.) parser.add_argument('-l', '--lim', metavar='flux_lim', nargs=2, required=False, help='limits for the intensity map', default=None) args = parser.parse_args() diff --git a/package/lib/reduction.py b/package/lib/reduction.py index 2553f13..36227f1 100755 --- a/package/lib/reduction.py +++ b/package/lib/reduction.py @@ -30,7 +30,7 @@ prototypes : Compute Stokes parameters I, Q and U and their respective correlated errors from data_array. - compute_pol(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers) -> P, debiased_P, s_P, s_P_P, PA, s_PA, s_PA_P - Compute polarisation degree (in %) and angle (in degree) and their respective errors. + Compute polarization degree (in %) and angle (in degree) and their respective errors. - rotate_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, ang, SNRi_cut) -> I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers Rotate I, Q, U given an angle in degrees using scipy functions. @@ -466,7 +466,7 @@ def get_error(data_array, headers, error_array=None, data_mask=None, sub_type=No mask = zeropad(mask_c, data[0].shape).astype(bool) background = np.zeros((data.shape[0])) - # wavelength dependence of the polariser filters + # wavelength dependence of the polarizer filters # estimated to less than 1% err_wav = data*0.01 # difference in PSFs through each polarizers @@ -963,7 +963,7 @@ def polarizer_avg(data_array, error_array, data_mask, headers, FWHM=None, scale= pol120, err120 = smooth_data(pol120_array, err120_array, data_mask, headers120, FWHM=FWHM, scale=scale, smoothing=smoothing) else: - # Sum on each polarisation filter. + # Sum on each polarization filter. pol0_t = np.sum([header['exptime'] for header in headers0]) pol60_t = np.sum([header['exptime'] for header in headers60]) pol120_t = np.sum([header['exptime'] for header in headers120]) @@ -1070,10 +1070,10 @@ def compute_Stokes(data_array, error_array, data_mask, headers, FWHM=None, scale total intensity Q_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - vertical/horizontal linear polarisation intensity + vertical/horizontal linear polarization intensity U_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - +45/-45deg linear polarisation intensity + +45/-45deg linear polarization intensity Stokes_cov : numpy.ndarray Covariance matrix of the Stokes parameters I, Q, U. """ @@ -1236,17 +1236,17 @@ def compute_Stokes(data_array, error_array, data_mask, headers, FWHM=None, scale 2 + Q_diluted**2*U_diluted_err**2 - 2.*Q_diluted*U_diluted*QU_diluted_err) for header in headers: - header['P_int'] = (P_diluted, 'Integrated polarisation degree') - header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarisation degree error') - header['PA_int'] = (PA_diluted, 'Integrated polarisation angle') - header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarisation angle error') + header['P_int'] = (P_diluted, 'Integrated polarization degree') + header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarization degree error') + header['PA_int'] = (PA_diluted, 'Integrated polarization angle') + header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarization angle error') return I_stokes, Q_stokes, U_stokes, Stokes_cov def compute_pol(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers): """ - Compute the polarisation degree (in %) and angle (in deg) and their + Compute the polarization degree (in %) and angle (in deg) and their respective errors from given Stokes parameters. ---------- Inputs: @@ -1255,10 +1255,10 @@ def compute_pol(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers): total intensity Q_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - vertical/horizontal linear polarisation intensity + vertical/horizontal linear polarization intensity U_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - +45/-45deg linear polarisation intensity + +45/-45deg linear polarization intensity Stokes_cov : numpy.ndarray Covariance matrix of the Stokes parameters I, Q, U. headers : header list @@ -1266,21 +1266,21 @@ def compute_pol(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers): ---------- Returns: P : numpy.ndarray - Image (2D floats) containing the polarisation degree (in %). + Image (2D floats) containing the polarization degree (in %). debiased_P : numpy.ndarray - Image (2D floats) containing the debiased polarisation degree (in %). + Image (2D floats) containing the debiased polarization degree (in %). s_P : numpy.ndarray - Image (2D floats) containing the error on the polarisation degree. + Image (2D floats) containing the error on the polarization degree. s_P_P : numpy.ndarray Image (2D floats) containing the Poisson noise error on the - polarisation degree. + polarization degree. PA : numpy.ndarray - Image (2D floats) containing the polarisation angle. + Image (2D floats) containing the polarization angle. s_PA : numpy.ndarray - Image (2D floats) containing the error on the polarisation angle. + Image (2D floats) containing the error on the polarization angle. s_PA_P : numpy.ndarray Image (2D floats) containing the Poisson noise error on the - polarisation angle. + polarization angle. new_headers : header list Updated list of headers corresponding to the reduced images accounting for the new orientation angle. @@ -1354,10 +1354,10 @@ def rotate_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, total intensity Q_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - vertical/horizontal linear polarisation intensity + vertical/horizontal linear polarization intensity U_stokes : numpy.ndarray Image (2D floats) containing the Stokes parameters accounting for - +45/-45deg linear polarisation intensity + +45/-45deg linear polarization intensity Stokes_cov : numpy.ndarray Covariance matrix of the Stokes parameters I, Q, U. data_mask : numpy.ndarray @@ -1379,10 +1379,10 @@ def rotate_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, accounting for total intensity new_Q_stokes : numpy.ndarray Rotated mage (2D floats) containing the rotated Stokes parameters - accounting for vertical/horizontal linear polarisation intensity + accounting for vertical/horizontal linear polarization intensity new_U_stokes : numpy.ndarray Rotated image (2D floats) containing the rotated Stokes parameters - accounting for +45/-45deg linear polarisation intensity. + accounting for +45/-45deg linear polarization intensity. new_Stokes_cov : numpy.ndarray Updated covariance matrix of the Stokes parameters I, Q, U. new_headers : header list @@ -1497,10 +1497,10 @@ def rotate_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, 2 + Q_diluted**2*U_diluted_err**2 - 2.*Q_diluted*U_diluted*QU_diluted_err) for header in new_headers: - header['P_int'] = (P_diluted, 'Integrated polarisation degree') - header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarisation degree error') - header['PA_int'] = (PA_diluted, 'Integrated polarisation angle') - header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarisation angle error') + header['P_int'] = (P_diluted, 'Integrated polarization degree') + header['P_int_err'] = (np.ceil(P_diluted_err*1000.)/1000., 'Integrated polarization degree error') + header['PA_int'] = (PA_diluted, 'Integrated polarization angle') + header['PA_int_err'] = (np.ceil(PA_diluted_err*10.)/10., 'Integrated polarization angle error') return new_I_stokes, new_Q_stokes, new_U_stokes, new_Stokes_cov, new_data_mask, new_headers diff --git a/package/overplot_MRK463E.py b/package/overplot_MRK463E.py index 79de0ca..fed7e2f 100755 --- a/package/overplot_MRK463E.py +++ b/package/overplot_MRK463E.py @@ -10,17 +10,11 @@ Stokes_Xr = fits.open("./data/MRK463E/Chandra/X_ray_crop.fits") levels = np.geomspace(1., 99., 7) -# A = overplot_chandra(Stokes_UV, Stokes_Xr) -# A.plot(levels=levels, SNRp_cut=3.0, SNRi_cut=20.0, zoom=1, savename='./plots/MRK463E/Chandra_overplot.pdf') - -B = overplot_chandra(Stokes_UV, Stokes_Xr, norm=LogNorm()) -B.plot(levels=levels, SNRp_cut=3.0, SNRi_cut=3.0, vec_scale=5, zoom=1, savename='./plots/MRK463E/Chandra_overplot_forced.pdf') -B.write_to(path1="./data/MRK463E/FOC_data_Chandra.fits", path2="./data/MRK463E/Chandra_data.fits", suffix="aligned") - -# C = overplot_pol(Stokes_UV, Stokes_IR) -# C.plot(SNRp_cut=3.0, SNRi_cut=20.0, savename='./plots/MRK463E/IR_overplot.pdf') +A = overplot_chandra(Stokes_UV, Stokes_Xr, norm=LogNorm()) +A.plot(levels=levels, SNRp_cut=3.0, SNRi_cut=3.0, vec_scale=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.*Stokes_UV[0].header['photflam'] -D = overplot_pol(Stokes_UV, Stokes_IR, norm=LogNorm()) -D.plot(levels=levels, SNRp_cut=3.0, SNRi_cut=3.0, vec_scale=5, norm=LogNorm(8.5e-18, 2.5e-15), savename='./plots/MRK463E/IR_overplot_forced.pdf') -D.write_to(path1="./data/MRK463E/FOC_data_WFPC.fits", path2="./data/MRK463E/WFPC_data.fits", suffix="aligned") +B = overplot_pol(Stokes_UV, Stokes_IR, norm=LogNorm()) +B.plot(levels=levels, SNRp_cut=3.0, SNRi_cut=3.0, vec_scale=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")