[GIS] How to pan-sharpen Landsat images for image-classification in GRASS

grasslandsatremote sensing

I'd like to learn how to do an unsupervised classification of a Landsat scene using i.cluster > i.maxlik in GRASS using pan-sharpened 15m-resolution images (the example given in its Wiki uses the ordinary 30m-resolutions).

I tried using i.pansharpen to create the pan-sharpened images first, but i.pansharpen can only output 3 channels that can be combined using d.rgb or r.composite. As far as I know, the image-classification process requires complete multispectral bands 1 to 7. How can I produce separate bands (band 1 to band 7), pan-sharpened using its 15m-resolution band 8 image, before inputting them to the image-classification process?

I've found a paper that demonstrates how they did this; basically they used some kind of Principal Component Analysis to somehow merge the 30m-multispectral bands with the 15m-panchromatic band. The exact quote would be:

"The method first resamples the 30-meter multispectral image to 15 meters. It then computes all six principal components of the multispectral image (We deleted the thermal band due to is courser resolution). Then the histogram of the panchromatic band (15-meter resolution) is rescaled to match the histogram
of the first principal component of the 30-meter image and the first component is replaced with the rescaled panchromatic band. This is justified because the first principal component represents overall brightness in a way that is similar to the broad spectral band of the panchromatic image. After replacement, then the six components are transformed back into the original data space, rendering an improvement in spatial resolution."

First off, the paper showed no algorithm/formulas whatsoever. I have no idea how to turn the above quote into a corresponding mathematical formula. I realized that I can use i.pansharpen with the PCA algorithm instead of the usual Brovey or IHS – but still – the output will only be 3 channels of red, blue and green – which I sadly have no clue on how to use them for image-classification..

So, before I even attempt to crack my head trying to write a new PCA algorithm manually, can anybody help me to point to an easier and better way to perform image-classification on pan-sharpened Landsat images? I mean – there should be an easier way, right? I feel like I'm missing something simple.

If the only route left is to write my own script, can you guys please point me to something out there that resembles an example of what I'm trying to do?

Any help is very much appreciated!

Best Answer

i.cluster needs at least two rasters. So, three outputs from i.pansharpen will be enough.

Overview

One working approach inside GRASS-GIS version 7 to get an acceptable color-balanced composite image after Pan-sharpening is

check if input data are 8-bit ranging inside [0, 255]
if the data are inside [0, 255] proceed then to pan-sharpening (i.pansharpen)
if the data are not inside [0, 255], rescale them to this range (r.rescale)
pan-sharpen with any of the featured methods (Brovey, IHS, PCA)
color-balance automatically by using the i.landsat.rgb module or manually adjusting the color tables of the bands of interest

Details and example instructions

Pan-Sharpening / Fusion

GRASS 7 holds a dedicated pan-sharpening module, i.pansharpen which features three techniques for sharpening, namely the Brovey transformation, the classical IHS method and one that is based on PCA.

i.pansharpen works fine with 8-bit raster maps as an input. If the data to be processed are out of this range, that is out of [0, 255], they can be rescaled to fit into this range by using GRASS' r.rescale module.

Given a set of 11-bit spectral bands (for example Blue, Green, Red, NIR and Pan) ranging between [0, 2047], querying the Blue band for example would return

r.info Blue_DNs -r
min=0
max=2047

Rescaling the Blue band to range between [0, 255]

r.rescale in=Blue_DNs out=Blue_DNs_255 from=0,2047 to=0,255

The same step applies to both the rest of the multi-spectral bands and the Panchromatic band of interest.

As usual when working with GRASS, it is required to set the region of interest, i.e. g.regionrast=Blue_DNs_255 to match the extent of the band(s) or else. The resolution itself is taken care in this particular case by the module and the resulting pan-sharpened raster maps will be of the same high(er) resolution as the Panchromatic band.

An example command for an IHS-based Pan-Sharpening action might look like

i.pansharpen pan=Pan_DNs_255 ms1=Blue_DNs_255 ms2=Green_DNs_255 ms3=Red_DNs_255 output=sharptest255 sharpen=ihs

Color Balancing

After the process completion, the module outputs

...
The following pan-sharpened output maps have been generated:
sharptest255_red
sharptest255_green
sharptest255_blue

To visualize output, run: g.region -p rast=sharptest255.red
d.rgb r=sharptest255_red g=sharptest255_green b=sharptest255_blue

Normally it should be enough to re-balance the colors after the pan-sharpening by using for example the i.landsat.rgb module or manual adjustment of each of the three bands that would compose an RGB image.

Screenshots

...to be added

[GIS] How to perform unsupervised classification of LANDSAT time series using Grass

Did you know this wiki? Maybe it helps you:

http://grasswiki.osgeo.org/wiki/Image_classification

Best Answer

Related Solutions

grass-gis – How to Color-Balance Pan-Sharpened Landsat Images in GRASS GIS

Overview

Details and example instructions

[GIS] How to perform unsupervised classification of LANDSAT time series using Grass

Related Question