- 积分
- 1868
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2015-1-1
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
求助各位大神,我绘制水汽通量散度的t检验时,大值中心没有通过显著性检验,反而在图片的边边角角(右下角网格处代表通过信度检验)通过了检验,脚本标红处为计算t检验的部分。
自己看脚本也没有找出错误,纠结一个星期了QAQ,忘各位大神指点迷津,不知究竟是我水汽通量散度没计算对,还是我的t处理部分有问题,忘大神给与指点Orz。
非常感谢!
begin
f1=addfile("E:/TPV/first/data/strong_summer.nc","r")
f2=addfile("E:/TPV/first/data/weak_summer.nc","r")
air= short2flt(f1->t(:,7,:,:))
printVarSummary(air)
u=short2flt(f1->u(:,7,:,:))
v=short2flt(f1->v(:,7,:,:))
rhum=short2flt(f1->r(:,7,:,:))
lat=f1->latitude
lon=f1->longitude
es = 6.112*exp(17.67*(air-273.16)/(air-29.65))
qs = 0.62197*es*1000/(700-0.378*es)
q = qs*rhum/100
q = smth9(q, 0.50, -0.25, False)
copy_VarCoords(air,es)
copy_VarCoords(air,qs)
copy_VarCoords(air,q)
qu = new((/12,61,121/),float) ;gen ju zi liao jing wei du xiu gai
qv = new((/12,61,121/),float)
do i = 0,11
do j=0,60
do k=0,120
qu(i,j,k) = q(i,j,k)*u(i,j,k)/9.8
qv(i,j,k) = q(i,j,k)*v(i,j,k)/9.8
end do
end do
end do
copy_VarCoords(air,qu)
copy_VarCoords(air,qv)
; qu = smth9(qu, 0.50, -0.25, False)
; qv = smth9(qv, 0.50, -0.25, False)
vapord = new((/12,61,121/),float)
vapord=uv2dv_cfd(qu,qv,lat,lon,2) ; u,v ==> divergence
vapord = vapord*100000
copy_VarCoords(air,vapord)
printVarSummary(vapord)
printVarSummary(qu)
;>=========================================================<
air2= short2flt(f2->t(:,7,:,:))
u2=short2flt(f2->u(:,7,:,:))
v2=short2flt(f2->v(:,7,:,:))
rhum2=short2flt(f2->r(:,7,:,:))
lat=f2->latitude
lon=f2->longitude
es2 = 6.112*exp(17.67*(air2-273.16)/(air2-29.65))
qs2 = 0.62197*es2*1000/(700-0.378*es2)
q2 = qs2*rhum2/100
q2 = smth9(q2, 0.50, -0.25, False)
copy_VarCoords(air2,es2)
copy_VarCoords(air2,qs2)
copy_VarCoords(air2,q2)
qu2 = new((/15,61,121/),float)
qv2 = new((/15,61,121/),float)
do i = 0,14
do j=0,60
do k=0,120
qu2(i,j,k) = q2(i,j,k)*u2(i,j,k)/9.8
qv2(i,j,k) = q2(i,j,k)*v2(i,j,k)/9.8
end do
end do
end do
copy_VarCoords(air2,qu2)
copy_VarCoords(air2,qv2)
; qu = smth9(qu, 0.50, -0.25, False)
; qv = smth9(qv, 0.50, -0.25, False)
vapord2 = new((/15,61,121/),float)
vapord2=uv2dv_cfd(qu2,qv2,lat,lon,2) ; u,v ==> divergence
vapord2 = vapord2*100000
copy_VarCoords(air2,vapord2)
;printVarSummary(vapord2)
quavg = dim_avg_n_Wrap(qu, 0)
printVarSummary(quavg)
quavg2 = dim_avg_n_Wrap(qu2, 0)
quchazhi = quavg - quavg2
copy_VarMeta(quavg, quchazhi)
qvavg = dim_avg_n_Wrap(qv, 0)
printVarSummary(qvavg)
qvavg2 = dim_avg_n_Wrap(qv2, 0)
qvchazhi = qvavg - qvavg2
copy_VarMeta(qvavg, qvchazhi)
vapordavg = dim_avg_n_Wrap(vapord, 0)
vapordavg2 = dim_avg_n_Wrap(vapord2, 0)
vapordchazhi = vapordavg - vapordavg2
copy_VarMeta(vapordavg, vapordchazhi)
;printVarSummary(vapordchazhi)
t1 = 4
t2 = 5
wsgd1_ave=dim_avg_n_Wrap(vapord,0)
wsgd2_ave=dim_avg_n_Wrap(vapord2,0)
wsgd1_var=dim_variance_n_Wrap(vapord,0)
wsgd2_var=dim_variance_n_Wrap(vapord2,0)
printMinMax(wsgd1_var, False)
printMinMax(wsgd2_var, False)
iflag = False
prob= ttest(wsgd1_ave,wsgd1_var,t1,wsgd2_ave,wsgd2_var,t2,iflag,False)
copy_VarMeta(vapordavg,prob)
printVarSummary(prob)
wks = gsn_open_wks("png","E:/TPV/first/picture/1shuiqi_ttest_try")
res = True
res@tiMainString =""
res@gsnMaximize = False
res@gsnDraw = False
res@gsnFrame = False
res@gsnAddCyclic = False
;>--------------------------------------------<
; set for the map
;>--------------------------------------------<
res@mpMinLatF = 0.
res@mpMaxLatF = 60.
res@mpMinLonF = 40.
res@mpMaxLonF = 160.
res@pmTickMarkDisplayMode = "Always"
;>--------------------------------------------<
res@gsnLeftString = ""
res@gsnContourZeroLineThicknessF= 2 ;加粗
res@gsnContourNegLineDashPattern= 2;5等值线负值虚线
res@gsnRightString =""
res@cnInfoLabelOn = False
res@cnFillOn = True
res@cnLinesOn = False
res@cnLineLabelsOn = False
res@cnLevelSelectionMode = "ManualLevels"
res@cnMinLevelValF = -0.3 ; set min contour level
res@cnMaxLevelValF= 0 ; set max contour level
res@cnLevelSpacingF= 0.05
res@cnFillPalette = "gsdtol"
res@lbLabelBarOn = True
res@lbBoxMinorExtentF = 0.15
res@tmXBLabelFontHeightF = 0.02
res@tmXBLabelFontThicknessF = 0.04
res@tmYLLabelFontHeightF = 0.02
res@tmYLLabelFontThicknessF = 0.04
vres = True
vres@gsnDraw = False
vres@gsnFrame = False
vres@gsnAddCyclic = False
vres@gsnMaximize = True
vres@gsnRightString = ""
vres@gsnLeftString = ""
vres@vcRefMagnitudeF = 0.5
vres@vcRefAnnoString1 = "0.5"
vres@vcRefAnnoString2 = ""
vres@vcRefAnnoFontHeightF = 0.015
vres@vcRefLengthF = 0.03
vres@vcGlyphStyle = "LineArrow" ;WindBarb CurlyVector LineArrow
vres@vcMinDistanceF = 0.03
vres@vcLineArrowThicknessF = 1.5
vres@vcRefAnnoOrthogonalPosF = -1.36 ; move ref vector
;==========================================================================
pres = True
pres@gsnMaximize = True
pres@gsnDraw = False ; don't draw
pres@gsnFrame = False ; don't advance frame
pres@gsnLeftString = ""
pres@gsnRightString= ""
pres@gsnAddCyclic = False
pres@tiMainString = ""
pres@cnInfoLabelOn = False
pres@cnFillOn = True
pres@cnLinesOn = False
pres@cnSmoothingOn = True
pres@cnMonoFillPattern = False
pres@cnLineLabelsOn = False
pres@cnLevelSelectionMode = "ExplicitLevels"
pres@cnLevels = (/0.1/) ;; set to significance level
pres@cnFillPatterns = (/15,-1/)
pres@cnFillColors = (/1,0/)
pres@lbLabelBarOn = False
;======================================
plot2 = gsn_csm_contour(wks,prob,pres)
vplot = gsn_csm_vector(wks,quchazhi,qvchazhi,vres)
plot=gsn_csm_contour_map(wks,vapordchazhi,res)
overlay(plot, vplot)
overlay(plot, plot2)
;>============================================================<
draw(plot)
frame(wks)
end
|
-
灰色填色代表水汽通量散度,网格处为通过90的信度检验
|