Chapter

1. Chapter Image Rectification and Restoration Analysis and applications of remote sensing imagery Instructor: Dr. Cheng-Chien Liu Department of Earth Sciences National Cheng Kung University Last updated: 12 November 2014

2. Introduction • Rectification 糾正  distortion 畸變 • Restoration 復原 degradation • Source • Digital image acquisition type • Platform • TFOV

3. Geometric correction • Geometric distortion 幾何畸變 • Altitude, attitude, velocity of sensor platform • Panoramic distortion, earth curvature, atmospheric refraction, relief displacement, nonlinearities in the sweep of a sensor’s IFOV

4. Geometric correction (cont.) • Two-step procedure • Systematic (predictable) • e.g. eastward rotation of the earth  skew distortion • Deskewing  offest each successive scan line slightly to the west  parallelogram image • Random (unpredictable) • e.g. random distortions and residual unknown systematic distortions • Ground control points (GCPs) • Highway intersections, distinct shoreline features,… • Two coordinate transformation equations • Distorted-image coordinate  Geometrically correct coordinate

5. Two coordinate transformation equations • Affine coordinate transform • Six factors • Transformation equation • x = a0 + a1X + a2Y • y = b0 + b1X + b2Y • (x, y): image coordinate • (X, Y): ground coordinate • Six parameters  six conditions  3 GCPs • If GCPs > 3  redundancy  LS solutions

6. Resampling • Resampling • Fig 7.1: Resampling process • Transform coordinate • Adjust DN value  perform after classification • Methods • Nearest neighbor • Bilinear interpolation • Bicubic convolution

7. Resampling (cont.) • Nearest neighbor • Fig 7.1: a  a΄ (shaded pixel) • Fig C.1: implement • Rounding the computed coordinates to the nearest whole row and column number • Advantage • Computational simplicity • Disadvantage • Disjointed appearance: feature offset spatially up to ½ pixel (Fig 7.2b)

8. Resampling (cont.) • Bilinear interpolation • Fig 7.1: a, b, b, b  a΄ (shaded pixel) • Takes a distance-weighted average of the DNs of the four nearest pixels • Fig C.2a: implement • Eq. C.2 • Eq. C.3 • Advantage • Smoother appearing (Fig 7.2c) • Disadvantage • Alter DN values • Performed after image classification procedures

9. Resampling (cont.) • Bicubic (cubic) interpolation • Fig 7.1: a, b, b, b, c, …  a΄ (shaded pixel) • Takes a distance-weighted average of the DNs of the four nearest pixels • Fig C.2b: implement • Eq. C.5 • Eq. C.6 • Eq. C.7 • Advantage (Fig 7.2d) • Smoother appearing • Provide a slightly sharper image than the bilinear interpolation image • Disadvantage • Alter DN values • Performed after image classification procedures

10. 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 characterstics

11. Radiometric correction (cont.) • 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

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

13. Radiometric correction (cont.) • Conversion of DNs to radiance values • Measure over time using different sensors • Different range of reflectance • e.g. land  water • Fig 7.4: radiometric response function • Linear • Wavelength-dependent • Characteristics are monitored using onboard calibration lamp • DN = GL + B • G: channel gain (slope) • B: channel offset (intercept) • Fig 7.5: inverse of radiometric response function • Equation • LMAX: saturated radiance • LMAX - LMIN: dynamic range for the channel

14. 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

15. 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

16. 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

17. Tutorial: image georeferencing and registration • Georeferenced Data and Image-Map • Image to Image Registration • Image to Map Registration • HSV Merge of Different Resolution Georeferenced Data Sets

18. Georeferenced Data and Image-Map • Georeferenced Data and Image-Map • Construct an image-map complete with map grids and annotation, and produce an output image • Start ENVI • Open and Display SPOT Data • bldr_reg subdirectory: bldr_sp.img • Edit Map Info in ENVI Header • Edit Map Information • The basic map information used by ENVI in georeferencing. • Click on the arrow next to the Projection/Datum field • Click on the active DMS or DDEG button • Cursor Location/Value

19. Georeferenced Data and Image-Map (cont.) • Overlay Map Grids • Overlay →Grid Lines. • File →Restore Setup • file bldr_sp.grd • Options →Edit Map Grid Attributes • Options →Edit Geographic Grid Attributes • Apply • Overlay Map Annotation • Overlay →Annotation • File →Restore Annotation • file bldr_sp.ann • Object • Output to Image or Postscript • Direct Printing

20. Image to Image Registration • Image to Image Registration • The georeferenced SPOT image will be used as the Base image, and a pixel-based Landsat TM image will be warped to match the SPOT. • Open and Display Landsat TM Image File • bldr_reg directory: file bldr_tm.img • Band 3 • Display the Cursor Location/Value • Start Image Registration and Load GCPs • Map → Registration → Select GCPs: • Base Image: Display #1 (the SPOT image) • Warp Image: Display #2 (the TM image). • SPOT image to 753, 826 • TM image to 331, 433 • Add Point • Show List • Try this for a few points to get the feel of selecting GCPs. Once you have at least 4 points, the RMS error is reported. • Options → Clear All Points to clear all of your points.

21. Image to Image Registration (cont.) • File → Restore GCPs from ASCII. • file name bldr_tm.pts • Working with GCPs • On/Off • Delete • Update • Predict • Warp Images • Options → Warp • Displayed Band. • Warp Method • RST • Resampling • Nearest Neighbor • filename bldr_tm1.wrp • repeat steps 1 and 2 still using RST warping but with both Bilinear, and Cubic Convolution resampling methods. • Output the results to bldr_tm2.wrp and bldr_tm3.wrp, respectively. • Repeat steps 1 and 2 twice more, this time performing a 1st degree Polynomial warp using Cubic Convolution resampling, and again using a Delaunay Triangulation warp with Cubic Convolution resampling. • Output the results to bldr_tm4.wrp and bldr_tm5.wrp, respectively.

22. Image to Image Registration (cont.) • Compare Warp Results • Tools → Link → Link Displays • Load bldr_tm2.wrp and bldr_tm3.wrp into new displays and use the image linking and dynamic overlays to compare the effect of the three different resampling methods: nearest neighbor, bilinear interpolation, and cubic convolution. • Note how jagged the pixels appear in the nearest neighbor resampled image. The bilinear interpolation image looks much smoother, but the cubic convolution image is the best result, smoother, but retaining fine detail. • Examine Map Coordinates • Tools → Cursor Location/Value • Close All Files

23. Image to Map Registration • Image to Map Registration • The map coordinates picked from the georeferenced SPOT image and a vector Digital Line Graph (DLG) will be used as the Base, and the pixel-based Landsat TM image will be warped to match the map data. • Open and Display Landsat TM Image File • File → Open Image File. • bldr_reg directory: file bldr_tm.img • Gray Scale • Band 3

24. Image to Map Registration (cont.) • Select Image-to-Map Registration and Restore GCPs • Map → Registration → Select GCPs: • Image to Map • UTM • enter 13 in the Zone text field. • Leave the pixel size at 30 m and click OK to start the registration. • Add Individual GCPs by moving the cursor position in the warp image to a ground location for which you know the map coordinate (either read from a map or ENVI vector file [see the next section]). • Enter the known map coordinates manually into the E (Easting) and N (Northing) text boxes and click Add Point to add the new GCP. • File → Restore GCPs from ASCII • file bldrtm_m.pts. • Show List

25. Image to Map Registration (cont.) • Select Image-to-Map Registration and Restore GCPs • Add Map GCPs Using Vector Display of DLGs • File → Open Vector File → USGS DLG. • bldr_rd.dlg • Memory • ROADS AND TRAILS: • BOULDER, CO file in the Available Vectors Layers • Load Selected • New Vector Window • Click and drag the left mouse button in the Vector Window #1 to activate a crosshair cursor. • Tools → Pixel Locator • 402, 418 • Apply. • In the Vector Window 477593.74, 4433240.0 • Select Export Map Location. The new map coordinates will appear in the Ground Control Points Selection dialog. • Add Point • observe the change in RMS error

26. Image to Map Registration (cont.) • RST and Cubic Convolution Warp • Options →Warp File • file name bldr_tm.img • select all 6 TM bands for warping. • Warp Method RST • Resampling Cubic Convolution • background value 255 • output file name bldrtm_m.img • Display Result and Evaluate • Close Selected Files

27. HSV Merge of Different Resolution Georeferenced Data Sets • HSV Merge of Different Resolution Georeferenced Data Sets • We will use the TM color-composite image registered above as the low-resolution color image and the georeferenced SPOT image as the high resolution image. The result is a color composite image with enhanced spatial resolution. • Display 30 m TM Color Composite • file bldrtm_m.img. • RGB load bands 4, 3, and 2 (R, G, and B) into a new display. • Display 10 m SPOT Data • file bldr_sp.img. • Gray Scale • New Display

28. HSV Merge of Different Resolution Georeferenced Data Sets (cont.) • Perform HSV Sharpening • Transform → Image Sharpening → HSV • Select Input Band SPOT image • HSV Sharpening Parameters dialog, enter the output file name: bldrtmsp.img • Display 10 m Color Image • Transforms → Image Sharpening → Color Normalized (Brovey), • Overlay Map Grid • Overlay → Grid Lines. • File → Restore Setup • bldrtmsp.grd • Overlay Annotation • Overlay → Annotation. • File → Restore Annotation • file bldrtmsp.ann • Output Image Map

29. Orthorectification • Orthorectification • Definition • The geometry of an image is made planimetric (map-accurate) by modeling the nature and magnitude of geometric distortions in the imagery • Steps • Build the interior orientation (aerial photograph only) • Build the exterior orientation • Orthorectify using a Digital Elevation Model (DEM)

30. Georeferencing Images Using Input Geometry • Georeferencing Images Using Input Geometry • Modern sensors → detailed acquisition (platform geometry) information → model-based geometric rectification and map registration • Users must have the IGM or GLT file as a minimum to conduct this form of geocorrection • The Input Geometry (IGM) file: the X and Y map coordinates for a specified map projection for each pixel in the uncorrected input image. • The Geometry Lookup (GLT) file: the sample and line that each pixel in the output image came from in the input image. • If the GLT value is positive, there was an exact pixel match. If the GLT value is negative, there was no exact match and the nearest neighboring pixel is used

31. Georeferencing Images Using Input Geometry (cont.) • Uncorrected HyMap Hyperspectral Data • HyMap • Aircraft-mounted commercial hyperspectral sensor • 126 spectral channels covering the 0.44 - 2.5 mm range with approximately 15nm spectral 162 resolution and 1000:1 SNR over a 512-pixel swath. Spatial resolution is 3-10 m • Gyro-stabilized platform • Open HyMap data • envidata/cup99hym directory • File: cup99hy_true.img • Examine Uncorrected Data • Cursor Location/Value • Examine IGM files • envidata/cup99hym directory • File: cup99hy_geo_igm • Available Bands List dialog • Gray Scale • IGM Input X Map • New Display • IGM Input Y Map • New Display

32. Georeferencing Images Using Input Geometry (cont.) • Uncorrected HyMap Hyperspectral Data (cont.) • Geocorrect Image Using IGM File • Map →Georeference from Input Geometry →Georeference from IGM • File: cup99hy.eff • Input Data File • File: cup99hy.eff • Spectral Subset • File Spectral Subset: band 109 • Input Data File • Input X Geometry Band: IGM Input X Map • Input Y Geometry Band: IGM Input Y Map • Geometry Projection Information • UTM, Zone 13, datum: North America 1927 • the same map projection as the input geometry. • Build Geometry Lookup File Parameters • background value of -9999, output filename • Display and Evaluate Correction Results • Available Bands List • Georef band • Cursor Location/Value • Examine GLT Files • GLT Sample Look-up • GLT Line Look-up

33. Georeferencing Images Using Input Geometry (cont.) • Geocorrect Image using GLT File • Map →Georeference from Input Geometry →Georeference from GLT • Input Geometry Lookup File: cup99hy_geo_glt • Input Data File: cup99hy.eff • Spectral Subset • File Spectral Subset: band 109 • Input Data File • Georeference from GLT Parameters -9999 • output filename • Display and Evaluate Correction Results • Available Bands List • Georef band. • Cursor Location/Value

34. Georeferencing Images Using Input Geometry (cont.) • Using Build GLT with Map Projection • File →Open Image • File: cup99hy_geo_igm • Input X Geometry Band • IGM Input X Map • Input Y Geometry Band • IGM Input Y Map • Geometry Projection Information • State Plane (NAD 27) • Set Zone • Nevada West (2703) • Build Geometry Lookup File Parameters • Overlay Map Grids

35. IKONOS and QuickBird Orthorectification • Orthorectification • Use the Rational Polynomial Coefficients (RPCs) provided by the data vendors with some products • Orthorectification 正射糾正 • Open files • File → Open Image File • ortho subdirectory • File: po_101515_pan_0000000.tif • File → Open External File → Digital Elevation → USGS DEM • File: CONUS_USGS.dem • USGS DEM Input Parameters dialog • output filename: ortho_dem.dat • New Display • Load Band

36. IKONOS and QuickBird Orthorectification (cont.) • Run the Orthorectification • Map → Orthorectification → Orthorectify IKONOS. • File: po_101515_pan_0000000.tif • Enter Orthorectification Parameters dialog • Image Resampling: Bilinear • Background 0 • Input Height • specifies whether a fixed elevation or a DEM (more accurate) value will be used for the entire image • ortho_dem.dat • DEM Resampling • Bilinear • Geoid Offset • The height of the geoid above mean sea level in the location of the image. • -35: means that the ellipsoid is about 35 meters above mean sea level in this area • Many institutions doing photogrammetric processing have their own software for geoid height determination, or you can obtain software from NOAA, NIMA, USGS, or other sources. A geoid height calculation can currently be found at the following URL: http://www.ngs.noaa.gov/cgi-bin/GEOID_STUFF/geoid99_prompt1.prl • Save Computed DEM • Orthorectified Image • File: ikonos_ortho.dat

37. IKONOS and QuickBird Orthorectification (cont.) • Examine the Orthorectification Results • Tools → Link Displays → Link • Notice the difference in geometry, especially in the upper right corner of the two images. That is the result of the orthorectification process