Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

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!

share|improve this question
First, it's not clear to me after reading the 1st paragraph in page 5, if the 2nd step of the presented method (PCA on all bands but the thermal one), uses the resampled to 15m bands or the original 30m bands. From the 3rd step (histogram matching of the Pan-band using as a reference the 1st PC wihch is of spatial resolution... ?), I guess the 2nd step (PCA) was applied on the original bands (30m). In the 4th step, the enhanced Pan-band replaces the 1st PC -- so, probably the 2nd step is applied on 15m-bands!-- and, finally, a backwards PCA derives an enhanced set of images. Is it so? – Nikos Alexandris Jan 15 '13 at 14:30
Do they apply twice a PCA, once for the original set of six bands (30m) and once for the resampled-to-15m set of bands? – Nikos Alexandris Jan 15 '13 at 14:32
up vote 3 down vote accepted

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

share|improve this answer
You mean I have to put all 3 outputs from i.pansharpen into Hmm, never thought of that (I thought requires rasters in separate bands).. The quality of the classification using only 3 rasters should be different from using the whole range though, right? Never mind, I'll try it out first. Thanks @Vladimir Naumov! – user8723 Oct 31 '12 at 14:11
I think that should do the trick for now! I knew I missed something simple :) – user8723 Oct 31 '12 at 16:16

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.