LearningENVI&IDL分享 http://blog.sciencenet.cn/u/dongyanqing Learning ENVI&IDL

博文

编写高效率的IDL程序

已有 3641 次阅读 2011-9-18 12:59 |系统分类:科研笔记| 程序, IDL, 高效率

程序的效率问题,在大数据或复杂运算的时候是不能忽略的。但在IDL程序的编写方式上,不能按照常规的循环for依次处理方式写,简单归纳下,提高效率的运行的写法注意下面两种方式。

1、 尽量避免或少用循环

2、 多用whereHistogram

说起来很容易,但实际写的时候一定要多斟酌斟酌。

举例1:对2000*2000的数组中大于100的值进行累加

代码:

PRO TEST_TIME 
  ;
  a =  DIST(2000,2000)
  sum = 0.
  sum1 = 0.
  start = SYSTIME(1)
  FOR i=0L,N_ELEMENTS(a)-1L DO BEGIN
    IF(a[i] GT 100.0) THEN  BEGIN
      sum = sum +a[i]      
    ENDIF
  ENDFOR
  fortime = SYSTIME(1)-start
  PRINT,'for time:',fortime
 
  i=0L
  start = SYSTIME(1)  
  WHILE i LT N_ELEMENTS(a)-1L DO BEGIN
    IF(a[i] GT 100.0)THEN sum = sum +a[i]    
    i++
  ENDWHILE
  whiletime = SYSTIME(1)-start
  PRINT,'while time:',whiletime

  start = SYSTIME(1)
  sum = TOTAL(a * (a GT 100.0))
  funtime = SYSTIME(1)-start
  PRINT,'function time:',funtime
  ;倍数
  print,'result for VS function:', fortime/funtime
  print,'result while VS function:',whiletime/funtime  
END

行后的输出:



到差别了吧,循环比函数直接运算慢至少一个数量级!

 

举例2:对一个图像中的特定值,若存在,则以该像素为中心,特定半径内的元素统一修改为某值。

IDL自带的一个图像为例,将数据值等于142的赋为0.

源码如下:

;≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌
;
;计算两个点的距离
;
function CalDistance, point1, point2
  compile_opt idl2    ;
  Return, SQRT((point1[0]-point2[0])^2+(point1[1]-point2[1])^2)
end
;≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌
;
;搜索当前坐标周围Distance内的下标,注意输入x和y方向的坐标范围xRange和yRange
;
function calIdxInDistance, curLoc, distance;,xRange,yRange
  ;初始化临时下标
  suitLoc = [0,0]
  ;循环一次,计算矩形范围内的符合要求下标
  for xLoc = curLoc[0]-distance, curLoc[0]+distance do begin
    for yLoc = curLoc[1]-distance, curLoc[1]+distance  do begin
      if calDistance(curLoc, [xLoc,yLoc]) LE distance then suitLoc = [[suitLoc],[xLoc,yLoc]]
    endfor
  end
  ;
  return, suitLoc[*,1:(N_Elements(suitLoc)/2-1)]
 
end
;≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌≌
;
;测试即调用主函数
pro test_process

  ;原数据
  file = FILEPATH('rbcells.jpg', $
    SUBDIRECTORY = ['examples', 'data'])
  READ_JPEG, file, data
  ;
  if size(data,/n_dimensions) ne 2 then return
  dims = size(data,/dimension)
  startTime  = systime(1)
 
  ;当值等于142时,半径5内的元素赋值为0
  eqValue = 142
  repValue = 0
  ;符合要求的坐标
  suitIdx = calIdxInDistance([0,0],5)
  ;新数据
  nData = data
  idxs = where(data eq eqValue,count)
  for curIdx =0,count-1 do begin
 
    nSuit = suitIdx
    ;转换为二维坐标
    suitLoc = ARRAY_INDICES(data, idxs[curIdx])
    nSuit[0,*]= nSuit[0,*]+suitLoc[0]
    nSuit[1,*]= nSuit[1,*]+suitLoc[1]
    ;下标要在数组自身范围内
    nSuit[0,*] = dims[0] < nSuit[0,*] > 0
    nSuit[1,*] = dims[1] < nSuit[1,*] > 0
    ;符合要求的位置赋值
    nData[nSuit[0,*],nSuit[1,*]] = repValue
  endfor
  ;输出花费时间
  print,systime(1) - startTime
  ;
  ns = dims[0]
  nl = dims[1]
  new = data
  new1 = data
  startTime = systime(1)
  for i=0,ns-1 do begin
    for j=0,nl-1 do begin
      if new[i,j] eq eqValue then begin
        for m=0,ns-1 do begin
          for n=0,nl-1 do begin
            if (m-i*1L)^2+(n-j*1L)^2 le 25 then new1[m,n]=0
          endfor
        endfor
      endif
    endfor
    
  endfor
  print,'common time:',systime(1) - startTime
 
  window,1,xSize = dims[0]*2,ySize = dims[1]
  tvscl,Data,0
  tvscl,ndata,1
  window,2,xSize = dims[0]*2,ySize = dims[1]
  tvscl,new,0
  tvscl,new1,1
 
end

运行后输出

输出图像如下:


结果一样,但循环所花费的时间多了1千多倍,所以方法或函数选择一定要注意下。

 



https://blog.sciencenet.cn/blog-344887-487680.html

上一篇:Envi下的特殊值替换扩展
下一篇:2011中国用户大会遥感专题共享之“IDL新特性与遥感GIS系统一体化
收藏 IP: 223.214.189.*| 热度|

2 陈小润 liguoshuai

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-4-25 09:59

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部