What does this code do?

Can anybody tell me what every line of this code do? I can't figure out. You can show me a numeric example? Thank you.

mean = conv2(in,ones(N)/N^2,'same');
variance = conv2(in.^2,ones(N)/N^2,'same')-mean.^2;
out = mean+max(variance-sigma_n^2,0)./variance.*(in-mean);

Could I write it in other simple mode?

Answer 1

[Image Analyst]

The first line has a very very badly chosen name for the output variable. mean is a built-in function that you Definitely do not want to destroy by assigning the output of conv2() to. Call it localMeans or blurredImage instead.

Then conv2 takes a kernel with is an N by N window where each value is 1/N^2. So that slides along the matrix one pixel at a time and multiplies the kernel by the matrix underneath it and adds up all the values. It has the effect of taking the mean value of the matrix at that location in an N-by-N neighborhood around it. For example if N = 3, and you're at element (12,34) in the output image, that value would be the mean of the input image in the square window from rows 11-13, and columns 33-35.

The variance is the same concept except since the matrix is squared, it's the variance. It's using the typical variance formula you see in your first statistics class. You could also do it with stdfilt() which computes the standard deviation image and then square that.

Not sure what the third line is doing - it's somewhat complicated.

max(variance-sigma_n^2,0)./variance

is basically making the variance in the range 0-1 - some kind of normalization, but it clips at 0 so no negatives are allowed. Then it multiplies it by the input image minus the mean. So if the input image mean was 80 and the values ranged from 60 - 140, then now out will be in the range -20 to +60. Finally it adds back in the mean (80) so now you have values of out between 60 to 140 again. So basically it looks like it non-linearly scaled the variance of the image to the same gray level range as the input image. I'm not sure why. It's non-linear because of the clipping to 0. This clipping will occur in areas of low variance - smooth regions of your image. Those smooth areas will all have the minimum gray level of the input image so you'll have large uniform regions and then regions where the values are the scaled variance. Again, not sure what the intent of all this is.

Answer 2

[Walter Roberson]

The mean = line is a way of implementing a sliding mean in an N x N window around each point.

The formula for mean is the sum of the values divided by the number of values, and sum of the values is the same as (sum of 1 times each value to be selected) divided by the number of values. Which is the same as sum of (1/number_of_values times each value to be selected). Selecting an N x N block is ones(N), number of such values is N^2, so the conv2() is taking the mean around each point.

Likewise from the definition of variance, you can see how the second line implements sliding variance.

However, the formula used has a limitation: more typically, you divide by (number of points minus 1) when computing the variance, reflecting a loss of 1 degree of freedom when you subtracted out the mean. The variance it is calculating is what would be calculated with var() with the w flag set to 1, not the more common w flag set to 0.

I do not know what the last line is for; you did not happen to show the definition of sigma_n