画图,subplot多子图

import matplotlib.pyplot as plt
import os
import numpy as np
from astropy.coordinates import SkyCoord
from astropy.table import Table
from astropy.io import fits
import warnings
warnings.filterwarnings('ignore')

T = Table.read('atoms_B3_12M_hotcore_candidate_modified.tbl',format='ascii')
#spws = ['spw23','spw25','spw27','spw29','spw31','spw33','spw35','spw37','spw39']
lines = ['H$^{13}$CN','H$^{13}$CO$^+$','CCH','SiO','HCN','HCO$^+$']
freqs = [86339.9,86754.3,87316.9,86846.96,88631.847,89188.526]

def seek_file(path,sourcename):
    files = []
    for filename in os.listdir(path):
        fp = os.path.join(filename)
        fmts = '.txt'
        if os.path.isfile(path+'/'+fp) and fmts in filename:
            files.append(fp)
    files.sort()
    return files

def show_spec(path,sourcename,ID,files,ra,dec):
    a = 0
    coord = SkyCoord(ra,dec,unit='deg',frame='fk5')
    hh,hm,hs = coord.ra.hms
    dd,dm,ds = coord.dec.dms
    dm,ds = abs(dm),abs(ds)
    for j in range(len(files)):
        filename = files[j]
        x,y = np.loadtxt(path+'/'+filename,usecols=(0,1),unpack=True)
        fitsfile = str(filename)[:-4]+'.fits'
        hdu = fits.open(path+'/'+fitsfile)
        windows = hdu[0].header['FILNAM04']
        plt.subplots_adjust(hspace=0.3)
        if len(x)>2000:
            a = a+1
            ax1 = plt.subplot(5,2,a)
            ax1.tick_params(labelsize=4,direction='in', length=1, width=0.2,pad=0.7)
            plt.title('%s_%s_%s' %(sourcename,ID,windows),loc='center',fontsize=4,pad=0.7)
            if a==9:
                plt.xlabel('Frequency (MHz)',fontsize=4)
            plt.ylabel('Temperature (K)',fontsize=4)
            a = a+1
            ax2 = plt.subplot(5,2,a)
            ax2.tick_params(labelsize=4,direction='in', length=1, width=0.2,pad=0.7)
            ax1.plot(x[:int(len(x)/2)],y[:int(len(x)/2)],c='k',linewidth=0.2)
            ax1.set_xlim(xmin=np.nanmin(x[:int(len(x)/2)]),xmax=np.nanmax(x[:int(len(x)/2)]))
            ax2.plot(x[int(len(x)/2):],y[int(len(x)/2):],c='k',linewidth=0.2)
            ax2.set_xlim(xmin=np.nanmin(x[int(len(x)/2):]),xmax=np.nanmax(x[int(len(x)/2):]))
            if a==10:
                plt.xlabel('Frequency (MHz)',fontsize=4)
            plt.title('%s_%s_%s' %(sourcename,ID,windows),loc='center',fontsize=4,pad=0.7)
            plt.ylabel('Temperature (K)',fontsize=4)
        else:
            a = a+1
            ax = plt.subplot(5,2,a)
            ax.tick_params(labelsize=4,direction='in', length=1, width=0.2,pad=0.7)
            ax.plot(x,y,c='k',linewidth=0.2)
            ax.set_xlim(xmin=np.nanmin(x),xmax=np.nanmax(x))
            #plt.xlabel('Frequency (MHz)',fontsize=2)
            plt.ylabel('Temperature (K)',fontsize=4)
            for n in range(len(lines)):
                if np.nanmin(x)<freqs[n]<np.nanmax(x):
                    plt.title('%s_%s_%s  %s' %(sourcename,ID,windows,lines[n]),loc='center',fontsize=4,pad=0.7)

    #plt.tight_layout() 
    plt.savefig('spec/figure/%s_%s.pdf' %(sourcename,ID))
    plt.close('all')


for i in range(len(T['sname'])):
    path = 'spec/'+T['sname'][i]+'/'+str(T['ID'][i])
    files = seek_file(path,T['sname'][i])
    show_spec(path,T['sname'][i],str(T['ID'][i]),files,T['ra'][i],T['dec'][i])

 

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值