;----------------------------------------------------------------- ; lee_ref.pro ; ; 06/2001 by Andreas Reigber ;------------------------------------------------------------------ ; Refined Lee filter for multiplicative image noise using local ; statistics. Better preservation of image structures than the ; original Lee filter by the use of directional windows. Very useful ; for speckle reduction of SAR amplitude images. More in: ; J.S.Lee: "Refined filtering of image noise using local statistics" ; Comp. Graph. Image Process. Vol. 15, pp. 380-389, 1981 ; ; Usage: ; res = lee_ref(image,size) ; ; image = floating point image amplitude ; size = size of the filter window ; ; Keywords: ; ; looks = Number of looks (default = 1) ; threshold = Threshold for switching to rectangular boxes ; ;------------------------------------------------------------------ ; This software has been released under the terms of the GNU Public ; license. See http://www.gnu.org/copyleft/gpl.html for details. ;------------------------------------------------------------------ function lee_ref,amp,smm,LOOKS=looks,THRESHOLD=threshold if not keyword_set(looks) then looks=1 if not keyword_set(threshold) then threshold=0.5 siz = size(amp) anz_rg = siz[1] anz_az = siz[2] delta = (smm-1)/2 ol = delta / 2 + 1 thr = threshold sig2 = 1.0/looks sfak = 1.0+sig2 out = amp-amp fbox = mbox(smm,anz_rg) sbox = intarr(9,delta^2) for i=0,2 do begin sbox[i*3,*] = mbox(delta,anz_rg)+ol*(i-1)*anz_rg-ol sbox[i*3+1,*] = mbox(delta,anz_rg)+ol*(i-1)*anz_rg sbox[i*3+2,*] = mbox(delta,anz_rg)+ol*(i-1)*anz_rg+ol endfor dummy = sbox[4,*] sbox[4,*] = sbox[8,*] sbox[8,*] = dummy sbox[0,*] = mbox(delta,anz_rg)-ol sbox[1,*] = mbox(delta,anz_rg)-ol*anz_rg-ol sbox[2,*] = mbox(delta,anz_rg)-ol*anz_rg sbox[3,*] = mbox(delta,anz_rg)-ol*anz_rg+ol sbox[4,*] = mbox(delta,anz_rg)+ol sbox[5,*] = mbox(delta,anz_rg)+ol*anz_rg-ol sbox[6,*] = mbox(delta,anz_rg)+ol*anz_rg sbox[7,*] = mbox(delta,anz_rg)+ol*anz_rg+ol sbox[8,*] = mbox(delta,anz_rg) cbox = intarr(8,smm*(delta+1)) for i=0,smm-1 do begin cbox[0,i*(delta+1):(i+1)*(delta+1)-1] = indgen(delta+1)+(i-delta)*anz_rg cbox[2,i*(delta+1):(i+1)*(delta+1)-1] = indgen(delta+1)*anz_rg+(i-delta) cbox[4,i*(delta+1):(i+1)*(delta+1)-1] = -indgen(delta+1)+(i-delta)*anz_rg cbox[6,i*(delta+1):(i+1)*(delta+1)-1] = -indgen(delta+1)*anz_rg+(i-delta) endfor pos = 0 for i=0,smm-1 do begin cbox[1,pos:pos+i] = (delta-indgen(i+1))*anz_rg + (i-delta) cbox[3,pos:pos+i] = (delta-indgen(i+1))*anz_rg - (i-delta) cbox[5,pos:pos+i] = (-delta+indgen(i+1))*anz_rg + (i-delta) cbox[7,pos:pos+i] = (-delta+indgen(i+1))*anz_rg - (i-delta) pos = pos + (i+1) endfor for i=delta,anz_rg-delta-1 do begin for j=delta,anz_az-delta-1 do begin pos = i+j*anz_rg stat = total(amp[pos+sbox],2)/9 grad = max(abs(stat/stat[8]-1.0),dir) if grad gt thr then win = cbox[dir,*] else win=fbox box = amp[pos+win] mean = avg(box) mean2= mean^2 z = total((box-mean)^2)/(n_elements(win)-1) z = ((z+mean2)/sfak)-mean2 > 0 out[pos] = mean + (amp[pos]-mean)*(z/(mean2*sig2+z)) endfor endfor return,out end