爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

搜索
查看: 11975|回复: 3

[作图] NCL求时间平均和距平相关问题

[复制链接]
发表于 2016-10-25 17:21:20 | 显示全部楼层 |阅读模式

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

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

x
大家好,我想做1990-2015年间的500hPa高度场的时间平均和距平,然后再对距平场做EOF分解,在官网上看了anom.ncl的和相关平均函数,
但是还是不太明白,一个下午都在调试这个错误,也没搞懂,来向大家请教这个问题:
dim_avg(x,0)函数是对变量x的最右边维进行平均,这里0代表什么意思,
dim_avg_n_Wrap(x,0)的解释是可以不用重置矩阵,可以对任一维度求平均,怎么使用,
我该如何设置我所想要的年份区间呢,是在输入资料的时候就设置好吗?

; eof_1.ncl
;
; Concepts illustrated:
;   - Calculating EOFs
;   - Using coordinate subscripting to read a specified geographical region
;   - Rearranging longitude data to span -180 to 180
;   - Calculating symmetric contour intervals
;   - Drawing filled bars above and below a given reference line
;   - Drawing subtitles at the top of a plot
;   - Reordering an array
; ==============================================================
; NCL V6.4.0 has new functions eofunc_n_Wrap and
; eofunc_ts_n_Wrap that allow you to calculate the EOFs without
; first having to first reorder the data. See eof_1_640.ncl.
; ==============================================================
; Calculate EOFs of the Sea Level Pressure over the North Atlantic.
; ==============================================================
; The slp.mon.mean file can be downloaded from:
; http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.surface.html
; ==============================================================
; These files are loaded by default in NCL V6.2.0 and newer
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
; ==============================================================
; User defined parameters that specify region of globe and
; ==============================================================
  latS   =  10.
  latN   =  40.
  lonL   =  80.
  lonR   =  160.
  yrStrt = 1990
  yrLast = 2015
  season = "JJA"    ; choose Dec-Jan-Feb seasonal mean
  neof   = 3        ; number of EOFs
  optEOF = True      
  optEOF@jopt = 0   ; This is the default; most commonly used; no need to specify.
;;optEOF@jopt = 1   ; **only** if the correlation EOF is desired
  optETS = False
; ==============================================================
; Open the file: Read only the user specified period
; ==============================================================
  f      = addfile ("/cygdrive/d/fnl/hgt.mon.mean.nc", "r")
  TIME   = f->time
  YYYY   = cd_calendar(TIME,-1)/100                 ; entire file
  iYYYY  = ind(YYYY.ge.yrStrt .and. YYYY.le.yrLast)
  slp    = f->hgt(510:824,12,:,:)  
;  printVarSummary(slp)                              ; variable overview
; ==============================================================
; dataset longitudes span 0=>357.5
; Because EOFs of the North Atlantic Oscillation are desired
; use the "lonFlip" (contributed.ncl) to reorder
; longitudes to span -180 to 177.5: facilitate coordinate subscripting
; ==============================================================
  slp    = lonFlip( slp )
  anom = dim_rmvmean_n_Wrap(slp,0)
  printVarSummary(slp)                              ; note the longitude coord
; ==============================================================
; compute desired global seasonal mean: month_to_season (contributed.ncl)
; ==============================================================


; SLP    = month_to_season (anom, season)
  
  nyrs   = dimsizes(anom&time)
; printVarSummary(SLP)
; =================================================================
; create weights:  sqrt(cos(lat))   [or sqrt(gw) ]
; =================================================================
  rad    = 4.*atan(1.)/180.
  clat   = f->lat           
  clat   = sqrt( cos(rad*clat) )                 ; gw for gaussian grid
; =================================================================
; weight all observations
; =================================================================
  wSLP   = anom                                   ; copy meta data
  wSLP   = anom*conform(anom, clat, 1)
  wSLP@long_name = "Wgt: "+wSLP@long_name
; =================================================================
; Reorder (lat,lon,time) the *weighted* input data
; Access the area of interest via coordinate subscripting
; =================================================================
  xw     = wSLP({lat|latS:latN},{lon|lonL:lonR},time|:)
  x      = wSLP(time|:,{lat|latS:latN},{lon|lonL:lonR})
  eof      = eofunc_Wrap(xw, neof, optEOF)      
  eof_ts   = eofunc_ts_Wrap (xw, eof, optETS)
  printVarSummary( eof )                         ; examine EOF variables
  printVarSummary( eof_ts )
; =================================================================
; Normalize time series: Sum spatial weights over the area of used
; =================================================================
  dimxw  = dimsizes( xw )
  mln    = dimxw(1)
  sumWgt = mln*sum( clat({lat|latS:latN}) )
  eof_ts = eof_ts/sumWgt
; =================================================================
; Extract the YYYYMM from the time coordinate
; associated with eof_ts [same as SLP&time]
; =================================================================
  yyyymm = cd_calendar(eof_ts&time,-2)/100  
;;yrfrac = yyyymm_to_yyyyfrac(yyyymm, 0.0); not used here
;============================================================
; PLOTS
;============================================================
  wks = gsn_open_wks("pdf","eof")         ; send graphics to PNG file
  plot = new(neof,graphic)                ; create graphic array
                                          ; only needed if paneling
; EOF patterns
  res                      = True         
  res@gsnDraw              = False        ; don't draw yet
  res@gsnFrame             = False        ; don't advance frame yet
  res@gsnAddCyclic         = False        ; plotted dataa are not cyclic

  res@mpFillOn             = False        ; turn off map fill
  res@mpMinLatF            = latS         ; zoom in on map
  res@mpMaxLatF            = latN
  res@mpMinLonF            = lonL
  res@mpMaxLonF            = lonR
  res@cnFillOn             = True         ; turn on color fill
  res@cnLinesOn            = False        ; True is default
;res@cnLineLabelsOn       = False        ; True is default
  res@cnFillPalette        = "BlWhRe"     ; set color map
  res@lbLabelBarOn         = False        ; turn off individual lb's
                                          ; set symmetric plot min/max
  symMinMaxPlt(eof, 16, False, res)       ; contributed.ncl
; panel plot only resources
  resP                     = True         ; modify the panel plot
  resP@gsnMaximize         = True         ; large format
  resP@gsnPanelLabelBar    = True         ; add common colorbar
  yStrt                    = yyyymm(0)/100
  yLast                    = yyyymm(nyrs-1)/100
  resP@txString            = "HGT_500hPa: "+season+": "+yStrt+"-"+yLast
;*******************************************
; first plot
;*******************************************
  do n=0,neof-1
     res@gsnLeftString  = "EOF "+(n+1)
     res@gsnRightString = sprintf("%5.1f", eof@pcvar(n)) +"%"
     plot(n)=gsn_csm_contour_map_ce(wks,eof(n,:,:),res)
  end do
  gsn_panel(wks,plot,(/neof,1/),resP)     ; now draw as one plot
;*******************************************
; second plot
;*******************************************
; EOF time series  [bar form]
  rts           = True
  rts@gsnDraw   = False       ; don't draw yet
  rts@gsnFrame  = False       ; don't advance frame yet
  rts@gsnScale  = True        ; force text scaling               
; these four rtsources allow the user to stretch the plot size, and
; decide exactly where on the page to draw it.
  rts@vpHeightF = 0.40        ; Changes the aspect ratio
  rts@vpWidthF  = 0.85
  rts@vpXF      = 0.10        ; change start locations
  rts@vpYF      = 0.75        ; the plot

  rts@tiYAxisString = "gpm"                    ; y-axis label      
  rts@gsnYRefLine           = 0.              ; reference line   
  rts@gsnXYBarChart         = True            ; create bar chart
  rts@gsnAboveYRefLineColor = "red"           ; above ref line fill red
  rts@gsnBelowYRefLineColor = "blue"          ; below ref line fill blue
; panel plot only resources
  rtsP                      = True            ; modify the panel plot
  rtsP@gsnMaximize          = True            ; large format
  rtsP@txString             = "HGT_500hPa: "+season+": "+yStrt+"-"+yLast
  year = yyyymm/100
; create individual plots
  do n=0,neof-1
     rts@gsnLeftString  = "EOF "+(n+1)
     rts@gsnRightString = sprintf("%5.1f", eof@pcvar(n)) +"%"
     plot(n) = gsn_csm_xy (wks,year,eof_ts(n,:),rts)
  end do
  gsn_panel(wks,plot,(/neof,1/),rtsP)     ; now draw as one plot
end
设置这样以后,提示在plot(n) = gsn_csm_xy (wks,year,eof_ts(n,:),rts)这一句出错:
(0)     gsn_csm_xy: Fatal: The x array is not monotonically increasing.
(0)     Cannot draw a bar chart.
fatal:Illegal right-hand side type for assignment
fatal:["Execute.c":8128]:Execute: Error occurred at or near line 207 in file /home/Administrator/eof_hgt500.ncl。
各路大神求解啊!
密码修改失败请联系微信:mofangbao
发表于 2017-1-6 16:34:10 | 显示全部楼层
dim_avg(x,0) 可能代表对变量的第一维求平均
密码修改失败请联系微信:mofangbao
发表于 2018-4-11 19:23:30 | 显示全部楼层
你好~我没有找到anom.ncl,请问您能给我发一份吗~
密码修改失败请联系微信:mofangbao
发表于 2020-12-3 21:38:55 | 显示全部楼层
百度搜ncl提取时间
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

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

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

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