2020-04-17

利用FNL資料畫(huà)Tlog-p圖 ,參考鏈接http://www.ncl.ucar.edu/Applications/skewt.shtml
代碼如下:

;**************************************************
; skewt_6.ncl
;
; Concepts illustrated:
;   - Reading RUC (Rapid Update Cycle)  GRIB data
;   - Using getind_latlon2d to determine grid locations
;   - Drawing Skew-T plots at nearest grid locations  
;                    to user specified locations
;**************************************************
;
; 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/contributed.ncl"
;
; This file still has to be loaded manually
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/skewt_func.ncl"
;***********************************************
; RUC Data downloaded from:
;     http://nomads.ncdc.noaa.gov/data.php?name=access#hires_weather_datasets
; RUC-ANL: http
;     Specifically: http://nomads.ncdc.noaa.gov/data/rucanl/201205/20120501/
;***********************************************
; The GRIB file's contents can be examined via:
;     ncl_filedump -itime ruc2anl_130_20120501_0000_000.grb2  | less
;***********************************************
 begin

; --- Read RUC GRIB file------------;
f = addfile("new.grib2", "r")
;printVarSummary(f)
time  = f->initial_time0_hours
YYYYMMDDHH = tostring(cd_calendar(time, -3))
dims = dimsizes(YYYYMMDDHH)
cst=new((dims),string)
cst1=new((dims),string)
do i = 0,dims-1
    utc = str_get_cols(YYYYMMDDHH(i), 0, 7)
   hr = str_get_cols(YYYYMMDDHH(i), 8, 9)
    cst(i) = systemfunc("date +'%Y%m%d%H' -d " + "'" + utc + " " + hr + " " + 8 + " hour'")
    cst1(i) = str_get_cols(cst(i), 2, 9)
end do
;  p     = f->lv_ISBL0               ;  ( lv_ISBL0 )
  p     = f->lv_ISBL5({10000:100000})
  print(p)
  ;time  = f->initial_time0_encoded  ; yyyymmddhh.hh_frac

                                    ; RUC grid point locations
  lat2d = f->lat_0              ; ( ygrid_0, xgrid_0 ) 
  lon2d = f->lon_0
;  xCopy = new ( (/dimsizes(lat)-1,dimsizes(lon)-1/), typeof(x))
;  do i=0,dimsizes(lat)-1
;;    do j=0,dimsizes(lon)-1
;;     lat2d(k,:,:) =lat(i)
;     lat2d(k,:,:) = lat()
;  end do
  print("lat2d: min="+min(lat2d)+"  ;  max="+max(lat2d))
  print("lon2d: min="+min(lon2d)+"  ;  max="+max(lon2d))
  printVarSummary(lat2d)
  p     = p*0.01                    ; change units
  p@units = "hPa"                   ; skewT, mixhum_ptrh use mb (hPa)
;--- Specify one or more locations
  printVarSummary(p)

nn=29.28
mm=113.5

;*************************
; create plot(s)
;*************************
do i = 0, dims-1
    skewtOpts                 = True
    skewtOpts@DrawColAreaFill = True    ; default is False
   
    dataOpts   = True
    dataOpts@PrintZ = True                          
    tk   = f->TMP_P0_L100_GLL0(i,{10000:100000},{nn},{mm})
    z    = f->HGT_P0_L100_GLL0(i,{10000:100000},{nn},{mm})
    rh   = f->RH_P0_L100_GLL0(i,{10000:100000},{nn},{mm})
    u    = f->UGRD_P0_L100_GLL0(i,{10000:100000},{nn},{mm})
    v    = f->VGRD_P0_L100_GLL0(i,{10000:100000},{nn},{mm})
    printVarSummary(tk)
    printVarSummary(rh)
    printVarSummary(u)
    printVarSummary(z)

         ; change units and calculate needed variables
    tc   = tk - 273.15
    tc@units= "degC"            

    q    = mixhum_ptrh (p, tk, rh, 2)
    q@units = "kg/kg"

    tdc  = dewtemp_trh(tk,rh) - 273.15 
    tdc@units = "degC"            

    wspd = sqrt((u*2.5)^2 + (v*2.5)^2)
    wdir = wind_direction(u,v,0)

     ;itime= toint(time)
    skewtOpts@tiMainString    = cst(i)+"("+nn+","+mm+")"
    skewtOpts@DrawFahrenheit  = False   ; default is True
          ; each location will have a different file name      
    wks  = gsn_open_wks ("png", "FNL_skewt_"+cst(i))
    skewt_bkgd = skewT_BackGround (wks, skewtOpts)
    skewt_data = skewT_PlotData   (wks, skewt_bkgd, p,tc,tdc,z \
                                       , wspd,wdir, dataOpts)
    draw (skewt_bkgd)
    draw (skewt_data)
    frame(wks)
end do
end  
;shox 肖沃特穩(wěn)定指數(shù) Showalter stability index 斜坡對(duì)流 sloping convection 
FNL_skewt_2020032120.png

(21日夜間湘北狂風(fēng)驟雨,多地出現(xiàn)冰雹天氣,其中臨湘還出現(xiàn)了超級(jí)單體)
值得注意的是FNL資料一直在更新,2020更新了lv_ISBL5垂直分布,共有34層,這樣Tk,z都是34層,而rh和 u、v是31層,比lv_ISBL0(31層)多了3層 ,在畫(huà)圖時(shí)要注意挑選層次{10000:100000}

最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點(diǎn),簡(jiǎn)書(shū)系信息發(fā)布平臺(tái),僅提供信息存儲(chǔ)服務(wù)。

友情鏈接更多精彩內(nèi)容