编写高效率的IDL程序

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

1、 尽量避免或少用循环

2、 多用where和Histogram；

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

举例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

运行后的输出

IDL> test_time % Compiled module: TEST_TIME. for time:       1.4510000 while time:       2.0430000 function time:     0.054000139 result for VS function:       26.870301 result while VS function:       37.833235

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

 

举例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

 

运行后输出：

% Compiled module: TEST_PROCESS.       0.10999990 common time:       116.80900

 

输出图像如下：

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