|
wolffd@0
|
1 % PLOTCOV2 - Plots a covariance ellipsoid with axes for a bivariate
|
|
wolffd@0
|
2 % Gaussian distribution.
|
|
wolffd@0
|
3 %
|
|
wolffd@0
|
4 % Usage:
|
|
wolffd@0
|
5 % [h, s] = plotcov2(mu, Sigma[, OPTIONS]);
|
|
wolffd@0
|
6 %
|
|
wolffd@0
|
7 % Inputs:
|
|
wolffd@0
|
8 % mu - a 2 x 1 vector giving the mean of the distribution.
|
|
wolffd@0
|
9 % Sigma - a 2 x 2 symmetric positive semi-definite matrix giving
|
|
wolffd@0
|
10 % the covariance of the distribution (or the zero matrix).
|
|
wolffd@0
|
11 %
|
|
wolffd@0
|
12 % Options:
|
|
wolffd@0
|
13 % 'conf' - a scalar between 0 and 1 giving the confidence
|
|
wolffd@0
|
14 % interval (i.e., the fraction of probability mass to
|
|
wolffd@0
|
15 % be enclosed by the ellipse); default is 0.9.
|
|
wolffd@0
|
16 % 'num-pts' - if the value supplied is n, then (n + 1)^2 points
|
|
wolffd@0
|
17 % to be used to plot the ellipse; default is 20.
|
|
wolffd@0
|
18 % 'label' - if non-empty, a string that will label the
|
|
wolffd@0
|
19 % ellipsoid (default: [])
|
|
wolffd@0
|
20 % 'plot-axes' - a 0/1 flag indicating if the ellipsoid's axes
|
|
wolffd@0
|
21 % should be plotted (default: 1)
|
|
wolffd@0
|
22 % 'plot-opts' - a cell vector of arguments to be handed to PLOT3
|
|
wolffd@0
|
23 % to contol the appearance of the axes, e.g.,
|
|
wolffd@0
|
24 % {'Color', 'g', 'LineWidth', 1}; the default is {}
|
|
wolffd@0
|
25 % 'fill-color' - a color specifier; is this is not [], the
|
|
wolffd@0
|
26 % covariance ellipse is filled with this color
|
|
wolffd@0
|
27 % (default: [])
|
|
wolffd@0
|
28 %
|
|
wolffd@0
|
29 % Outputs:
|
|
wolffd@0
|
30 % h - a vector of handles on the axis lines
|
|
wolffd@0
|
31 %
|
|
wolffd@0
|
32 % See also: PLOTCOV3
|
|
wolffd@0
|
33
|
|
wolffd@0
|
34 % Copyright (C) 2002 Mark A. Paskin
|
|
wolffd@0
|
35 %
|
|
wolffd@0
|
36 % This program is free software; you can redistribute it and/or modify
|
|
wolffd@0
|
37 % it under the terms of the GNU General Public License as published by
|
|
wolffd@0
|
38 % the Free Software Foundation; either version 2 of the License, or
|
|
wolffd@0
|
39 % (at your option) any later version.
|
|
wolffd@0
|
40 %
|
|
wolffd@0
|
41 % This program is distributed in the hope that it will be useful, but
|
|
wolffd@0
|
42 % WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
wolffd@0
|
43 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
wolffd@0
|
44 % General Public License for more details.
|
|
wolffd@0
|
45 %
|
|
wolffd@0
|
46 % You should have received a copy of the GNU General Public License
|
|
wolffd@0
|
47 % along with this program; if not, write to the Free Software
|
|
wolffd@0
|
48 % Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
|
|
wolffd@0
|
49 % USA.
|
|
wolffd@0
|
50 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
|
wolffd@0
|
51
|
|
wolffd@0
|
52 function [h, s] = plotcov2New(mu, Sigma, varargin)
|
|
wolffd@0
|
53
|
|
wolffd@0
|
54 h = [];
|
|
wolffd@0
|
55 s = [];
|
|
wolffd@0
|
56
|
|
wolffd@0
|
57 if size(Sigma) ~= [2 2], error('Sigma must be a 2 by 2 matrix'); end
|
|
wolffd@0
|
58 if length(mu) ~= 2, error('mu must be a 2 by 1 vector'); end
|
|
wolffd@0
|
59
|
|
wolffd@0
|
60 Sigma = checkpsd(Sigma);
|
|
wolffd@0
|
61
|
|
wolffd@0
|
62 [p, ...
|
|
wolffd@0
|
63 n, ...
|
|
wolffd@0
|
64 label, ...
|
|
wolffd@0
|
65 plot_axes, ...
|
|
wolffd@0
|
66 plot_opts, ...
|
|
wolffd@0
|
67 fill_color] = process_options(varargin, 'conf', 0.9, ...
|
|
wolffd@0
|
68 'num-pts', 20, ...
|
|
wolffd@0
|
69 'label', [], ...
|
|
wolffd@0
|
70 'plot-axes', 1, ...
|
|
wolffd@0
|
71 'plot-opts', {}, ...
|
|
wolffd@0
|
72 'fill-color', []);
|
|
wolffd@0
|
73 holding = ishold;
|
|
wolffd@0
|
74 % Compute the Mahalanobis radius of the ellipsoid that encloses
|
|
wolffd@0
|
75 % the desired probability mass.
|
|
wolffd@0
|
76 k = conf2mahal(p, 2);
|
|
wolffd@0
|
77 % Scale the covariance matrix so the confidence region has unit
|
|
wolffd@0
|
78 % Mahalanobis distance.
|
|
wolffd@0
|
79 Sigma = Sigma * k;
|
|
wolffd@0
|
80 % The axes of the covariance ellipse are given by the eigenvectors of
|
|
wolffd@0
|
81 % the covariance matrix. Their lengths (for the ellipse with unit
|
|
wolffd@0
|
82 % Mahalanobis radius) are given by the square roots of the
|
|
wolffd@0
|
83 % corresponding eigenvalues.
|
|
wolffd@0
|
84 [V, D] = eig(full(Sigma));
|
|
wolffd@0
|
85 V = real(V);
|
|
wolffd@0
|
86 D = real(D);
|
|
wolffd@0
|
87 D = abs(D);
|
|
wolffd@0
|
88
|
|
wolffd@0
|
89 % Compute the points on the boundary of the ellipsoid.
|
|
wolffd@0
|
90 t = linspace(0, 2*pi, n);
|
|
wolffd@0
|
91 u = [cos(t(:))'; sin(t(:))'];
|
|
wolffd@0
|
92 w = (V * sqrt(D)) * u;
|
|
wolffd@0
|
93 z = repmat(mu(:), [1 n]) + w;
|
|
wolffd@0
|
94 h = [h; plot(z(1, :), z(2, :), plot_opts{:})];
|
|
wolffd@0
|
95 if (~isempty(fill_color))
|
|
wolffd@0
|
96 s = patch(z(1, :), z(2, :), fill_color);
|
|
wolffd@0
|
97 end
|
|
wolffd@0
|
98
|
|
wolffd@0
|
99 % Plot the axes.
|
|
wolffd@0
|
100 if (plot_axes)
|
|
wolffd@0
|
101 hold on;
|
|
wolffd@0
|
102 L = sqrt(diag(D));
|
|
wolffd@0
|
103 h = plot([mu(1); mu(1) + L(1) * V(1, 1)], ...
|
|
wolffd@0
|
104 [mu(2); mu(2) + L(1) * V(2, 1)], plot_opts{:});
|
|
wolffd@0
|
105 h = [h; plot([mu(1); mu(1) + L(2) * V(1, 2)], ...
|
|
wolffd@0
|
106 [mu(2); mu(2) + L(2) * V(2, 2)], plot_opts{:})];
|
|
wolffd@0
|
107 end
|
|
wolffd@0
|
108
|
|
wolffd@0
|
109
|
|
wolffd@0
|
110 if (~isempty(label))
|
|
wolffd@0
|
111 th = text(mu(1), mu(2), label);
|
|
wolffd@0
|
112 set(th, 'FontSize', 18);
|
|
wolffd@0
|
113 set(th, 'FontName', 'Times');
|
|
wolffd@0
|
114 set(th, 'FontWeight', 'bold');
|
|
wolffd@0
|
115 set(th, 'FontAngle', 'italic');
|
|
wolffd@0
|
116 set(th, 'HorizontalAlignment', 'center');
|
|
wolffd@0
|
117 end
|
|
wolffd@0
|
118
|
|
wolffd@0
|
119 if (~holding & plot_axes) hold off; end
|
