求助，MI画假相当位温

lixiliang2008

下边这个计算假相当位温有错误，请指教！数据是ERA5的。
f = addfile('E:/20220211-13.nc')
levels = f['level'][::-1]
t = f['t'][::-1,::-1,'40.95','115.27']
q = f['q'][::-1,::-1,'40.95','115.27']
udata = f['u'][::-1,::-1,'40.95','115.27']
vdata = f['v'][::-1,::-1,'40.95','115.27']

prs = levels
es = 6.1078*exp(17.2693882*(t-273.16)/(t-35.86))
e = prs*q/(0.62197+q)+1e-10
tlcl = 55.0+2840.0/(3.5*log(t)-log(e)-4.805)
theta = t*pow((1000/prs),(0.2854*(1.0-0.28*q)))
eqt = theta*exp(((3376./tlcl)-2.54)*q*(1.0+0.81*q))
eqt = eqt.T

lev1 = eqt.dimvalue(0)
lev1=lev1[::-1]
lev2 = meteo.pressure_to_height_std(lev1)
lev2 = lev2[:]/1000
eqt.setdimvalue(0, lev2)

u = udata.T
v = vdata.T
lev3 = u.dimvalue(0)
lev3 = lev3[::-1]
lev4 = meteo.pressure_to_height_std(lev3)
lev4 = lev4[:]/1000
u.setdimvalue(0, lev4)
lev5 = v.dimvalue(0)
lev5 = lev5[::-1]
lev6 = meteo.pressure_to_height_std(lev5)
lev6 = lev6[:]/1000
v.setdimvalue(0, lev6)

layer1 = barbs(u , v, color = 'k')
layer2 = contour(eqt, 38, color='r')
clabel(layer2)

xaxis(axistype='normal')
xlim(f.gettime(10), f.gettime(52))
ylim(lev4.min(),6)
xticks([f.gettime(10), f.gettime(13), f.gettime(16), f.gettime(19), f.gettime(22), f.gettime(25), f.gettime(28), f.gettime(31), f.gettime(34), f.gettime(37), f.gettime(40), f.gettime(43), f.gettime(46), f.gettime(49), f.gettime(52)], ['13-21','13-18','13-15','13-12','13-09','13-06','13-03','13-00','12-21','12-18','12-15','12-12','12-09','12-06','12-03'])
levy = array([850, 825, 800, 775, 750, 700, 650, 600, 550, 500, 450, 400, 350, 300, 250, 225, 200])
yticks(eqt.dimvalue(0), levy)
xaxis(tickin=False,tickfontsize=20)
yaxis(tickin=False,tickfontsize=20)
xlabel('Time(dd-hh)',fontsize=20)
ylabel('Pressure(hPa)',fontsize=20)

lixiliang2008

需要在levels = f['level'][::-1]加一个0.0就可以了
改为levels = f['level'][::-1]+0.0
研究了半天

MeteoInfo

lixiliang2008 发表于 2022-7-14 16:10
需要在levels = f['level'][::-1]加一个0.0就可以了
改为levels = f['level'][::-1]+0.0
研究了半天

为了节省内存，数组的slice操作通常只是生成一个新的index，并不会对原始数据进行拷贝，这样的话在数组计算的函数中需要很小心处理，否则可能会从内存中取错数，也许是哪个函数的bug造成你之前的问题，+0.0会生成新的数组（会进行内存拷贝），这样从内存中取数就很容易，其实用 .copy() 函数就可以，比如  a = a.copy() 。当然如果你能定位出哪个函数有bug就更好了，我这里好进行修正和更新。