Matlab - Scripts

From XennisWiki
Jump to: navigation, search

A collection of short Matlab scripts

Plots

Normal (Gaussian) distribution (normpdf)

Matlab - Normal distribution.png

x = [-6:0.01:8];
N_0_2 = normpdf(x,0,sqrt(2));  % Normal distribution N(0,2) 
N_3_1 = normpdf(x,3,1);        % Normal distribution N(3,1)
 
figure
    plot(x, N_0_2, x, N_3_1)
    grid on
    legend('N(0,2)', 'N(3,1)')
    xlabel('x')
    ylabel('\Phi_{\mu,\sigma^2}(x)')

Lagrange multiplier (contour)

Matlab - Lagrange multiplier.png

f = @(x,y) x.^2 + y.^2 + 1;
g = @(x,y) x + y + 2*sqrt(2);

x = linspace(-3,3);
y = linspace(-3,3);
[X,Y] = meshgrid(x, y);
Z = f(X,Y);

figure
    hold on
        contour(X, Y, Z, 1:10)
        % g(x,y) = 0 -> x = - y -2*sqrt(2)
        x = - y - 2*sqrt(2);
        plot(x, y, 'LineWidth', 2)
    hold off
    grid on
    xlabel('x')
    ylabel('y')
    legend('f(x,y)', 'g(x,y)')

Two dimensional data points (marker)

Matlab - Data points.png

% Given data points
w1 = [ 4  2  1;
      -3 -2  3];
w2 = [-1 -4 -2;
      -3  3  2];

figure
    hold on
        plot(w1(1,:), w1(2,:), 'bo', 'MarkerSize', 15)
        plot(w2(1,:), w2(2,:), 'rx', 'MarkerSize', 15)
        plot(4, 3, 'd', 'MarkerSize', 15)
    hold off
    title('Data points');
    legend('Class \omega _1', 'Class \omega _2', '(4,3)^T', 'Location', 'southwest')
    axis([-6 6 -6 6]);
    xlabel('x');
    ylabel('y');
    grid on

Norms

Own implementation of the 1-norm, 2-norm and the infinite-norm and plot of these three norms using surfc.

Matlab Normen.png

A = meshgrid(-2:0.2:2); % Untere Grenze: -2, Obere Grenze: 2, Intervall: 0.2
B = A;
C = A;

%{
Funktion l1
    Die Funktion erhält als Eingabeparameter die zwei Einträge eines
    Vektor aus dem R^2 und gibt die Summe der Beträge (Befehl: abs) der
    Vektoreinträge zurück.

Funktion l2
    Die Funktion erhält als Eingabeparameter die zwei Einträge eines
    Vektors aus dem R^2. Jeder Vektoreintrag wird quadiert und die Summe
    aus diesen gebildet. Aus der Summe wird die Wurzel (Befehl sqrt)
    gezogen und zurückgegeben.

Funktion lunendlich
    Die Funktion erhält als Eingabeparameter die zwei Einträge eines
    Vektors aus dem R^2 und gibt den maximalen Wert (Befehl max()) der
    Beträge (Befehl abs) der Vektoreinträge zurück.
%}
l1 = @(x, y) abs(x) + abs(y);
l2 = @(x, y) sqrt(x^2 + y^2);
lunendlich = @(x, y) max(abs(x), abs(y));

% For-Schleife durchläuft die Zeilen der Matrizen
for n = 1:length(A)
    
    % For-Schleife durchläuft die Spalten der Matrizen
    for m = 1:length(A)
        % Falls sich der Eintrag (n, m) auf der Hauptdiagonalen oder
        % oberhalb von dieser befindet
        if n <= m
            % Die entsprechende Norm wird aus den Einträgen (m, n) und
            % (n, m) der jeweiligen Matrix berechnet uns in dieser
            % gespeichert
            A(n, m) = l1( A(m, n), A(n, m) );
            B(n, m) = l2( B(m, n), B(n, m) );
            C(n, m) = lunendlich( C(m, n), C(n, m) );
        % Eintrag (n, m) befindet sich unterhalb der Hauptdiagonalen
        else
            % Die Symmetrie wird ausgenutzt
            A(n, m) = A(m, n);
            B(n, m) = B(m, n);
            C(n, m) = C(m, n);
        end
    end
end


figure

% 1-Norm zeichnen
subplot(1, 3, 1)
surfc(A)
title('1-Norm')

% 2-Norm zeichnen
subplot(1, 3, 2)
surfc(B)
title('2-Norm')

% Unendlich-Norm zeichnen
subplot(1, 3, 3)
surfc(C)
title('Unendlich-Norm')

Bode plot (bode)

Matlab - Bode plot.png

s = tf('s');
G = 100 / ( s*(1+0.1*s)*(1+0.5*s) );
bode(G)

Nyquist plot (nyquist)

Matlab - Nyquist plot.png

s = tf('s');
G = 100 / ( s*(1+0.1*s)*(1+0.5*s) );
nyquist(G)

Root locus (rlocus)

Matlab - Root locus.png

s = tf('s');
G = 100 / ( s*(1+0.1*s)*(1+0.5*s) );
rlocus(G)

Parallelepiped (plot3d)

Plot a parallelepiped for a given set of points. In the extended version the parallelepiped gets rotated.

Matlab Parallelepiped.png

% The given data points
E1 = [0 0 0]';
E2 = [3 0 0]';
E3 = [0 1 0]';
E4 = [3 1 0]';
E5 = [0 0 1]';
E6 = [3 0 1]';
E7 = [0 1 1]';
E8 = [3 1 1]';

% Assemble points to matrix E
E  = [E1 E2 E6 E2 E4 E8 E4 E3 E7 E3 E1 E5 E6 E8 E7 E5];

xE = E(1,:);
yE = E(2,:);
zE = E(3,:);

% Plot
figure()
    plot3(xE, yE, zE, '.-')
    axis([0,3,0,3,0,3])
    title('Parallelepiped')
    xlabel('x')
    ylabel('y')
    zlabel('z')
    grid on

Code extension

Matlab Parallelepiped extension.png

% The given data points
E1 = [0 0 0]';
E2 = [3 0 0]';
E3 = [0 1 0]';
E4 = [3 1 0]';
E5 = [0 0 1]';
E6 = [3 0 1]';
E7 = [0 1 1]';
E8 = [3 1 1]';

% Assemble points to matrix E
E  = [E1 E2 E6 E2 E4 E8 E4 E3 E7 E3 E1 E5 E6 E8 E7 E5];
xE = E(1,:);
yE = E(2,:);
zE = E(3,:);

%
A = [0 0.5 0; 0.4 0 0; 1 -1 0];
EA = A * E;
xEA = EA(1,:);
yEA = EA(2,:);
zEA = EA(3,:);

w1 = 1/sqrt(2) * [1 1 0]';
w2 = 1/sqrt(2) * [1 -1 0]';
w3 = [0 0 1]';
B = [w1 w2 w3];
EB = B * E;
xEB = EB(1,:);
yEB = EB(2,:);
zEB = EB(3,:);

% Plot
figure()
   p = plot3(xE, yE, zE, '.-b', xEA, yEA, zEA, '.-r', xEB, yEB, zEB, '.-g');
   set(p, 'LineWidth', 2)
   axis([-1,3,-1,3,-1,3])
   title('Parallelepiped')
   xlabel('x')
   ylabel('y')
   zlabel('z')
   legend('P', 'S(P)', 'T(P)')
   grid on

Scripts only

Bayer matrix

The RGB components Ib, Ir and Ig of the Bayer matrix are given and they should be assembled to the picture I.

load('Bayer.mat') % Contains Ib, Ir and Ig

% Filter mask rot and blue
h_rb = [0.25	0.5 	0.25;
		0.5 	1 	0.5;
		0.25 	0.5 	0.25];
% Filter mask green
h_g = [	0	0.25	0;
	0.25 	1 	0.25;
	0 	0.25 	0];

b = conv2(Ib, h_rb, 'valid');
r = conv2(Ir, h_rb, 'valid');
g = conv2(Ig, h_g,  'valid');

I(:,:,1)=r;
I(:,:,2)=g;
I(:,:,3)=b;

imagesc(uint8(I))

Save figure as PDF

function [] = saveFigureAsPdf( figure_handle, pdf_file_name )
% saveFigureAsPdf saves a figure in a pdf.
%
%   Example: saveFigureAsPdf('gcf', 'myFile.pdf');
%
%	figure_handle: a figure handle, e.g. 'gcf' for current figure handle

    if nargin < 2
        error('Not enough input arguments.');
    elseif isempty(pdf_file_name)
        error('Argument is empty.');
    end
    
    %
    % Expand axes to fill figure
    %
    style = hgexport('factorystyle');
    style.Bounds = 'tight';
    hgexport(figure_handle, '-clipboard', style, 'applystyle', true);

    % Set the orientation to landscape
    set(figure_handle,'PaperOrientation','landscape');
    %
    set(figure_handle,'PaperPosition', [1 1 28 19]);
    % Save PDF
    print(gcf, '-dpdf', pdf_file_name);

end

Other

Image gradient

Matlab - Image gradient.png

I = double(imread('image.jpg'));
[F_X F_Y] = gradient(I);
F_XY = sqrt(F_X.^2 + F_Y.^2);

figure
	subplot(2,2,1)
		imshow(I, [])
		title('Original photo')
	subplot(2,2,2)
		imshow(F_XY, [])
		title('Gradient')
	subplot(2,2,3)
		imshow(F_X, [])
		title('Gradient in x direction')
	subplot(2,2,4)
		imshow(F_Y, [])
		title('Gradient in y direction')

Key point detection

Matlab - Key point detection.png

keyPointDetection.m

function [K, H, locationsLocalMaxima] = keyPointDetection(I, threshold)
% Key point detector
%
% Input
%	I: Image
%	threshold: Threshold
%
% Output:
%	K:
%	H: 
%	locationsLocalMaxima:
%
[imageHeight, imageWidth] = size(I);
K = zeros(imageHeight, imageWidth);
H = zeros(imageHeight, imageWidth);

gauss = fspecial('gaussian', 15);

% 1. Compute the partial derivatives I_x and I_y of the image itensity I(x,y)
[I_x, I_y] = gradient(I);
% 2. Compute the product images I_xx, I_yy, I_xy --------------------------
% 3. Convolve the product image with a Gaussian kernel (the width of which
% is a parameter)
I_x2 = imfilter(I_x .^ 2, gauss);
I_y2 = imfilter(I_y .^ 2, gauss);
I_xy = imfilter(I_x .* I_y, gauss);

for x = 1:imageWidth
	for y = 1:imageHeight
		J = [I_x2(y, x)   I_xy(y, x);
			 I_xy(y, x) I_y2(y, x)];
		H(y,x) = trace(J);
		K(y,x) = det(J);
	end
end

% 5. find the locations of local maxima of the interest measure; maybe use
% non-maxima supression for sparser key points
locationsLocalMaxima = K > threshold + imdilate(K, [1 1 1; 1 0 1; 1 1 1]);

end

main.m

I = double(imread('image.jpg'));
[K, H, locationsLocalMaxima] = keyPointDetection(I, 900);

figure
	imshow(I, [])
	hold on
	[col, row]=find(locationsLocalMaxima);
	plot(row, col, 'r.')
	hold off
	title(['K, threshold=', num2str(threshold)])

See also