爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 10706|回复: 2

相关显著性检验

[复制链接]

新浪微博达人勋

发表于 2021-4-6 17:46:28 | 显示全部楼层 |阅读模式
30金钱
做秋雨指数与850hpa温度的相关分析,秋雨指数处理成每年一个数,共计60年,一个是风场的显著性,一个是温度的显著性,显著性区域都成方块状,有人知道是为什么吗?
begin

  latS     = 0              
  latN     = 90
  lonL     = 0
  lonR     = 150
  siglvl   = 0.05
  ;;-----------------------------------------------------------------------
   f1      = addfile("pre2.nc", "r") ;降水

   pre     = f1->pre(:,{latS:latN},{lonL:lonR})
   index   = dim_avg_n_Wrap(pre, (/1,2/));指数,只剩时间维度,后面才能做相关; [time | 60]
   ;printVarSummary(index)

   ;;---------------------------------------------------------------------

   f2      = addfile("air.mean.nc", "r") ;(都是每年一个数)
   airtmp     = f2->air(:,{850},{latS:latN},{lonL:lonR})  ; [time | 60] x [lat | 21] x [lon | 37]
   printVarSummary(airtmp)
   dimx = dimsizes(airtmp)
   ntim = dimx(0)               

   ;;-----------------------------------------------------------------------
   mxlag = 0 ;同时相关
   rairtmp  = esccr(index, airtmp({lat|:},{lon|:},{time|:}), mxlag); [21] x [37] x [1]
   printVarSummary(rairtmp)
   rairtmp1 :=(/rairtmp(:,:,0)/)
   copy_VarMeta(airtmp(0,:,:), rairtmp1);把经纬度元数据赋值给rv1
   ;;---------------------------------------------------------------------------------
   n = ntim
   df = n-2
   t = rairtmp1*sqrt((n-2)/(1-rairtmp1^2)) ;t检验
   p = student_t(t, df) ;两侧student-t分布概率值
   printVarSummary(p)   ; [21] x [37]
   psig = 0.05
   copy_VarMeta(airtmp(0,:,:), p)


   ;;----------------------------------------------------------------------------
  wks = gsn_open_wks("png", "corr.airtmp")
  gsn_define_colormap(wks, "GMT_panoply")

  res = True
  res@gsnDraw        = False
  res@gsnFrame       = False
  res@gsnMaximize    = True
  res@gsnLeftString  = ""
  res@gsnRightString = ""
  res@pmTitleDisplayMode = "Always"
  res@pmTickMarkDisplayMode = "Always"

;;set map;;
  mpres                             = res
  mpres@mpDataSetName               = "Earth..4"
  mpres@mpDataBaseVersion           = "MediumRes"
  mpres@mpOutlineOn                 = True
  mpres@mpOutlineSpecifiers         = (/"China:states","Taiwan"/)
  mpres@mpGeophysicalLineThicknessF = 2
  mpres@mpNationalLineThicknessF    = 2
  mpres@mpFillDrawOrder             = "PostDraw"
  mpres@mpFillOn                    = False


;;set area;;
  mpres@mpMinLatF                   = 0
  mpres@mpMaxLatF                   = 90
  mpres@mpMinLonF                   = 0
  mpres@mpMaxLonF                   = 150
  map     = gsn_csm_map(wks,mpres)
;;--------------------------------------------------------------------------相关系数线
  cnres                             = res
  cnres@cnFillDrawOrder             = "PreDraw"
  cnres@cnFillOn                    = False
  cnres@cnLinesOn                   = True
  cnres@cnLineLabelsOn              = True
  cnres@cnLineThicknesses           = 3.     ;等值线粗细
  cnres@cnSmoothingOn               = True
  cnres@cnSmoothingDistanceF        = 0.001
  cnres@cnSmoothingTensionF         = -2
  cnres@cnInfoLabelOn               = False  ;右下方等值线信息标签
  cnres@gsnContourNegLineDashPattern = 0
  cnres@gsnContourPosLineDashPattern = 2

  contour = gsn_csm_contour(wks,rairtmp1,cnres)
  contour = ColorNegDashZeroPosContour(contour,"black","black","black")
;;-------------------------------------------------------------------------------------
  sres                             = True                                                           ; set up a second resource list
  sres@gsnDraw            = False                                                          ; do not draw the plot
  sres@gsnFrame           = False                                                          ; do not advance the frame
  sres@cnLineLabelsOn  = False                                                          ; do not use line labels
  sres@cnLinesOn          = False                                                          ; do not draw contour lines
  sres@cnInfoLabelOn    = False
  sres@cnFillOn              = False                                                           ; color fill
  sres@lbLabelBarOn      = False
  sres@cnInfoLabelOn    = False
  sres@gsnLeftString = ""
  sres@gsnRightString = ""

  sres@cnFillDotSizeF  = 0.005
;;---------------------------------------------------------------
  opt = True
  opt@gsnShadeFillType = "color"
  opt@gsnShadeMid     = "gray"

  Rxt_test = where(p.le.0.05, rairtmp1, rairtmp1@_FillValue)
  copy_VarMeta(rairtmp1(:,:), Rxt_test(:,:))
  printVarSummary(Rxt_test)

  Rxt_max = max(Rxt_test)
  Rxt_min = min(Rxt_test)

  plot1                = gsn_csm_contour(wks,Rxt_test,sres)
  plot1                = gsn_contour_shade(plot1,Rxt_min,Rxt_max,opt)

  overlay(map,plot1)
  overlay(map,contour)
  draw(map)
  frame(wks)
end


corr.rainfall & H850uv.png
corr.airtmp.png
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2021-4-12 00:12:09 来自手机 | 显示全部楼层
可能因为你的数据是格点的,建议用打点的形式会好看一些么,但会和风场重叠,所以最好打点用其他颜色
密码修改失败请联系微信:mofangbao
回复

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2021-4-16 21:54:16 | 显示全部楼层
杨松 发表于 2021-4-12 00:12
可能因为你的数据是格点的,建议用打点的形式会好看一些么,但会和风场重叠,所以最好打点用其他颜色

后面解决了,是因为gsn_contour_shade这个函数,用等值线方法画就好了
密码修改失败请联系微信:mofangbao
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

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

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

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