请选择 进入手机版 | 继续访问电脑版
爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 7013|回复: 8

[作图] 关于EOF画海温异常的疑问

[复制链接]

新浪微博达人勋

发表于 2022-3-3 22:02:48 | 显示全部楼层 |阅读模式

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

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

x

热带太平洋冬季海温EOF

热带太平洋冬季海温EOF

我画了海温的EOF分解,但是明显图a应该海温是正的,问题出在哪里?
下面是我的代码:
begin
f = addfile("sst.DJF.nc", "r")
sst = f->sst
lat = f->latitude
sst = lonFlip(sst)


sst = dtrend_msg_n(sst&time, sst, False, False, 0)   ;--去趋势
ssta = dim_rmvmean_n_Wrap(sst, 0)
ny = dimsizes(sst&time)

; ;--------面积权重-------------
; rad = get_d2r("f")
; wgt = sqrt(cos(lat*rad))

; ;nlat = dimsizes(sst&latitude)
; sstw = sst
; sstw = sst*conform(sst, wgt, 1)
; printVarSummary(sstw)

TPsst = ssta(:,{-20:25},{100:280})

;===================
;EOF
;===================
neof    = 2        ; number of EOFs

optEOF    = True
optEOF@jopt = 0
eof     = eofunc_n_Wrap(TPsst,neof,optEOF,0)
;eof = smth9_Wrap(eof, 0.5, 0.5, True)


optETS    = False
eof_ts  = eofunc_ts_n_Wrap(TPsst, eof, optETS,0)
printVarSummary(eof_ts)

do ne=0,neof-1
     eof(ne,:,:) = eof(ne,:,:)*sqrt(eof@eval(ne))
end do
printVarSummary(eof)
eof_ts = dim_standardize_Wrap(eof_ts, 0)

pcvar = eof@pcvar
prinfo = True
sig_ev = eofunc_north(eof@eval, ny, prinfo)

wks = gsn_open_wks("eps","eof")      
plot = new(4,graphic)               

res                      = True         
res@gsnDraw              = False      
res@gsnFrame             = False        
res@gsnStringFont               = "times-roman"

res@cnFillPalette        = "temp_19lev"   ; choose colormap
res@cnFillOn                    = True
res@cnLinesOn                   = False             ;--不画等值线
res@cnLineLabelsOn              = False             ;--不显示等值线数值

res@lbLabelBarOn                = False             ;--关闭lbar box

res@cnLevelSelectionMode        = "ManualLevels"
res@cnMaxLevelValF              = 1
res@cnMinLevelValF              = -1
res@cnLevelSpacingF             = 0.1
;------------------------------设置刻度线内容-----------------------------------------
res@tmXBLabelFontHeightF        = 0.02
res@tmYLLabelFontHeightF        = 0.02
res@tmBorderThicknessF          = 1                 ;--设置边框厚度
res@tmXBMajorThicknessF         = 1                 ;--主刻度线厚度
res@tmYLMajorThicknessF         = 1
res@tmXBMinorThicknessF         = 1     
res@tmYLMinorThicknessF         = 1
res@tmXTOn                      = False             ;--不显示x上刻度
res@tmYROn                      = False             ;--不显示y右刻度
;res@tmYLLabels                  = (/"20~S~o~N~N","30~S~o~N~N","40~S~o~N~N","50~S~o~N~N"/)  
;res@tmXBLabels                  = (/"100~S~o~N~E","110~S~o~N~E","120~S~o~N~E","130~S~o~N~E"/)  
;res@pmTickMarkDisplayMode       = "Always"
res@vpHeightF                   = 0.45
res@vpWidthF                    = 0.6


res@mpMaxLonF                   = 280
res@mpMinLonF                   = 100
res@mpMaxLatF                   = 25
res@mpMinLatF                   = -20
res@mpCenterLonF                = 190

res@gsnAddCyclic                = False              ;--是否要添加周期点使纬圈完整, 不是全球数据设为 False

res@mpOutlineOn                 = True
res@mpDataBaseVersion           = "Ncarg4_1"
res@mpDataSetName               = "Earth..4"

res@mpGeophysicalLineThicknessF = 1.0               ;--这两行是为了加粗边界和国界线
res@mpNationalLineThicknessF    = 1.0

;res@mpOceanFillColor        = 0
;res@mpInlandWaterFillColor  = 0
res@mpLandFillColor         = 1
res@gsnLeftStringFontHeightF  = 0.02

;symMinMaxPlt(eof, 16, False, res)
lbal = (/"(a)","(b)","(c)","(d)"/)
do n=0,neof-1
   res@gsnLeftString  = lbal(n)+" EOF"+(n+1)+"("+sprintf("%5.1f", eof@pcvar(n)) +"%)"
   plot(n) = gsn_csm_contour_map(wks,eof(n,:,:),res)
end do

;==========================
;时间序列
;==========================

yy = ispan(1951, 2010, 1)
rts             = True
rts@gsnDraw     = False
rts@gsnFrame    = False
rts@gsnScale    = True
rts@tmBorderThicknessF          = 2                 ;--设置边框厚度
rts@tmXBMajorThicknessF         = 2                 ;--主刻度线厚度
rts@tmYLMajorThicknessF         = 2
rts@tmXBMinorThicknessF         = 2     
rts@tmYLMinorThicknessF         = 2
rts@tmXTOn                      = False             ;--不显示x上刻度
rts@tmYROn                      = False             ;--不显示y右刻度
rts@gsnStringFont               = "times-roman"


rts@vpHeightF = 0.35        ; Changes the aspect ratio
rts@vpWidthF  = 0.55
rts@tiYAxisString               = ""
rts@trXMinF                     = min(yy)-1
rts@trXMaxF                     = max(yy)+1
; rts@trYMinF                     = min(yy)-1
; rts@trYMaxF                     = max(yy)+1

rts@gsnYRefLine           = 0.              ; reference line   
rts@gsnXYBarChart         = True
rts@gsnAboveYRefLineColor = "red"           ; above ref line fill red
rts@gsnBelowYRefLineColor = "blue"          ; below ref line fill blue
rts@gsnLeftStringFontHeightF  = 0.02
rts@tiXAxisString             = "Time"

do i = 0, 1
  rts@gsnLeftString = lbal(i+2)+ " PC"+(i+1)
  plot(i+2) = gsn_csm_xy(wks, yy, eof_ts(i,:), rts)
end do

gres    = True
gres@gsLineThicknessF               = 3
gres@gsLineColor                    = "black"
; lnid1 = gsn_add_polyline(wks, plot(0), (/80,100/), (/-15,-15/), gres)
; lnid2 = gsn_add_polyline(wks, plot(0), (/80,100/), (/-25,-25/), gres)
; lnid3 = gsn_add_polyline(wks, plot(0), (/80,80/), (/-25,-15/), gres)
; lnid4 = gsn_add_polyline(wks, plot(0), (/100,100/), (/-25,-15/), gres)

; lnid5 = gsn_add_polyline(wks, plot(0), (/45,75/), (/-30,-30/), gres)
; lnid6 = gsn_add_polyline(wks, plot(0), (/45,75/), (/-45,-45/), gres)
; lnid7 = gsn_add_polyline(wks, plot(0), (/45,45/), (/-45,-30/), gres)
; lnid8 = gsn_add_polyline(wks, plot(0), (/75,75/), (/-45,-30/), gres)


pnlres = True
pnlres@amJust = "TopLeft"
pnlres@gsnPanelLabelBar           = True
pnlres@pmLabelBarWidthF           = 0.5

gsn_panel(wks,plot,(/2,2/),pnlres)        ; draw all 'neof' as one plot
; system("rm -f "+"eof_ts.nc")
; outfile = addfile("eof_ts.nc", "c")
; outfile->eof_ts = eof_ts

end

密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2022-3-26 21:47:56 | 显示全部楼层
本人暂时给出的解决方案是根据Nino3.4指数判断海温距平是否为负值,若为负,则eof和eof_ts都×-1
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

发表于 2022-3-4 10:29:03 | 显示全部楼层
本帖最后由 井中月 于 2022-3-4 10:33 编辑

空间型和时间序列都乘-1不就行了。好好学习学习EOF的基本原理吧
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2022-3-4 13:11:53 | 显示全部楼层
井中月 发表于 2022-3-4 10:29
空间型和时间序列都乘-1不就行了。好好学习学习EOF的基本原理吧

单个模式可以这么弄,但是几十个模式要做EOF相关还要一个一个确定有点不现实,我如何可以通过不用画图方式来保持都是正的呢
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2022-3-26 15:09:58 | 显示全部楼层
循环,判断你想要的关键区里eof那个模态是正是负,如果是负就把eof和pc都乘-1
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2023-3-20 17:16:00 | 显示全部楼层
谢谢分享代码
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2023-3-31 15:30:49 | 显示全部楼层
给nino3.4区域的EOF算一个区域平均,如果小于0,eof和pc都*-1
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2023-3-31 16:03:38 | 显示全部楼层
ncl的玄学
密码修改失败请联系微信:mofangbao
回复

使用道具 举报

新浪微博达人勋

发表于 2023-4-4 23:38:03 来自手机 | 显示全部楼层
EOF得到的模态结果是对称的,得到相反的模态乘-1就行了
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

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

本版积分规则

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

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

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