笔者前言
本指南一共9 章,由 Valentina Plekhanova 博士撰写,链接地址http://osiris.sunderland.ac.uk/%7Ecs0vpl/GIP-VP%20Tutorials.htm。笔者花费了大量时间来找寻关于 Matlab 图像处理方面的资料,很遗憾的是中文资料实在太少,可以说根本没有!全国众多高等院校、科研院所,竟然如此缺乏有奉献精神的知识分子,笔者深感遗憾!
本指南描述性的句子不多,文笔流畅,所以笔者放弃了翻译的念头。算是转载一下吧,为国内关注图像处理的研发人员尽绵薄之力。指南中用到的图像可在
推荐比较好的论坛:仿真论坛 http://www.simwe.com/forum,如果你还在为 C++ 与 Matlab 混合编程以及调用发愁,请参考《深入浅出MATLAB 7.x混合编程》一书。
1-Intoduction to Matlab & Curve Drawing
I. Entering and Loading Data
Examples:
Enter a 4 x 1 vector
>> x = [1; 2; 3; 4]
Enter a 4 x 1 vector
>> x= [1, 2, 3, 4]
Enter a 4 x 4 matrix
>> A = [1,2,3,4; 11,22,33,44; 111,222,333,444; 1111,2222,3333,4444]
Note that the rows of a matrix are separated by semicolons, while the entries on a row are separated by spaces (or commas).
4. Enter a vector:
>> V=[0 0.1 0.2 0.3]
Note that a semicolon after a written statement presents the echo on the screen; and it may be useful if long vectors are entered.
The vector V can be defined as follows:
>> V= 0 : 0.1 : 0.3
Creat a 4 x1 vector of ones
>> iota=ones(4,1)
Get the diagonal elements of a matrix
>> DA= diag(A)
Display the names of all defined variable and thier types
>> whos
Deleting Rows and Columns
You can delete rows and columns from a matrix using just a pair of square brackets.
>> X=A;
>> X(:, 2)=[]
% to delete the second column of X
Arithmetic Operators
+
Addition
*
Multiplication
-
Subtraction
/
Division
^
Power
>=
Equal or greater than
./
Element-by-element division
.*
Element-by-element multiplication
Enter the following sequence of commands (press Enter after each command):
1
5*5+2*2
2
9*(1/(12-3)-(1/3^2))
3
s=2; h=3; g=s+h
4
pr=s*h
5
u=[1,2,3]
6
u=[1 2 3]
7
v=[-1, 0, -3]
8
w=u-2*v
9
range =
1:12
10
odd=1:2:12
11
down=20:-0.5:0
12
even=odd+1
13
u'
14
v'
15
w'
16
pi
17
xgrid=0:.05:1
18
x=xgrid*pi
19
y=sin(x)
20
a=2.; b=a^2
21
sqrt(9)
22
Z= zeros(2.5)
The
Colon
Operator
It occurs in different forms. Try the following exercises:
>> 1:10
(it is a row vector containing the integers from 1 to 10)
>> 100 : -7 : 50
>> 0 : pi/4 : pi
>>sum(A(1 : 4,4))
>>sum (A(:, end))
(it computes the sum of the elements in the last column of A)
>> sum(
1:16
)/4
II. Basic Plotting
Creating a Plot
The
plot
function has different forms and depends on the input arguments:
plot(y)
, where y is a vector
plot (x,y)
, where x and y are vectors.
Plotting in Polar Coordinates
polar
is the function that is used for polar coordinate plot: e.g.
polar ([0 2*pi], [0 1])
Controlling the Axis
* Matlab selects axis limits on the range of the plotted data. To specify the limits manually you need to use the axis command:
axis ([xmin, xmax, ymin, ymax])
.
* The command
axis('equal')
makes the x- and y-axes equal in length. To make the x- and y-data units equal use the command:
axis('square')
.
* To add a title to the plot use:
title ('Title of the Plot')
* To add x-, y- labels use:
xlabel('x-Axis')
and
ylabel('y-Axis')
Grid Lines
The
grid
command sets grid lines.
III. Examples
Make m-files that plot the following functions:
III.1
Explain how the following function works. Find/define a problem.
function drawline1
x = [1 2];
y = [1 4];
y= 2.*x+3.;
plot(x,y);
xlabel('x-Axis');
ylabel('y-Axis');
title('Plot of the Function','FontSize',12);
III.2
function drawline2
x =[1,2];
y= 2.*x+3;
plot(x,y, 'g');
xlabel('x-Axis');
ylabel('y-Axis');
title('Plot of the Function','FontSize', 12);
axis([-10 5 -5 20]);
III.3
function drawline3
x=[0 15];
y= x+1.;
plot(y,'r');
xlabel('x-Axis');
ylabel('y-Axis');
title('Plot of the Function','FontSize', 12);
axis([-10 5 -5 20]);
grid
III.4
function drawsin1
xlabel('x=0:2\pi')
ylabel('Sine of x')
title('Plot of the Sine Function','FontSize',12)
x = 0:pi/100:2*pi;
y=sin(x);
plot(x,y)
III.5
function drawsin2
plot(sin(0:.01:10));
III.6
function drawsin3
x=0:.01:10;
plot(sin(x));
III.7
% generate a spiral in polar coordinates
theta= 0:0.2:5*pi;
rho=theta.^2;
polar(theta, rho, 'go')
III.8
function Drawsquare1;
%To draw a square using nodes
%First of all define the nodes
Node = [-2 -2; -2 2; 2 2; 2 -2;-2 -2];
%Draw the square
plot(Node(:,1),Node(:,2))
axis([-5,5,-5,5]);
axis square;
III.9
function Drawsquare2;
%To draw a square using nodes
%First of all define the nodes
Node = [-2 -2; -2 2; 2 2; 2 -2;-2 -2];
%Draw the square
fill(Node(:,1),Node(:,2),'red')
axis([-5,5,-5,5]);
axis square;
2-Intoduction to Matlab: Drawing the Curves and 2D Objects
I-1 Specifying Line Styles and Colours. Plotting Lines and Markers.
If we specify a marker style but not a line style, Matlab draws only the marker.
To specify colour use command
plot (x,y, 'colour_style_marker')
, where colour strings are
'c'
cyan
'r'
red
'w'
white
'g'
green
'y'
yellow
'b'
blue
'k'
black
'm'
magenta
Line style
strings are
' - '
solid
' - - '
dashed
' : '
dotted
' - . '
dash-dot
For example:
plot (x, y, 'r : ')
plots a red dotted line.
Themarker typesare
' + '
+
' s '
square
' x '
x
' p '
pentagram
' 0 '
0
' v '
down triangle
' * '
*
' d '
diamond
For example:
plot (x, y, 'r :s ')
plots a red dotted line and places square markers at each data point.
I-2 Knowing where you are:
Use the
pwd
command to know where you are
Use the
dir
command to list the files in your current directory
Use the
cd
command to change directory
I-3 Reminder: Read "Introduction to Matlab" and do Computer Graphics Exercises - pages 11-13.
II-1. Examples
Make m-files that plot the following graphs:
II.1
function drawParabola;
% To draw parabola y=(x+3)^2+25 in a figure window.
% Define the start and end points on the x-axis
x= -
18:15
;
y=(x+3).^2+25;
plot(x,y, 'g');
Explain the program and results.
II.2
function drawEllipse;
% to draw the ellipse (x^2)/16 +(y^2)/4 = 1
%you need to get parametric representation of the ellipse, i.e. x=a cos(t); y=b sin(t)
t=0:0.01:2*pi;
a=16.;
b=4.;
x=a*cos(t);
y=b*sin(t);
plot(x,y, 'g'), grid
Explain the program and results. Is there any problem in the program?
II.3
function drawHyperbola;
% Draw hyperbola y=3/x in a figure window.
x=0.01:0.01:1;
y=3./x;
plot(x,y,'r');
Explain the program and results.
II-2 Examples
II.4
Make
NewObject.m
file
function NewObject;
%To draw an object using fill
%First of all define the nodes
X1 = [-2 -2 0];
Y1 = [-1 1 0];
X2 = [-1 1 0];
Y2 = [2 2 0];
X3 = [2 2 0];
Y3 = [1 -1 0];
X4 = [1 -1 0];
Y4 = [-2 -2 0];
%Draw the object
fill(X1,Y1,'green',X2,Y2,'green',X3,Y3,'green',X4,Y4,'green');
axis([-5,5,-5,5]);
axis square;
Explain the program and results.
II.5
function drawFig (N);
% Hello I'm a function drawFig
R = 10;
t = 0:2*pi/N:2*pi;
x = R*sin(t);
y = R*cos(t);
hold on
plot(x,y,'linewidth',2)
for
J = 1:N
for
K = 1:J
XP = [x(K) x(J)];
YP = [y(K) y(J)];
plot(XP,YP,'linewidth',2)
end
end
;
axis square
1) Save this program with
drawFig.m
2) In Matlab Command Window type
help drawFig
3) Now type, e.g.
drawFig(4)
Explain the program and results.
II.6
function CircleSimple(R);
%Draws a circle centred on the origin of radius R
K = 0;
for
t = 0:pi/36:2*pi;
K = K+1;
X(K) = R*cos(t);
Y(K) = R*sin(t);
end
fill(X,Y,'w');
axis square
Explain the program and results.
3-2D Transforms & Data Structures
I. Matlab Section
Use the Matlab
help
system to read about the following commands:
tic, toc, drawnow , clock
and
etime.
II-1. Examples
II.1 Programming the rotation
First step
: Make the following files (see programs below):
testT1.m
,
Rotatesquare.m
Second step
: Run
testT1.m
Third step
: Explain the program and results
Make
testT1.m
file
function
testT1
for
theta = 0:360
Rotatesquare(theta);
drawnow;
t = clock;
while
etime(clock,t)<0.15
end
end
Make
Rotatesquare.m
file
function
Rotatesquare(theta);
%First of all define the nodes in the homogeneous format
% see Topic 2 & "Mathematical Pages"
Node = [-2 -2 2 2 -2; -2 2 2 -2 -2;1 1 1 1 1]
%Now convert the rotation angle to radians
%Remember 180 degs. = pi radians
thetarad = pi*theta/180;
%Now let us define the rotation matrix
Ctheta = cos(thetarad);
Stheta = sin(thetarad);
Rot = [Ctheta -Stheta 0; Stheta Ctheta 0; 0 0 1]
%Now rotate the square
NewNode = Rot*Node;
%Draw the square
plot(NewNode(1,:),NewNode(2,:));
%Set the limits to the axes and ensure they are of equal dimension
axis([-5,5,-5,5]);
axis square;
Explain the program and results.
II.2
First step:
Make the following files (see programs below):
fig1.
txt
,
Tut3.m
, rotate.m, scale.m, Readobject.m
Second step
: Run
Tut3.m
Third step
: Explain the program and results
Make
fig1.txt
file
0,6
6,0
0,-6
-6,0
Make
Tut3.m
file
function
Tut3
N = readobject('
fig1.txt
');
NP = N;
hold off
patch(NP(:,1),NP(:,2),'green','edgecolor','green');
hold on
axis([-10,10,-10,10]);
axis square
tic
timeinterval = 0;
while
timeinterval<1
timeinterval = toc
end
NR = rotate(N,45);
NR = scale(NR',0.8,0.8);
NP = NR';
patch(NP(:,1),NP(:,2),'blue','edgecolor','blue');
%
NR = rotate(N,150);
NR = scale(NR',0.5,0.5);
NP = NR';
patch(NP(:,1),NP(:,2),'red','edgecolor','red');
%
hold off
Make
rotate.m
file
function
NR = rotate(N,theta)
%To rotate the node array N by the angle theta
thetarad = pi*theta/180;
Ctheta = cos(thetarad);
Stheta = sin(thetarad);
R = zeros(3,3);
R(1,1) = Ctheta;
R(2,2) = Ctheta;
R(1,2) = -Stheta;
R(2,1) = Stheta;
R(3,3) = 1;
NR = R*N';
Make
scale.m
file
function
NS = scale(N,SX,SY)
%To scale the node array N by SX and SY
S = zeros(3,3);
S(1,1) = SX;
S(2,2) = SY;
S(3,3) = 1;
NS = S*N';
Make
Readobject.m
file
function
N = Readobject(FN)
%To load data from a text file 'FN' into a node array N
temp = dlmread(FN);
[a,b] = size(temp);
N = ones(a,3);
N(:,1)=temp(:,1);
N(:,2)=temp(:,2);
Explain the program and results.
II-3
function
Spiral
%To draw a 3D spiral
%For efficiency we calculate X,Y and Z for a single turn
R = 5;
C = 2;
thetastep = 10*pi/180;
theta = [0:thetastep:2*pi];
X = R*cos(theta);
Y = R*sin(theta);
Z = C*theta;
%Now extend to the second term of the spiral
[dummy numpoints] = size(X);
X(numpoints+1:2*numpoints)= X;
Y(numpoints+1:2*numpoints)= Y;
Z(numpoints+1:2*numpoints)= Z + C*2*pi;
plot3(X,Y,Z,'linewidth',5);
Explain the program and results.
4-3D Transforms & Data Structures
I. Matlab Section
Use the Matlab
help
system to read about the following commands:
dlmread, rotate3d on, patch,
Nan
, set
II-1. Examples
II.1
First step
: Make the following file (see program below):
CubeProject.m
Second step
: Run
CubeProject.m
Third step
: Explain the program and results
Make
CubeProject.m
file
function
CubeProject
%To illustrate normal perspective projection
%Currently this progam only rotates and displays the cube
%First define the basic cube
V = [2 -2 -2 1; 2 -2 2 1; 2 2 2 1; 2 2 -2 1; -2 -2 -2 1; -2 -2 2 1; -2 2 2 1; -2 2 -2 1];
F = [1 2 3 4; 3 7 8 4; 7 6 5 8; 6 2 1 5; 2 6 7 3; 4 8 5 1];
VT = V';
%Define sutable rotation angles to orientate the cube theta
thetaY = 60*pi/180;
thetaX = 30*pi/180;
%Set up a rotation matrix for the Y axis
YROT = [cos(thetaY) 0 -sin(thetaY) 0 ; 0 1 0 0; sin(thetaY) 0 cos(thetaY) 0; 0 0 0 1];
XROT = [1 0 0 0 ; 0 cos(thetaX) -sin(thetaX) 0; 0 sin(thetaX) cos(thetaX) 0; 0 0 0 1];
%Now do the rotations
VT = XROT*YROT*VT;
%Strip off the homogeneous co_ordinate to use the patch command
VP = zeros(3,8);
VP=VT(1:3,:);
colourset = [0 1 0;1 0 0;0 0 1;1 1 0;0 1 1;1 1 1];
patch('Vertices',VP','Faces',F,'FaceVertexCData',colourset,'facecolor','flat')
axisarray = [-5 5 -5 5 -5 5 0 1];
axis(axisarray)
II.2
Step 1:
Make the following files (see programs below):
DisplayBuilding.m; Nhouse.
txt
; Phouse.
txt
Step 2
: Run
DisplayBuilding.m
Step 3
: Explain the program and results
Make
Nhouse.txt
file
-10,-20,0
-10,-20,20
0,-20,30
10,-20,20
10,-20,0
-10,20,0
-10,20,20
0,20,30
10,20,20
10,20,0
Make
Phouse.txt
file
1,2,3,4,5
5,4,9,10,
NaN
6,7,2,1,
NaN
4,3,8,9,
NaN
7,8,3,2,
NaN
10,9,8,7,6
1,5,10,6,
NaN
Make
DisplayBuilding.m
file
function
DisplayBuilding
%To display a basic Building
Node = dlmread('Nhouse.txt');
Face = dlmread('Phouse.txt');
patch('vertices',Node,'faces',Face,'facecolor','b');
%Depending on the size of the building you may need to change the axes range
axis([-100 100 -100 100 -100 100 0 1]);
grid;
axis square
set(gcf,'renderer','zbuffer')
II-3
Step 1:
Make the following file:
HouseDraw1.m;
and check that you already have files:
NHouse.txt
&
PHouse.txt
Step 2
: Run
HouseDraw1.m
Step 3
: Explain the program and results
function
HouseDraw1
%Draws a simple house
%Data on nodes is read from a text file NHouse.txt
%Data on patches is read from a text file PHouse.txt
Node = dlmread('NHouse.txt');
Face = dlmread('PHouse.txt');
%Draw the house
colourset = [0.2 0.8 0.2];
patch('vertices',Node,'faces',Face,'FaceVertexCData',colourset,'FaceColor','flat');
axis_data = [-40 40 -40 40 -40 40 0 1];
grid;
axis(axis_data);
axis square;
rotate3d on
II-4.
Step 1:
Make the following file:
HouseDraw2.m;
and check that you already have files:
NHouse.txt
&
PHouse.txt
Step 2
: Run
HouseDraw2.m
Step 3
: Explain the program and results
Step 4
: Compare results with Example II-3
Step 5:
Delete
NaN
s from
PHouse.txt
file
Step 6
: Explain and Compare results with Example II-3
Make
HouseDraw2.m
file:
function
HouseDraw2
%Draws a simple house
%Data on nodes is read from a text file NHouse.txt
%Data on patches is read from a text file PHouse.txt
Node = dlmread('NHouse.txt');
%Now read in the uncorrected patch data
UPatch = dlmread('PHouse.txt');
%First find how many patches we need
%We don't care about the other dimension X
[numpatches,X]=size(UPatch);
%Now lets sort into 4 node & 5 node faces
num4 = 0;
for
k = 1:numpatches
if
UPatch(k,5)==0
num4 = num4 + 1;
end
end
Face4 = zeros(num4,4);
Face5 = zeros(numpatches-num4,5);
j = 1;
i=1;
for
k = 1:numpatches
if
UPatch(k,5)==0
Face4(j,:)=UPatch(k,1:4);
j=j+1;
else
Face5(i,:)=UPatch(k,:);
i = i+1;
end
end
colourset = [0];
%
patch('Vertices',Node,'Faces',Face4,'FaceVertexCData',colourset,'FaceColor','r')
patch('Vertices',Node,'Faces',Face5,'FaceVertexCData',colourset,'FaceColor','flat')
axis_data = [-40 40 -40 40 -40 40 0 1];
grid;
axis(axis_data);
axis square;
rotate3d on
II-5.
First step:
Make the following files:
testCone.m
&
Cone.m
Second step
: Run
testCone.m
Third step
: Explain the program and results
Make
testCone.m
file:
function
testCone;
% How to colour a face
R = 4;
h = 2;
L = 10;
[X,Y,Z] = Cone;
surf(X,Y,Z,'facecolor',[0.2 0.8 0.2])
axis equal
axis square
rotate3D on
Make
Cone.m
file:
function
[X,Y,Z] = Cone(R,h,L)
% To create X,Y,Z data to draw a truncated cone.
% R = base radius of the cone.
% h = height to which the cone is drawn
% L = apex of the cone
% Note setting h = L will draw the full cone
if
nargin < 3, L = 1; end
if
nargin < 2, h = L; end
if
nargin == 0, R = 1; end
stepsize = h/10;
t = 0:stepsize:h;
[X,Y,Z] = cylinder((R*(1-t/L)));
5-3D Transforms; Shading,Light & Colour
I. Matlab Section
Use the Matlab
help
system to read about the following commands:
light, lighting gouraud, cylinder, sphere
II-1. Examples
II.1
First step
: Make the following file (see program below):
LightDemo.m
Second step
: Run
LightDemo.m
Third step
: Explain the program and results
Make
LightDemo.m
file
function
LightDemo
%To demonstrate some of the properties of lighting
%Create a sphere with a red colour and draw it with the surf command
[X,Y,Z] = sphere(30);
surf(X,Y,Z,'FaceColor','red','EdgeColor','none');
%Try changing the material to dull, shiny, metal and see the effect material shiny;
%material metal;
%material dull;
%Light is incident along the X axis and is located an infinite distance away
%Try changing the colour of the light (e.g. 'white', 'green', 'blue' etc.)
light('Position',[1 0 0],'Style','infinite','color','white');
%Try the effects of flat, gouraud or phong shading
lighting gouraud
axis square
rotate3D on
II.2
Step 1:
Make the following files:
goblet.
txt
&
gobletprofile.m
&
RevObject.m
Step
2: Run
gobletprofile.m
Step 3:
Run
RevObject(6)
Step 4
: Explain the program and results
Make
gobletprofile.m
file:
function
gobletprofile
%To display the outline profile of a the goblet
U = dlmread('goblet.txt');
[a,b] = size(U);
plot(U(:,1),U(:,2),'-o','linewidth',2)
axis([0 12 0 12])
axis square
Make
goblet.
txt
file (i.e. the digitised profile of the object is presented as follows):
5,10
5,6
2,4
2,2
5,0
Make
RevObject.m
file:
function
RevObject(n)
%To generate a 3D object from the digitised profile in 'goblet.txt'
%NOTE this version only works with 5 nodes in the profile
%First read in the data, set the number of nodes and faces
U = dlmread('goblet.txt');
[a,b] = size(U);
numnodes = a*n;
numpatch = (a-1)*n;
Node = zeros(numnodes,3);
PSet = zeros(numpatch,4);
%Now calculate all the node values
for
k = 0:n-1
for
L = 1:a
theta = 2*k*pi/n;
Node(L+a*k,1) = U(L,1)*cos(theta);
Node(L+a*k,2) = U(L,1)*sin(theta);
Node(L+a*k,3) = U(L,2);
end
;
end
;
%Uncomment the following line if you want to see the node list
%Node
%Now assign nodes to faces (or patches)
for
k = 1:n term = 5*(k-1);
for
L = 1:4
Pset(4*(k-1)+L,:) = [(term+L) (term+L+1) (term+L+a+1) (term+L+a)];
end
;
end
;
%Finally ensure that the last but one profile connects to the first
%Pset
for
k = numpatch-3:numpatch
Pset(k,3)=Pset(k,3)-numnodes;
Pset(k,4) = Pset(k,4) - numnodes;
end
%Uncomment the following line if you want to see the patch array
%Pset
patch('Vertices',Node,'Faces',Pset,'facecolor',[0.2 0.8 0.2],'edgecolor',[0 0 0])
%alpha(0.3)
%light('Position',[1 0 0],'Style','infinite');
axis square
rotate3d on
6-Image Analysis - Basics
I. Matlab Section
I-1. Loading & Displaying an Images
We can load an image with the ‘
imread
’ function and display an image with the ‘
imshow
’ function.
A = imread(filename,fmt)
[X,map] = imread(filename,fmt)
In the first form the image is read into the array
A
and the format of the image ‘
fmt
’ is optional.
In the second form the image is read into the array
X
and the associated indexed map into ‘
map
’ scaled 0to 1as double.
imshow(A)
imshow(X,map)
imshow
has many other formats
I-2. Image Type Conversions
MatLab allows easy conversion between image types using :
rgb2hsv
to convert to a HSV image
rgb2gray
to convert to a grey scale image (note American spelling).
rgb2ind
to convert to an indexed image
You should check the details in the MatLab help files.
With the addition of the following line (and a variable name change) the previous example can be used to display the HSV components.
B = rgb2hsv(A);
I-3. Inspecting and Recording Image Colours
improfile
: which will reveal the colour intensity profile along a line defined on the image.
Simply placing
improfile
on a line will allow you to interactively explore the colour profile alternatively using :
c = improfile(image,xi,yi,n);
will record in the array ‘
c
’ the RGB values from image
image
with line segment end points defined by
xi
&
yi
using ‘
n
’ equally spaced points
pixval
: which will interactively record the location and colour components of each pixel as you move a cursor across the image.
impixel
: returns the RGB values for a specified pixel.
I-4. Histogram Manipulation
The histogram of an image is a simple way to see how the gray level (i.e. the information) is spread over the available range of grey levels. This in turn can be used to highlight problems during the acquisition. The histogram can also be modified to make better use of the available range for further processing and display.
Matlab has several functions available to calculate and process histograms:
imhist
: Display the histogram of an image. Synthax:
hist = imhist(image,nb_boxes)
imhist
works with uint8 types images.
histeq
: histogram equalisation and specification. Synthax:
ima1 = histeq(image,hgram)
where
hgram
is a specified histogram. If
hgram
is not present, histogram equalisation is performed.
II-1. Examples
Make m-files
II.1
function
colourdisplay
%Displaying the individual colours
%To demonstrate the splitting of an image into its primary colours
% Save an image, for example
Egik.jpg
A = imread('Egik.jpg');
subplot(2,2,1);
imshow(A);
title('RGB image');
Redimage = A(:,:,1);
subplot(2,2,2);
imshow(Redimage);
title('Red image');
Greenimage = A(:,:,2);
subplot(2,2,3);
imshow(Greenimage);
title('Green image');
Blueimage = A(:,:,3);
subplot(2,2,4);
imshow(Blueimage);
title('Blue image');
Explain the program and results.
Select color which gives the best contrast.
II-2.
function
ImageTest
A = imread('Egik.jpg');
%Check that it was indeed loaded
whos
% Display the image
imshow(A)
% Convert the variable into double
A=im2double(A);
% Check that the variable indeed was converted into double
whos
% The next procedure cuts out the upper left corner
% (i.e. the leaf) of the image
% and stores the reduced image as Ared
for
i=1:145
for
j=1:180
Ared(i,j)=A(i,j);
end
end
%Check what variables you now have stored
whos
% Display the reduced image
imshow (Ared)
Explain the program and results.
II.3
function
AnalyseImage
%To demonstrate using arithmetic on image data
% Save an image
Egik.jpg
myfilename = 'Egik.jpg';
A = imread(myfilename);
R = A(:,:,1);
G = A(:,:,2);
B = A(:,:,3);
subplot(2,2,1);imshow(A);
title('Original');
subplot(2,2,2);imshow(R);
title('Red');
subplot(2,2,3);imshow(G);
Title('Green');
Rconverted = double(R) + 1;
Gconverted = double(G) + 1;
Dconverted = Rconverted - Gconverted;
D = uint8(Dconverted - 1);
subplot(2,2,4);imshow(D);
title('Difference red - green');
Explain the program and results.
II.4
function
TestTry
%First load an image and display
filename = ['Egik.jpg'];
J = imread(filename);
figure;
subplot(2,2,1),imshow(J);
Title('Original');
% In order to convert the indexed image into an intensity
% (gray scale) image use the
ind2gray
command
J=rgb2gray(J);
whos % now the size indicates that our image is a regular matrix.
subplot(2,2,2),imhist(J);
Title('imhist');
%Use only one of the following 3 lines at a time
%I = (J>35) & (J<60);
%I = (J>70) & (J<110);
%I = (J>110);
%figure, imshow(I);
lowin = 35/255;
highin = 60/255;
K = imadjust(J,[lowin highin],[ ]);
subplot(2,2,3),imshow(K);
Title('imadjust');
Explain the program and results.
II-5.
function
ThreshHold
%To demonstrate thresholding an image
J = imread('Egik.jpg');
subplot(2,1,1);imshow(J);
BW = im2bw(J,0.6);
subplot(2,1,2);imshow(BW);
II-6.
function
Equalise
%To show the effects of histogram equalisation
J = imread('Egik.jpg');
J=rgb2gray(J);
subplot(2,1,1);imshow(J);
K = histeq(J);
subplot(2,1,2);imshow(K);
Explain the program and results.
7-Image Segmentation
I. Matlab Section
I-1.
im2bwconverts image to binary image by thresholding.
I-2.
graythreshis used to determine a threshhold for converting the image to binary.
I-3.
bwlabelsearches for connected components and label them with unique numbers.bwlabeltakes a binary input image and a value specifying the connectivity of objects.
I-4.
STATS = regionprops (L,PROPERTIES)measures a set of properties for each labeled region in the label matrix L.
Positive integer elements of L correspond to different regions. For example, the set of elements of L equal to 1 corresponds to region 1; the set of elements of L equal to 2 corresponds to region 2; and so on.
STATSis a structure array of lengthmax(L(:)).
The fields of the structure array denote different properties for each region, as specified byPROPERTIES.
Use the Matlabhelpsystem to read more about this command.
I-5.
You can usefindfunction in conjunction withbwlabelto return vectors of indices for the pixels that make up a specific object.
I-6.
Use the Matlabhelpsystem to read about the following commands:max, min, find
II-1. Examples
Make m-files
II-1.functionHist1
% Basic Global Thresholding
% Select a Threshold Value from the Histogram
A=imread('Egik.jpg');
% imshow(A);
figure;
image(A);
A=rgb2gray(A);
imhist(A);
% figure;
% to define normalised gray level, seeimhiste.g.k1=200/255
k1=200/255;
BW1=im2bw(A,k1);
% imshow(BW1);
% figure;
k2=20/255;
BW2=im2bw(A,k2);
% imshow(BW2);
% figure;
k3=150/255;
BW3=im2bw(A,k3);
% imshow(BW3);
figure;
subplot(2,2,1), imhist(A);
Title('imhist of Image');
subplot(2,2,2), imshow(BW1);
Title(['BW1 Image with k1=', num2str(k1)]);
subplot(2,2,3), imshow(BW2);
Title(['BW2 Image with k2=', num2str(k2)]);
subplot(2,2,4), imshow(BW3);
Title(['BW3 Image with k3=', num2str(k3)]);
Explain the program and results.
II-2.functionGlobalThresh
% Compute global image threshold usingOtsu's method
A = imread('Egik.jpg');
A=rgb2gray(A);
level=graythresh(A);
level
BW6 = im2bw(A,level);
figure, imshow(BW6);
Title(['BW6 Image with global threshold level =', num2str(level)]);
Explain the program and results.
II-3.functionEdgeEgik
%To edge detect an Egik
A = imread('Egik.jpg');
A=rgb2gray(A);
[BW,thresh] = edge(A,'sobel');
imshow(BW);
thresh
Explain the program and results.
II-4.functionConnectedObjects
A = imread('bean.jpg');
A=rgb2gray(A);
[BW,thresh] = edge(A,'sobel');
imshow(BW);
Title('BW Image');
thresh
% usebwlabelto return innumthe number of connected objects found inBW
[L,num] = bwlabel(BW,8);
% convert a label matrix into anRGBimage for the purpose of visualising
% the labeled regions
% use the following two lines to view the regions found and
% to decide on the selection criteria
RGB=label2rgb(L);
figure, imshow(RGB);
Title('RGB Image');
num
Explain the program and results.
II-5.functionSelectRegion
% savebean.jpg
A = imread('bean.jpg');
A=rgb2gray(A);
BW = edge(A,'sobel',0.04);
imshow(BW);
Title('BW Image');
[L,num] = bwlabel(BW,8);
RGB=label2rgb(L);
figure, imshow(RGB);
Title('RGB Image');
num
% Measure properties of image regions.
% Consider the following approach to selecting the image profile
STATS = regionprops(L,'BoundingBox','MajorAxisLength');
forj=1:num
ifSTATS(j).MajorAxisLength>100
maxobject=j;
end
end
maxobject
boxsize=STATS(maxobject).BoundingBox;
% The bounding box represents the smallest rectangle that can contain a region.
%The four element vector returned by the BoundingBox field.
% Two first elements show the upper left corner of the bounding box
% and two last ones represent a width and a hight of the box.
boxsize
xb1 = round(boxsize(1));
xb2= round(boxsize(1)+boxsize(3));
yb1=round(boxsize(2));
yb2= round(boxsize(2)+boxsize(4));
BWW=BW(yb1:yb2, xb1:xb2);
figure, imshow(BWW)Explain the program and results.
8、9-Morphological Image Processing
I. Matlab Section
I-1.
Use the Matlabhelpsystem to read about the following commands:strel, imopen, imclose
I-2.
Reminder:impixelcan return values and/or coordinates (in ROW/COL)
Use the Matlabhelpsystem to read more about this command.
II-1. Examples
Make m-files
II-1.functionMorDisk
% Create morphological structuring element with a disk of the given radious
clear, close all,
A = imread('bean.jpg');
A=rgb2gray(A);
BW = edge(A,'sobel',0.04);
imshow(BW);
Title('BW Image');
[L,num] = bwlabel(BW,8); % Label components
RGB=label2rgb(L);
figure, imshow(RGB);
Title('RGB Image');
num
figure,
subplot(2,2,1), imshow(A);
Title('Original');
% Perform a morphological opening operation by calling imopen with a disk-shaped
% structuring element with radiouses of10, 25, 30, 40.
% The structuring element is created by thestrelfunction.
% The morphological opening has the effect of removing objects that cannot
% completely contain a disk of the given radious.
% tryRL=imopen(A,strel('disk', 40));
%RL=imopen(A,strel('disk', 25));
RL=imopen(A,strel('disk', 3));
subplot(2,2,2), imshow(RL)
Title('RL Image');
Explain the program and results.
II-2.functionMorSquare
% Create morphological structuring element with a sguare
clear, close all,
A = imread('bean.jpg');
A=rgb2gray(A);
BW = edge(A,'sobel',0.04);
imshow(BW);
Title('BW Image');
[L,num] = bwlabel(BW,8); % Label components
RGB=label2rgb(L);
figure, imshow(RGB);
Title('RGB Image');
num
figure,
subplot(2,2,1), imshow(A);
Title('Original');
% Perform a morphological opening operation by calling imopen with a sguare structuring element
%SE = strel('square',W) creates a square structuring element whose width isWpixels.Wmust be a nonnegative integer scalar.
% tryRL = imopen(A,strel('square',100));
%RL = imopen(A,strel('square',10));
RL = imopen(A,strel('square',40));
subplot(2,2,2), imshow(RL)
Title('RL Image');
Explain the program and results.
II-3.functionTestMor
% Step 1: Threshold the image
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
% tryBW = im2bw(A,graythresh(A));
imshow(A), title('Original')
figure, imshow(BW);
Title('Step 1: Thresholded Image')
%Step 2: Create morphological structuring element with a disk-shaped
% structuring element with a given radius
MR = strel('disk',6);
% Step 3: Close with a disk of radius6to merge
% together small features that are close together.
BW2 = imclose(BW,MR);
figure, imshow(BW2);
Title('Step 3: Closing')
% Step 4: Follow with an opening to remove the isolated white pixels.
BW3 = imopen(BW2,MR);
figure, imshow(BW3);
Title('Step 4: Opening')
Explain the program and results.
II-4.function MorNum
% Step 1: Threshold the image
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
imshow(A), title('Original')
figure, imshow(BW);
Title('Step 1: Thresholded Image')
%Step 2: Create morphological structuring element with a disk-shaped
% structuring element with a given radius
MR = strel('disk',6);
% Step 3: Close with a disk of radius 6 to merge
% together small features that are close together.
BW2 = imclose(BW,MR);
figure, imshow(BW2);
Title('Step 3: Closing')
% Step 4: Follow with an opening to remove the isolated white pixels.
BW3 = imopen(BW2,MR);
figure, imshow(BW3);
Title('Step 4: Opening')
% Determine the Number of Objects in the Image
[L,num] = bwlabel(BW2,8);
% compare the number of beans on the image with num that you have received
% after opening process
num
STATS = regionprops(L,'BoundingBox','MajorAxisLength');
forj=1:num
ifSTATS(j).MajorAxisLength>100
maxobject=j;
end
end
maxobject
boxsize=STATS(maxobject).BoundingBox;
boxsize
xb1 = round(boxsize(1));
xb2= round(boxsize(1)+boxsize(3));
yb1=round(boxsize(2));
yb2= round(boxsize(2)+boxsize(4));
BWW=BW(yb1:yb2, xb1:xb2);
figure, imshow(BWW)
Explain the program and results.
II-5.functionMorColor
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
imshow(A), title('Original')
figure, imshow(BW);
Title('Step 1: Thresholded Image')
% Step 2: Create morphological structuring element
MR = strel('disk',6);
BW2 = imclose(BW,MR);
figure, imshow(BW2);
Title('Step 3: Closing')
BW3 = imopen(BW2,MR);
figure, imshow(BW3);
Title('Step 4: Opening')
[L,num] = bwlabel(BW2,8);
num
% Now draw theoutline profile(note that an outline profile is a graphic summary of the object presented by the line by which the object is defined or bounded; contour)
imwrite(BW2,'BW.jpg'); % saves the image
D = imread('BW.jpg');
ED = edge(D,'sobel'); % creates a binary image using the Sobel approximation
imshow(ED);% displays image
Title('Outline profiles - Image');
%
RGB=label2rgb(L);
figure, imshow(RGB);
Title('RGB Image');
Explain the program and results.
II-6.Run the following program with the processes: Closing-Opening, i.e.functionMorCO
A = imread('Egik.jpg');
BW = ~im2bw(A,graythresh(A));
MR = strel('disk',6);
BW2 = imclose(BW,MR);
BW3 = imopen(BW2,MR);
figure,
subplot (2,2,1), imshow(A);
Title('Original')
subplot(2,2,2), imshow(BW);
Title('Step 1: Thresholded Image')
% Note that Step 2: Create morphological structuring element
subplot(2,2,3), imshow(BW2);
Title('Step 3: Closing')
subplot(2,2,4), imshow(BW3);
Title('Step 4: Opening')
[L,num] = bwlabel(BW2,4);
num
Explain the program and results.
II-7.Now we modify and run the program II-6with the processes: Opening-Closing, i.e.functionMorOC
A = imread('Egik.jpg');
BW = ~im2bw(A,graythresh(A));
MR = strel('disk',6);
BW2 = imopen(BW,MR);
BW3 = imclose(BW2,MR);
figure,
subplot (2,2,1), imshow(A);
Title('Original')
subplot(2,2,2),imshow(BW);
Title('Step 1: Thresholded Image')
% Note that Step 2: Create morphological structuring element
subplot(2,2,3),imshow(BW2);
Title('Step 3: Opening')
subplot(2,2,4), imshow(BW3);
Title('Step 4: Closing')
[L,num] = bwlabel(BW2,4);
num
Explain the program and results.
II-8.Run the following program with the processes: Closing-Opening, i.e.function MorCOBean
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
MR = strel('disk',6);
BW2 = imclose(BW,MR);
BW3 = imopen(BW2,MR);
[L,num] = bwlabel(BW2,4);
num
% find All Beans
beandata=regionprops(L,'basic');
allbeans=[beandata.Area];
% find max & min Beans
maxbean=max(allbeans);
maxbean
minbean=min(allbeans);
minbean
% find Big Bean
bigbean=find(allbeans==maxbean);
bigbean
% find Small Bean
smallbean=find(allbeans==minbean);
smallbean
Explain the program and results (find in Matlab Command Window).
II-9.Now we modify and run the program II-8 with the processes: Opening-Closing, i.e.functionMorOCBean
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
MR = strel('disk',6);
BW2 = imopen(BW,MR);
BW3 = imclose(BW2,MR);
[L,num] = bwlabel(BW2,4);
num
% find All Beans
beandata=regionprops(L,'basic');
allbeans=[beandata.Area];
% find max & min Beans
maxbean=max(allbeans);
maxbean
minbean=min(allbeans);
minbean
% find Big Bean
bigbean=find(allbeans==maxbean);
bigbean
% find Small Bean
smallbean=find(allbeans==minbean);
smallbean
Explain the program and results (find in Matlab Command Window).
II-10.functionCoordinateGenerator
% generate a Vector Model of the object profile
clear, close all,
A = imread('bean.jpg');
BW = ~im2bw(A,graythresh(A));
MR = strel('disk',6);
BW2 = imclose(BW,MR);
[L,num] = bwlabel(BW2,8);
num
STATS = regionprops(L,'BoundingBox','MajorAxisLength');
forj=1:num
ifSTATS(j).MajorAxisLength>100
maxobject=j;
end
end
maxobject
boxsize=STATS(maxobject).BoundingBox;
boxsize
xb1 = round(boxsize(1));
xb2 = round(boxsize(1)+boxsize(3));
yb1 = round(boxsize(2));
yb2= round(boxsize(2)+boxsize(4));
BWW=BW2(yb1:yb2, xb1:xb2);
figure,
subplot(2,2,1), imshow(BWW)
% to digitise an object (e.g. image, profile) means to convert this object into numbers
% to digitise the outline profile generate Coordinates for a Vector Model as follows:
% define/set a number of steps, e.g. 100
ystepsnum= 100;
% to define the step size use a hight of the box,
% i.e. boxsize(4)-see Tutorial 7, Example II-5.
step=round((boxsize(4))/ystepsnum);
% to generate coordinates of the bean useimpixel
% Note: a Vector Model of the bean profile is defined by a set of coordinates
forI=1:ystepsnum
ystep=1+(I-1)*step;
forK=1:boxsize(3)
Px=impixel(BWW,K,ystep);
ifPx>0
X(I)= K- boxsize(3)/2;
Y(I)= boxsize(4)- ystep;
break
end
end
end
%use theVector Model to plot the profile of the bean
subplot (2,2,2), plot(X,Y)
axis equal
请注意!引用、转贴本文应注明原收集者:Rosen Jiang 以及出处:http://www.blogjava.net/rosen