麻辣GIS微信平台

更多 GIS 干货

微信关注不错过

IDL遥感应用入门(14):IDL读取netcdf数据(.nc)

netcdf数据也是遥感应用中经常见到的一种图像格式,关于这种影像格式的一些介绍,你可以看这篇文章:关于netcdf文件(.nc)的一点概述,在这篇文章中,有一个pdf文件,专门介绍IDL读取netcdf的,本文还是以实际实例为基础,对IDL读取netcdf做一个简单的介绍。

一些说明

本文所用到的netcdf数据,你可以从这里下载的到。

此处内容已被隐藏
为了防止资源恶意爬取被和谐,关注微信公众号即可免费解锁
麻辣GIS公众号 微信扫一扫关注

1. 微信扫码关注公众号,回复 资源下载 获取验证码

2. 在下方输入验证码并提交查看隐藏内容:

也可在微信搜索公众号“麻辣GIS”关注获取

详细代码

;Author:Sailor
;2013-12-25
PRO Course_13

;定义文件路径
MyRootDir='D:\10\'

;遍历文件夹
filearr = file_search(MyRootDir,'*.nc',count=num);
  FOR fileindex=0,num-1,1 DO BEGIN
    nid = ncdf_open(filearr[fileindex], /nowrite )

    ; inquire about this file; returns structure
    file_info = ncdf_inquire( nid )

    ; print out the dimensions of this file
    FOR dimid=0, file_info.ndims -1 DO BEGIN
      ncdf_diminq, nid, dimid, name, size
      print, ' ---> dimension ' + name  + ' is: ', size 
    ENDFOR

    FOR varid=0, file_info.nvars-1 DO BEGIN
      ; inquire about the variable; returns structure
      var = ncdf_varinq( nid, varid )
      print,var
      print,'========================'
      ;read all attributes
      FOR var_att_id=0,var.natts -1 DO BEGIN
        att_name = ncdf_attname( nid, varid, var_att_id )
        print,att_name
        ncdf_attget, nid, varid, att_name, tematt
        print,string(tematt)
      ENDFOR
    ENDFOR

    ;read sst
    sstid = ncdf_varid(nid, 'analysed_sst')
    ncdf_varget, nid, sstid, sst
    sst=sst*0.01+273.150
    landindex=where(sst eq min(sst))
    sst[landindex]=max(sst)+1
    Img=image(sst,rgb_table=39,title=filearr(fileindex),grid_units=1,POSITION=[0.1,0.2,0.9,0.9])
    xaxis=axis('X',LOCATION=[0,0],AXIS_RANGE=[-180,180],MINOR=0, MAJOR=19,COORD_TRANSFORM=[-180,360.0/4096.0],title='Longitude(°)')
    yaxis=axis('Y',LOCATION=[0,0],AXIS_RANGE=[-90,90],MINOR=0, MAJOR=7,COORD_TRANSFORM=[-90,180.0/2048.0],title='Latitude(°)')
    c1 = COLORBAR(TARGET=Img, ORIENTATION=0,TITLE='(K) ',POSITION=[0.1,0.1,0.9,0.15])

  ENDFOR
END

效果图

代码解释

;定义文件路径
MyRootDir='D:\10\'

这里定义了文件路径,当然你要改成你自己的。

filearr = file_search(MyRootDir,'*.nc',count=num);
FOR fileindex=0,num-1,1 DO BEGIN

这句代码段是对上面的路径做一个遍历,找出所有的.nc文件,之所以这么写是为了方便批处理,因为不可能只处理一个文件。不过如果仅仅是做测试,这里可以改为num-num只处理一个文件测试效果,效果满意再做批处理。关于file_search函数的详细使用,你可以看这篇文章:IDL使用file_search函数遍历某个文件夹下的所有文件

; inquire about this file; returns structure
file_info = ncdf_inquire( nid )

这里是读取了整个nc文件的文件信息,返回的是一个结构体。详细参数有:

Ndims:文件的维数 Nvars:netcdf文件定义的变量 Ngatts:netcdf文件的全局属性 RecDim:未限制维数的ID,如果没有则为-1

; print out the dimensions of this file
FOR dimid=0, file_info.ndims -1 DO BEGIN
  ncdf_diminq, nid, dimid, name, size
  print, ' ---> dimension ' + name  + ' is: ', size 
ENDFOR

这段代码是读取文件的的所有维的信息,包括名字和大小。

FOR varid=0, file_info.nvars-1 DO BEGIN
  ; inquire about the variable; returns structure
  var = ncdf_varinq( nid, varid )
  print,var
  print,'========================'
  ;read all attributes
  FOR var_att_id=0,var.natts -1 DO BEGIN
    att_name = ncdf_attname( nid, varid, var_att_id )
    print,att_name
    ncdf_attget, nid, varid, att_name, tematt
    print,string(tematt)
  ENDFOR
ENDFOR

读取所有变量和变量的属性信息。

;read sst
sstid = ncdf_varid(nid, 'analysed_sst')
ncdf_varget, nid, sstid, sst
sst=sst*0.01+273.150
landindex=where(sst eq min(sst))
sst[landindex]=max(sst)+1
Img=image(sst,rgb_table=39,title=filearr(fileindex),grid_units=1,POSITION=[0.1,0.2,0.9,0.9])
xaxis=axis('X',LOCATION=[0,0],AXIS_RANGE=[-180,180],MINOR=0, MAJOR=19,COORD_TRANSFORM=[-180,360.0/4096.0],title='Longitude(°)')
yaxis=axis('Y',LOCATION=[0,0],AXIS_RANGE=[-90,90],MINOR=0, MAJOR=7,COORD_TRANSFORM=[-90,180.0/2048.0],title='Latitude(°)')
c1 = COLORBAR(TARGET=Img, ORIENTATION=0,TITLE='(K) ',POSITION=[0.1,0.1,0.9,0.15])

最后出图,其中这里数据处理所做的处理来自文件本身,也就是在上面读取所有变量和变量的属性信息所读取出的数据。

一些说明

这里只是对IDL读取netcdf数据的读取做一个示例,如有问题,欢迎交流。

相关阅读

麻辣GIS-麻辣GIS

作者:

GIS爱好者,学GIS,更爱玩GIS。

声明

1.本文所分享的所有需要用户下载使用的内容(包括但不限于软件、数据、图片)来自于网络或者麻辣GIS粉丝自行分享,版权归该下载资源的合法拥有者所有,如有侵权请第一时间联系本站删除。

2.下载内容仅限个人学习使用,请切勿用作商用等其他用途,否则后果自负。

手机阅读
公众号关注
知识星球
手机阅读
麻辣GIS微信公众号关注
最新GIS干货
关注麻辣GIS知识星球
私享圈子

已有 38 条评论

  1. 亖木子
    很有帮助。我想知道,如果我想把每一层单独存储成tif,后面怎么写。我写出来的没有坐标信息,而且纬度是颠倒的。
      1. 桑榆非晚
        桑榆非晚
        同学,你的问题解决了吗?我现在也是在为这个问题烦恼。
  2. hellosyndy
    hellosyndy 2#
    请问,在给nid赋值时nid = ncdf_open(filearr[fileindex], /nowrite ),运行会报Attempt to subscript FILEARR with FILEINDEX is out of range错误怎么解决?
  3. hikari
    hikari 3#
    请问 sst[landindex]=max(sst)+1 这部分有什么意义
      1. 麻辣GIS
        麻辣GIS 作者
        landindex表示陆地,这里表示将陆地的值统一成所以数据的最大值+1,可以将陆地对比显示出来。
  4. Linda
    Linda 4#
    我在运行这个程序时,出现Attempt to call undefined procedure/function: 'READ_WSPD'.想问一下是怎么回事
  5. 花花
    花花 5#
    我也出现了和楼上一样的错误,请问是怎么解决的?
      1. 麻辣GIS
        麻辣GIS 作者
        IDL版本是多少?
          1. 花花
            花花 回复 @麻辣GIS
            5.1
              1. 麻辣GIS
                麻辣GIS 作者 回复 @花花
                我的程序是使用IDL8.2.3写的,可能与此有关。
                  1. 花花
                    花花 回复 @麻辣GIS
                    嗯嗯 谢谢
  6. ShirleyYuan
    ShirleyYuan 6#
    你好,我的运行的时候出现下面是怎么回事呀?
    IDL Version 8.3, Microsoft Windows (Win32 x86_64 m64). (c) 2013, Exelis Visual Information Solutions, Inc.
    Installation number: 6046848432.
    Licensed for use by:

    % Restored file: ENVI.
    % Loaded DLM: PNG.
    % Loaded DLM: HPGRAPHICS.
    % Loaded DLM: XML.
    The specified language is not available in the language catalog.
    The available languages are:

    eng
    ENVI> nc文件读取
    % Compiled module: COURSE_13.

    nc文件读取
    ^
    % Syntax error.
    ENVI>
      1. 麻辣GIS
        麻辣GIS 作者
        目测是不支持中文
          1. 时光
            时光 回复 @麻辣GIS
            您好,我想问一下,sstid = ncdf_varid(nid, 'analysed_sst'),如果是自己的程序,analysed_sst应该怎么换?谢谢你
              1. zww
                zww 回复 @时光
                你好,请问这个问题解决了吗
              2. eyes
                eyes 回复 @时光
                相同的问题啊,请问如何解决呢?
                  1. 麻辣GIS
                    麻辣GIS 作者 回复 @eyes
                    不同的数据,换成数据对应的不就行了。
                  2. 6idle
                    6idle 回复 @eyes
                    您好请问在哪里换?
              3. 不用谢
                不用谢 回复 @时光
                你好,解决了吗,分享一下呗
              4. 6idle
                6idle 回复 @时光
                “ncdf_varget变量查询失败”怎么办?换成对应数据的话,这个需要在哪里查看呀求解?谢谢!
  7. 宝恋枫
    宝恋枫 7#
    请问我用IDL 读出来的数据是颠倒的,这个该如何解决
      1. 麻辣GIS
        麻辣GIS 作者
        把数组转一下。
  8. Susanna
    Susanna 8#
    请问您用的什么气象数据?
      1. 麻辣GIS
        麻辣GIS 作者
        这个测试数据是老师给的,很久之前,忘记是什么的了。
  9. happy
    happy 9#
    您好,我运行的结果总提示:NCDF_VARGET: Variable Inquiry failed, -1 is not a valid variable id.为什么?
      1. 麻辣GIS
        麻辣GIS 作者
        IDL版本是?
      2. naonao
        naonao
        你好,请问你的这个问题解决了吗
  10. 夏天里游泳
    夏天里游泳 10#
    我和一楼有同样的问题。我想知道,如果我想把每一层单独存储成tif,后面怎么写。我写出来的没有坐标信息,而且纬度是颠倒的。
    还有,sstid = ncdf_varid(nid, 'analysed_sst'),如果是自己的程序,analysed_sst应该怎么换?

    还有,我想加入你那个群的,但是付款需要用财付通,我没有开通,所以加不了。
  11. XJ
    XJ 11#
    强烈建议出一个IDL 读取OMI he5 数据的教程
  12. 初学者
    初学者 12#
    编译运行的时候错误提示指在了1,DO 的位置,for fileindex=0,num-1, 1, DO BEGIN
      1. 不用谢
        不用谢
        sstid = ncdf_varid(nid, 'analysed_sst'),如果是自己的程序,analysed_sst应该怎么换?
          1. 麻辣GIS
            麻辣GIS 作者 回复 @不用谢
            你先看看ncdf_varid这个函数是做什么的,你就知道应该怎么换了
  13. 陈小白
    陈小白 13#
    今天才看到博主发的帖,我想读取在cryosat下载的nc数据,在sstid = ncdf_varid(nid, 'analysed_sst')上有点不懂,怎么看自己用的数据的变量名呢?比如在cryosat上下载的nc数据他的变量名是什么呢?望解答;谢谢啦
  14. VVVV
    VVVV 14#
    你好,我想请问一下可以实现.nc数据转成hdf5数据吗
  15. LucyChen
    LucyChen 15#
    您好,
    通过您的教程我读取了一个自己的.nc数据并成功抓取里面需要的数值。其中一个WindSpeed数值我绘图成功了,然而另一个nrcs数值没有成功。我发现这两个dimension不一样,其中Windspeed 是 float[254,194], 而NRCS 是 float [2,254,194]:

    IDL> size(NRCS)
    3 2 254 194 4 98552
    IDL> size(WindSpeed)
    2 254 194 4 49276


    有什么办法对NRCS的数值也进行绘图呢,或者降维?
      1. 麻辣GIS
        麻辣GIS 作者
        我理解降维画就可以了

留言板(小编看到第一时间回复)