nc文件的读取与处理

1、.nc 文件的查看和处理2013-08-04 12:45阅读：2,168 最近在学习冰后回弹模型（http:/www.atmosp.physics.utoronto.ca/peltier/data.php），里面给出的数据为nc（netcdf）格式的，读取是个问题，在言深深的帮助下，成功读取，在这里将其总结出来，以备后用也希望可以方便有用之人。有两种方式可以进行读取：1、Excel的一个插件，NetCDF4Excel可以完成相关的查看工作，需要的可以进行下载( . 416233&third=15 （华为网盘太不靠谱，不知道百度靠谱到合适）安装之后，我们可以得到，一个快捷方式，

2、如图：打开之后，我们可以看到下面的情形：然后打开nc文件，我们就可以进行查看操作：2、通过Grads进行相关的读取，主要是在言深深同学的帮助下完成，深深的感谢！第一步，安装相关的grads软件，这个这里不做介绍，相关可以论坛网站；第二步，打开grads软件，如图，回车之后的那个窗口千万不能关，那个窗口是图像显示的，关闭之后输入命令之后直接就退出了，刚开始我老犯这个低级错误（对菜鸟而言）。第三步，可以通过两个命令简单的看到数据的头文件包含信息：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程序将1.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：这里的参数要与自己的数据一致，比如mno，以及最后的输出，有的是1°，有的是2.5°。然后我们就可以得到相关的数据文件如下：如何用ncdump查看.nc文件由于用grads打开.nc文件时遇到问题，提示SDF file has no discernable X coordinate，所以在网上找了好久解决的办法，几乎所有办法都指向一个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文件？谢谢 复制链接 |关注本帖取消最新回复取消置顶回复取消最新编辑duan2008 duan2008 当前离线 积分1023贡献精华在线时间小时注册时间2012-5-7最后登录1970-1-1窥视卡 雷达卡 电梯直达 楼主 发表于 2012-11-30 19:12:24 | 只看该作者 |倒

7、序浏览 |阅读模式 |取消关注该作者的回复登录后查看更多精彩内容您需要 登录 才可以下载或查看，没有帐号？立即注册 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空间 新浪微博 腾讯微博 人人网 更多 收藏2 分享 淘帖0 赞0

9、踩0 给自己设置一个个性的签名吧点击前往GrADS官网新地址回复 使用道具 举报 提升卡置顶卡沉默卡喧嚣卡变色卡抢沙发千斤顶显身卡river river 当前离线 积分24285贡献精华在线时间小时注册时间2012-8-28最后登录1970-1-1窥视卡 雷达卡 沙发 发表于 2012-11-30 20:11:19 | 只看该作者 |取消关注该作者的回复'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'还不行就把错误提示发上来啊，这样才能对症下药啊给自己设置一个个性的签名吧点击前往GrADS官网新地址回复 支持 反对 使用道具 举报 显身卡duan2008 duan2008 当前离线 积分1023贡

11、献精华在线时间小时注册时间2012-5-7最后登录1970-1-1窥视卡 雷达卡 板凳  楼主| 发表于 2012-11-30 20:30:14 | 只看该作者 |取消关注该作者的回复river 发表于 2012-11-30 20:11                                 登录/注册后可看大图'reinit''sdfopen D:TDDOWNLOADhgt.1981.nc&#

12、39;'sdfopen D:TDDOWNLOADhgt.1982.nc'谢谢啊，还是不行。错误如下 登录/注册后可看大图 1.JPG (23.67 KB, 下载次数: 0)下载附件  保存到相册 2012-11-30 20:27 上传登录/注册后可看大图 2.JPG (26.74 KB, 下载次数: 0)下载附件  保存到相册 2012-11-30 20:29 上传我在想是时间的设置问题吗？因为第一个T=1是1981年1月1日，而第二个文件是1982年1月1日，而且有可能是365天和366天，在求和的时候有数据缺失？ 谢谢啊。给自己设置一个个性的签名吧点击前

13、往GrADS官网新地址回复 支持 反对 使用道具 举报 显身卡river river 当前离线 积分24285贡献精华在线时间小时注册时间2012-8-28最后登录1970-1-1窥视卡 雷达卡 地板 发表于 2012-11-30 21:18:42 | 只看该作者 |取消关注该作者的回复duan2008 发表于 2012-11-30 20:30                                 登录/注册后可看大图谢谢啊，

14、还是不行。错误如下把那句set t 1去掉，把define那句改成define hgt12=hgt.1(t=1)+hgt.2(t=1)试试给自己设置一个个性的签名吧点击前往GrADS官网新地址回复 支持 反对 使用道具 举报 显身卡duan2008 duan2008 当前离线 积分1023贡献精华在线时间小时注册时间2012-5-7最后登录1970-1-1窥视卡 雷达卡 5#  楼主| 发表于 2012-11-30 21:57:08 | 只看该作者 |取消关注该作者的回复river 发表于 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天的，谢谢。给自己设置一个个性的签名吧点击前往GrADS官网新地址回复 支持 反对 使用道具 举报 显身卡river river 当前离线 积分24285贡献精华在线时间小时注册时间2012-8-2

16、8最后登录1970-1-1窥视卡 雷达卡 6# 发表于 2012-11-30 22:31:41 | 只看该作者 |取消关注该作者的回复duan2008 发表于 2012-11-30 21:57                                 登录/注册后可看大图这样可以计算出一个值，但是对于366天的值怎么计算了，把hgt12存入一个NC文件，但是T是366天的，谢谢。用循环可以吧，加几句'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'但是多出来一天就不好办了啊给自己设置一个个性的签名吧点击前往GrADS官网新地址回复 支持 反对 使用道具 举报 显身卡xjliutao xjliutao 当前离线 积分56贡献精华在线时间小时注册时间2013-3-21最后登录1970-1-1窥视卡 雷达卡 7# 发表于 2013-3-21 20:18:20 | 只看该作者 |取消关注该作者的回复学习了，我一个都没有处理好 你都可以处理多个了论坛里

19、常常有人对nc这种内置ctl的文件数据处理 提出问题，这里总结下常见的一些问题及对策。nc数据是一类内置ctl的文件，也就说没有独立的ctl文件，那么打开命令也不是常见的open，而是sdfopen命令，那么要处理这类数据就得知道变量名、变量层次分布等信息，用q ctlinfo命令就可以查看全部的ctl了。为什么提查看ctl 呢，很多人在处理这类数据的时候会想当然的认定一个变量名，缺省值、甚至精度都会不注意，那么无论怎么提取，出来的图都不会对的。下面总结下容易出的问题：1、从nc文件提取出二进制数据（fwrite）设定需要注意经纬度，层次，时间是否是你需要的部分，尤其是用循环提取数据的朋友，尤

20、其要注意。有的时候设定经纬度（set lat 0 359.8）出了问题，可是尝试定义x y格点（set x 1 180）的方法试试，实在不行，可以用循环来输出不同的经纬度的数值。2、注意提取数据的grads命令'set gxout fwrite''set fwrite c:ss00.grd'.'disablefwrite'最后一行的这个命令，如果忘了加就容易数据不正常，还有的朋友，一个gs文件提取了两个数据文件，两组这类命令中间那个disable fwrite丢掉了，数据就彻底乱套了。3、提取的数据，写ctl，需要注意与nc的ctl对应，包括缺测

21、值，精度，这个问题只要照着nc的ctl写，层次格点什么的写成你数据对应的就好了。4、数据提取了，ctl写好了，画画图看是不是与nc的一致，设定的时间、层次，图是不是不一样，因为用循环出数据的朋友容易犯 写了好几个时间点，数据是一样的错误。5、还有的朋友需要进一步提取这个grd或者dat文件，这就涉及到二进制数据的读取与操作问题了。二进制数据是无格式的数据，读取与写入都不需要格式说明，里面的数据无法直接看到，其内部无非是0和1，而且数据之间没有间隔，读取写入只要read（文件号） write（文件号）就好读取种要注意，grads写入的数据循环从内到外依次为x 、y、z、变量、t，如果你需要确定某

22、个数据的位置，那么，x是从西到东，y是从南到北，z是从地面到高空，t是从前到后的顺序。另外此类文件读取目前论坛里主要有两种方法，一种是直接声明form=binary的，这种形式简单，如果数据是常规的顺序读取，选择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另一种是无格式直接读取的方式，无格式直接读取主要是用于大量的数据需要处理，你只需要其中部分，可以用记录的跳跃选取你要的数据，如下面这段：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最后总结下检查错误的思路：1、数据没提取出来，先检查gs是否没有写disable fwrite，或者文件名.grd写错了；2、数据出来后，ctl写好，画图与nc对比，如果发现，出来的图上只有个别区域有大

25、条粗线，整幅图都是黑的，说明你的ctl里缺测值或者精度不够对应，q ctlinfo 查询下；3、图对应上了，改变设定时间或层次，没有变化，说明，你提取数据的gs错了，层次或者时间都是输出的同一个点；4、前三步都没有问题，恭喜你，数据提取成功了。Fortran编程的时候如果报错说endof file during read，很可能是你文件名错了，或者你的你读取数据的顺序不是按照grads的数据顺序读的，检查你的各个循环的次数，及循环的嵌套顺序是否正确。5、如果你是顺序的读取数据的话，建议你选择binary的方式读取，犯错几率要小；6、如果以上，还是报错，剩下的原因就可能是你马虎导致微小的错误了，

26、你要认真的从头到尾检查程序，不行就让别人给你检查，因为有些问题太细微了，自己已经习惯了。求问，我想计算两个NDVI数据集中相同一个地区2001年到2006年12月-2月平均值的差值图，结果老是出错，不知大家可否帮忙看一下程序，先谢谢大家了！这是程序：'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'结果老是说“Error in gagchk: axis sizes are not the sameOperation err

30、or:  Incompatable grids ”但我的两个NDVI数据集都是0.5X0.5的，不知道是什么原因？不过gimms数据集是从lat -61.75  lon -179.75开始; modis数据集是从lat 0  lon -20 开始，不知是否是这个缘故？如果是的话，应该怎么解决呢？由于本人新手，都没有积分的，先大大的谢谢大家！=07/29/2014解决方法因为两个NDVI数据的起始网格点不同，需要转换到同一的网格点上，运用函数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上显示没有此函数，后来重新装了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