# Binary processing

## Moments of Regions

### Moments

The moment of a region is computed thanks to next equation:
${m}_{pq}=\sum _{x=0}^{M}\sum _{y=0}^{N}I\left(x,y\right)\cdot {x}^{p}{y}^{q}$

### Aera of regions

The aera of a regions can be compute with the moment of order p=0 et q=0:
$A\left(R\right)=\sum _{x=0}^{M}\sum _{y+0}^{N}I\left(x,y\right)\cdot {x}^{0}{y}^{0}={m}_{00}\left(R\right)$
where A(R) is the aera of the region.
Here is the program:

```clear all; close all; %x = imread('objetinv.PNG'); x = imread('objetinv2.PNG'); %x = imread('forme.PNG'); x=x(:,:,1); [M N] = size(x); %binarisation level = graythresh(x); img = im2bw(x,0.1); %labellisation CC = bwconncomp(img); L = labelmatrix(CC); %affihcage figure subplot(221) imshow(img) subplot(222) imshow(label2rgb(L)); %nombre d'objet, number of objects nbobj = max(max(L)); %parametres des objets, object's parameters param=zeros(nbobj,6); for i = 1 : nbobj    param(i,1)=i; end %param = [ numero , nb pixel , boite_gauche , boite_droite , boite_haut , %boite_bas] for i = 1 : M     for j = 1 : N         if L(i,j)~=0             numero = L(i,j);             param(numero,2)= param(numero,2) + 1 ;             if param(numero,3)==0 || param(numero,3)>j                 param(numero,3)=j;             end             if param(numero,4)==0 || param(numero,4)<j                 param(numero,4)=j;             end             if param(numero,5)==0 || param(numero,5)>i                 param(numero,5)=i;             end             if param(numero,6)==0 || param(numero,6)<i                 param(numero,6)=i;             end         end     end end % supression petit objet var =1; for i = 1 : nbobj     if param(i,2)>10         parametre(var,:)=param(i,:);         var=var+1;     end end nbobj = size(parametre,1);   %########### %########### %%%%%%moment moment_0_0=zeros(nbobj,1); p=0; q=0; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_0_0(i,1)=moment_0_0(i,1)+ (m^p*n^q);             end         end                 end end % affichage centre de gravité %subplot(223) figure('color',[0.3137 0.3137 0.5098]) imshow(img) hold on for i = 1 : nbobj     subplot(3,4,i)     imshow(img(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4)))     title(strcat('A(R)=',num2str(moment_0_0(i)))) end```

An example

### Gravity center

The gravity center is computed with moments of order p=1 and q=0 for the coordinate on x and p=0 and q=1 for y. The moments found are then normalized by the moment of order 0 of the region (A(R)). we have:
${x}_{centre}=\genfrac{}{}{0.1ex}{}{1}{A\left(R\right)}\cdot \sum _{x=0}^{M}\sum _{y+0}^{N}I\left(x,y\right)\cdot {x}^{1}{y}^{0}=\genfrac{}{}{0.1ex}{}{{m}_{10}\left(R\right)}{{m}_{00}\left(R\right)}$
${y}_{centre}=\genfrac{}{}{0.1ex}{}{1}{A\left(R\right)}\cdot \sum _{x=0}^{M}\sum _{y+0}^{N}I\left(x,y\right)\cdot {x}^{0}{y}^{1}=\genfrac{}{}{0.1ex}{}{{m}_{01}\left(R\right)}{{m}_{00}\left(R\right)}$
Here is the program:

```%centroide moment_1_0=zeros(nbobj,1); p=1; q=0; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_1_0(i,1)=moment_1_0(i,1)+ (m^p*n^q);             end         end                 end end moment_1_0=round(moment_1_0./moment_0_0); moment_1_0= moment_1_0+parametre(:,5); % moment_0_1=zeros(nbobj,1); p=0; q=1; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_0_1(i,1)=moment_0_1(i,1)+ (m^p*n^q);             end         end                 end   end moment_0_1=round(moment_0_1./moment_0_0); moment_0_1= moment_0_1+parametre(:,3); % affichage centre de gravité %subplot(223) figure('color',[0.3137 0.3137 0.5098]) imshow(img) hold on for i = 1 : nbobj     plot(moment_0_1(i),moment_1_0(i),'r+','LineWidth',2) end title('center of gravity')```

An example

### Central Moments

They are computed for coordinates taken with respect to the coordinates of the center of gravity:
${\mu }_{pq}=\sum _{x=0}^{M}\sum _{y+0}^{N}I\left(x,y\right)\cdot {\left(x-{x}_{centre}\right)}^{p}{\left(y-{y}_{centre}\right)}^{q}$

### Normalized Central Moments

To normalize a moment, we multiply the inverse of the moment of order 0 (the aera) raised to the power (p+q+2)/2. We find:
${\stackrel{‾}{\mu }}_{pq}\left(R\right)={\mu }_{pq}\cdot {\left(\genfrac{}{}{0.1ex}{}{1}{{\mu }_{00}\left(R\right)}\right)}^{\left(\genfrac{}{}{0.1ex}{}{p+q+2}{2}\right)}$

### Orientation of regions

The orientation of a region shows the axe passing trough its centroïde and by the major part of the region. The angle $\theta$ represent the angle of rotation of the axe with respect to the horizontal. We have:
$\theta \left(R\right)=\genfrac{}{}{0.1ex}{}{1}{2}\cdot ta{n}^{-1}\left(\genfrac{}{}{0.1ex}{}{2\cdot {\mu }_{11}\left(R\right)}{{\mu }_{20}\left(R\right)-{\mu }_{02}\left(R\right)}\right)$
Thanks to this angle , we can plot the line segment passing through the center of gravity of the region. The factor $\lambda$ allows to defined the length of the line segment.
$\left(\genfrac{}{}{0}{}{x\text{'}}{y\text{'}}\right)=\left(\genfrac{}{}{0}{}{{x}_{c}entre}{{y}_{c}entre}\right)+\lambda \left(\genfrac{}{}{0}{}{\mathrm{cos}\left(\theta \left(R\right)\right)}{\mathrm{sin}\left(\theta \left(R\right)\right)}\right)$
Here is the program:

```% orientation theta=zeros(nbobj,1); moment_0_0=zeros(nbobj,1); moment_1_0=zeros(nbobj,1); moment_0_1=zeros(nbobj,1); moment_1_1=zeros(nbobj,1); moment_2_0=zeros(nbobj,1); moment_0_2=zeros(nbobj,1); for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6)+1,parametre(i,3):parametre(i,4)+1);     label = parametre(i,1);         for m = 1 : size(boite,1)         for n = 1 : size(boite,2)                         if boite(m,n)==label                 p=0;                 q=0;                 moment_0_0(i,1)=moment_0_0(i,1)+ (m^p)*(n^q);                 p=1;                 q=0;                 moment_1_0(i,1)=moment_1_0(i,1)+ (m^p)*(n^q);                 p=0;                 q=1;                 moment_0_1(i,1)=moment_0_1(i,1)+ (m^p)*(n^q);             end         end       end             moment_1_0(i,1)=moment_1_0(i,1)/moment_0_0(i,1);     moment_0_1(i,1)=moment_0_1(i,1)/moment_0_0(i,1);         for m = 1 : size(boite,1)         for n = 1 : size(boite,2)                         if boite(m,n)==label                 p=1;                 q=1;                 moment_1_1(i,1)=moment_1_1(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);                 p=2;                 q=0;                 moment_2_0(i,1)=moment_2_0(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);                 p=0;                 q=2;                 moment_0_2(i,1)=moment_0_2(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);             end         end       end     moment_1_1(i,1)=moment_1_1(i,1)*(1/moment_0_0(i,1))^2;         moment_2_0(i,1)=moment_2_0(i,1)*(1/moment_0_0(i,1))^2;         moment_0_2(i,1)=moment_0_2(i,1)*(1/moment_0_0(i,1))^2;             a=2*moment_1_1(i,1);     b=moment_2_0(i,1)-moment_0_2(i,1);     theta(i)= atan(a/b)/2;     if moment_2_0(i,1)<moment_0_2(i,1)          theta(i)=theta(i)+pi/2;     end     end % affichage orientation des momentsd %subplot(223) figure('color',[0.3137 0.3137 0.5098]) a=axes x=0; y=0; lambda=50; imshow(img) hold on for i = 1 : nbobj         x(1)=moment_1_0(i)+parametre(i,5);     x(2)=(x(1) + cos(theta(i))*lambda);     y(1)=moment_0_1(i)+parametre(i,3);     y(2)=(y(1) + sin(theta(i))*lambda);     x=round(x);     y=round(y);     plot(y(1),x(1),'r+','LineWidth',2)     plot([y(1) y(2)],[x(1) x(2)] ,'r') end title('moment orientation')  ```

An example

### Eccentricity of regions

$Ecc\left(R\right)=\genfrac{}{}{0.1ex}{}{{\mu }_{20}+{\mu }_{02}+\sqrt{{\left({\mu }_{20}-{\mu }_{02}\right)}^{2}+4\cdot {\mu }_{11}^{2}}}{{\mu }_{20}+{\mu }_{02}-\sqrt{{\left({\mu }_{20}-{\mu }_{02}\right)}^{2}+4\cdot {\mu }_{11}^{2}}}=\genfrac{}{}{0.1ex}{}{{a}_{1}}{{a}_{2}},$
$ra={\left(\genfrac{}{}{0.1ex}{}{2{a}_{1}}{A\left(R\right)}\right)}^{\genfrac{}{}{0.1ex}{}{1}{2}}$
$rb={\left(\genfrac{}{}{0.1ex}{}{2{a}_{2}}{A\left(R\right)}\right)}^{\genfrac{}{}{0.1ex}{}{1}{2}}$
where ra and rb are the ellipse parameters of the region. We can plot the ellipse thanks to ellipse equations:
$\left(\genfrac{}{}{0}{}{x\left(t\right)}{y\left(t\right)}\right)=\left(\genfrac{}{}{0}{}{{x}_{c}entre}{{y}_{c}entre}\right)+\left(\begin{array}{cc}\mathrm{cos}\left(\theta \right)& -\mathrm{sin}\left(\theta \right)\\ \mathrm{sin}\left(\theta \right)& \mathrm{cos}\left(\theta \right)\end{array}\right)\cdot \left(\genfrac{}{}{0}{}{ra\cdot \mathrm{cos}\left(t\right)}{rb\cdot \mathrm{sin}\left(t\right)}\right)$
where $0\le t\le 2\pi$
the program:

```% excentricité for i = 1 : nbobj     a1(i)=moment_2_0(i)+moment_0_2(i)+sqrt( (moment_2_0(i)-moment_0_2(i))^2 + 4*moment_1_1(i)^2);     a2(i)=moment_2_0(i)+moment_0_2(i)-sqrt( (moment_2_0(i)-moment_0_2(i))^2 + 4*moment_1_1(i)^2);     ra(i)=(2*a1(i)/moment_0_0(i))^0.5;     rb(i)=(2*a2(i)/moment_0_0(i))^0.5; end   figure('color',[0.3137 0.3137 0.5098]) a=axes x=0; y=0; imshow(img) hold on for i = 1 : nbobj     a(1)=moment_1_0(i)+parametre(i,5);       b(1)=moment_0_1(i)+parametre(i,3);     a=round(a);     b=round(b);     ra(i)=ra(i)*moment_0_0(i);     rb(i)=rb(i)*moment_0_0(i);     t=0:1/1000:2*pi;     for j = 1 : length(t)         [x(1:2,j)]= [a;b]+[cos(theta(i)) -sin(theta(i));sin(theta(i)) cos(theta(i))]*[ra(i)*cos(t(j));rb(i)*sin(t(j))];     end     plot(b,a,'r+','LineWidth',2)     plot(x(2,:),x(1,:) ,'r') end title('moment orientation')  ```

An example

### Complete Program

```clear all; close all; x = imread('objetinv2.PNG'); x=x(:,:,1); [M N] = size(x); %binarisation level = graythresh(x); img = im2bw(x,0.1); %labellisation CC = bwconncomp(img); L = labelmatrix(CC); %affihcage figure subplot(221) imshow(img) subplot(222) imshow(label2rgb(L)); %nombre d'objet, number of objects nbobj = max(max(L)); %parametres des objets, object's parameters param=zeros(nbobj,6); for i = 1 : nbobj    param(i,1)=i; end %param = [ numero , nb pixel , boite_gauche , boite_droite , boite_haut , %boite_bas] for i = 1 : M     for j = 1 : N         if L(i,j)~=0             numero = L(i,j);             param(numero,2)= param(numero,2) + 1 ;             if param(numero,3)==0 || param(numero,3)>j                 param(numero,3)=j;             end             if param(numero,4)==0 || param(numero,4)<j                 param(numero,4)=j;             end             if param(numero,5)==0 || param(numero,5)>i                 param(numero,5)=i;             end             if param(numero,6)==0 || param(numero,6)<i                 param(numero,6)=i;             end         end     end end % supression petit objet var =1; for i = 1 : nbobj     if param(i,2)>10         parametre(var,:)=param(i,:);         var=var+1;     end end nbobj = size(parametre,1);   %%%%%%moment moment_0_0=zeros(nbobj,1); p=0; q=0; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_0_0(i,1)=moment_0_0(i,1)+ (m^p*n^q);             end         end                 end end % affichage centre de gravité %subplot(223) figure('color',[0.3137 0.3137 0.5098]) imshow(img) hold on for i = 1 : nbobj     subplot(3,4,i)     imshow(img(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4)))     title(strcat('A(R)=',num2str(moment_0_0(i)))) end %centroide moment_1_0=zeros(nbobj,1); p=1; q=0; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_1_0(i,1)=moment_1_0(i,1)+ (m^p*n^q);             end         end                 end end moment_1_0=round(moment_1_0./moment_0_0); moment_1_0= moment_1_0+parametre(:,5); % moment_0_1=zeros(nbobj,1); p=0; q=1; for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6),parametre(i,3):parametre(i,4));     label = parametre(i,1);     for m = 1 : size(boite,1)         for n = 1 : size(boite,2)             if boite(m,n)==label                 moment_0_1(i,1)=moment_0_1(i,1)+ (m^p*n^q);             end         end                 end   end moment_0_1=round(moment_0_1./moment_0_0); moment_0_1= moment_0_1+parametre(:,3); % affichage centre de gravité %subplot(223) figure('color',[0.3137 0.3137 0.5098]) imshow(img) hold on for i = 1 : nbobj     plot(moment_0_1(i),moment_1_0(i),'r+','LineWidth',2) end title('center of gravity')   % orientation theta=zeros(nbobj,1); moment_0_0=zeros(nbobj,1); moment_1_0=zeros(nbobj,1); moment_0_1=zeros(nbobj,1); moment_1_1=zeros(nbobj,1); moment_2_0=zeros(nbobj,1); moment_0_2=zeros(nbobj,1); for i = 1 : nbobj     boite = L(parametre(i,5):parametre(i,6)+1,parametre(i,3):parametre(i,4)+1);     label = parametre(i,1);         for m = 1 : size(boite,1)         for n = 1 : size(boite,2)                         if boite(m,n)==label                 p=0;                 q=0;                 moment_0_0(i,1)=moment_0_0(i,1)+ (m^p)*(n^q);                 p=1;                 q=0;                 moment_1_0(i,1)=moment_1_0(i,1)+ (m^p)*(n^q);                 p=0;                 q=1;                 moment_0_1(i,1)=moment_0_1(i,1)+ (m^p)*(n^q);             end         end       end             moment_1_0(i,1)=moment_1_0(i,1)/moment_0_0(i,1);     moment_0_1(i,1)=moment_0_1(i,1)/moment_0_0(i,1);         for m = 1 : size(boite,1)         for n = 1 : size(boite,2)                         if boite(m,n)==label                 p=1;                 q=1;                 moment_1_1(i,1)=moment_1_1(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);                 p=2;                 q=0;                 moment_2_0(i,1)=moment_2_0(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);                 p=0;                 q=2;                 moment_0_2(i,1)=moment_0_2(i,1)+ ((m-moment_1_0(i,1))^p)*((n-moment_0_1(i,1))^q);             end         end       end     moment_1_1(i,1)=moment_1_1(i,1)*(1/moment_0_0(i,1))^2;         moment_2_0(i,1)=moment_2_0(i,1)*(1/moment_0_0(i,1))^2;         moment_0_2(i,1)=moment_0_2(i,1)*(1/moment_0_0(i,1))^2;             a=2*moment_1_1(i,1);     b=moment_2_0(i,1)-moment_0_2(i,1);     theta(i)= atan(a/b)/2;     if moment_2_0(i,1)<moment_0_2(i,1)          theta(i)=theta(i)+pi/2;     end     end % affichage orientation des momentsd %subplot(223) figure('color',[0.3137 0.3137 0.5098]) a=axes x=0; y=0; lambda=50; imshow(img) hold on for i = 1 : nbobj         x(1)=moment_1_0(i)+parametre(i,5);     x(2)=(x(1) + cos(theta(i))*lambda);     y(1)=moment_0_1(i)+parametre(i,3);     y(2)=(y(1) + sin(theta(i))*lambda);     x=round(x);     y=round(y);     plot(y(1),x(1),'r+','LineWidth',2)     plot([y(1) y(2)],[x(1) x(2)] ,'r') end title('moment orientation') % excentricité for i = 1 : nbobj     a1(i)=moment_2_0(i)+moment_0_2(i)+sqrt( (moment_2_0(i)-moment_0_2(i))^2 + 4*moment_1_1(i)^2);     a2(i)=moment_2_0(i)+moment_0_2(i)-sqrt( (moment_2_0(i)-moment_0_2(i))^2 + 4*moment_1_1(i)^2);     ra(i)=(2*a1(i)/moment_0_0(i))^0.5;     rb(i)=(2*a2(i)/moment_0_0(i))^0.5; end   figure('color',[0.3137 0.3137 0.5098]) x=0; y=0; imshow(img) hold on for i = 1 : nbobj     a(1)=moment_1_0(i)+parametre(i,5);       b(1)=moment_0_1(i)+parametre(i,3);     a=round(a);     b=round(b);     ra(i)=ra(i)*moment_0_0(i);     rb(i)=rb(i)*moment_0_0(i);     t=0:1/1000:2*pi;     for j = 1 : length(t)         [x(1:2,j)]= [a;b]+[cos(theta(i)) -sin(theta(i));sin(theta(i)) cos(theta(i))]*[ra(i)*cos(t(j));rb(i)*sin(t(j))];     end     plot(b,a,'r+','LineWidth',2)     plot(x(2,:),x(1,:) ,'r') end title('moment orientation')```