Re: STANDARD DEVIATON

[Alex Schuster <alex(at)pet.mpin-koeln.mpg.de>]

Ben Marriage wrote:

> I want to calculate the standard deviation of a 3x3 pixel area for each
> element of a satellite image (5000x2000 pixels). At the moment I use for
> loops to scan through the image and then use the built-in IDL standard
> deviation routines. I've tried to speed this up with the use of
> convolution and whatnot but never seemed to solve the problem. Does
> anyone know of a simple and elegant (read "quick") method of doing this
> without using loops?

Without loops yes, elegant, well yes of course, but quick I don't know?
Our fastest SUN takes 8.5 seconds for a dist( 5000, 2000 ). How fast was
your routine using convolution, and how fast should it be?

For some 3x3 pixels x_i (and a mean of x_mean), the standard deviation
is 
s = 1/9 * sqrt( sum( (x_i-x_mean)^2 ) ). img=img-smooth(img,3) does the
inner subtraction of the x_means, img=img^2 squares each element,
convol( img, k ) sums them all up Then qrt(img)/9 gives the result.

function foo, img, width
return, sqrt( convol( ( img - smooth( img, width ) )^2, $
                      fltarr( width, width ) + 1.0 ) ) / width^2
end

        Alex
-- 
  Alex Schuster     Wonko@weird.cologne.de          PGP Key available
                    alex@pet.mpin-koeln.mpg.de