[GRASS-dev] Implementation of the High Pass Filter Additive Fusion technique (i.fusion.hpf)

Thu Nov 21 23:04:25 PST 2013

Because I am doing different pan sharpening algorithms, I think I got the equation backward previously. Now I am using

r.mapcalc "histmatched= ( pan - pan_mean ) / pan_sd * colorband_sd + colorband_mean

where pan=high res panchromatic band and colorband=lower resolution restricted frequency bands (e.g., red, green, blue)

This gives me a sharpened image. But I don't seem much difference from matching by the more standard cumulative distribution function except perhaps the mean/stdev method is a bit faster (and maybe especially so for bands with many grey levels).

Here is the new code. You can just drop this into the scripts directory for $GRASS_ADDON_BASE to test.

Here are links for comparing the results using the Landsat 7 images from the sc_08 demo data. I've zoomed in to shows the sharpened resolution and colors. The file names reference 3 different sharpening methods (brovey, IHS, and PCA), and 2 different histogram matching methods (CDF=cdf and mean/stdev=msd)

https://www.dropbox.com/s/m3me8yyqrjqihmv/brovey_msdhist.png
https://www.dropbox.com/s/3rqvrzdhn6g8w75/brovey_cdfhist.png

https://www.dropbox.com/s/2rvjd44dya8nrq8/ihs_msdhist.png
https://www.dropbox.com/s/hja31kvudb0yerz/ihs_cdfhist.png

https://www.dropbox.com/s/5yxp75or1ybj510/pca_msdhist.png
https://www.dropbox.com/s/7bfe71i7sbrze39/pca_cdfhist.png

Michael
____________________
C. Michael Barton
Director, Center for Social Dynamics & Complexity 
Professor of Anthropology, School of Human Evolution & Social Change
Arizona State University

voice: 	480-965-6262 (SHESC), 480-965-8130/727-9746 (CSDC)
fax:          480-965-7671 (SHESC),  480-727-0709 (CSDC)
www: http://www.public.asu.edu/~cmbarton, http://csdc.asu.edu

On Nov 21, 2013, at 1:20 PM, Michael Barton <michael.barton at asu.edu> wrote:

> I tried to implement this as an alternate way of histogram matching in i.pansharpen (attached as i.pansharpen3). It works and gives different colors, but it no longer sharpens. I’m attaching the code for you to look at (new method matchhist_mean_sd). I’m tied up the rest of the afternoon. So maybe someone else can puzzle this out.
> 
> Michael
> 
> 
> ____________________
> C. Michael Barton
> Director, Center for Social Dynamics & Complexity 
> Professor of Anthropology, School of Human Evolution & Social Change
> Arizona State University
> 
> voice:   480-965-6262 (SHESC), 480-965-8130/727-9746 (CSDC)
> fax:          480-965-7671 (SHESC),  480-727-0709 (CSDC)
> www: http://www.public.asu.edu/~cmbarton, http://csdc.asu.edu
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> On Nov 21, 2013, at 12:38 AM, Nikos Alexandris <nik at nikosalexandris.net> wrote:
> 
> > Nick Ves wrote:
> > 
> >>> The HPF algorithm about the histrogram matching states "Stretch the
> >>> new multispectral image to match the mean and standard deviation of
> >>> the original (input) multispectral image"
> >>> In that context why it is wrong to do:
> > 
> >>> Ouput <- output/sddev(output)*sddev(input)
> >>> Output <- Output - mean(output) + mean(input)
> > 
> > Moritz Lennert:
> > 
> >> To give the 'new' image the same mean and stddev of the 'input' image:
> >> (new - mean(new)) / stddev(new) * stddev(input) + mean(input)
> > 
> > Not entirely sure, I think it works. Taken a fragment from the publicly 
> > available QuickBird2 image over Sri Lanka (the extent shown at 
> > <http://grasswiki.osgeo.org/wiki/File:RGB_04APR05050541-M2AS-000000186011_01_P001.jpg>),
> > 
> > 
> > get mean and sd for r, g and b
> > 
> > # red
> > mean_red=321.375
> > sd_red=172.756
> > 
> > # hpf_red
> > mean_red_hpf=321.364
> > sd_red_hpf=238.391
> > 
> > # green
> > mean_green=434.033
> > sd_green=159.9
> > 
> > # hpf_green
> > mean_green_hpf=434.02
> > sd_green=258.249
> > 
> > # blue
> > mean_blue=285.168
> > sd_blue=75.0218
> > 
> > # hpf_blue
> > mean_blue_hpf=285.161
> > sd_blue_hpf=145.395
> > 
> > 
> > and then
> > 
> > r.mapcalc "hpf_red_histomatched = ( hpf_red - $mean_red_hpf ) / $sd_red_hpf * 
> > $sd_red + $mean_red" --o
> > 
> > r.mapcalc "hpf_blue_histomatched = ( hpf_blue - $mean_blue_hpf ) / 
> > $sd_blue_hpf * $sd_blue + $mean_blue" --o
> > 
> > r.mapcalc "hpf_green_histomatched = ( hpf_green - $mean_green_hpf ) / 
> > $sd_green_hpf * $sd_green + $mean_green" --o
> > 
> > 
> > check output stats
> > 
> > # r.univar hpf_red_histomatched
> > mean: 321.375
> > standard deviation: 172.756
> > 
> > # r.univar hpf_green_histomatched
> > mean: 434.033
> > standard deviation: 159.9
> > 
> > # r.univar hpf_blue_histomatched
> > mean: 285.168
> > standard deviation: 75.0217
> > 
> > # reset colors to match originals doesn't play exactly nice, so
> > r.colors hpf_red_histomatched color=grey -e
> > # repeat for green, blue
> > 
> > # draw
> > d.rgb ...
> > 
> > 
> > Looks nice! Is it that?
> > 
> > Nikos
> 
> <i.pansharpen3>

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20131122/986cef8b/attachment-0002.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: i.pansharpen3
Type: application/octet-stream
Size: 19069 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20131122/986cef8b/attachment-0001.obj>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20131122/986cef8b/attachment-0003.html>