Average of the surrounding pixels for each 3x3 matrix, excluding the "outliers" from the count.

Dear all,

I ve got several 134x134 double class temperature data matrices.

For each pixel, I need to calculate the average of the surrounding 8 pixels (excluding the central pixel from the calculation).

I do this simply by:

MyMtx % My 134x134 matrix
kernel = ones(3, 3) / 8;
kernel(2, 2) = 0;
bg= nanconv(MyMtx, kernel, 'same');

What I need to do is to compute the average excluding (set to NaN?), the values in each sub matrix falling outside the upper and lower limits, namely those grater than the mean+standard deviation of the 3x3 matrix and those smaller than the mean-standard deviation of the 3x3 matrix, respectively.

For instance, my sub 3x3 matrix is:

MySubMtx=[284.536443397838,285.257928640519,284.887393876757;...
    281.380348448367,285.774025332031,284.514456427377;...
    274.226936908716,278.793284550825,283.871352790957] % Note: the central pixel is excluded from the average (its weight is 0)

MySubMtxStd=nanstd(MySubMtx(:))     % 3.8919 - Note: the central pixel has to be excluded from the std 
MySubMtxAvg=nanmean(MySubMtx(:))    % 282.1835 - Note: the central pixel has to be excluded from the average
UpperLimit=MySubMtxAvg+MySubMtxStd  % 286.0754
LowerLimit=MySubMtxAvg-MySubMtxStd  % 278.2916

% This would leave me with the average (as resulting from nanconv) computed on this "residual" matrix
MySubMtx=[284.536443397838,285.257928640519,284.887393876757;...
    281.380348448367,NaN,284.514456427377;...
    NaN,278.793284550825,283.871352790957]

Any help is grately appreciated!

*EDIT:

I managed to write a code that does the job, but it is extremely slow. Would you be able to help me amend the following code to enhance its performance reducing its computational time:

A=magic(134); % Assuming this is my 134x134 matrix
background=NaN(size(A));
Mtx_Siz=134*134;
xmin=1;
xmax=3;
ymin=1;
ymax=3;
for i=1:Mtx_Siz
    if ymax<135;
        Mask=zeros(size(A));
        Mask(ymin:ymax,xmin:xmax)=1;
        Mask(ymax-1,xmax-1)=0;
        jmask=find(Mask==1);
        B=A(jmask);
        B_Mean=nanmean(B);
        B_std=nanstd(B);
        B_Upper=B_Mean+B_std;
        B_Lower=B_Mean-B_std;
        jout=find(B<B_Lower|B>B_Upper);
        B(jout)=NaN;
        B_Mean2=nanmean(B);
        B_std2=nanstd(B);
        background(ymax-1,xmax-1)=B_Mean2;
        xmin=xmin+1;
        xmax=xmax+1;
        if xmax==135;
            xmin=1;
            xmax=3;
            ymin=ymin+1;
            ymax=ymax+1;
        else
        end
    else
    end
    keep A background Mtx_Siz ymin ymax xmin xmax
end

Answer 1

[Antoni Garcia-Herreros]

Hello,

You can try something like this:

Your code takes arround 0.9 seconds to run on my laptop. This one takes 0.06 seconds

clear all
tic
A=magic(134);
[sx,sy]=size(A); 
matMEAN=zeros(size(A)); % Matrix containing the mean of the neighboring values
for i=1:sx
    for j=1:sy
        if i>1 && i<sx && j>1 && j<sy % If index in the central region
            subMat=A(i-1:i+1,j-1:j+1); % 3x3 neighbor matrix
            uplim=nanmean(subMat(:))+nanstd(subMat(:));
            lowlim=nanmean(subMat(:))-nanstd(subMat(:));
            B=subMat>lowlim & subMat<uplim;
            B(2,2)=0;
            matMEAN(i,j)=mean(subMat(B));
        else % If in the boundary
            if i==1 && j==1 % Top left corner
                subMat=A(i:i+1,j:j+1);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(1,1)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif i==1 && j>1 && j<sy %Top boundary
                subMat=A(i:i+1,j-1:j+1);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(1,2)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif j==1 && i>1 && i<sx    %Left boundary
                subMat=A(i-1:i+1,j:j+1);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(2,1)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif j==sy && i>1 && i<sx    %Right boundary
                subMat=A(i-1:i+1,j-1:j);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(2,2)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif i==sx && j>1 && j<sy    %Bottom boundary
                subMat=A(i-1:i,j-1:j+1);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(2,2)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif i==sx && j==1 % Bottom right corner
                subMat=A(i-1:i,j:j+1);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(2,1)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif i==sx && j==sy % Bottom left corner
                subMat=A(i-1:i,j-1:j);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(2,2)=0;
                matMEAN(i,j)=mean(subMat(B));
            elseif i==1 && j==sy % Top left corner
                subMat=A(i:i+1,j-1:j);
                uplim=nanmean(subMat(:))+nanstd(subMat(:));
                lowlim=nanmean(subMat(:))-nanstd(subMat(:));
                B=subMat>lowlim & subMat<uplim;
                B(1,2)=0;
                matMEAN(i,j)=mean(subMat(B));
            
            end

            
        end
    end
end
toc

Hope this helps

[Simone A.]

Hi @Antoni Garcia-Herreros, thanks for taking the time to help me. I ran the code you kindly wrote but I get a strange output. the first and the last columns as well as the first and the last rows have numbers, whilst all the rest of the output matrix (matMEAN) has zeros.

Any suggestions on what might cause this?

Thanks once again

matMEAN (ex. for the first 10 rows/columns):

matMEAN_First_Ten_Rows_Columns=[15815.5000000000,15861.3333333333,15930.3333333333,15999.3333333333,16068.3333333333,16137.3333333333,16206.3333333333,16275.3333333333,16344.3333333333,16413.3333333333;...
    15861,0,0,0,0,0,0,0,0,0;...
    15929,0,0,0,0,0,0,0,0,0;...
    15997,0,0,0,0,0,0,0,0,0;...
    16065,0,0,0,0,0,0,0,0,0;...
    16133,0,0,0,0,0,0,0,0,0;...
    16201,0,0,0,0,0,0,0,0,0;...
    16269,0,0,0,0,0,0,0,0,0;...
    16337,0,0,0,0,0,0,0,0,0;...
    16405,0,0,0,0,0,0,0,0,0]

[Antoni Garcia-Herreros]

Yes, you are right! There was a mistake in the first " if " statement, I just eddited the code.