西南科技大学oj答案:IDL直接读取ENVI格式数据代码[修订版]

来源:百度文库 编辑:偶看新闻 时间:2024/04/30 00:11:15

[资源共享] IDL直接读取ENVI格式数据代码[修订版]

在帖子http://bbs.esrichina-bj.cn/ESRI/thread-44373-1-1.html中楼主提供了一个非常好的读取ENVI格式数据的函数read_envi_image.pro,该函数在读取无扩展名ENVI数据时报错,为方便大家使用现进行了完善并添加了部分中文注释。
      ;
; Copyright (c) 2003,Institute of Photogrammetry and Remote Sensing, ;(IPF),
; Technical University of Vienna. Unauthorised reproduction prohibited.
;
;+
; NAME:
; read_envi_file
;
; PURPOSE:
; IDL program, which reads standard ENVI image files (*.img).
;
; CALLING SEQUENCE:
; read_envi_file, infile, img, xs, ys, type,offset
;
; INPUTS:
; infile - 传入数据文件名
;
; OPTIONAL INPUTS:
; None
;
; KEYWORD PARAMETERS:
; None
;
; OUTPUTS:--这些均是返回值,注意是位置参数
; img - ENVI的图像文件;
; xs - 列号;
; ys - 行号;
; type - 数据类型代码
; offset - 头文件偏移量;
; mapinfo - map及地理坐标信息;
;
; EXAMPLE:  read_envi_image, 'C:\Program Files\ITT\IDL708\products\envi46\data\ddd', data;;;data即读取的数据值
;
; MODIFICATION HISTORY:
; Written by: Carsten Pathe, cp@ipf.tuwien.ac.at
; Date: 25.08.2003
;
;   Modified By DYQ→2009-8-18
;   修正:数据无扩展名会读取错误
;         简化数组创建代码
;-
PRO read_envi_image, infile, img, xs, ys, type, offset, mapinfo
  COMPILE_OPT idl2
  image = infile
  ;Add by DYQ-for数据文件无扩展名的bug
  pointExist = STRPOS(infile,'.')
  IF pointExist[0] NE -1 THEN BEGIN
    header = strsplit(infile,'.',/extract)
    header = header[N_ELEMENTS(header)-2]+'.hdr'
  END ELSE header = infile + '.hdr'
  ;打开解析头文件
  OPENR, unit, header, /get_lun
  header_line = ''
  ;文件未读取完之前一直循环
  WHILE NOT EOF(unit) DO BEGIN
    READF, unit, header_line
    tmp = strsplit(header_line[0], '=', /extract)
    header_keyword = strsplit(tmp[0], ' ', /extract)
    ;解析头文件中的信息
    IF header_keyword[0] EQ 'samples' THEN xs = LONG(tmp[1])
    IF header_keyword[0] EQ 'lines' THEN ys = LONG(tmp[1])
    IF header_keyword[0] EQ 'header' THEN offset = LONG(tmp[1])
    IF header_keyword[0] EQ 'data' THEN type = LONG(tmp[1])
    ;如读取map信息则解析map信息
    IF header_keyword[0] EQ 'map' THEN BEGIN
      mapinfo_tmp=strsplit(tmp[1],'{',/extract)
      mapinfo_tmp=strsplit(mapinfo_tmp[1],',',/extract)
      mapinfo={ulx:0.,uly:0.,spacing:0.}
      mapinfo.ulx=mapinfo_tmp[3]
      mapinfo.uly=mapinfo_tmp[4]
      mapinfo.spacing=mapinfo_tmp[5]
    ENDIF
  ENDWHILE
  ;关闭头文件
  CLOSE,unit & FREE_LUN, unit
  ;modified by dyq
  img = MAKE_ARRAY(xs,ys,type = type)
  ;读取数据文件中的数据
  OPENR, unit,image, /get_lun
  POINT_LUN, unit, offset
  READU, unit, img
  FREE_LUN, unit
END