diff --git a/plots/NGC1068_x274020/NGC1068_FOC_bin10px.jpg b/plots/NGC1068_x274020/NGC1068_FOC_bin10px.jpg new file mode 100644 index 0000000..8442ea8 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_FOC_bin10px.jpg differ diff --git a/plots/NGC1068_x274020/NGC1068_FOC_bin10px.png b/plots/NGC1068_x274020/NGC1068_FOC_bin10px.png new file mode 100644 index 0000000..9407274 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_FOC_bin10px.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px.png new file mode 100644 index 0000000..a0ce617 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I.png new file mode 100644 index 0000000..14e7988 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_IQU.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_IQU.png new file mode 100644 index 0000000..d9151cf Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_IQU.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I_err.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I_err.png new file mode 100644 index 0000000..6431e9c Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_I_err.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P.png new file mode 100644 index 0000000..47e9a49 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_err.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_err.png new file mode 100644 index 0000000..4cdfab0 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_err.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_flux.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_flux.png new file mode 100644 index 0000000..d874c3e Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_P_flux.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRi.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRi.png new file mode 100644 index 0000000..9c76b46 Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRi.png differ diff --git a/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRp.png b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRp.png new file mode 100644 index 0000000..eb2c5db Binary files /dev/null and b/plots/NGC1068_x274020/NGC1068_K_FOC_bin10px_SNRp.png differ diff --git a/src/FOC_reduction.py b/src/FOC_reduction.py index 1b4365a..6e60d54 100755 --- a/src/FOC_reduction.py +++ b/src/FOC_reduction.py @@ -18,12 +18,13 @@ from astropy.wcs import WCS ##### User inputs ## Input and output locations -#globals()['data_folder'] = "../data/NGC1068_x274020/" -#globals()['infiles'] = ['x274020at_c0f.fits','x274020bt_c0f.fits','x274020ct_c0f.fits', -# 'x274020dt_c0f.fits','x274020et_c0f.fits','x274020ft_c0f.fits', -# 'x274020gt_c0f.fits','x274020ht_c0f.fits','x274020it_c0f.fits'] -##psf_file = 'NGC1068_f253m00.fits' -#globals()['plots_folder'] = "../plots/NGC1068_x274020/" +globals()['data_folder'] = "../data/NGC1068_x274020/" +#globals()['infiles'] = ['xn1c400.fits','xn2c400.fits','xn3c400.fits'] +globals()['infiles'] = ['x274020at_c0f.fits','x274020bt_c0f.fits','x274020ct_c0f.fits', + 'x274020dt_c0f.fits','x274020et_c0f.fits','x274020ft_c0f.fits', + 'x274020gt_c0f.fits','x274020ht_c0f.fits','x274020it_c0f.fits'] +#psf_file = 'NGC1068_f253m00.fits' +globals()['plots_folder'] = "../plots/NGC1068_x274020/" #globals()['data_folder'] = "../data/IC5063_x3nl030/" #globals()['infiles'] = ['x3nl0301r_c0f.fits','x3nl0302r_c0f.fits','x3nl0303r_c0f.fits'] @@ -90,13 +91,13 @@ from astropy.wcs import WCS #globals()['plots_folder'] = "../plots/3C273_x0u20/" #BEWARE: 5 observations separated by 1 year each (1995, 1996, 1997, 1998, 1999) -globals()['data_folder'] = "../data/M87/POS1/" +#globals()['data_folder'] = "../data/M87/POS1/" #globals()['infiles'] = ['x2py010ct_c0f.fits','x2py010dt_c0f.fits','x2py010et_c0f.fits','x2py010ft_c0f.fits'] #1995 #globals()['infiles'] = ['x3be010ct_c0f.fits','x3be010dt_c0f.fits','x3be010et_c0f.fits','x3be010ft_c0f.fits'] #1996 #globals()['infiles'] = ['x43r010km_c0f.fits','x43r010mm_c0f.fits','x43r010om_c0f.fits','x43r010rm_c0f.fits'] #1997 #globals()['infiles'] = ['x43r110kr_c0f.fits','x43r110mr_c0f.fits','x43r110or_c0f.fits','x43r110rr_c0f.fits'] #1998 -globals()['infiles'] = ['x43r210kr_c0f.fits','x43r210mr_c0f.fits','x43r210or_c0f.fits','x43r210rr_c0f.fits'] #1999 -globals()['plots_folder'] = "../plots/M87/POS1/" +#globals()['infiles'] = ['x43r210kr_c0f.fits','x43r210mr_c0f.fits','x43r210or_c0f.fits','x43r210rr_c0f.fits'] #1999 +#globals()['plots_folder'] = "../plots/M87/POS1/" #BEWARE: 5 observations separated by 1 year each (1995, 1996, 1997, 1998, 1999) #globals()['data_folder'] = "../data/M87/POS3/" @@ -128,28 +129,28 @@ def main(): # Data binning rebin = True if rebin: - pxsize = 0.05 - px_scale = 'arcsec' #pixel, arcsec or full + pxsize = 10 + px_scale = 'pixel' #pixel, arcsec or full rebin_operation = 'sum' #sum or average # Alignement align_center = 'image' #If None will align image to image center display_data = False # Smoothing smoothing_function = 'combine' #gaussian_after, weighted_gaussian_after, gaussian, weighted_gaussian or combine - smoothing_FWHM = 0.10 #If None, no smoothing is done + smoothing_FWHM = None #If None, no smoothing is done smoothing_scale = 'arcsec' #pixel or arcsec # Rotation rotate_stokes = True #rotation to North convention can give erroneous results rotate_data = False #rotation to North convention can give erroneous results # Final crop crop = False #Crop to desired ROI - final_display = True + final_display = False # Polarization map output - figname = 'M87_POS1_1999_FOC' #target/intrument name - figtype = '_combine_FWHM010' #additionnal informations + figname = 'NGC1068_K_FOC' #target/intrument name + figtype = '_bin10px' #additionnal informations SNRp_cut = 3. #P measurments with SNR>3 SNRi_cut = 30. #I measurments with SNR>30, which implies an uncertainty in P of 4.7%. - step_vec = 0 #plot all vectors in the array. if step_vec = 2, then every other vector will be plotted + step_vec = 1 #plot all vectors in the array. if step_vec = 2, then every other vector will be plotted # if step_vec = 0 then all vectors are displayed at full length ##### Pipeline start diff --git a/src/lib/plots.py b/src/lib/plots.py index 5935c60..6d46f63 100755 --- a/src/lib/plots.py +++ b/src/lib/plots.py @@ -1703,12 +1703,14 @@ class pol_map(object): ax.add_artist(self.north_dir) def display(self, fig=None, ax=None): + norm = None if self.display_selection is None: self.display_selection = "total_flux" if self.display_selection.lower() in ['total_flux']: self.data = self.I*self.convert_flux - vmin, vmax = 0., np.max(self.data[self.data > 0.]) + vmin, vmax = np.min(self.data[self.cut])/10., np.max(self.data[self.data > 0.]) + norm = LogNorm(vmin, vmax) label = r"$F_{\lambda}$ [$ergs \cdot cm^{-2} \cdot s^{-1} \cdot \AA^{-1}$]" elif self.display_selection.lower() in ['pol_flux']: self.data = self.I*self.convert_flux*self.P @@ -1738,12 +1740,18 @@ class pol_map(object): ax = self.ax if hasattr(self, 'im'): self.im.remove() - self.im = ax.imshow(self.data, vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno') + if not norm is None: + self.im = ax.imshow(self.data, norm=norm, aspect='equal', cmap='inferno') + else: + self.im = ax.imshow(self.data, vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno') self.cbar = plt.colorbar(self.im, cax=self.cbar_ax, label=label) fig.canvas.draw_idle() return self.im else: - im = ax.imshow(self.data, vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno') + if not norm is None: + im = ax.imshow(self.data, norm=norm, aspect='equal', cmap='inferno') + else: + im = ax.imshow(self.data, vmin=vmin, vmax=vmax, aspect='equal', cmap='inferno') ax.set_xlim(0,self.data.shape[1]) ax.set_ylim(0,self.data.shape[0]) plt.colorbar(im, pad=0.025, aspect=80, label=label) diff --git a/src/lib/reduction.py b/src/lib/reduction.py index 955dca8..2bd5d10 100755 --- a/src/lib/reduction.py +++ b/src/lib/reduction.py @@ -1475,24 +1475,23 @@ def rotate_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, data_mask, headers, new_U_stokes = np.zeros(shape) new_Stokes_cov = np.zeros((*Stokes_cov.shape[:-2],*shape)) - for i in range(shape[0]): - for j in range(shape[1]): - new_I_stokes[i,j], new_Q_stokes[i,j], new_U_stokes[i,j] = np.dot(mrot, np.array([I_stokes[i,j], Q_stokes[i,j], U_stokes[i,j]])).T - new_Stokes_cov[:,:,i,j] = np.dot(mrot, np.dot(Stokes_cov[:,:,i,j], mrot.T)) - #Rotate original images using scipy.ndimage.rotate - new_I_stokes = sc_rotate(new_I_stokes, ang, order=1, reshape=False, cval=0.) - new_Q_stokes = sc_rotate(new_Q_stokes, ang, order=1, reshape=False, cval=0.) - new_U_stokes = sc_rotate(new_U_stokes, ang, order=1, reshape=False, cval=0.) + new_I_stokes = sc_rotate(I_stokes, ang, order=1, reshape=False, cval=0.) + new_Q_stokes = sc_rotate(Q_stokes, ang, order=1, reshape=False, cval=0.) + new_U_stokes = sc_rotate(U_stokes, ang, order=1, reshape=False, cval=0.) new_data_mask = sc_rotate(data_mask.astype(float)*10., ang, order=1, reshape=False, cval=0.) new_data_mask[new_data_mask < 2] = 0. new_data_mask = new_data_mask.astype(bool) for i in range(3): for j in range(3): - new_Stokes_cov[i,j] = sc_rotate(new_Stokes_cov[i,j], ang, order=1, - reshape=False, cval=0.) + new_Stokes_cov[i,j] = sc_rotate(Stokes_cov[i,j], ang, order=1, reshape=False, cval=0.) new_Stokes_cov[i,i] = np.abs(new_Stokes_cov[i,i]) + for i in range(shape[0]): + for j in range(shape[1]): + new_I_stokes[i,j], new_Q_stokes[i,j], new_U_stokes[i,j] = np.dot(mrot, np.array([new_I_stokes[i,j], new_Q_stokes[i,j], new_U_stokes[i,j]])).T + new_Stokes_cov[:,:,i,j] = np.dot(mrot, np.dot(new_Stokes_cov[:,:,i,j], mrot.T)) + #Update headers to new angle new_headers = [] mrot = np.array([[np.cos(-alpha), -np.sin(-alpha)],