1 / 36

Chapter 2

Chapter 2. Image Restoration Analysis and applications of remote sensing imagery Instructor: Dr. Cheng-Chien Liu Department of Earth Sciences National Cheng Kung University Last updated: 10 April 2005. Introduction. Why need restoration ( 復原 ) ? Degradation  Restoration

august
Download Presentation

Chapter 2

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Chapter 2 Image Restoration Analysis and applications of remote sensing imagery Instructor: Dr. Cheng-Chien Liu Department of Earth Sciences National Cheng Kung University Last updated: 10 April 2005

  2. Introduction • Why need restoration (復原)? • Degradation  Restoration • Source of degradation • Digital image acquisition type • Platform • TFOV

  3. Radiometric correction • Radiometric correction 輻射校正 • Varies with sensors • Mosaics of images taken at different times  require radiometric correction • Influence factors • Scene illumination • Atmospheric correction • Viewing geometry • Instrument response characteristics

  4. Scene illumination • Sun elevation correction • Fig 7.3: seasonal variation • Normalize by calculating pixel brightness values assuming the sun was at the zenith on each date of sensing • Multiply by cosq0 • Earth-Sun distance correction • Decrease as the square of the Earth-Sun distance • Divided by d2 • Combined influence

  5. Atmospheric correction • Atmospheric correction • Atmospheric effects • Attenuate (reduce) the illuminating energy • Scatter and add path radiance • Combination • Haze compensation  minimize Lp • Band of zero Lp (e.g. NIR for clear water) • Path length compensation • Off-nadir pixel values are normalized to their nadir equivalents • Atmospheric correction in ENVI

  6. Atmospheric correction in ENVI • Flat Field Correction • Normalize images to an area of known “flat” reflectance • ROI: Locate a large, spectrally flat, spectrally uniform area in the AVIRIS data • Assumed to be composed of primarily atmospheric effects and the solar spectrum • The average AVIRIS radiance spectrum from the ROI is used as the reference spectrum, which is then divided into the spectrum at each pixel of the image. • The result is apparent reflectance data that can be compared with laboratory spectra

  7. Exercise 1 • Working with Cuprite Radiance Data • File → Open Image File • Subdirectory: c95avsub • File: cup95_rd.int • Display a Gray Image • Band 193 • Display a Color Image • Band 183, Band 193, and Band 207 • Link Two Displays • Extract Radiance Spectra • Tools → Pixel Locator • sample 590 and line 570 • Tools → Profiles → Z-Profile (Spectrum)

  8. Exercise 1 (cont.) • Working with Cuprite Radiance Data (cont.) • Options → Collect Spectra • Stonewall Playa 590 570 • Varnished Tuff 凝灰岩 435 555 • Silica Cap 494 514 • Opalite Zone with Alunite 531 541 • Strongly Argillized 硬頁岩質 Zone with Kaolinite 502 589 • Buddingtonite 水銨長石 Zone 448 505 • Calcite 260 613 • Options → Stack Plots • Plot Key • Edit → Data Parameters

  9. Exercise 2 • Flat Field Correction • Define an ROI over a spectrally flat area • Select Basic Tools → Preprocessing → Calibration Utilities → Flat Field. • When the Flat Field Calibration Parameters dialog appears, select the average spectrum to use for the flat field calibration by clicking on the desired ROI name in the column labeled Select Calibration from Regions. The selected ROI name appears at the bottom of the window in the text box labeled Selected Item. • Memory output. • Click OK to start the calibration. • Compare to the precalculated file: cup95_ff.int • Link display • Show z-profile • Edit plot parameter to fix the tick label of y-axis

  10. Atmospheric correction in ENVI (cont.) • Internal Average Relative Reflectance (IARR) • Normalize images to a scene average spectrum • Effective in an area where no ground measurements exist and little is known about the scene • Works best for arid areas with no vegetation. • Calculating an average spectrum for the entire AVIRIS scene and using this as the reference spectrum • Apparent reflectance is calculated for each pixel of the image by dividing the reference spectrum into the spectrum for each pixel

  11. Exercise 3 • Internal Average Relative Reflectance (IARR) • Select Basic Tools → Preprocessing → Calibration Utilities → IAR Reflectance. • When the Calibration Input File dialog appears, select the input file and perform any subsetting. • When the IARR Calibration Parameters dialog appears, select output to File or Memory. • Compare to the precalculated file: cup95_ia.int • Link display • Show z-profile • Edit plot parameter to fix the tick label of y-axis

  12. Atmospheric correction in ENVI (cont.) • Empirical Line Calibration • Force image data to match selected field reflectance spectra • Requires ground measurements and/or knowledge • Two or more ground targets are identified and reflectance is measured in the field (one light and one dark area) • The same two targets are identified in the AVIRIS images and average spectra are extracted for Regions of Interest • A linear regression is calculated between the field reflectance spectra and the image radiance spectra to determine a linear transform from radiance to reflectance for each band of the AVIRIS data set • Gains and offsets calculated in the regression are applied to the radiance spectra for each pixel to produce apparent reflectance on a pixel-by-pixel basis.

  13. Exercise 4 • Empirical Line Calibration • Select Basic Tools → Preprocessing → Calibration Utilities → Empirical Line → Compute Factors and Calibrate. • When the Empirical Line Input File dialog appears, select the input file and perform any subsetting. • Click OK. The Empirical Line Spectra dialog appears. • Collect image spectra and reference spectra, and pair spectra using the procedures described in the following sections. • Collecting Data (Image) Spectra • Collecting Field Spectra • Pairing Spectra • Completing the Calibration • Compare to the precalculated file: cup95_el.int • Link display • Show z-profile • Edit plot parameter to fix the tick label of y-axis • Calibrating Using Existing Factors

  14. Self test • Empirical Line Calibration • File → Open Image File • Subdirectory: c95avsub • File: cup95_rd.int • Use the spectral library file: cup95GroundTruth.sli to do the Empirical Line Calibration • Playa → Illite5.spc • Opalite Zone with Alunite → Allunite5.spc • Compare to the precalculated file: cup95_el.int • Display a Color Image (40%) • R → Band 183, G → Band 193, and B → Band 207 • Link Two DisplaysOverlay the existing ROI file: cup95_av.roi (20%) • Compare and collect the spectra at each ROI (30%) • Match the color of each ROI and spectrum pair (10%)

  15. Rigorous atmospheric correction

  16. Rigorous atmospheric correction (cont.) • Compare Radiance and ATREM • Background: ATREM Calibration • ATmospheric REMoval Program (ATREM) • A radiative transfer model-based technique for deriving scaled surface reflectance from AVIRIS data without a priori knowledge of surface characteristics • 0.94 and 1.1 micrometer water vapor bands → calculate water vapor on a pixel by pixel basis • From • the AVIRIS data • the solar irradiance curve above the atmosphere • transmittance spectra for each of the atmospheric gases CO2, O3, N2O, CO, CH4, and O2 • Unavailable for distribution • ATREM is not included as part of ENVI • The other calibration methods

  17. Exercise 5 • Compare Radiance and ATREM (cont.) • Load ATREM Apparent Reflectance Data and Start the Z Profiler • File → Open Image File • File: cup95_at.int • Gray Scale • Band 193 • Tools → Profiles → Z Profile (Spectrum) • Link Images and Compare Spectra • Tools → Link → Link Displays • Tools → Link → Dynamic Overlay Off • Tools → Profiles → Additional Z Profile • Options → Stack Plots • Stonewall Playa 590 570 • Varnished Tuff 435 555 • Silica Cap 494 514 • Opalite Zone with Alunite 531 541 • Strongly Argillized Zone with Kaolinite 502 589 • Buddingtonite Zone 448 505 • Calcite 260 613 • Close the Windows

  18. Rigorous atmospheric correction (cont.) • FLAASH

  19. Noise removal • Noise • Definition • Sources • Periodic drift, malfunction of a detector, electronic interference, intermittent hiccups in the data transmission and recording sequence • Influence • Degrade or mask the information content

  20. Noise removal (cont.) • Systematic noise • Striping or banding • e.g. Landsat MSS six detectors drift • Destriping (Fig 7.6) • Compile a set of histograms • Compare their mean and median values  identify the problematic detectors • Gray-scale adjustment factors • Line drop • Line drop correction (Fig 7.7) • Replace with values averaged from the above and below

  21. Noise removal (cont.) • Random noise • Bit error  spikey  salt and pepper or snowy appearance • Moving windows • Fig 7.8: moving window • Fig 7.9: an example of noise suppression algorithm • Fig 7.10: application to a real imagey

  22. Noise removal using IDL • Removing salt-and-pepper noise • Salt and pepper noise • Generate salt-and-pepper noise using randomu • Removing salt-and-pepper noise using median • Write a median kernel by yourself

  23. Exercise 6 • Removing salt-and-pepper noise • Subset a 300 x 300 image • Read in the subset image • ImageSize = [300, 300] • file = FILEPATH('Nagasaki', ROOT_DIR='D:',SUBDIRECTORY = ['examples','data']) • Image = READ_BINARY(file,DATA_DIMS=ImageSize) • Setup the display • DEVICE, DECOMPOSED=0 • LOADCT, 0 • WINDOW, 0, XSIZE = 300*3, YSIZE = 300*1, TITLE = 'Chapter 2: Exercise 6' • Generate salt-and-pepper noise using randomu • NoisyImage = Image • points = randomu(seed, 1800) * 300 * 300 • NoisyImage[points] = 255 • points = randomu(seed, 1800) * 300 * 300 • NoisyImage[points] = 0

  24. Exercise 6 (cont.) • Removing salt-and-pepper noise • Removing salt-and-pepper noise using median • MedianImage = MEDIAN(NoisyImage, 3) • Display the images for comparison • TVSCL, Image, /ORDER, 0 • TVSCL, NoisyImage, /ORDER, 1 • TVSCL, MedianImage, /ORDER, 2 • Save all images for further comparison • file = FILEPATH('OriginalImage', ROOT_DIR='D:',SUBDIRECTORY = ['Temp']) • WRITE_IMAGE,file,'BMP',Image • file = FILEPATH('NoisyNagasaki', ROOT_DIR='D:',SUBDIRECTORY = ['Temp']) • WRITE_IMAGE,file,'BMP',NoisyImage • file = FILEPATH('MedianImage', ROOT_DIR='D:',SUBDIRECTORY =['Temp']) • WRITE_IMAGE,file,'BMP', MedianImage

  25. Exercise 7 • Removing salt-and-pepper noise with your own kernel • Create a median kernel to convolve the NoisyImage • ProcessedImage = NoisyImage • KernelSize = 3 • KernelHalfSize = (KernelSize - 1) / 2 • Kernel = BYTARR(KernelSize, KernelSize) • SampleStart = 0 + KernelHalfSize • SampleEnd = ImageSample - KernelHalfSize - 1 • LineStart = 0 + KernelHalfSize • LineEnd = ImageLine - KernelHalfSize - 1 • FOR i = SampleStart, SampleEnd, 1 do begin • FOR j = LineStart, LineEnd, 1 do begin • KernelSampleStart = i - KernelHalfSize • KernelLineStart = j - KernelHalfSize • FOR k = 0, KernelSize-1, 1 do begin • FOR l = 0, KernelSize-1, 1 do begin • Kernel[k, l] = NoisyImage[KernelSampleStart + k, KernelLineStart + l] • ENDFOR • ENDFOR • MedianValue = MEDIAN(Kernel) • ProcessedImage[i, j] = MedianValue • ENDFOR • ENDFOR

  26. Self test • Write an IDL program to remove the salt-and-pepper noise • NoisyImage: NoisyNagasaki.img • Dimension: 500 x 500 • Use the Median function to filter the salt-and pepper noise • Show the pre- and post-processing images on your screen

  27. Noise removal using IDL • Fast Fourier Transform (FFT) • Spatial domain  Frequency domain • Decompose an image into sines and cosines of varying amplitudes and phases • The values of the resulting transform represent the amplitudes of particular horizontal and vertical frequencies • Low frequencies  gradual variations in the image  overall shape or pattern • High frequencies  abrupt variations in the image  details • Rationale of noise removal using FFT • FFT of a normal image  low frequencies  brightest pixels expressed in image or a spiky peak expressed in surface • Forward FFT  identify and mask the unwanted frequencies  Inverse FFT

  28. Exercise 8 • Refer to the source code on page 362 of book “Image Processing in IDL” • Example Code: Removing Noise with the FFT • Forward FFT • Import the image from the file C:\RSI\IDL60\examples\data\abnorm.dat • ImageSize = [64, 64] • File = FILEPATH('abnorm.dat',SUBDIRECTORY = ['examples','data']) • Image = READ_BINARY(File, DATA_DIMS = ImageSize) • Initialize a display size parameter to resize the image when displaying it • DisplaySize = ImageSize * 5 • Initialize the display • DEVICE, DECOMPOSED = 0 • LOADCT, 0 • Create a window and display the image • WINDOW, 0, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'Original Image' • TVSCL, CONGRID(Image, DisplaySize[0], DisplaySize[1]) • Transform the image into the frequency domain • FFTImage = FFT(Image)

  29. Exercise 8 (cont.) • Forward FFT (cont.) • Shift the zero frequency location from (0, 0) to the center of the display • Center = ImageSize / 2 + 1 • ShiftFFTImage = SHIFT(FFTImage, Center) • Calculate the horizontal and vertical frequency values, which will be used as the values for the axes of the display (refer to the FFT definition for storage scheme) • interval = 1. • hFrequency = INDGEN(ImageSize[0]) • hFrequency[Center[0]] = Center[0] - ImageSize[0] + FINDGEN(Center[0] - 2) • hFrequency = hFrequency / (ImageSize[0] / interval) • hFrequencyShifted = SHIFT(hFrequency, -Center[0]) • vFrequency = INDGEN(ImageSize[0]) • vFrequency[Center[1]] = Center[1] - ImageSize[1] + FINDGEN(Center[1] - 2) • vFrequency = vFrequency / (ImageSize[1] / interval) • vFrequencyShifted = SHIFT(vFrequency, -Center[1])

  30. Exercise 8 (cont.) • Forward FFT (cont.) • Create another window and display the frequency transform • WINDOW, 1, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'FFT Power Spectrum' • SHADE_SURF, ShiftFFTImage, hFrequencyShifted, vFrequencyShifted, $ • /XSTYLE, /YSTYLE, /ZSTYLE, TITLE = 'Transform of Image', $ • XTITLE = 'Horizontal Frequency', YTITLE = 'Vertical Frequency', ZTITLE = 'Amplitude', CHARSIZE = 1.5 • Create another window and display the frequency transform within the data (z) range of 0 to 5 • WINDOW, 1, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'FFT Power Spectrum' • SHADE_SURF, ShiftFFTImage, hFrequencyShifted, vFrequencyShifted, $ • /XSTYLE, /YSTYLE, /ZSTYLE, TITLE = 'Transform of Image', $ • XTITLE = 'Horizontal Frequency', YTITLE = 'Vertical Frequency', ZTITLE = 'Amplitude', CHARSIZE = 1.5, ZRANGE = [0.,5.] • Compute the power spectrum of the transform • PowerSpectrum = ABS(ShiftFFTImage)^2.0

  31. Exercise 8 (cont.) • Forward FFT (cont.) • Apply a logarithmic scale to values of the power spectrum • ScaledPowerSpectrum = ALOG10(PowerSpectrum) • Create another window and display the power spectrum as a surface • WINDOW, 2, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'FFT Power Spectrum: Logarithmic Scale (surface)' • SHADE_SURF, ScaledPowerSpectrum, hFrequencyShifted, vFrequencyShifted, $ • /XSTYLE, /YSTYLE, /ZSTYLE, TITLE = 'Transform of Image', XTITLE = 'Horizontal Frequency', YTITLE = 'Vertical Frequency', ZTITLE = 'Amplitude', CHARSIZE = 1.5 • Create another window and display the log-scaled transform as an image • WINDOW, 3, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'FFT Power Spectrum: Logarithmic Scale (image)' • TVSCL, CONGRID(ScaledPowerSpectrum, DisplaySize[0], DisplaySize[1])

  32. Exercise 8 (cont.) • Thresholding and masking • Scale the power spectrum to a zero maximum • ScaledPS0 = ScaledPowerSpectrum - MAX(ScaledPowerSpectrum) • Create another window and display the scaled transform as a surface • WINDOW, 4, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'Power Spectrum Scaled to a Zero Maximum' • SHADE_SURF, ScaledPS0, hFrequencyShifted, vFrequencyShifted, $ • /XSTYLE, /YSTYLE, /ZSTYLE, TITLE = 'Zero Maximum Power Spectrum', XTITLE = 'Horizontal Frequency', YTITLE = 'Vertical Frequency', ZTITLE = 'Max-Scaled(Log(Power Spectrum))', CHARSIZE = 1.5 • Threshold the image using -5.25, which is just below the peak of the transform, to remove the noise • Mask = REAL_PART(ScaledPS0) GT -5.25 • Mask the transform to exclude the noise • MaskedTransform = ShiftFFTImage * Mask • Create another window and display the power spectrum of the masked transform. • WINDOW, 5, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'Power Spectrum of Masked Transform and Results' • TVSCL, CONGRID(ALOG10(ABS(MaskedTransform^2)), DisplaySize[0], DisplaySize[1]), 0, /NAN

  33. Exercise 8 (cont.) • Inverse FFT • Shift the masked transform to the position of the original transform • ShiftedMaskedTransform = SHIFT(MaskedTransform, -Center) • Apply the inverse transformation to masked transform • InverseTransform = REAL_PART(FFT(ShiftedMaskedTransform, /INVERSE)) • Display results of inverse transformation • WINDOW, 6, XSIZE = DisplaySize[0], YSIZE = DisplaySize[1], TITLE = 'Noise Removed Image' • TVSCL, CONGRID(InverseTransform, DisplaySize[0], DisplaySize[1])

More Related