标准化画图【cmap】【imshow】【EllipseSkyRegion】【contour】【set_ticklabels】【coords】【colorbar】

标准化画图，包括colorbar，imshow数据，坐标投影到世界坐标系，画分辨率beam，contour轮廓。由以下代码提供，文件来源MPIFR网址

from astropy.io import fits as pf
from reproject import reproject_interp
from astropy.wcs.wcs import WCS
from matplotlib.pyplot import plot, sca
import numpy as np
from astropy import wcs
from matplotlib import pylab as plt
import matplotlib as mpl
from mpl_toolkits.axes_grid1 import make_axes_locatable, axes_size
from copy import deepcopy
from regions import PixCoord
from regions import EllipseSkyRegion, EllipsePixelRegion
from astropy import units as u
from astropy.coordinates import SkyCoord
import matplotlib.colors as colors
import math
import sys

#数据导入
data = pf.getdata('m31_6cm.fits')[0]
hdr0 = pf.getheader('m31_6cm.fits')

#修改头文件，其实也可以不修改，因为原来的fits文件的beam信息没有给出所以这边修改成自己想要的。header各个信息参考以下文章
#https://blog.csdn.net/weixin_42842694/article/details/118723782

hdr = pf.Header({ 'SIMPLE':   True,
                  'BITPIX':    -32,
                   'NAXIS':       2,
                  'NAXIS1':  hdr0['NAXIS1'],
                  'NAXIS2':  hdr0['NAXIS2'],
                  'CDELT1':  hdr0['CDELT1'],
                  'CRVAL1':  hdr0['CRVAL1'],
                  'CRPIX1':  hdr0['CRPIX1'],
                  'CTYPE1':  hdr0['CTYPE1'],
                  'CDELT2':  hdr0['CDELT2'],
                  'CRPIX2':  hdr0['CRPIX2'],
                  'CRVAL2':  hdr0['CRVAL2'],
                  'CTYPE2':  hdr0['CTYPE2'],
                  'BMAJ'  :  1.5000E-01,
                  'BMIN'  :  1.5000E-01,
                  'BPA'   :  0.00,
                  'EQUINOX': 2000.00, 
                  'EXTEND':True
                             })


#设置画布的尺寸大小
plt.rcParams.update({'font.size': 30})
fig = plt.figure(figsize=(12,12))

wcs = WCS(hdr)
ax = plt.subplot(1,1,1, projection=wcs)

#colorbar样式选择
cm1 = plt.get_cmap('cubehelix_r')
#norm = mpl.colors.LogNorm(vmax=vmax_v, vmin=vmin_v)

#imshow参考
#https://blog.csdn.net/weixin_42842694/article/details/119249749
vmax_v = np.percentile(data,95)     
vmin_v = np.percentile(data,4)
im = ax.imshow(data, origin='lower', cmap=cm1, vmax=vmax_v, vmin=vmin_v)

ax.coords['ra'].set_axislabel('Right Ascension')
ax.coords['dec'].set_axislabel('Declination')

#绘制分辨率beam。单口径天线是一个圆，如果是干涉阵是一个有倾角的椭圆
word = wcs.wcs_pix2world([[20,20]],0)  #选择椭圆要绘制的位置，以图中pix为准
center_sky = SkyCoord(word[0][0], word[0][1], unit='deg', frame='fk5')
region_sky = EllipseSkyRegion(center=center_sky,
                              height=hdr['BMAJ'] * u.deg, width=hdr['BMIN'] * u.deg,
                              angle=hdr['BPA'] * u.deg)
pixel_region = region_sky.to_pixel(wcs)
pixel_region.plot(color='b')

#画contour，level可以设置需要每个contour的间隔
level =np.linspace(vmin_v, vmax_v, 5)
ax.contour(data,levels=level,colors=['black'])
 
#设置colorbar的位置和大小
position=fig.add_axes([0.77, 0.11, 0.04, 0.77])#位置[左,下,右,上]
#位置左和下相当于确定位置（下是从下边开始算0），右和上是在确定了左，下的增量
#说的不是很好，参考我文章的【colorbar位置参数详解】
cb = fig.colorbar(im, cax = position)

#或者可以考虑用这个来设置colorbar的宽和长，区别在自动帮你设置好colorbar的位置
cb = fig.colorbar(im, fraction =0.04, pad = 0.03)

#设置colorbar自定义刻度和label
colors_range = np.linspace(vmin_v, vmax_v, 5)
cb.set_ticks(colors_range)
cb.set_ticklabels([1,2,3,4,5])
cb.set_label(r'mJy beam$^{-1}$')

plt.show()