Lab exercise 6, due right before the next class
GEO5083: Remote Sensing Image Processing and Analysis, UTSA

 

Student name: ______________

 

 

Spatial filtering  

Purpose

 

            This lab is to practice one type of image enhancement techniques: spatial filtering. 

           

 

Step 1. Preparation

 

    the data you will used for step 3 can be downloaded from \\129.115.25.240\XIE_misc\ees5083, under lab6 subdirectory. this image  (p33r36sep14_2000) is a subset of Path33/Row 36, acquired on September 14, 2000 in central New Mexico.

 

           

 Step 2. Spatial filtering: convolution based

           

        Open and load any Landsat image you used before. from the ENVI main menu, click Filter -> Convolution and Morphology. a new window looks like below.

Figure 1

 

    Actually, we even do not talk about the Morphology in the lecture, so we only use the Convolutions. If you put you mouse to the Convolutions of the above window, you will see a lot of methods we did talked about in the lecture such as: High pass, low pass, Laplacian, directional,  Sobel, and Roberts. Click any of the first 4, you will see a Editable Kernel (3 x 3), which is also called convolution mask in our lecture or moving window sometime. This kernel, mask, or window is editable (values can be changed) and size can be changed to 5 x 5, 7 x 7, 9 x 9, 11 x 11, ... through the Kernel Size. In our lecture, we use the 3 x 3 as an example, but you will practice other sizes in this lab. You can first use the default kernels in the above window. Click the Apply to File, it will lead you to select your image. In this exercise, you can only select one of the 6 bands image using the Spectral Subset (from the input file window as below Figure 2), remembering you even can select your Spatial Subset from this window. This two buttons (Spectral Subset and Spatial Subset) are very useful, I am sure most of you already used them in previous labs.

Figure 2

 

    Click OK in the window, and output result to memory. then you will have your first high pass image (or high frequency image). Open it in a new window. Link and compare with the original image (the same band) to see the high frequency areas have been captured in the frequency image. Now change your Kernel size to 5 x 5, click the Quick Apply in figure 1. you will get a Quick-Apply Input Band window (figure 3). In this image, you can directly select the original band you did for the 3 x 3 high pass. You also have the choice to do a Spatial Subset. Click OK, you will have a new high pass new window opened. This window even does not output to a memory file (I do not know where it is stored). Change Kernel Size to 7 x 7, or 9 x 9, ... click the Quick Apply, what happen? you can keep increase the Kernel size to see how the high frequency will change compare with the 3 x 3 window.

Figure 3

 

    You can practice the low pass, Laplacian, Directional in the same way and even replace the default Kernel values using the ones we talked in the lecture. Please refer to the Help to learn more about how to use them and refer to the book or lecture slides for the principles and algorithms of them.

 

    If you click Sobel and Roberts edge detector filters. they does not allow you to see the Kernel and edit them since they used the standard Kernels as you learned from the Lecture. I want you to run these two methods to the same band that you used for the previous works. do a quick comparison and evaluation of these two methods for edge detection.

 

 

Step 3. Spatial filtering: FFT transform

   

    Before you do this step, I suggest you to close all opened windows and delete all memory files. Under the FFT filtering of the Filter menu, you will see Forward FFT, Filter Definition, and Inverse FFT. the Forward FFT converts the spatial domain image to frequency domain image. Filter Definition is used to define a filter to remove some noises you can see in the frequency image created from the Forward FFT. Some noises are difficult to be identified or removed from spatial domain but easy to do so in the frequency domain. Inverse FFT will convert the frequency domain image back to spatial domain after some noises removed through Filter Definition step. Please refer to my image registration paper using FFT and inverse FFT for the same way but for different purpose.

    Now, open the image you downloaded in step 1, make a 321 as RGB image. You will see the image looks like figure 4 below: obviously these noises are from different bands. I want to you open each band individually to check which bands contribute noises. then you will do a Forward FFT, and then you will check to see how the frequency images response these noises on band by band basis.  You may do an analysis on how and where these noises may come from based on the Lecture.

Figure 4

 

    Now, if you click the Filter Definition, you will open the Filter Definition window (Figure 5). Click the Filter_Type, you will see several filters to be used to define a filter to remove the noises. "Pass" is passing the low-frequency, removing the high frequency, "Cut" is cutting the low-frequency and keeping the high frequency. In this case, I am not sure which one will work, it maybe none of them will work except the User Defined Cut might work. Once you defined your filter to remove the noise from the frequency domain, then you will use the Inverse FFT to convert back to spatial domain, then the noises should be removed. However, this work needs to make a annotation file required for the User Defined Cut. To develop a algorithm using IDL and ENVI to easily remove it could be a good class project.

        

Figure 5

 

 

Write a report about the procedure, answer some questions brought out in the lab, and do necessary analysis to your results. for the step 3, you need to show some details as well as the final noise-removed combination  image (such as 321 combination).

 

(two IDL codes for those who have extra energy and are interested in learning more about IDL programming)

 

PRO Lab6_1
; illustrating convolution in the spatial and frequency domains

envi_select, title='Choose multispectral image', fid=fid, dims=dims, pos=pos
if (fid EQ -1) THEN BEGIN
PRINT, 'cancelled'
RETURN
ENDIF

num_cols=(c=dims[2]-dims[1]+1)
num_rows=dims[4]-dims[3]+1
;num_bands=n_elements(pos)
;num_pixels=num_cols*num_rows

; pick out the center row of pixels
image=envi_get_data(fid=fid, dims=dims, pos=pos[0]) ;pos[0] is the first band, pos[1] the second band, ...
g=float(image[*, num_rows/2])

; define a FIR kernel of length m=5
h=[1,2,3,2,1]

; convolve in the spatial domain
window, 11
PLOT, convol(g,h,center=0)

; pad the arrays to c+m-1
g=[g, [0,0,0,0]]
hp=g*0
hp[0:4]=h

; convolove in the frequency domain
window, 12
PLOT, fft(c*fft(g)*fft(hp),1)

END

 

 

PRO Lab6_2
; displaying the power spectrum of an image band

envi_select, title='Choose multispectral image', fid=fid, dims=dims, pos=pos
if (fid EQ -1) THEN BEGIN
PRINT, 'cancelled'
RETURN
ENDIF

num_cols=(c=dims[2]-dims[1]+1)
num_rows=dims[4]-dims[3]+1
;num_bands=n_elements(pos)
;num_pixels=num_cols*num_rows

; get the image band from ENVI
g=envi_get_data(fid=fid, dims=dims, pos=pos[0]) ;pos[0] is the first band, pos[1] the second band, ...

; transform it
g_hat=fft(g)

; create a Gauss filter in the frequency domain
sigma=50
d=dist(num_cols, num_rows)
h_hat=exp(-d^2/sigma^2)

; output surface plot of filter as EPS file
thisDevice=!D.Name
set_plot, 'PS'
Device, Filename='fig4_1.eps', xsize=3, ysize=3, /inches, /encapsulated
shade_surf, shift(h_hat, num_cols/2, num_rows/2)
device, /close_file
set_plot, thisDevice

; multiply, do inverse FFT and returen result to ENVI
envi_enter_data, fft(g_hat*h_hat,1) ; low-pass
envi_enter_data, fft(g_hat*(1-h_hat), 1) ; high-pass



END