nc文件的讀取與處理

1、.nc 文件的查看和處理2013-08-04 12:45閱讀：2,168 最近在學(xué)習(xí)冰后回彈模型（http:/www.atmosp.physics.utoronto.ca/peltier/data.php），里面給出的數(shù)據(jù)為nc（netcdf）格式的，讀取是個問題，在言深深的幫助下，成功讀取，在這里將其總結(jié)出來，以備后用也希望可以方便有用之人。有兩種方式可以進(jìn)行讀?。?、Excel的一個插件，NetCDF4Excel可以完成相關(guān)的查看工作，需要的可以進(jìn)行下載( . 416233&third=15 （華為網(wǎng)盤太不靠譜，不知道百度靠譜到合適）安裝之后，我們可以得到，一個快捷方式，

2、如圖：打開之后，我們可以看到下面的情形：然后打開nc文件，我們就可以進(jìn)行查看操作：2、通過Grads進(jìn)行相關(guān)的讀取，主要是在言深深同學(xué)的幫助下完成，深深的感謝！第一步，安裝相關(guān)的grads軟件，這個這里不做介紹，相關(guān)可以論壇網(wǎng)站；第二步，打開grads軟件，如圖，回車之后的那個窗口千萬不能關(guān)，那個窗口是圖像顯示的，關(guān)閉之后輸入命令之后直接就退出了，剛開始我老犯這個低級錯誤（對菜鳥而言）。第三步，可以通過兩個命令簡單的看到數(shù)據(jù)的頭文件包含信息：sdfopen d:1.ncq ctlinfo然后編寫gs文件，如下：'reinit''sdfopen D:1.nc'&#

3、39;set gxout fwrite''se t fwrite D:1.dat''set x 1 360''set y 1 179''set z 1''set t 1 1''d dsea_250''disable fwrite'ps1：路徑自己可以按需更改，不要有空格。ps2：這里的dsea_250可不是隨便取的，看上面的ctlinfo心里里最后一行，必須要一致！然后就可以生成一個1.dat的文件。第四步，通過fortran程序?qū)?.dat讀為txt文檔parameter(

4、m=179,n=360,o=1)dimension x(m,n,o)open(1,file='kaka.dat',access='direct',form='binary',recl=n*m*o*4)read(1,rec=1) (x(i,j,k),j=1,n),i=1,m),k=1,o)close(1)open(2,file='data.txt')do i=1,mprint*,ido j=1,nwrite(2,'(2f10.2,2x,f6.2)')(i)*1.-90,(j-1)*1.,x(i,j,1)enddoen

5、ddoclose(2)endps：這里的參數(shù)要與自己的數(shù)據(jù)一致，比如mno，以及最后的輸出，有的是1°，有的是2.5°。然后我們就可以得到相關(guān)的數(shù)據(jù)文件如下：如何用ncdump查看.nc文件由于用grads打開.nc文件時遇到問題，提示SDF file has no discernable X coordinate，所以在網(wǎng)上找了好久解決的辦法，幾乎所有辦法都指向一個ncdump的使用，于是又開始尋找使用ncdump的辦法，最后終于得到一個可行的方法，如下：需要以下兩個文件ncdump.exe和netcdf.dll讀nc文件的步驟：1、將ncdump.exe和netcdf.

6、dll放在一個文件夾里，然后將要讀取的nc文件也放到相同的文件夾中2、將此文件夾放到一個目錄下，例如放到d：netcdf目錄下3、在dos命令行窗口中輸入：d：cd netcdfncdump wrf.nc或者ncdump wrf.nc|more通過以上三步即可查看nc文分享資料 GRADS怎么處理多個NC文件？謝謝 復(fù)制鏈接 |關(guān)注本帖取消最新回復(fù)取消置頂回復(fù)取消最新編輯duan2008 duan2008 當(dāng)前離線 積分1023貢獻(xiàn)精華在線時間小時注冊時間2012-5-7最后登錄1970-1-1窺視卡 雷達(dá)卡 電梯直達(dá) 樓主 發(fā)表于 2012-11-30 19:12:24 | 只看該作者 |倒

7、序瀏覽 |閱讀模式 |取消關(guān)注該作者的回復(fù)登錄后查看更多精彩內(nèi)容您需要 登錄 才可以下載或查看，沒有帳號？立即注冊 x 'reinit''sdfopen D:TDDOWNLOADhgt.1981.nc''set dfile 1''set lat 0 60''set lon 70 180''set lev 500''set t 1''sdfopen D:TDDOWNLOADhgt.1982.nc''set dfile 2''set lat 0 6

8、0''set lon 70 180''set lev 500''set t 1''define hgt12=hgt.1+hgt.2''set t 1''d hgt12''enable print D:TDDOWNLOAD1981-1982HGT.gmf'請問高手  我的程序那里錯誤了？我想合并這兩年的HGT。程序無法求取兩年的和，有可能是第一年的時間抽和第二年的不一樣，謝謝 分享到： QQ空間 新浪微博 騰訊微博 人人網(wǎng) 更多 收藏2 分享 淘帖0 贊0

9、踩0 給自己設(shè)置一個個性的簽名吧點擊前往GrADS官網(wǎng)新地址回復(fù) 使用道具 舉報 提升卡置頂卡沉默卡喧囂卡變色卡搶沙發(fā)千斤頂顯身卡river river 當(dāng)前離線 積分24285貢獻(xiàn)精華在線時間小時注冊時間2012-8-28最后登錄1970-1-1窺視卡 雷達(dá)卡 沙發(fā) 發(fā)表于 2012-11-30 20:11:19 | 只看該作者 |取消關(guān)注該作者的回復(fù)'reinit''sdfopen D:TDDOWNLOADhgt.1981.nc''sdfopen D:TDDOWNLOADhgt.1982.nc''set lat 0 60'&#

10、39;set lon 70 180''set lev 500''set t 1''define hgt12=hgt.1+hgt.2''d hgt12''enable print D:TDDOWNLOAD1981-1982HGT.gmf''print''disable print'還不行就把錯誤提示發(fā)上來啊，這樣才能對癥下藥啊給自己設(shè)置一個個性的簽名吧點擊前往GrADS官網(wǎng)新地址回復(fù) 支持 反對 使用道具 舉報 顯身卡duan2008 duan2008 當(dāng)前離線 積分1023貢

11、獻(xiàn)精華在線時間小時注冊時間2012-5-7最后登錄1970-1-1窺視卡 雷達(dá)卡 板凳  樓主| 發(fā)表于 2012-11-30 20:30:14 | 只看該作者 |取消關(guān)注該作者的回復(fù)river 發(fā)表于 2012-11-30 20:11                                 登錄/注冊后可看大圖'reinit''sdfopen D:TDDOWNLOADhgt.1981.nc&#

12、39;'sdfopen D:TDDOWNLOADhgt.1982.nc'謝謝啊，還是不行。錯誤如下 登錄/注冊后可看大圖 1.JPG (23.67 KB, 下載次數(shù): 0)下載附件  保存到相冊 2012-11-30 20:27 上傳登錄/注冊后可看大圖 2.JPG (26.74 KB, 下載次數(shù): 0)下載附件  保存到相冊 2012-11-30 20:29 上傳我在想是時間的設(shè)置問題嗎？因為第一個T=1是1981年1月1日，而第二個文件是1982年1月1日，而且有可能是365天和366天，在求和的時候有數(shù)據(jù)缺失？ 謝謝啊。給自己設(shè)置一個個性的簽名吧點擊前

13、往GrADS官網(wǎng)新地址回復(fù) 支持 反對 使用道具 舉報 顯身卡river river 當(dāng)前離線 積分24285貢獻(xiàn)精華在線時間小時注冊時間2012-8-28最后登錄1970-1-1窺視卡 雷達(dá)卡 地板 發(fā)表于 2012-11-30 21:18:42 | 只看該作者 |取消關(guān)注該作者的回復(fù)duan2008 發(fā)表于 2012-11-30 20:30                                 登錄/注冊后可看大圖謝謝啊，

14、還是不行。錯誤如下把那句set t 1去掉，把define那句改成define hgt12=hgt.1(t=1)+hgt.2(t=1)試試給自己設(shè)置一個個性的簽名吧點擊前往GrADS官網(wǎng)新地址回復(fù) 支持 反對 使用道具 舉報 顯身卡duan2008 duan2008 當(dāng)前離線 積分1023貢獻(xiàn)精華在線時間小時注冊時間2012-5-7最后登錄1970-1-1窺視卡 雷達(dá)卡 5#  樓主| 發(fā)表于 2012-11-30 21:57:08 | 只看該作者 |取消關(guān)注該作者的回復(fù)river 發(fā)表于 2012-11-30 21:18         &#

15、160;                       登錄/注冊后可看大圖把那句set t 1去掉，把define那句改成define hgt12=hgt.1(t=1)+hgt.2(t=1)試試這樣可以計算出一個值，但是對于366天的值怎么計算了，把hgt12存入一個NC文件，但是T是366天的，謝謝。給自己設(shè)置一個個性的簽名吧點擊前往GrADS官網(wǎng)新地址回復(fù) 支持 反對 使用道具 舉報 顯身卡river river 當(dāng)前離線 積分24285貢獻(xiàn)精華在線時間小時注冊時間2012-8-2

16、8最后登錄1970-1-1窺視卡 雷達(dá)卡 6# 發(fā)表于 2012-11-30 22:31:41 | 只看該作者 |取消關(guān)注該作者的回復(fù)duan2008 發(fā)表于 2012-11-30 21:57                                 登錄/注冊后可看大圖這樣可以計算出一個值，但是對于366天的值怎么計算了，把hgt12存入一個NC文件，但是T是366天的，謝謝。用循環(huán)可以吧，加幾句'reinit'

17、'sdfopen D:TDDOWNLOADhgt.1981.nc''sdfopen D:TDDOWNLOADhgt.1982.nc''enable print D:TDDOWNLOAD1981-1982HGT.gmf''set lat 0 60''set lon 70 180''set lev 500'it=1while（it<=365）'set t 'it'''define hgt12=hgt.1(t='it')+hgt.2(t=

18、9;it')'draw title 'it'''print''c'it=it+1'disable print'但是多出來一天就不好辦了啊給自己設(shè)置一個個性的簽名吧點擊前往GrADS官網(wǎng)新地址回復(fù) 支持 反對 使用道具 舉報 顯身卡xjliutao xjliutao 當(dāng)前離線 積分56貢獻(xiàn)精華在線時間小時注冊時間2013-3-21最后登錄1970-1-1窺視卡 雷達(dá)卡 7# 發(fā)表于 2013-3-21 20:18:20 | 只看該作者 |取消關(guān)注該作者的回復(fù)學(xué)習(xí)了，我一個都沒有處理好 你都可以處理多個了論壇里

19、常常有人對nc這種內(nèi)置ctl的文件數(shù)據(jù)處理 提出問題，這里總結(jié)下常見的一些問題及對策。nc數(shù)據(jù)是一類內(nèi)置ctl的文件，也就說沒有獨立的ctl文件，那么打開命令也不是常見的open，而是sdfopen命令，那么要處理這類數(shù)據(jù)就得知道變量名、變量層次分布等信息，用q ctlinfo命令就可以查看全部的ctl了。為什么提查看ctl 呢，很多人在處理這類數(shù)據(jù)的時候會想當(dāng)然的認(rèn)定一個變量名，缺省值、甚至精度都會不注意，那么無論怎么提取，出來的圖都不會對的。下面總結(jié)下容易出的問題：1、從nc文件提取出二進(jìn)制數(shù)據(jù)（fwrite）設(shè)定需要注意經(jīng)緯度，層次，時間是否是你需要的部分，尤其是用循環(huán)提取數(shù)據(jù)的朋友，尤

20、其要注意。有的時候設(shè)定經(jīng)緯度（set lat 0 359.8）出了問題，可是嘗試定義x y格點（set x 1 180）的方法試試，實在不行，可以用循環(huán)來輸出不同的經(jīng)緯度的數(shù)值。2、注意提取數(shù)據(jù)的grads命令'set gxout fwrite''set fwrite c:ss00.grd'.'disablefwrite'最后一行的這個命令，如果忘了加就容易數(shù)據(jù)不正常，還有的朋友，一個gs文件提取了兩個數(shù)據(jù)文件，兩組這類命令中間那個disable fwrite丟掉了，數(shù)據(jù)就徹底亂套了。3、提取的數(shù)據(jù)，寫ctl，需要注意與nc的ctl對應(yīng)，包括缺測

21、值，精度，這個問題只要照著nc的ctl寫，層次格點什么的寫成你數(shù)據(jù)對應(yīng)的就好了。4、數(shù)據(jù)提取了，ctl寫好了，畫畫圖看是不是與nc的一致，設(shè)定的時間、層次，圖是不是不一樣，因為用循環(huán)出數(shù)據(jù)的朋友容易犯 寫了好幾個時間點，數(shù)據(jù)是一樣的錯誤。5、還有的朋友需要進(jìn)一步提取這個grd或者dat文件，這就涉及到二進(jìn)制數(shù)據(jù)的讀取與操作問題了。二進(jìn)制數(shù)據(jù)是無格式的數(shù)據(jù)，讀取與寫入都不需要格式說明，里面的數(shù)據(jù)無法直接看到，其內(nèi)部無非是0和1，而且數(shù)據(jù)之間沒有間隔，讀取寫入只要read（文件號） write（文件號）就好讀取種要注意，grads寫入的數(shù)據(jù)循環(huán)從內(nèi)到外依次為x 、y、z、變量、t，如果你需要確定某

22、個數(shù)據(jù)的位置，那么，x是從西到東，y是從南到北，z是從地面到高空，t是從前到后的順序。另外此類文件讀取目前論壇里主要有兩種方法，一種是直接聲明form=binary的，這種形式簡單，如果數(shù)據(jù)是常規(guī)的順序讀取，選擇binary讀取方法簡單的多，也不容易出錯。如下面這段：program ex01real ss(360,180)integeri,jopen(100,file='s01.dat',form='binary')open(200,file='data0.txt')doj=1,180       read

23、(100) (ss(i,j),i=1,360)        write(200,*)(ss(i,j),i=1,360)        write(*,*) j,"  ",ss(1,j)enddoend另一種是無格式直接讀取的方式，無格式直接讀取主要是用于大量的數(shù)據(jù)需要處理，你只需要其中部分，可以用記錄的跳躍選取你要的數(shù)據(jù)，如下面這段：program readsstinteger,parameter: ix=360,iy=180,iz=

24、7real sst(ix,iy,iz)open(9,file='f:samareacchadleycentersst_update.dat',form='unformatted',access='direct',recl=ix*iy*4)read(9,rec=1) (sst(i,j),i=1,ix),j=1,iy)write(*,*) sstend最后總結(jié)下檢查錯誤的思路：1、數(shù)據(jù)沒提取出來，先檢查gs是否沒有寫disable fwrite，或者文件名.grd寫錯了；2、數(shù)據(jù)出來后，ctl寫好，畫圖與nc對比，如果發(fā)現(xiàn)，出來的圖上只有個別區(qū)域有大

25、條粗線，整幅圖都是黑的，說明你的ctl里缺測值或者精度不夠?qū)?yīng)，q ctlinfo 查詢下；3、圖對應(yīng)上了，改變設(shè)定時間或?qū)哟危瑳]有變化，說明，你提取數(shù)據(jù)的gs錯了，層次或者時間都是輸出的同一個點；4、前三步都沒有問題，恭喜你，數(shù)據(jù)提取成功了。Fortran編程的時候如果報錯說endof file during read，很可能是你文件名錯了，或者你的你讀取數(shù)據(jù)的順序不是按照grads的數(shù)據(jù)順序讀的，檢查你的各個循環(huán)的次數(shù)，及循環(huán)的嵌套順序是否正確。5、如果你是順序的讀取數(shù)據(jù)的話，建議你選擇binary的方式讀取，犯錯幾率要小；6、如果以上，還是報錯，剩下的原因就可能是你馬虎導(dǎo)致微小的錯誤了，

26、你要認(rèn)真的從頭到尾檢查程序，不行就讓別人給你檢查，因為有些問題太細(xì)微了，自己已經(jīng)習(xí)慣了。求問，我想計算兩個NDVI數(shù)據(jù)集中相同一個地區(qū)2001年到2006年12月-2月平均值的差值圖，結(jié)果老是出錯，不知大家可否幫忙看一下程序，先謝謝大家了！這是程序：'reinit'*gimms'sdfopen /Users/geog493/ndvi/ndvi0.5deg.2006.nc''sdfopen /Users/geog493/ndvi/ndvi0.5deg_modis_westafrica.nc''set lat 13 20.0''

27、;set lon -15.0 20.0' * DJF'define gy01=ave(ndvi.1,t=467,t=472)' 'define gy02=ave(ndvi.1,t=491,t=496)''define gy03=ave(ndvi.1,t=515,t=520)''define gy04=ave(ndvi.1,t=539,t=544)'  'define gy05=ave(ndvi.1,t=563,t=568)''define gy06=ave(ndvi.1,t=58

28、7,t=592)' 'define gyclim=(gy01+gy02+gy03+gy04+gy05+gy06)/6.'*'modis ndvi''set dfile 2''set t 1 307'* DJF'define my01=ave(ndvi.2,t=15,t=20)' 'define my02=ave(ndvi.2,t=38,t=43)''define my03=ave(ndvi.2,t=61,t=66)''define my04=ave(ndvi.2,t=

29、84,t=89)'  'define my05=ave(ndvi.2,t=107,t=112)''define my06=ave(ndvi.2,t=130,t=135)''define myclim=(my01+my02+my03+my04+my05+my06)/6.''define difyclim=gyclim-myclim''d gyclim-myclim'結(jié)果老是說“Error in gagchk: axis sizes are not the sameOperation err

30、or:  Incompatable grids ”但我的兩個NDVI數(shù)據(jù)集都是0.5X0.5的，不知道是什么原因？不過gimms數(shù)據(jù)集是從lat -61.75  lon -179.75開始; modis數(shù)據(jù)集是從lat 0  lon -20 開始，不知是否是這個緣故？如果是的話，應(yīng)該怎么解決呢？由于本人新手，都沒有積分的，先大大的謝謝大家！=07/29/2014解決方法因為兩個NDVI數(shù)據(jù)的起始網(wǎng)格點不同，需要轉(zhuǎn)換到同一的網(wǎng)格點上，運用函數(shù)re可做到這步：define gyclim1=re(gyclim, 70, linear, -1

31、5, 0.5, 14, linear, 13, 0.5, ba)define myclim1=re(myclim, 70, linear, -15, 0.5, 14, linear, 13, 0.5, ba)但不知什么原因，re在我最初的grads上顯示沒有此函數(shù)，后來重新裝了openGrADS才可行。openGrADS的安裝方法如下：Installing openGrADS /wiki/index.php?title=Installing_the_OpenGrADS_Bundle1: mv Contents $HOME/opengrads2: export PATH=$HOME/opengrads:$PATHnote: in order to e