import xlrd
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.gridspec import GridSpec
import scipy.stats as stats
from decimal import Decimal
import matplotlib as mpl
import matplotlib
import matplotlib.cm as cm
from matplotlib.colors import ListedColormap, LinearSegmentedColormap
import math

wb = xlrd.open_workbook(filename='realdata1.xlsx')  # 读取Excel观测数据
realdata = wb.sheet_by_index(0)  # 通过索引获取表格

# read winter pol obs
TSK_obs_01_22_24 = TSK_obsx = realdata.col_values(22)[148:291:2]  # 第23列温度 p
TSK_obs_01_08_14 = TSK_obsx = realdata.col_values(22)[17793:18128:2]
TSK_obs_01_08_10 = TSK_obsx = realdata.col_values(22)[17793:17936:2]#18128:2] #17936 c
TSK_obs_01_11 = TSK_obsx = realdata.col_values(22)[17937:17984:2] #17936 p
TSK_obs_01_12 = TSK_obsx = realdata.col_values(22)[17985:18032:2] #17936  c
TSK_obs_01_13_14 = TSK_obsx = realdata.col_values(22)[18033:18128:2] #17936 p
TSK_obs_01_20_21 = TSK_obsx = realdata.col_values(22)[52:147:2]  # c
TSK_obs_01_15 = TSK_obsx = realdata.col_values(22)[18129:18176:2] #17936 p
T_obs_01_22_24 = T_obsx = realdata.col_values(7)[148:291:2]  # 第8列温度 p
rh_obs_01_22_24 = rh_obsx = realdata.col_values(8)[148:291:2]  # 9湿度 p
T_obs_01_08_14 = T_obsx = realdata.col_values(7)[17793:18128:2]
T_obs_01_08_10 = T_obsx = realdata.col_values(7)[17793:17936:2]#18128:2] #17936 c
T_obs_01_11 = T_obsx = realdata.col_values(7)[17937:17984:2] #17936 p
T_obs_01_12 = T_obsx = realdata.col_values(7)[17985:18032:2] #17936  c
T_obs_01_13_14 = T_obsx = realdata.col_values(7)[18033:18128:2] #17936 p
T_obs_01_20_21 = T_obsx = realdata.col_values(7)[52:147:2]  # c
rh_obs_01_08_14 = rh_obsx = realdata.col_values(8)[17793:18128:2]
rh_obs_01_20_21 = rh_obsx = realdata.col_values(8)[52:147:2]
rh_obs_01_08_10 = rh_obsx = realdata.col_values(8)[17793:17936:2]#18128:2] #17936 c
rh_obs_01_11 = rh_obsx = realdata.col_values(8)[17937:17984:2] #17936 p
rh_obs_01_12 = rh_obsx = realdata.col_values(8)[17985:18032:2] #17936  c
rh_obs_01_13_14 = rh_obsx = realdata.col_values(8)[18033:18128:2] #17936 p
T_obs_01_15 = T_obsx = realdata.col_values(7)[18129:18176:2] #17936 p
rh_obs_01_15 = rh_obsx = realdata.col_values(8)[18129:18176:2]  # p

T_obs_jan = T_obs_01_08_14 + T_obs_01_15
T_obs_jan2 = T_obs_01_20_21 + T_obs_01_22_24
#T_obs_janc = T_obs_01_08_14 + T_obs_01_20_21
#T_obs_janp = T_obs_01_15 + T_obs_01_22_24
T_obs_janc = T_obs_01_08_10 + T_obs_01_12 + T_obs_01_20_21  # 6 days
T_obs_janp = T_obs_01_11 + T_obs_01_13_14 + T_obs_01_15 + T_obs_01_22_24  # 7 days
rh_obs_jan = rh_obs_01_08_14 + rh_obs_01_15
rh_obs_jan2 = rh_obs_01_20_21 + rh_obs_01_22_24
rh_obs_janc = rh_obs_01_08_10 + rh_obs_01_12 + rh_obs_01_20_21
rh_obs_janp = rh_obs_01_11 + rh_obs_01_13_14 + rh_obs_01_15 + rh_obs_01_22_24
TSK_obs_jan = TSK_obs_01_08_14 + TSK_obs_01_15
TSK_obs_jan2 = TSK_obs_01_20_21 + TSK_obs_01_22_24
#T_obs_janc = T_obs_01_08_14 + T_obs_01_20_21
#T_obs_janp = T_obs_01_15 + T_obs_01_22_24
TSK_obs_janc = TSK_obs_01_08_10 + TSK_obs_01_12 + TSK_obs_01_20_21  # 6 days
TSK_obs_janp = TSK_obs_01_11 + TSK_obs_01_13_14 + TSK_obs_01_15 + TSK_obs_01_22_24  # 7 days



ddate1 = ['01-07', '01-08', '01-09', '01-10', '01-11', '01-12', '01-13', '01-14',
          '01-15', '01-19', '01-20', '01-21', '01-22', '01-23']  # data date
model1 = ['AC_0.01_2014', 'ACM2_MM5_0.2_2014', 'ACM2_MM5_0.5_2014','AC_0.8_2014','AC_1.0_2014', 'BL_0.1_2014',
          'BL_0.2_2014', 'BL_0.5_2014', 'ACM2_kf_2014', 'ACM2_0.55_2014']
model1_pf = ['AC_0.01', 'ACM2_MM5_0.2', 'ACM2_MM5_0.5','AC_0.8','AC_1.0', 'BL_0.1',
             'BL_0.2', 'BL_0.5', 'ACM2_kf', 'ACM2_0.55']
model2 = ['ACM2_0.01', 'ACM2_0.2', 'ACM2_0.5', 'ACM2_0.8', 'ACM2_1.0', 'BL_0.1', 'BL_0.2', 'BL_0.5', 'ACM2_CMAQ', 'AC_0.55']
model2diff = ['ACM2_0.01', 'ACM2_0.2 - ACM2_0.01', 'ACM2_0.5 - ACM2_0.01',
              'ACM2_0.8 - ACM2_0.01', 'ACM2_1.0 - ACM2_0.01', 'BL_0.1', 'BL_0.2', 'BL_0.5', 'ACM2_CMAQ', 'AC_0.55']
linetype = ['-', '-', '-','-', '-', '--', '--', '--', '-']
ptype = ['x', 'x', 'x','x', 'x', '^', '^', '^', 'o']
colorline = ['r', 'b', 'g', 'c', 'm', 'r', 'b', 'g', 'y']
eta_levels = [1.000, 0.997, 0.994, 0.991,
                0.988, 0.985, 0.980, 0.975, 0.970, 0.960,
                0.950, 0.940, 0.930, 0.920, 0.910, 0.895,
                0.880, 0.865, 0.850, 0.825, 0.800, 0.775,
                0.750, 0.720, 0.690, 0.660, 0.630, 0.600,
                0.570, 0.540, 0.510, 0.475, 0.440,
                0.405, 0.370, 0.330, 0.290, 0.250,
                0.210, 0.175, 0.145, 0.115, 0.090,
                0.060, 0.045, 0.025, 0.010, 0.000,]
p_top = 5000  # hpa
SIGMA = 5.67e-8  # SIGMA sb
emiss = 0.98 # emiss
KZM = [0.01, 0.2, 0.5, 0.1, 0.2, 0.5]

def rh2(modeldata):
    psfc = modeldata[:, 7]
    q2 = modeldata[:, 1]
    T = modeldata[:, 2]
    e = list(map(lambda x, y: x * y / 0.622 / 100, psfc, q2))
    e0 = list(map(lambda x: 6.11 * 10 ** ((7.45 * (x - 273.15)) / (x - 273.15 + 235)), T))
    rh2_sc = list(map(lambda x, y: x / y * 100, e, e0))
    return rh2_sc


def rmse(predictions, targets):
    return np.sqrt(((predictions - targets) ** 2).mean())


def ioa(p, o):
    ioa = 1 - np.sum((p - o) ** 2)/np.sum((np.fabs(p - o.mean())+np.fabs(o-o.mean())) ** 2)
    return ioa

def wind(uu, vv):
    u = np.array(uu)
    v = np.array(vv)
    ws = np.sqrt(u ** 2 + v ** 2)
    theta = u
    n = len(theta)
    for i in range(n):
        if u[i] > 0 and v[i] > 0:
            theta[i] = 3/2 * math.pi - np.arctan(v[i]/u[i])
        elif u[i] < 0 and v[i] > 0:
            theta[i] = 1/2 * math.pi - np.arctan(v[i]/u[i])
        elif u[i] > 0 and v[i] < 0:
            theta[i] = 3/2 * math.pi - np.arctan(v[i]/u[i])
        elif u[i] < 0 and v[i] < 0:
            theta[i] = 1/2 * math.pi - np.arctan(v[i]/u[i])
        elif u[i] == 0 and v[i] > 0:
            theta[i] = math.pi
        elif u[i] == 0 and v[i] < 0:
            theta[i] = 0
        elif u[i] > 0 and v[i] == 0:
            theta[i] = 3/2 * math.pi
        elif u[i] < 0 and v[i] == 0:
            theta[i] = 1/2 * math.pi
        elif u[i] == 0 and v[i] == 0:
            theta[i] = -999
    theta1 = theta/(2*math.pi)*360

    return ws, theta1


def s_plo_sizeset(yname, ys, ye):
    ax1.set_xticks(np.arange(0, 192, 24))
    ax1.set_xticklabels(['jan08', '09', '10', '11', '12', '13', '14', '15', '16'])
    ax1.set_ylabel(yname, fontsize=25)
    ax2.set_xticks(np.arange(0, 121, 24))
    ax2.set_xticklabels(['20', '21', '22', '23', '24'])
    ax1.set_xlim(0, 192)
    ax2.set_xlim(0, 120)  #ax2.set_xlim(72, 119)
    ax1.set_ylim(ys, ye)
    ax2.set_ylim(ys, ye)
    ax1.spines['right'].set_visible(False)
    #  ax5.spines['right'].set_visible(True)
    ax2.spines['left'].set_visible(False)
    ax2.tick_params(labelleft=False)  # don't put tick labels at the top
    ax1.tick_params(labelsize=25)
    ax2.tick_params(labelsize=25)
    bwith = 2
    ax1.spines['top'].set_linewidth(bwith)
    ax1.spines['left'].set_linewidth(bwith)
    ax1.spines['bottom'].set_linewidth(bwith)
    ax2.spines['top'].set_linewidth(bwith)
    ax2.spines['right'].set_linewidth(bwith)
    ax2.spines['bottom'].set_linewidth(bwith)
    #  ax5.yaxis.tick_right()
    ax2.yaxis.set_visible(False)
    d = .005  # how big to make the diagonal lines in axes coordinates
    # arguments to pass to plot, just so we don't keep repeating them
    kwargs = dict(transform=ax1.transAxes, color='k', clip_on=False)
    ax1.plot((1 - d, 1 + d), (1 - d, 1 + d), **kwargs)        # B-R diagonal
    ax1.plot((1 - d, 1 + d), (-d, +d), **kwargs)  # top-right diagonal
    kwargs = dict(transform=ax2.transAxes, color='k', clip_on=False)
    ax2.plot((- d, + d), (1 - d, 1 + d), **kwargs)        # B-R diagonal
    ax2.plot((- d, + d), (-d, +d), **kwargs)  # top-right diagonal
#    ax1.annotate('s-pol', xy=(0.077, 0.883), xycoords='figure fraction', fontsize=15)


def box(xdata, ydata, xlabel, ylabel, s, e, col, labeln, axzanshi, themarker, size):
    # r1 = np.corrcoef(xdata, ydata)
    # print(r1)
    r2, p = stats.pearsonr(xdata, ydata)
    print(r2, p)
    r2_3 = Decimal(r2).quantize(Decimal('0.000'))  # 四舍五入保留小数点3位
    #  r2_3 = round(r2, 3)
    rr = str(r2_3)  # 转字符串
    # print(r2_3, rr)
    rm = rmse(np.array(xdata), np.array(ydata))
    rm_3 = Decimal(rm).quantize(Decimal('0.000'))
    #  rm_3 = round(rm, 3)
    rrm = str(rm_3)
    io = ioa(np.array(xdata), np.array(ydata))
    io_3 = Decimal(io).quantize(Decimal('0.000'))
    iio = str(io_3)
    ax = axzanshi
    ax.plot([-1000, 1000], [-1000, 1000], '--', color='k', linewidth=2.0)
    ax.scatter(xdata, ydata, color=col, s=size, marker=themarker, linewidths=0.00001,
               label=labeln+'  R='+rr+'  RMSE='+rrm+'  IOA='+iio, edgecolors='none')
    plt.xlim(s, e)
    plt.ylim(s, e)
    plt.yticks(fontsize=15)
    plt.xticks(fontsize=15)
    plt.xlabel(xlabel, fontsize=15)
    plt.ylabel(ylabel, fontsize=15)
    leg = ax.legend(loc='upper left', frameon=False)  # legend设置
    ltext = leg.get_texts()
    plt.setp(ltext, fontsize=15)
    # end legend
    fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
    fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
    fig.subplots_adjust(top=0.97)  # 设置窗口右内边距


def lk(modeldata0,psta1,pend1,psta2,pend2,psta3,pend3,col):
    t2p0 = np.concatenate((modeldata0[psta1:pend1, col], modeldata0[psta2:pend2, col]), axis=0)
    t2p0 = np.concatenate((t2p0, modeldata0[psta3:pend3, col]), axis=0)
    return t2p0


def lkrh(modeldata0,psta1,pend1,psta2,pend2,psta3,pend3):
    t2p0 = np.concatenate((modeldata0[psta1:pend1], modeldata0[psta2:pend2]), axis=0)
    t2p0 = np.concatenate((t2p0, modeldata0[psta3:pend3]), axis=0)
    return t2p0



def meanvar(T):

    T2 = T.reshape(int(len(T)/24),24)
    T2m = T2.mean(axis=0)
    return T2m

def avelevel(vardata,modeldata,p_top,eta_level,level):  #vardata= porfile ,eta_level[i] 第i层 大于1  位温转温度
    vardata = vardata.transpose(1, 0)
    psfc = modeldata[:, 7]
    theta = np.concatenate((vardata[level-1,0:192],vardata[level-1,216:336]), axis=0)
    p0=eta_level*(psfc-p_top)+p_top
    p=np.concatenate((p0[0:192],p0[216:336]), axis=0)
    p1=((100000/p)**0.286)
    T = theta/p1
    return T

def tstheta(modeldata):  #vardata= porfile ,eta_level=第几层 大于1
    psfc = modeldata[:, 7]
    TSK = modeldata[:, 16]
    TStheta = TSK*(100000/psfc)**0.286
    TStheta1 = np.concatenate((TStheta[0:192],TStheta[216:336]), axis=0)
    return TStheta1

var = ['EXCH_H', 'T', 'Q']
t0 = 300
for k in range(3):  # var
    for j in range(10):  # model
        vars()[model1[j]] = []
        for i in range(len(ddate1)):  # date
            datazanshi = []
            datazanshi = np.loadtxt(open(model1_pf[j] + '-' + var[k] + '-porfile-' + ddate1[i], 'rb'))
            if i < 1:
                vars()[model1[j] + '_' + var[k]] = datazanshi
            else:
                vars()[model1[j] + '_' + var[k]] = np.concatenate((vars()[model1[j] + '_' + var[k]], datazanshi), axis=0)
        if k == 1:
            vars()[model1[j] + '_' + 'theta'] = vars()[model1[j] + '_' + var[k]] + t0
varname = ['EXCH_H', 'theta', 'Q']

#  model data read
for j in range(10):
    vars()[model1[j]] = []
    for i in range(len(ddate1)):
        datazanshi = []
        datazanshi = np.loadtxt(open(model1[j] + '-' + ddate1[i]+'_datachem.csv', 'rb'),
                                delimiter=',', skiprows=0)
        if i < 1:
            vars()[model1[j]] = datazanshi
        else:
            vars()[model1[j]] = np.concatenate((vars()[model1[j]], datazanshi), axis=0)

#  model data size -> si
modeldata = vars()[model1[0]]
si = modeldata.shape
leftspace = 0.00    # 设置窗口左内边距
rightspace = 0.997  # 设置窗口右内边距
topspace = 0.99

a=[0,1,2,3,4]

# plot mean T2
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((T_obs_jan[0:192], T_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)


for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 2], modeldata[216:336, 2]),axis=0)
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm - 273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
    print(Tm - 273.15)
    print(rmse(Tm-273.15,Tmobs))
    print(ioa(Tm-273.15,Tmobs))
    r2, p = stats.pearsonr(Tm-273.15, Tmobs)
    print(r2, p)
ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('Tem (℃)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

# plot mean TSK
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 16], modeldata[216:336, 16]),axis=0)  # 16 TSK
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm - 273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('TSK (℃)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距


# plot mean diff HFX
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
ad=[1,2,3,4]
var=5
for i in np.array(ad):
    modeldata0=vars()[model1[0]]
    T0 = np.concatenate((modeldata0[0:192, var], modeldata0[216:336, var]),axis=0)  # 5 HFX
    Tm0 = meanvar(T0)
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, var], modeldata[216:336, var]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm-Tm0, linetype[i], color=colorline[i], label=model2diff[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(np.arange(-25,26,5),fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('HFX-diff (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

# plot mean diff HFX
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#d=[1,2]
var=13
for i in np.array(ad):
    modeldata0=vars()[model1[0]]
    T0 = np.concatenate((modeldata0[0:192, var], modeldata0[216:336, var]),axis=0)  # 5 HFX
    Tm0 = meanvar(T0)
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, var], modeldata[216:336, var]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm-Tm0, linetype[i], color=colorline[i], label=model2diff[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(np.arange(-25,26,5),fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('LH-diff (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

'''''
# plot mean diff SD
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#ad=[1,2]
var=14
for i in np.array(ad):
    modeldata0=vars()[model1[0]]
    T0 = np.concatenate((modeldata0[0:192, var], modeldata0[216:336, var]),axis=0)  # 5 HFX
    Tm0 = meanvar(T0)
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, var], modeldata[216:336, var]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm-Tm0, linetype[i], color=colorline[i], label=model2diff[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(np.arange(-25,26,5),fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('SD-diff (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
'''''

# plot mean diff LD
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#ad = [1, 2]
var = 17
for i in np.array(ad):
    modeldata0=vars()[model1[0]]
    TTSK0 = np.concatenate((modeldata0[0:192, 16], modeldata0[216:336, 16]),axis=0)  # 16 TSK
    LU0 = emiss*SIGMA*TTSK0**4
    LUm0 = meanvar(LU0)
    T0 = np.concatenate((modeldata0[0:192, var], modeldata0[216:336, var]),axis=0)  # 5 HFX
    Tm0 = meanvar(T0)
    to0 = LUm0 - Tm0
    modeldata=vars()[model1[i]]
    TTSK = np.concatenate((modeldata[0:192, 16], modeldata[216:336, 16]),axis=0)  # 16 TSK
    LU = emiss*SIGMA*TTSK**4
    LUm = meanvar(LU)
    T = np.concatenate((modeldata[0:192, var], modeldata[216:336, var]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    to = LUm - Tm
    ax.plot(list(range(0,24)), to0-to, linetype[i], color=colorline[i], label=model2diff[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
plt.ylim(-25, 25)
plt.yticks(np.arange(-25, 26, 5), fontsize=15)
plt.xticks(np.arange(0, 25, 6), fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('NL-diff (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
'''''
# plot mean LU
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 16], modeldata[216:336, 16]),axis=0)  # 16 TSK
    LU = emiss*SIGMA*T**4
    LUm = meanvar(LU)
    ax.plot(list(range(0,24)), LUm, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('LU (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
'''
# plot mean diff HFX
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#ad=[1,2]
var=19
for i in np.array(ad):
    modeldata0=vars()[model1[0]]
    T0 = np.concatenate((modeldata0[0:192, var], modeldata0[216:336, var]),axis=0)  # 5 HFX
    Tm0 = meanvar(T0)
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, var], modeldata[216:336, var]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm-Tm0, linetype[i], color=colorline[i], label=model2diff[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(np.arange(-25,26,5),fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('G-diff (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
# meanheatdiff(5,linetype,colorline,model1,model2,'HFX-diff (W/m^2)')
# meanheatdiff(13,linetype,colorline,model1,model2,'LH-diff (W/m^2)')
# meanheatdiff(14,linetype,colorline,model1,model2,'SD-diff (W/m^2)')
# meanheatdiff(17,linetype,colorline,model1,model2,'LD-diff (W/m^2)')
# meanheatdiff(19,linetype,colorline,model1,model2,'G-diff (W/m^2)')
'''''
# plot mean HFX
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 5], modeldata[216:336, 5]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('HFX (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
'''''
'''''
# plot mean reduation downward
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 14], modeldata[216:336, 14]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    #T1 = np.concatenate((modeldata[0:192, 17], modeldata[216:336, 17]),axis=0)  # 5 HFX
    #Tm1 = meanvar(T1)
    ax.plot(list(range(0,24)), Tm, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
#ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('Downward short wave radiation (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

# plot mean reduation downward
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 17], modeldata[216:336, 17]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
#ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('Downward long wave radiation (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

# plot mean LH
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
Tobs = np.concatenate((TSK_obs_jan[0:192], TSK_obs_jan2[0:120]),axis=0)
Tmobs = meanvar(Tobs)
#a=[0,6]
for i in np.array(a):
    modeldata=vars()[model1[i]]
    T = np.concatenate((modeldata[0:192, 13], modeldata[216:336, 13]),axis=0)  # 5 HFX
    Tm = meanvar(T)
    ax.plot(list(range(0,24)), Tm, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
#ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
#ax.scatter(list(range(0,24)), Tmobs, color='k', s=20, label='T_obs')
plt.xlim(0, 24)
#plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('LH (W/m^2)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距
'''''

# plot mean T_level1
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
for i in np.array(a):
    modeldata=vars()[model1[i]]
    vardata = vars()[model1[i] + '_' + varname[1]]
    T_1 = avelevel(vardata, modeldata,p_top,eta_levels[1],1)
    Tm = meanvar(T_1)
    ax.plot(list(range(0,24)), Tm-273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
plt.xlim(0, 24)
plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('T_level1 (℃)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

# calculate temperature profile
for i in np.array(a):
    modeldata = vars()[model1[i]]
    vardata = vars()[model1[i] + '_' + varname[1]]
    vars()[model1[i] + '_t_pf'] = np.zeros((47,312))  # chu shi hua
    for j in np.arange(47):
        vars()[model1[i] + '_t_pf'][j, :] = avelevel(vardata, modeldata, p_top, eta_levels[j+1], j+1)

# plot mean theta_level1
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
for i in np.array(a):
    modeldata=vars()[model1[i]]
    vardata = vars()[model1[i] + '_' + varname[1]]
    T_1 = np.concatenate((vardata[0:192, 0], vardata[216:336, 0]),axis=0)
    Tm = meanvar(T_1)
    ax.plot(list(range(0,24)), Tm-273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
plt.xlim(0, 24)
plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('theta_level1 (℃)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距

'''''
# plot mean TStheta
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
for i in np.array(a):
    modeldata=vars()[model1[i]]
    vardata = vars()[model1[i] + '_' + varname[1]]
    T_0 = np.concatenate((vardata[0:192, 0], vardata[216:336, 0]),axis=0)
    Tm0 = meanvar(T_0)
    T_1 = tstheta(modeldata)
    Tm = meanvar(T_1)
    ax.plot(list(range(0,24)), Tm-Tm0, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)
ax.plot(list(range(0,24)), Tm-Tm, '--', color='k', linewidth=2.5)
plt.xlim(0, 24)
plt.ylim(-10, 10)
plt.yticks(fontsize=15)
plt.xticks(np.arange(0,25,6),fontsize=15)
plt.xlabel('Local Standard Time', fontsize=15)
plt.ylabel('TStheta (℃)', fontsize=15)
leg = ax.legend(loc='upper left', frameon=False)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
fig.subplots_adjust(left=0.13)  # 设置窗口左内边距
fig.subplots_adjust(right=0.97)  # 设置窗口右内边距
fig.subplots_adjust(top=0.97)  # 设置窗口右内边距


#  plot t2
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
a = [0,6]#[0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 2] - 273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 2] - 273.15, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2

ax1.scatter(x1, T_obs_jan[0:len(x1)], color='k', s=20, label='T_obs')  # plot real data
ax2.scatter(x2, T_obs_jan2[0:len(x2)], color='k', s=20, label='T_obs')  # plot real data

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('T2 (℃)', -10, 11)

psta1 = 72  # 11
pend1 = 96   # 11 end
psta2 = 120   # 13
pend2 = 192   # 15 end
psta3 = 264   # 22
pend3 = 336   # 24 end
###################################################################
#### kzmin/kz #####################################################
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata0 = []
modeldata3 = []
modeldata4 = []
modeldata6 = []
vardata3 = vars()[model1[3] + '_' + var[0]]  # eval()字符串转变量名
vardata3 = vardata3.transpose(1, 0)
vardata4 = vars()[model1[4] + '_' + var[0]]  # eval()字符串转变量名
vardata4 = vardata4.transpose(1, 0)
modeldata0 = vars()[model1[0]]
modeldata3 = vars()[model1[4]]
modeldata4 = vars()[model1[5]]
modeldata6 = vars()[model1[6]]
t2p0 = lk(modeldata0,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p3 = lk(modeldata3,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p4 = lk(modeldata4,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p6 = lk(modeldata6,psta1,pend1,psta2,pend2,psta3,pend3,2)
box((KZM[4]/vardata4[1, 0:192]), modeldata4[0:192, 2]-modeldata3[0:192, 2], 'KZMIN/KZ', 'T2(℃)_Simulation_DIFF', -3, 3, 'b', 'AC0.01', ax, 'x', 50)
kzd=KZM[4]/vardata4[1, 0:192]
t2diff = modeldata4[0:192, 2]-modeldata3[0:192, 2]
print(kzd[np.where(kzd<0.95)])
fig, ax = plt.subplots(figsize=(6, 6))
ax.boxplot([t2diff[np.where(kzd<0.95)], t2diff[np.where(kzd>0.95)]], 0, 'rs', 1, whis=(0, 100), showmeans=True,
           meanline=True, labels=['kzmin/kz<0.95','kzmin/kz>0.95'])
ax.yaxis.grid(True)
ax.plot()
plt.ylim(-3, 3)
# box(T_obs_janp, t2p1 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC0.2', ax, 'x', 50)
# box(T_obs_janp, t2p2 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'g', 'AC0.5', ax, '^', 50)
#box(T_obs_janp, t2p6 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC_new', ax, 'o', 50)
# box(T_obs_janp, np.concatenate((modeldata6[psta1:pend1, 2], modeldata6[psta2:pend2, 2]), axis=0)-273.15, 'T2-OBS-Pollution', 'SIM', -10, 10, 'yellow', 'ACf', ax, 'o', 50)

################################################################


fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata0 = []
modeldata1 = []
modeldata2 = []
modeldata6 = []
modeldata0 = vars()[model1[0]]
modeldata1 = vars()[model1[1]]
modeldata2 = vars()[model1[2]]
modeldata6 = vars()[model1[6]]
t2p0 = lk(modeldata0,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p1 = lk(modeldata1,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p2 = lk(modeldata2,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p6 = lk(modeldata6,psta1,pend1,psta2,pend2,psta3,pend3,2)
box(T_obs_janp, t2p0 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'b', 'AC0.01', ax, 'x', 50)
# box(T_obs_janp, t2p1 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC0.2', ax, 'x', 50)
# box(T_obs_janp, t2p2 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'g', 'AC0.5', ax, '^', 50)
box(T_obs_janp, t2p6 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC_new', ax, 'o', 50)
# box(T_obs_janp, np.concatenate((modeldata6[psta1:pend1, 2], modeldata6[psta2:pend2, 2]), axis=0)-273.15, 'T2-OBS-Pollution', 'SIM', -10, 10, 'yellow', 'ACf', ax, 'o', 50)


csta1 = 0  # 8
cend1 = 72  # 10 end
csta2 = 96  # 12
cend2 = 120  # 12 end
csta3 = 216  # 20
cend3 = 264  # 21 end
fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata0 = []
modeldata1 = []
modeldata2 = []
modeldata0 = vars()[model1[0]]
modeldata1 = vars()[model1[1]]
modeldata2 = vars()[model1[2]]
modeldata6 = vars()[model1[6]]
t2c0 = lk(modeldata0,csta1,cend1,csta2,cend2,csta3,cend3,2)
t2c1 = lk(modeldata1,csta1,cend1,csta2,cend2,csta3,cend3,2)
t2c2 = lk(modeldata2,csta1,cend1,csta2,cend2,csta3,cend3,2)
t2c6 = lk(modeldata6,csta1,cend1,csta2,cend2,csta3,cend3,2)
box(T_obs_janc, t2c0 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'b', 'AC0.01', ax, 'x', 50)
# box(T_obs_janc, t2c1 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC0.2', ax, 'x', 50)
# box(T_obs_janc, t2c2 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'g', 'AC0.5', ax, '^', 50)
box(T_obs_janc, t2c6 - 273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'AC_new', ax, 'o', 50)
# box(T_obs_janc, np.concatenate((modeldata6[csta1:cend1, 2], modeldata6[csta2:cend2, 2]), axis=0)-273.15, 'T2-OBS-Clean', 'SIM', -10, 10, 'yellow', 'ACf', ax, 'o', 50)


fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata3 = []
modeldata4 = []
modeldata5 = []
modeldata6 = []
modeldata3 = vars()[model1[3]]
modeldata4 = vars()[model1[4]]
modeldata5 = vars()[model1[5]]
modeldata6 = vars()[model1[6]]
t2p3 = lk(modeldata3,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p4 = lk(modeldata4,psta1,pend1,psta2,pend2,psta3,pend3,2)
t2p5 = lk(modeldata5,psta1,pend1,psta2,pend2,psta3,pend3,2)
box(T_obs_janp, t2p3-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'b', 'BL0.1', ax, 'x', 50)
box(T_obs_janp, t2p4-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'BL0.2', ax, 'x', 50)
box(T_obs_janp, t2p5-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'g', 'BL0.5', ax, '^', 50)
# box(T_obs_janp, modeldata6[120:192, 2]-273.15, 'T-OBSP', 'SIM', -10, 10, 'yellow', 'ACf', ax, 'o', 50)

fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata3 = []
modeldata4 = []
modeldata5 = []
modeldata6 = []
modeldata3 = vars()[model1[3]]
modeldata4 = vars()[model1[4]]
modeldata5 = vars()[model1[5]]
modeldata6 = vars()[model1[6]]
t2c3 = lk(modeldata3,csta1,cend1,csta2,cend2,csta3,cend3,2)
t2c4 = lk(modeldata4,csta1,cend1,csta2,cend2,csta3,cend3,2)
t2c5 = lk(modeldata5,csta1,cend1,csta2,cend2,csta3,cend3,2)
box(T_obs_janc, t2c3-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'b', 'BL0.1', ax, 'x', 50)
box(T_obs_janc, t2c4-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'r', 'BL0.2', ax, 'x', 50)
box(T_obs_janc, t2c5-273.15, 'T2(℃)_Observation', 'T2(℃)_Simulation', -10, 10, 'g', 'BL0.5', ax, '^', 50)
# box(T_obs_janc, modeldata6[0:120, 2]-273.15, 'T-OBSC', 'SIM', -10, 10, 'yellow', 'ACf', ax, 'o', 50)

# print(rmse(T_obs_janp, modeldata0[120:192, 2]-273.5))   # test
# aaaa = Decimal(12.2225).quantize(Decimal('0.000'))
# print(aaaa)
# bbbb = round(12.2225, 3)
# print(bbbb)

#  plot RH
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距

for i in range(0, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    rh2_sc = rh2(modeldata)
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, rh2_sc[0:192], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, rh2_sc[216:336], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2

ax1.scatter(x1, rh_obs_jan[0:len(x1)], color='k', s=20, label='T_obs')  # plot real data
ax2.scatter(x2, rh_obs_jan2[0:len(x2)], color='k', s=20, label='T_obs')  # plot real data

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)
s_plo_sizeset('RH2 (%)', 0, 100)

fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata0 = []
modeldata1 = []
modeldata2 = []
modeldata6 = []
modeldata0 = rh2(vars()[model1[0]])
modeldata1 = rh2(vars()[model1[1]])
modeldata2 = rh2(vars()[model1[2]])
modeldata6 = rh2(vars()[model1[6]])
rh2p0 = lkrh(modeldata0,psta1,pend1,psta2,pend2,psta3,pend3)
rh2p1 = lkrh(modeldata1,psta1,pend1,psta2,pend2,psta3,pend3)
rh2p2 = lkrh(modeldata2,psta1,pend1,psta2,pend2,psta3,pend3)
rh2p6 = lkrh(modeldata6,psta1,pend1,psta2,pend2,psta3,pend3)
box(rh_obs_janp, rh2p0, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'b', 'AC0.01', ax, 'x', 50)
# box(rh_obs_janp, rh2p1, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'AC0.2', ax, 'x', 50)
# box(rh_obs_janp, rh2p2, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'g', 'AC0.5', ax, '^', 50)
box(rh_obs_janp, rh2p6, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'AC_new', ax, 'o', 50)
#box(rh_obs_janp, np.concatenate((modeldata6[psta1:pend1], modeldata6[psta2:pend2]), axis=0), 'RH-OBS-Pollution', 'SIM', 0, 100, 'yellow', 'ACf', ax, 'o', 50)
print(rh_obs_janp)

fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata0 = []
modeldata1 = []
modeldata2 = []
modeldata6 = []
modeldata0 = rh2(vars()[model1[0]])
modeldata1 = rh2(vars()[model1[1]])
modeldata2 = rh2(vars()[model1[2]])
modeldata6 = rh2(vars()[model1[6]])
rh2c0 = lkrh(modeldata0,csta1,cend1,csta2,cend2,csta3,cend3)
rh2c1 = lkrh(modeldata1,csta1,cend1,csta2,cend2,csta3,cend3)
rh2c2 = lkrh(modeldata2,csta1,cend1,csta2,cend2,csta3,cend3)
rh2c6 = lkrh(modeldata6,csta1,cend1,csta2,cend2,csta3,cend3)
box(rh_obs_janc, rh2c0, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'b', 'AC0.01', ax, 'x', 50)
# box(rh_obs_janc, rh2c1, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'AC0.2', ax, 'x', 50)
# box(rh_obs_janc, rh2c2, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'g', 'AC0.5', ax, '^', 50)
box(rh_obs_janc, rh2c6, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'AC_new', ax, 'o', 50)
#box(rh_obs_janc, np.concatenate((modeldata6[csta1:cend1], modeldata6[csta2:cend2]), axis=0), 'RH-OBSC', 'SIM', 0, 100, 'yellow', 'ACf', ax, 'o', 50)

fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata3 = []
modeldata4 = []
modeldata5 = []
modeldata6 = []
modeldata3 = rh2(vars()[model1[3]])
modeldata4 = rh2(vars()[model1[4]])
modeldata5 = rh2(vars()[model1[5]])
modeldata6 = rh2(vars()[model1[6]])
rh2c3 = lkrh(modeldata3,csta1,cend1,csta2,cend2,csta3,cend3)
rh2c4 = lkrh(modeldata4,csta1,cend1,csta2,cend2,csta3,cend3)
rh2c5 = lkrh(modeldata5,csta1,cend1,csta2,cend2,csta3,cend3)
box(rh_obs_janc, rh2c3, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'b', 'BL0.01', ax, 'x', 50)
box(rh_obs_janc, rh2c4, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'BL0.2', ax, 'x', 50)
box(rh_obs_janc, rh2c5, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'g', 'BL0.5', ax, '^', 50)
# box(rh_obs_janc, modeldata6[0:120], 'RH-OBSC', 'SIM', 0, 100, 'yellow', 'ACf', ax, 'o', 50)


fig = plt.figure(figsize=(8, 8))
ax = plt.subplot()
modeldata3 = []
modeldata4 = []
modeldata5 = []
modeldata6 = []
modeldata3 = rh2(vars()[model1[3]])
modeldata4 = rh2(vars()[model1[4]])
modeldata5 = rh2(vars()[model1[5]])
modeldata6 = rh2(vars()[model1[6]])
rh2p3 = lkrh(modeldata3,psta1,pend1,psta2,pend2,psta3,pend3)
rh2p4 = lkrh(modeldata4,psta1,pend1,psta2,pend2,psta3,pend3)
rh2p5 = lkrh(modeldata5,psta1,pend1,psta2,pend2,psta3,pend3)
box(rh_obs_janp, rh2p3, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'b', 'BL0.01', ax, 'x', 50)
box(rh_obs_janp, rh2p4, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'r', 'BL0.2', ax, 'x', 50)
box(rh_obs_janp, rh2p5, 'RH2(%)_Observation', 'RH2(%)_Simulation', 0, 100, 'g', 'BL0.5', ax, '^', 50)
# box(rh_obs_janp, modeldata6[120:192], 'RH-OBSP', 'SIM', 0, 100, 'yellow', 'ACf', ax, 'o', 50)


#  plot ws
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口右内边距
for i in range(6):
    modeldata = []
    modeldata = vars()[model1[i]]
    ws, wd = wind(modeldata[:, 3], modeldata[:, 4])
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, ws[0:192], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, ws[216:336], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2

# ax1.scatter(x1, T_obs_jan[0:len(x1)], color='k', s=7, label='T_obs')  # plot real data
# ax2.scatter(x2, T_obs_jan2[0:len(x2)], color='k', s=7, label='T_obs')  # plot real data

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)

leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)
s_plo_sizeset('WS (m s$^{-1}$)', 0, 15)

#  plot wd
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口右内边距
for i in range(6):
    modeldata = []
    modeldata = vars()[model1[i]]
    ws, wd = wind(modeldata[:, 3], modeldata[:, 4])
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.scatter(x1, wd[0:192], marker=ptype[i], color=colorline[i], label=model2[i], s=50)  # [:,2]=T2
    ax2.scatter(x2, wd[216:336], marker=ptype[i], color=colorline[i], label=model2[i], s=50)  # [:,2]=T2

# ax1.scatter(x1, T_obs_jan[0:len(x1)], color='k', s=7, label='T_obs')  # plot real data
# ax2.scatter(x2, T_obs_jan2[0:len(x2)], color='k', s=7, label='T_obs')  # plot real data

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='center left', frameon=False, bbox_to_anchor=(0.0, 0.8), ncol=6)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)
s_plo_sizeset('WD', 0, 360)

#  plot HFX
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 5], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 5], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
ax1.plot(x1, modeldata[0:192, 16]-modeldata[0:192, 16], '--', color='k', linewidth=1)  # [:,2]=T2
ax2.plot(x2, modeldata[216:336, 16]-modeldata[216:336, 16], '--', color='k', linewidth=1)  # [:,2]=T2
# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('HFX (W/m2)', -100, 300)
'''
''''
#  plot SWDOWN
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 17], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 17], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('lWDOWN (W/m2)', -100, 700)
'''''
''''
#  plot RDOWN
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 14]+modeldata[0:192, 17], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 14]+modeldata[216:336, 17], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('RDOWN (W/m2)', -100, 700)

#  plot TSK-T2
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 16]-modeldata[0:192, 2], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 16]-modeldata[216:336, 2], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax1.plot(x1, modeldata[0:192, 16]-modeldata[0:192, 16], '--', color='k', linewidth=1)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 16]-modeldata[216:336, 16], '--', color='k', linewidth=1)  # [:,2]=T2

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('TSK-T2', -10, 10)

#  plot TSK
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 16]-273, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 16]-273, linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('TSK', -15, 15)

#  plot GRDFLX
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
fig.subplots_adjust(top=topspace)  # 设置窗口左内边距
#a = [0,1,2,3,4,5]
for i in np.array(a):  #range(2, 6, 1):
    modeldata = []
    modeldata = vars()[model1[i]]
    x1 = list(range(0, 192))
    x2 = list(range(0, 120))
    ax1.plot(x1, modeldata[0:192, 19], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
    ax2.plot(x2, modeldata[216:336, 19], linetype[i], color=colorline[i], label=model2[i], linewidth=3.5)  # [:,2]=T2
ax1.plot(x1, modeldata[0:192, 16]-modeldata[0:192, 16], '--', color='k', linewidth=1)  # [:,2]=T2
ax2.plot(x2, modeldata[216:336, 16]-modeldata[216:336, 16], '--', color='k', linewidth=1)  # [:,2]=T2

# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
leg = ax1.legend(loc='upper left', frameon=False, ncol=2)  # legend设置
ltext = leg.get_texts()
plt.setp(ltext, fontsize=25)   # legend字体大小
s_plo_sizeset('GRDFLX (W/m2)', -300, 300)


def sens(kz1, v1, kz2, v2):
    v1 = np.array(v1)
    v2 = np.array(v2)
    sen = (v1-v2)/(kz1-kz2)*kz2/v2
    return sen


# sensitive - T2
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
k1 = 0.5
k2 = 0.55
modeldata1 = []
modeldata1 = vars()[model1[2]]
modeldata2 = []
modeldata2 = vars()[model1[7]]
sst = sens(k1, modeldata1[:, 2], k2, modeldata2[:, 2])

x1 = list(range(0, 192))
x2 = list(range(0, 120))
i = 0
ax1.plot(x1, sst[0:192]*1000, ptype[i], color=colorline[i], label=model2[i], linewidth=2.0)  # [:,2]=T2
ax2.plot(x2, sst[216:336]*1000, ptype[i], color=colorline[i], label=model2[i], linewidth=2.0)  # [:,2]=T2
ax1.plot(x1, np.zeros(192), '--', color='k', linewidth=2.0)  # [:,2]=T2
ax2.plot(x2, np.zeros(120), '--', color='k', linewidth=2.0)  # [:,2]=T2
sstp = lkrh(sst,psta1,pend1,psta2,pend2,psta3,pend3)
sstc = lkrh(sst,csta1,cend1,csta2,cend2,csta3,cend3)
meanc = Decimal((sstc.mean())*1000).quantize(Decimal('0.000'))
meanp = Decimal((sstp.mean())*1000).quantize(Decimal('0.000'))
ax1.annotate('clean = ' + str(meanc), xy=(0.12, 0.84), xycoords='figure fraction', fontsize=12)
ax1.annotate('pollution = ' + str(meanp), xy=(0.12, 0.81), xycoords='figure fraction', fontsize=12)
# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
# leg = ax1.legend(loc='upper left', frameon=False)  # legend设置
# ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
s_plo_sizeset('T2-sensitivity x1000', -10, 10)

# sensitive - RH2
fig = plt.figure(constrained_layout=False, figsize=(21, 5.5))
gs = GridSpec(1, 14, figure=fig)
ax1 = fig.add_subplot(gs[0, 1:9])
# identical to ax1 = plt.subplot(gs.new_subplotspec((0, 0), colspan=3))
ax2 = fig.add_subplot(gs[0, 9:14])
fig.subplots_adjust(left=leftspace)  # 设置窗口左内边距
fig.subplots_adjust(right=rightspace)  # 设置窗口右内边距
k1 = 0.5
k2 = 0.55
modeldata1 = []
modeldata1 = vars()[model1[2]]
modeldata2 = []
modeldata2 = vars()[model1[7]]
ssrh = sens(k1, rh2(modeldata1), k2, rh2(modeldata2))
x1 = list(range(0, 192))
x2 = list(range(0, 120))
i = 0
ax1.plot(x1, ssrh[0:192]*100, ptype[i], color=colorline[i], label=model2[i], linewidth=2.0)  # [:,2]=T2
ax2.plot(x2, ssrh[216:336]*100, ptype[i], color=colorline[i], label=model2[i], linewidth=2.0)  # [:,2]=T2
ax1.plot(x1, np.zeros(192), '--', color='k', linewidth=2.0)  # [:,2]=T2
ax2.plot(x2, np.zeros(120), '--', color='k', linewidth=2.0)  # [:,2]=T2
ssrhp = lkrh(ssrh,psta1,pend1,psta2,pend2,psta3,pend3)
ssrhc = lkrh(ssrh,csta1,cend1,csta2,cend2,csta3,cend3)
meanc = Decimal((ssrhc.mean())*100).quantize(Decimal('0.000'))
meanp = Decimal((ssrhp.mean())*100).quantize(Decimal('0.000'))
ax1.annotate('clean = ' + str(meanc), xy=(0.12, 0.84), xycoords='figure fraction', fontsize=12)
ax1.annotate('pollution = ' + str(meanp), xy=(0.12, 0.81), xycoords='figure fraction', fontsize=12)
# plotline('T2', 'T2', 'T2', -15, 10, 'upper left', 0)
# leg = ax1.legend(loc='upper left', frameon=False)  # legend设置
# ltext = leg.get_texts()
plt.setp(ltext, fontsize=15)
s_plo_sizeset('RH-sensitivity x100', -50, 50)
'''
plt.show()
