開始一個(gè)簡(jiǎn)單的NCL

為了更好地對(duì)NCL有一個(gè)整體的認(rèn)知葵硕,首先竞穷,舉一個(gè)最簡(jiǎn)單的例子梳庆。
畫出1948年以來(lái)500hPa高度上的平均等位勢(shì)線暖途。
數(shù)據(jù)可從以下網(wǎng)址中獲取:https://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.pressure.html

Fig.1.hgt.mon.mean.nc數(shù)據(jù)下載.png

有了數(shù)據(jù)膏执,知道了要畫什么驻售,首先你要做的是:知道數(shù)據(jù)的存儲(chǔ)格式,也就是知道數(shù)據(jù)長(zhǎng)什么樣子更米。
像我們剛剛下載的數(shù)據(jù)芋浮,屬于比較正規(guī)的數(shù)據(jù),其網(wǎng)站上面有一些說(shuō)明:

Fig.2.數(shù)據(jù)格式說(shuō)明.png

Brief Description: NCEP/NCAR再分析資料
Temporal Coverage: 從1948年1月1日起壳快,一天4次/逐日/逐月數(shù)據(jù)
-----> 我們下載的是monthly也就是逐月數(shù)據(jù)
Spatial Coverage: 分辨率為2.5×2.5°的格點(diǎn)纸巷,0.0E~357.5E, 90.0N~90.0S
Levels: 垂直分為17層,分別是1000hPa, 925hPa......以此類推

那除了這些信息眶痰,更多的數(shù)據(jù)信息是存儲(chǔ)在.nc文件本身里的瘤旨,對(duì)于這些數(shù)據(jù),我們可以直接通過ncl進(jìn)行查看竖伯。

數(shù)據(jù)信息的查看

1存哲、ncl_filedump命令(未進(jìn)入ncl):ncl_filedump文檔
打開終端,輸入ncl_filedump hgt.mon.mean.nc七婴,注意數(shù)據(jù)文件的路徑K钔怠!要么先cd到當(dāng)前目錄下(如下圖)打厘,要么在.nc文件名前加入相對(duì)路徑修肠!

Fig.3.png

2、ncl導(dǎo)入數(shù)據(jù)后再查看:
打開終端户盯,輸入ncl按回車進(jìn)入ncl命令行模式中嵌施。輸入:

f = addfile("..../hgt.mon.mean.nc", "r")  ; 將數(shù)據(jù)讀取并寫入f這個(gè)變量中
printVarSummary(f)  ; 輸出f變量的主要信息
print(f)  ; 輸出f變量的所有信息

可以對(duì)比一下兩種查看方式的異同點(diǎn)~


Fig.4.對(duì)比.png

數(shù)據(jù)信息含義

由以上兩種方式輸出的數(shù)據(jù)信息如下(以注釋形式在代碼中標(biāo)出,每行;之后的內(nèi)容為注釋):

Variable: f
Type: file
filename:   hgt.mon.mean
path:   hgt.mon.mean.nc
   file global attributes:
      description :  Data from NCEP initialized reanalysis (4x/day).  These are interpolated to pressure surfaces from model (sigma) surfaces.
      platform : Model
      Conventions : COARDS
      NCO : 20121012
      history : Created by NOAA-CIRES Climate Diagnostics Center (SAC) from the NCEP
reanalysis data set on 07/07/97 by calc.mon.mean.year.f using
/Datasets/nmc.reanalysis.derived/pressure/hgt.mon.mean.nc
from /Datasets/nmc.reanalysis/pressure/hgt.79.nc to hgt.95.nc
Converted to chunked, deflated non-packed NetCDF4 2014/09
      title : monthly mean hgt from the NCEP Reanalysis
      References : http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.derived.html
      dataset_title : NCEP-NCAR Reanalysis 1
   dimensions:
      level = 17        ; level這個(gè)維度有17層
      lat = 73            ; 緯度73個(gè)
      lon = 144        ; 精度144個(gè)
      time = 832  // unlimited
   variables:
      float level ( level )        ; 變量名:level
         units :    millibar    ; 單位:hPa
         long_name :    Level
         positive : down
         GRIB_id :  100
         GRIB_name :    hPa
         actual_range : ( 1000, 10 )
         axis : Z

      float lat ( lat )            ; 變量名:lat
         units :    degrees_north    ; 單位: 北緯
         actual_range : ( 90, -90 )    ; 值的范圍
         long_name :    Latitude
         standard_name :    latitude
         axis : Y

      float lon ( lon )
         units :    degrees_east
         long_name :    Longitude
         actual_range : (  0, 357.5 )
         standard_name :    longitude
         axis : X

      double time ( time )
         long_name :    Time
         delta_t :  0000-01-00 00:00:00
         avg_period :   0000-01-00 00:00:00
         prev_avg_period :  0000-00-01 00:00:00
         standard_name :    time
         axis : T
         units :    hours since 1800-01-01 00:00:0.0
         actual_range : ( 1297320, 1904352 )

      float hgt ( time, level, lat, lon )      ; 變量名:hgt
         long_name :    Monthly mean geopotential height
         valid_range :  ( -700, 35000 )
         units :    m
         add_offset :    0
         scale_factor :  1
         missing_value :    -9.96921e+36      ; 缺失值填補(bǔ)為這個(gè)值
         precision :    0
         least_significant_digit :  0
         GRIB_id :  7
         GRIB_name :    HGT
         var_desc : Geopotential height
         level_desc :   Multiple levels
         statistic :    Mean
         parent_stat :  Other
         dataset :  NCEP Reanalysis Derived Products
         actual_range : ( -354.4583, 32321.1 )
         _FillValue :   -9.96921e+36

可以看到莽鸭,其中hgt這個(gè)變量就是我們要用的吗伤,在某等壓面下的位勢(shì)數(shù)據(jù)啦。
float hgt ( time, level, lat, lon )說(shuō)明hgt這個(gè)變量是個(gè)四維(832*17*73*144)的變量硫眨。
首先將這個(gè)變量從f中抽取出來(lái)足淆,即輸入:

ncl 0> f = addfile("hgt.mon.mean.nc", "r")
ncl 1> hgt = f->hgt
ncl 2>

如果不放心還可以使用printVarSummary()函數(shù)輸出一下變量

ncl 0> f = addfile("hgt.mon.mean.nc", "r")
ncl 1> hgt = f->hgt
ncl 2> printVarSummary(hgt)

Variable: hgt
Type: float
Total Size: 594726912 bytes
            148681728 values
Number of Dimensions: 4
Dimensions and sizes:   [time | 832] x [level | 17] x [lat | 73] x [lon | 144]
Coordinates:
            time: [1297320..1904352]
            level: [1000..10]
            lat: [90..-90]
            lon: [ 0..357.5]
Number Of Attributes: 17
  long_name :   Monthly mean geopotential height
  valid_range : ( -700, 35000 )
  units :   m
  add_offset :   0
  scale_factor :     1
  missing_value :   -9.96921e+36
  precision :   0
  least_significant_digit : 0
  GRIB_id : 7
  GRIB_name :   HGT
  var_desc :    Geopotential height
  level_desc :  Multiple levels
  statistic :   Mean
  parent_stat : Other
  dataset : NCEP Reanalysis Derived Products
  actual_range :    ( -354.4583, 32321.1 )
  _FillValue :  -9.96921e+36
ncl 3>

數(shù)據(jù)的切片和整合

命題為:【畫出1948年以來(lái)500hPa高度上的平均等位勢(shì)線】
目前hgt這個(gè)變量time有832維,level也是包括17個(gè)等壓面的數(shù)據(jù)礁阁,我們需要:
1巧号、將時(shí)間維度做平均
2、切片只留500hPa這個(gè)level的數(shù)據(jù)

1氮兵、時(shí)間維度平均

ncl 3> mean_hgt = dim_avg_n_Wrap(hgt, 0)    ; 0表示對(duì)hgt的第一維做平均

輸出的mean_hgt如下:

ncl 4> printVarSummary(mean_hgt)

Variable: mean_hgt
Type: float
Total Size: 714816 bytes
            178704 values
Number of Dimensions: 3
Dimensions and sizes:   [level | 17] x [lat | 73] x [lon | 144]
Coordinates:
            level: [1000..10]
            lat: [90..-90]
            lon: [ 0..357.5]
Number Of Attributes: 18
  _FillValue :  -9.96921e+36
  long_name :   Monthly mean geopotential height
  valid_range : ( -700, 35000 )
  units :   m
  add_offset :   0
  scale_factor :     1
  precision :   0
  least_significant_digit : 0
  GRIB_id : 7
  GRIB_name :   HGT
  var_desc :    Geopotential height
  level_desc :  Multiple levels
  statistic :   Mean
  parent_stat : Other
  dataset : NCEP Reanalysis Derived Products
  actual_range :    ( -354.4583, 32321.1 )
  missing_value :   -9.96921e+36
  average_op_ncl :  dim_avg_n over dimension(s): time

2裂逐、切片只留level為500hPa的數(shù)據(jù)
切片操作我想大家都會(huì),那么問題來(lái)了泣栈,500hPa是第幾層卜高?
level是從高到底存的,還是從低到高南片?
float hgt ( time, level, lat, lon )這行代碼同時(shí)也意味著掺涛,level這個(gè)維度,對(duì)應(yīng)的坐標(biāo)是存在level這個(gè)變量中的疼进,其順序和level變量中的存儲(chǔ)順序是相同的薪缆。也就是說(shuō),我們看level這個(gè)變量里是什么樣的順序就好了伞广。

ncl 5> print(f->level)


Variable: level (file variable)
Type: float
Total Size: 68 bytes
            17 values
Number of Dimensions: 1
Dimensions and sizes:   [level | 17]
Coordinates:
            level: [1000..10]
Number Of Attributes: 7
  units :       millibar
  long_name :   Level
  positive :    down
  GRIB_id :     100
  GRIB_name :   hPa
  actual_range :        ( 1000, 10 )
  axis :        Z
(0)     1000
(1)     925
(2)     850
(3)     700
(4)     600
(5)     500            ; 看這里看這里看這里<鹈薄L鄣纭!
(6)     400
(7)     300
(8)     250
(9)     200
(10)    150
(11)    100
(12)    70
(13)    50
(14)    30
(15)    20
(16)    10

看到了吧减拭,500hPa是第6個(gè)值蔽豺,但是由于序號(hào)從0開始標(biāo),所以index=5.

ncl 6> mean_500_hgt = mean_hgt(5, :, :)
ncl 7> print(mean_500_hgt)


Variable: mean_500_hgt
Type: float
Total Size: 42048 bytes
            10512 values
Number of Dimensions: 2
Dimensions and sizes:   [lat | 73] x [lon | 144]
Coordinates:
            lat: [90..-90]
            lon: [ 0..357.5]
Number Of Attributes: 19
  level :       500
  average_op_ncl :      dim_avg_n over dimension(s): time
  missing_value :       -9.96921e+36
  actual_range :        ( -354.4583, 32321.1 )
  dataset :     NCEP Reanalysis Derived Products
  parent_stat : Other
  statistic :   Mean
  level_desc :  Multiple levels
  var_desc :    Geopotential height
  GRIB_name :   HGT
  GRIB_id :     7
  least_significant_digit :     0
  precision :   0
  scale_factor :         1
  add_offset :   0
  units :       m
  valid_range : ( -700, 35000 )
  long_name :   Monthly mean geopotential height
  _FillValue :  -9.96921e+36

現(xiàn)在這個(gè)變量是一個(gè)二維變量了~ 我們已經(jīng)可以將其畫在地圖上

畫圖設(shè)置

1拧粪、設(shè)置工作臺(tái)

    wks = gsn_open_wks("png", "500hPa mean hgt")

這是最簡(jiǎn)單的工作臺(tái)設(shè)置方式修陡,
第一個(gè)參數(shù)表示生成的格式,"png"表示png格式的圖片
第二個(gè)參數(shù)是生成的文件名稱
這樣設(shè)置出來(lái)以后意味著會(huì)生成一個(gè)500hPa mean hgt.png的文件可霎,一會(huì)兒畫出來(lái)圖就去相應(yīng)的路徑下找這個(gè)文件就好啦~

2魄鸦、設(shè)置圖片屬性
這一部分對(duì)于強(qiáng)迫癥患者來(lái)說(shuō)真的是最最最惡心的部分了,NCL語(yǔ)言設(shè)計(jì)中最最最惡心的地方就是在這里癣朗。圖片屬性中的每一個(gè)值都有對(duì)應(yīng)的參數(shù)名稱拾因,你想修改任意一個(gè)點(diǎn)都必須找到相應(yīng)的參數(shù),然后給予設(shè)置斯棒,不然就都是默認(rèn)值6苤隆!
以下是最簡(jiǎn)單的圖片屬性設(shè)置:

    res = True      ; resource變量為True
                      
    res@cnFillOn = True    ; 是否填充顏色
    res@cnLinesOn = False    ; 是否繪制等值線
    res@cnLineLabelsOn = False  ; 是否顯示等值線上的值標(biāo)簽
      ; cn這一溜表示這些參數(shù)是contour圖的屬性

畫圖

設(shè)置好了工作臺(tái)和圖片屬性荣暮,畫圖其實(shí)就是一行代碼的事庭惜。

plot = gsn_csm_contour_map(wks, mean_500_hgt, res)

gsn_csm_contour_map是一個(gè)繪圖函數(shù),_map表示繪制在地圖上穗酥,_contour表示畫圖的種類护赊,除了contour還有xy圖等等等好多,以后會(huì)講砾跃。
第一個(gè)參數(shù)骏啰,填入工作臺(tái)的變量。
第二個(gè)參數(shù)抽高,是你要繪制的變量判耕。
第三個(gè)參數(shù),是存儲(chǔ)圖片屬性的變量翘骂。

繪制出來(lái)的圖片:

Fig.5.500hPa平均位勢(shì)高度.png

....看吧壁熄,位勢(shì)高度場(chǎng)的年平均值……畫出來(lái)真的是毫無(wú)意義啊……有興趣的同學(xué)可以畫逐月的值,或者是夏季平均值~

結(jié)束語(yǔ)

只是為了讓大家快速地對(duì)NCL整體繪圖流程有一個(gè)模糊的概念碳竟,我相信中間有很多地方大家有疑惑的草丧,我會(huì)努力盡快整理出來(lái)~
大家對(duì)文章有哪些建議或者意見,還請(qǐng)直接留言評(píng)論~
有哪些迫切想要了解的點(diǎn)莹桅,也可以告訴我昌执,我可以提前梳理出來(lái)。

附:完整代碼
; ************************************
;   Chapter.2.開始一個(gè)簡(jiǎn)單的NCL畫圖程序
; ************************************
; 直接打開一個(gè)編輯器,另存為.ncl的文件就行了
; 自己回顧一下完整的.ncl代碼要怎么跑懂拾!
; *************************************
begin
    f = addfile("hgt.mon.mean.nc","r")
    ; printVarSummary(f)

    hgt = f->hgt

    mean_hgt = dim_avg_n_Wrap(hgt, 0)
    ; printVarSummary(mean_hgt)

    ;print(f->level)
    mean_500_hgt = mean_hgt(5, :, :)

    ; ------ set the work station ------
    wks = gsn_open_wks("png", "500hPa mean hgt")

    ; ------ set the resources ------
    res = True

    res@cnFillOn = True
    res@cnLinesOn = False
    res@cnLineLabelsOn = False

    plot= gsn_csm_contour_map(wks, mean_500_hgt, res)
end
最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
  • 序言:七十年代末煤禽,一起剝皮案震驚了整個(gè)濱河市,隨后出現(xiàn)的幾起案子委粉,更是在濱河造成了極大的恐慌呜师,老刑警劉巖,帶你破解...
    沈念sama閱讀 222,183評(píng)論 6 516
  • 序言:濱河連續(xù)發(fā)生了三起死亡事件贾节,死亡現(xiàn)場(chǎng)離奇詭異,居然都是意外死亡衷畦,警方通過查閱死者的電腦和手機(jī)栗涂,發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 94,850評(píng)論 3 399
  • 文/潘曉璐 我一進(jìn)店門,熙熙樓的掌柜王于貴愁眉苦臉地迎上來(lái)祈争,“玉大人斤程,你說(shuō)我怎么就攤上這事∑谢欤” “怎么了忿墅?”我有些...
    開封第一講書人閱讀 168,766評(píng)論 0 361
  • 文/不壞的土叔 我叫張陵,是天一觀的道長(zhǎng)沮峡。 經(jīng)常有香客問我疚脐,道長(zhǎng),這世上最難降的妖魔是什么邢疙? 我笑而不...
    開封第一講書人閱讀 59,854評(píng)論 1 299
  • 正文 為了忘掉前任棍弄,我火速辦了婚禮,結(jié)果婚禮上疟游,老公的妹妹穿的比我還像新娘呼畸。我一直安慰自己,他們只是感情好颁虐,可當(dāng)我...
    茶點(diǎn)故事閱讀 68,871評(píng)論 6 398
  • 文/花漫 我一把揭開白布蛮原。 她就那樣靜靜地躺著,像睡著了一般另绩。 火紅的嫁衣襯著肌膚如雪儒陨。 梳的紋絲不亂的頭發(fā)上,一...
    開封第一講書人閱讀 52,457評(píng)論 1 311
  • 那天板熊,我揣著相機(jī)與錄音框全,去河邊找鬼。 笑死干签,一個(gè)胖子當(dāng)著我的面吹牛津辩,可吹牛的內(nèi)容都是我干的。 我是一名探鬼主播,決...
    沈念sama閱讀 40,999評(píng)論 3 422
  • 文/蒼蘭香墨 我猛地睜開眼喘沿,長(zhǎng)吁一口氣:“原來(lái)是場(chǎng)噩夢(mèng)啊……” “哼闸度!你這毒婦竟也來(lái)了?” 一聲冷哼從身側(cè)響起蚜印,我...
    開封第一講書人閱讀 39,914評(píng)論 0 277
  • 序言:老撾萬(wàn)榮一對(duì)情侶失蹤莺禁,失蹤者是張志新(化名)和其女友劉穎,沒想到半個(gè)月后窄赋,有當(dāng)?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體哟冬,經(jīng)...
    沈念sama閱讀 46,465評(píng)論 1 319
  • 正文 獨(dú)居荒郊野嶺守林人離奇死亡,尸身上長(zhǎng)有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點(diǎn)故事閱讀 38,543評(píng)論 3 342
  • 正文 我和宋清朗相戀三年忆绰,在試婚紗的時(shí)候發(fā)現(xiàn)自己被綠了浩峡。 大學(xué)時(shí)的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片。...
    茶點(diǎn)故事閱讀 40,675評(píng)論 1 353
  • 序言:一個(gè)原本活蹦亂跳的男人離奇死亡错敢,死狀恐怖翰灾,靈堂內(nèi)的尸體忽然破棺而出,到底是詐尸還是另有隱情稚茅,我是刑警寧澤纸淮,帶...
    沈念sama閱讀 36,354評(píng)論 5 351
  • 正文 年R本政府宣布,位于F島的核電站亚享,受9級(jí)特大地震影響咽块,放射性物質(zhì)發(fā)生泄漏。R本人自食惡果不足惜虹蒋,卻給世界環(huán)境...
    茶點(diǎn)故事閱讀 42,029評(píng)論 3 335
  • 文/蒙蒙 一糜芳、第九天 我趴在偏房一處隱蔽的房頂上張望。 院中可真熱鬧魄衅,春花似錦峭竣、人聲如沸。這莊子的主人今日做“春日...
    開封第一講書人閱讀 32,514評(píng)論 0 25
  • 文/蒼蘭香墨 我抬頭看了看天上的太陽(yáng)。三九已至哲银,卻和暖如春扛吞,著一層夾襖步出監(jiān)牢的瞬間,已是汗流浹背荆责。 一陣腳步聲響...
    開封第一講書人閱讀 33,616評(píng)論 1 274
  • 我被黑心中介騙來(lái)泰國(guó)打工滥比, 沒想到剛下飛機(jī)就差點(diǎn)兒被人妖公主榨干…… 1. 我叫王不留,地道東北人做院。 一個(gè)月前我還...
    沈念sama閱讀 49,091評(píng)論 3 378
  • 正文 我出身青樓盲泛,卻偏偏與公主長(zhǎng)得像濒持,于是被迫代替她去往敵國(guó)和親。 傳聞我的和親對(duì)象是個(gè)殘疾皇子寺滚,可洞房花燭夜當(dāng)晚...
    茶點(diǎn)故事閱讀 45,685評(píng)論 2 360

推薦閱讀更多精彩內(nèi)容

  • Android 自定義View的各種姿勢(shì)1 Activity的顯示之ViewRootImpl詳解 Activity...
    passiontim閱讀 172,304評(píng)論 25 707
  • linux資料總章2.1 1.0寫的不好抱歉 但是2.0已經(jīng)改了很多 但是錯(cuò)誤還是無(wú)法避免 以后資料會(huì)慢慢更新 大...
    數(shù)據(jù)革命閱讀 12,175評(píng)論 2 33
  • Spring Cloud為開發(fā)人員提供了快速構(gòu)建分布式系統(tǒng)中一些常見模式的工具(例如配置管理柑营,服務(wù)發(fā)現(xiàn),斷路器村视,智...
    卡卡羅2017閱讀 134,708評(píng)論 18 139
  • 孟圓的筆記閱讀 304評(píng)論 0 0
  • 上周天氣開始徹底熱了起來(lái)官套,溫度灼燒,幾乎又開啟了周末白天泡在家里蚁孔,傍晚才出門的模式奶赔。 連著兩天去游泳池,一個(gè)選在4...
    Silly妞一枚閱讀 603評(píng)論 0 1