Lee filter using MATLAB

Lee filter for speckle reduction has been realized in MATLAB. For explanation check:
https://www.imageeprocessing.com/2014/08/lee-filter.html

MATLAB code:

%Read an Image

I = imread('coins.png');

%Add multiplicative noise to the image

J = imnoise(I,'speckle',0.01);

%Apply Lee filter

K = Lee_filter(I,J,[5 5]);

figure,subplot(121),imshow(J);title('Noisy Image');
subplot(122),imshow(uint8(K));title('Despeckled Image');

MATLAB Function:

function Y = Lee_filter(R,E,sz)

%R is the Reference Image

%E is the Error or Noisy Image

%K is the Kernel or Window

%Y is the Output Image

% Y = mean(K)+W*(C-mean(K);

% W = variance(K)/(variance(K)+variance(R))

%Define the type

R = double(R);

E = double(E);

%Preallocate the Output Matrix

Y = zeros(size(R));

mn = round((sz-1)/2);

Tot = sz(1,1)*sz(1,2);

EImg = padarray(E,mn);

%Variance of the reference Image

Rvar = var(R(:));

Indx = floor(median(1:Tot));

for i = 1:size(R,1)

    for j = 1:size(R,2)

        K = EImg(i:i+sz(1,1)-1,j:j+sz(1,2)-1);

        varK = var(K(:));

        meanK = mean(K(:));

        W = varK./(varK+Rvar);

       

        Y(i,j) = meanK + W * (K(Indx)-meanK);

      

    end

end

end

 EXPLANATION:

Save the function Lee-filter.m in separate file

Texture Measures from GLCM – MATLAB CODE

We have seen how to generate the GLCM matrix and as a continuation, we will focus on the different measures that can be obtained from the GLCM matrix. These measures are broadly classified into three different groups.

Contrast Group:

• Contrast (Con)
• Dissimilarity (Diss)
• Homogeneity (Hom)

Orderliness Group:

• Angular Second Moment (Asm)
• Maximum Probability (Max)
• Entropy (Entropy)

Statistics Group:

• GLCM Mean
• GLCM Variance
• GLCM Correlation

GLCM matrix contains the number of occurrences of two grey levels and let’s perform normalization to measure the probability of occurrence of these two pixels.

NORMALIZATION:

Divide the GLCM matrix by total counts of the pixel occurrences. Sum of the values in GLCM matrix = 20

MATLAB CODE:

GProb = glcm./sum(glcm(:));

MATLAB CODE:

%INITIALIZE THE VARIABLES

Con=0;

Diss=0;

Hom=0;

Asm=0;

Meanx=0;

Meany=0;

Varx=0;

Vary=0;

GLCMCorrelation=0;

%NORMALIZATION

GProb = glcm./sum(glcm(:));

for i = 1:size(glcm,1)

    for j = 1:size(glcm,2)

        ij = (i-j);

        Con = Con +(GProb(i,j).*(ij).^2);

        Diss = Diss +(GProb(i,j).*abs(ij));

        Hom = Hom +(GProb(i,j)/(1 + ij * ij));

        Asm = Asm+(GProb(i,j).^2);

        Meanx = Meanx + (GProb(i,j)*i);

        Meany = Meany + (GProb(i,j)*j);

       

    end

end

Energy = sqrt(Asm);

MaxGLCM = max(GProb(:));

Entropy = sum(-GProb .* log10(GProb+0.00001));

CG=sprintf("Contrast Group:\nContrast:%f\nDissimilarity:%f\nHomogeneity:%f\n",Con,Diss,Hom);

OG=sprintf("Orderliness Group:\nAngular Second Moment :%f\nMax Probability:%f\nEntropy:%f\nEnergy:%f\n",Asm,MaxGLCM,Entropy,Energy);

display(CG);

display(OG);

MATLAB CODE:

GLCMMean = (Meanx+Meany)/2;

%COMPUTE VARIANCE

for i = 1:size(glcm,1)

    for j = 1:size(glcm,2)

        Varx = Varx+(GProb(i,j)*(i-Meanx).^2);

        Vary = Vary+(GProb(i,j)*(j-Meany).^2);

       

    end

end

GLCMVar = (Varx+Vary)/2;

%COMPUTE CORRELATION

for i = 1:size(glcm,1)

    for j = 1:size(glcm,2)

        GLCMCorrelation = GLCMCorrelation +(GProb(i,j)*(i-Meanx)*(j-Meany)/sqrt(Varx*Vary));

    end

end

SG=sprintf("Statistics Group:\nGLCM Mean:%f\nGLCM Variance:%f\nGLCM Correlation:%f\n",GLCMMean,GLCMVar,GLCMCorrelation);

display(SG); 

EXPLANATION:

The above listed measures are computed on the normalized GLCM matrix.  

MATLAB CODE:

A = imread('cameraman.tif');

%OFFSET AND DIRECTION

Displacementx = -1;

Displacementy = 1;

  

As = (A-min(A(:)))+1;

NumQuantLevels = max(As(:));

glcm = zeros([NumQuantLevels,NumQuantLevels]);

if ( Displacementx < 0 ) 

    sx=abs(Displacementx)+1;

    ex=size(As,1);

else

    sx=1;

    ex=size(As,1)-Displacementx;

end

if ( Displacementy < 0 ) 

    sy=abs(Displacementy)+1;

    ey=size(As,2);

else

    sy=1;

    ey =size(As,2)-Displacementy;

end

for i=sx:ex

    for j=sy:ey

        glcm(As(i,j),As(i+(1*Displacementx),j+(1*Displacementy)))=glcm(As(i,j),As(i+(1*Displacementx),j+(1*Displacementy)))+1;

    end

end

GLCMProb = glcm./sum(glcm(:));

[jj,ii]=meshgrid(1:size(glcm,1),1:size(glcm,2));

ij=ii-jj;

Con = sum(sum(GLCMProb.*(ij).^2));

Diss = sum(sum(GLCMProb.*abs(ij)));

Hom = sum(sum(GLCMProb./(1+(ij).^2)));

Asm = sum(sum(GLCMProb.^2));     

Meanx = sum(sum(GLCMProb.*ii)); 

Meany = sum(sum(GLCMProb.*jj));        

Energy = sqrt(Asm);

MaxGLCM = max(GLCMProb(:));

Entropy = sum(sum(-GLCMProb .* log10(GLCMProb+0.00001)));

GLCMMean = (Meanx+Meany)/2;

Varx = sum(sum(GLCMProb.*(ii-Meanx).^2));

Vary = sum(sum(GLCMProb.*(jj-Meany).^2));

GLCMVar = (Varx+Vary)/2;

GLCMCorrelation = sum(sum(GLCMProb.*(ii-Meanx).*(jj-Meany)/sqrt(Varx*Vary)));

CG=sprintf("Contrast Group:\nContrast:%f\nDissimilarity:%f\nHomogeneity:%f\n",Con,Diss,Hom);

OG=sprintf("Orderliness Group:\nAngular Second Moment :%f\nMax Probability:%f\nEntropy:%f\nEnergy:%f\n",Asm,MaxGLCM,Entropy,Energy);

SG=sprintf("Statistics Group:\nGLCM Mean:%f\nGLCM Variance:%f\nGLCM Correlation:%f\n",GLCMMean,GLCMVar,GLCMCorrelation);

display(CG);

display(OG);

display(SG);

EXPLANATION:

The above code works also for the diagonal directions.  In this example, the diagonal direction 45 is considered with displacement of 1 pixel. The minimum value of the chosen example image, ‘cameraman.tif’ is 7 and maximum is 253 and it is scaled between 1 and 247. The quantization level is chosen as 247 so the GLCM matrix will be of size 247x247.  

NOTE:

·       Some of the measures obtained will be different from the computing using ‘graycoprops’ in MATLAB due to some difference in the equations adopted. For more reference, cross check with the MATLAB documentation

·       The input image or matrix is scaled between 1 and the maximum pixel value of the image. And the Quantization level is set based on it.

Grey Level Co-occurrence Matrix (GLCM) in MATLAB

Let’s see in these series of posts on how to extract the texture features from Grey Level Co-occurrence Matrix (GLCM) in MATLAB.  In simple terms, GLCM gives the spatial relationship between adjacent or neighbouring pixels. And from this GLCM Matrix, we will measure some texture features.

Let’s consider a simple example and start coding the steps in terms of algorithmic as well as programming in MATLAB.

Procedure to Generate GLCM Matrix:

MATLAB CODE:

clear all

%GLCM MATRIX

A=[1 2 3 2 3; 1 2 3 4 5; 4 4 1 1 4; 6 7 6 7 6; 1 8 1 8 8];

Displacement = 1; %NUMBER OF STEPS BETWEEN TWO PIXELS

NumQuantLevels = 8; %QUANTIZATION LEVEL

glcm = zeros([NumQuantLevels,NumQuantLevels]); %PREALLOCATE GLCM MATRIX

for i = 1:size(A,1)

    for j = 1 :size(A,2)-1

        glcm(A(i,j),A(i,j+1))=glcm(A(i,j),A(i,j+1))+1; %INCREMENT BY 1

    end

end

EXPLANATION:

The number of occurrences of a pixel and its neighboring pixel is calculated. The offset is 1 and it can be incremented based on the requirement. Here the direction of the pixel movement is towards the right. The maximum value of the pixel in the matrix A is 8 so the quantization level is fixed to 8. For instance, in the example given above, the pixel 2 and its neighbor 3 occurs 3 times in matrix A and it is updated in the GLCM matrix at the position (2,3) as 3.


Continue Reading: Texture Measures from GLCM - MATLAB CODE

RGB Image clustering using Mahalanobis distance

The Mahalanobis distance classification is widely used in clustering. The equation has a covariance matrix that works on the variation of the classes to create similarity.

In Matlab, we have the function ‘mahal’ that can calculate the distance between a point and a sample subset.

Let’s use the Mahal() function to cluster a RGB image,

Let’s make four clusters, for the image ‘flower8.jpg’ based on the colors(ie. RGB components).

Take a small subset of the yellow color flower as shown in the picture.

We will use the Red, Green and Blue components of this subset as the sample. The mean and the covariance of this subset is calculated.

Now the distance of each pixels in the image is computed from the sample subset.

The figure below shows the distance of each pixels in the image computed from the subset(yellow color).

Note that the distance is minimal for the yellow color and it is high for other colors.

Likewise, the subset of the blue and pink color is also considered as samples and the computation is done accordingly.

The minimum distance with respect to the subsets is generated and the corresponding pixel position is updated with a unique cluster number.

MATLAB CODE:

clear all

clc

%READ AN RGB IMAGE

A = imread('flower8.jpg');

R = A(:,:,1); %RED COMPONENT

G = A(:,:,2); %GREEN COMPONENT

B = A(:,:,3); %BLUE COMPONENT

%PREALLOCATE THE VECTOR

mahal_dist  = zeros([numel(B),4]);

%CROP THE SMALL SUBSET OF THE YELLOW FLOWER

A1 = imcrop(A,[824.5 203.5 33 22]);

subset_R = A1(:,:,1);

subset_G = A1(:,:,2);

subset_B = A1(:,:,3);

%CALCULATE MAHALANOBIS DISTANCE FOR THE RED, BLUE AND GREEN COMPONENTS OF THE WHOLE IMAGE AND ALSO THE

%RED,GREEN AND BLUE COMPONENTS OF THE SMALL SUBSET

mahal_dist(:,1) = mahal(double([R(:),G(:),B(:)]),double([subset_R(:),subset_G(:),subset_B(:)])); %Yellow

%CROP THE SMALL SUBSET OF THE BLUE FLOWER

A1 = imcrop(A,[449.5 44.5 49 33]);

subset_R = A1(:,:,1);

subset_G = A1(:,:,2);

subset_B = A1(:,:,3);

mahal_dist(:,2) = mahal(double([R(:),G(:),B(:)]),double([subset_R(:),subset_G(:),subset_B(:)])); %Blue

%CROP THE SMALL PORTION OF THE PINK FLOWER

A1 = imcrop(A,[765.5 243.5 60 50]);

subset_R = A1(:,:,1);

subset_G = A1(:,:,2);

subset_B = A1(:,:,3);

%CALCULATE MAHALANOBIS DISTANCE FOR THE RED, BLUE AND GREEN COMPONENTS OF THE WHOLE IMAGE AND ALSO THE

%RED,GREEN AND BLUE COMPONENTS OF THE SMALL SUBSET

mahal_dist(:,3) = mahal(double([R(:),G(:),B(:)]),double([subset_R(:),subset_G(:),subset_B(:)])); %Pink

%CROP THE SMALL SUBSET OF THE OTHER PARTS

A1 = imcrop(A,[505.5 389.5 49 33]);

subset_R = A1(:,:,1);

subset_G = A1(:,:,2);

subset_B = A1(:,:,3);

%CALCULATE MAHALANOBIS DISTANCE FOR THE RED, BLUE AND GREEN COMPONENTS OF THE WHOLE IMAGE AND ALSO THE

%RED,GREEN AND BLUE COMPONENTS OF THE SMALL SUBSET

mahal_dist(:,4) = mahal(double([R(:),G(:),B(:)]),double([subset_R(:),subset_G(:),subset_B(:)]));

%FIND THE MINIMUM VALUE AND GIVE A CLUSTER NUMBER

[val,ind]=min(mahal_dist,[],2);

mask = zeros(size(B));

for k=1:size(mahal_dist,2)

    mask(ind(:)==k)=k;

end

figure, subplot(121),imagesc(A); title('RGB IMAGE'); subplot(122),imagesc(mask);colormap(jet);colorbar;title('CLUSTERS');

figure,

for k=1:size(mahal_dist,2)

    subplot(2,2,k);

    imshow(uint8(double(A).*logical(mask==k)));

   

end

 You can also revisit ‘K_means Clustering for RGB image’ and compare both the methods.