reproduce contours from Morganti2018

src/lib/plots.py

@@ -580,7 +580,7 @@ class align_maps(object):
     
     def on_close_align(self, event):
         self.aligned = True
-        print(self.get_aligned_wcs())
+        #print(self.get_aligned_wcs())
     
     def align(self):
         self.fig.canvas.draw()
@@ -595,7 +595,7 @@ class overplot_maps(align_maps):
     Class to overplot maps from different observations.
     Inherit from class align_maps in order to get the same WCS on both maps.
     """
-    def overplot(self, other_vmin, other_vmax, other_num, SNRp_cut=3., SNRi_cut=30.):
+    def overplot(self, other_levels, SNRp_cut=3., SNRi_cut=30., savename=None):
         #Get Data
         obj = self.Stokes_UV[0].header['targname']
         stkI = self.Stokes_UV[np.argmax([self.Stokes_UV[i].header['datatype']=='I_stokes' for i in range(len(self.Stokes_UV))])]
@@ -606,6 +606,10 @@ class overplot_maps(align_maps):
         
         other_data = self.other_map[0].data
         other_unit = self.other_map[0].header['bunit']
+        other_convert = 1.
+        if other_unit.lower() == 'jy/beam':
+            other_unit = r"mJy/Beam"
+            other_convert = 1e3
         other_freq = self.other_map[0].header['crval3']
         
         convert_flux = self.Stokes_UV[0].header['photflam']
@@ -639,15 +643,18 @@ class overplot_maps(align_maps):
         self.ax.autoscale(False)
 
         #Display other map as contours
-        other_cont = self.ax.contour(other_data, transform=self.ax.get_transform(self.wcs_other), levels=np.linspace(other_vmin, other_vmax, other_num), colors='grey')
+        other_cont = self.ax.contour(other_data*other_convert, transform=self.ax.get_transform(self.wcs_other), levels=other_levels*other_convert, colors='grey')
         self.ax.clabel(other_cont, inline=True, fontsize=8)
 
-        self.ax.set(xlabel="Right Ascension (J2000)", ylabel="Declination (J2000)", title="HST/FOC UV polarization map of {0:s} overplotted with {1:.2e}Hz map in {2:s}.".format(obj, other_freq, other_unit))
+        self.ax.set(xlabel="Right Ascension (J2000)", ylabel="Declination (J2000)", title="HST/FOC UV polarization map of {0:s} overplotted with {1:.2f}GHz map in {2:s}.".format(obj, other_freq*1e-9, other_unit))
+
+        if not(savename is None):
+            self.fig2.savefig(savename,bbox_inches='tight',dpi=200)
 
         self.fig2.canvas.draw()
     
-    def plot(self, other_vmin, other_vmax, other_num, SNRp_cut=3., SNRi_cut=30.) -> None:
+    def plot(self, levels, SNRp_cut=3., SNRi_cut=30., savename=None) -> None:
         self.align()
         if self.aligned:
-            self.overplot(other_vmin, other_vmax, other_num, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut)
+            self.overplot(other_levels=levels, SNRp_cut=SNRp_cut, SNRi_cut=SNRi_cut, savename=savename)
             plt.show(block=True)

src/lib/test_overplot.py

@@ -1,14 +1,20 @@
+#!/usr/bin/python3
 from astropy.io import fits
+import numpy as np
 from plots import overplot_maps
 
 Stokes_UV = fits.open("../../data/IC5063_x3nl030/IC5063_FOC_combine_FWHM020.fits")
 Stokes_18GHz = fits.open("../../data/IC5063_x3nl030/radio/IC5063.18GHz.fits")
 Stokes_24GHz = fits.open("../../data/IC5063_x3nl030/radio/IC5063.24GHz.fits")
 
-A = overplot_maps(Stokes_UV, Stokes_18GHz)
-A.plot(5e-3, 5e-2, 8, SNRp_cut=6.0, SNRi_cut=180.0)
-A.fig2.savefig('../../plots/IC5063_x3nl030/1.8GHz_overplot.png',bbox_inches='tight',overwrite=True)
+levelsMorganti = np.array([1.,2.,3.,8.,16.,32.,64.,128.])
 
+#levels18GHz = np.array([0.6, 1.5, 3, 6, 12, 24, 48, 96])/100.*Stokes_18GHz[0].data.max()
+levels18GHz = levelsMorganti*0.28*1e-3
+A = overplot_maps(Stokes_UV, Stokes_18GHz)
+A.plot(levels=levels18GHz, SNRp_cut=6.0, SNRi_cut=180.0, savename='../../plots/IC5063_x3nl030/18GHz_overplot.png')
+
+#levels24GHz = np.array([1.,1.5, 3, 6, 12, 24, 48, 96])/100.*Stokes_24GHz[0].data.max()
+levels24GHz = levelsMorganti*0.46*1e-3
 B = overplot_maps(Stokes_UV, Stokes_24GHz)
-B.plot(1e-3, 3e-2, 8, SNRp_cut=6.0, SNRi_cut=180.0)
-B.fig2.savefig('../../plots/IC5063_x3nl030/2.4GHz_overplot.png',bbox_inches='tight',overwrite=True)
+B.plot(levels=levels24GHz, SNRp_cut=6.0, SNRi_cut=180.0, savename='../../plots/IC5063_x3nl030/24GHz_overplot.png')