爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 4730|回复: 1

[作图] 求各位大神指点迷津,ncl作图,大值中心未通过t检验,求助Orz

[复制链接]

新浪微博达人勋

发表于 2018-11-26 21:12:26 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

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的信度检验

灰色填色代表水汽通量散度,网格处为通过90的信度检验
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2018-11-26 22:59:10 | 显示全部楼层
没有看代码,仅发表一点儿意见供参考


算出来多大的水汽通量散度  和  通不通过检验 好像本来就没关系

各种检验是对统计模型适用性的检验,“假设变量分布满足XX分布时,XXX的分布是XXX的,这时XXX是XXX”
所以可以一致很低,甚至就是常量场,绝对能通过均值检验的那种

至于很高的地方……可能真的变化和波动太大了吧……

建议你拉出来两个代表点看下时间序列,画一下频率分布histogram看看
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表