Skip to content

Commit

Permalink
Update to use more common option name
Browse files Browse the repository at this point in the history
Switch the showplot option to output.
  • Loading branch information
Daniel Skoog committed Sep 20, 2017
1 parent bf6ca81 commit 266e1b4
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 18 deletions.
2 changes: 1 addition & 1 deletion src/ClusterPlot.mm
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ elif hastype(X, {'rtable', 'DataFrame'}) and cluster = NULL then
ldatasplit := Array(1..numdf);
for i to numdf do
ldatasplit(i) := ldata[i];
p1(i) := MVEE(ldatasplit[i][..,1..2], 'showplot', 'color' = plots:-setcolors()[i], _rest);
p1(i) := MVEE(ldatasplit[i][..,1..2], 'output' = 'plot', 'color' = plots:-setcolors()[i], _rest);
end do;
else
p1 := [NULL];
Expand Down
38 changes: 21 additions & 17 deletions src/MVEE.mm
Original file line number Diff line number Diff line change
Expand Up @@ -2,25 +2,30 @@
MVEE := proc(XY,
{tolerance::positive:= 1e-4}, #Convergence Criterion
{maxiterations::posint := 100},
{showplot::truefalse := false},
{output::{identical(data,plot),list(identical(data,plot))} := data},
{filled::truefalse := false}
)

local alpha, count, evalues, evectors, i, l_error, ldata, ldataext, M, maxvalindex, n, ncols, nrows, semiaxes, stepsize, U, U1, x, X, y;
local A, center; #Output
local alpha, evalues, evectors, i, l_error, ldata, ldataext, M, maxvalindex, n, ncols, nrows, p1, semiaxes, stepsize, U, U1, x, X, y;
local A, center, l_output; #Output

if hastype(output, 'list') then
l_output := output;
else
l_output := [output];
end if;

kernelopts(opaquemodules=false):

ldata := Statistics:-PreProcessData(XY, 2, 'copy');

nrows, ncols := upperbound(ldata);
ldataext := Matrix([ldata, Vector[column](nrows, 'fill' = 1)], 'datatype = float');
ldataext := Matrix([ldata, Vector[column](nrows, ':-fill' = 1)], 'datatype = float');

if ncols <> 2 then
error "expected 2 columns of data, got %1", ncols;
end if;

count := 1;
l_error := 1;

U := Vector[column](1..nrows, 'fill' = 1/nrows);
Expand All @@ -35,29 +40,28 @@
U1 := (1 - stepsize) * U;
U1[maxvalindex] := U1[maxvalindex] + stepsize;
l_error := LinearAlgebra:-Norm(LinearAlgebra:-DiagonalMatrix(U1 - U));
count := count + 1;
U := U1;

end do;

A := (1/ncols) * LinearAlgebra:-MatrixInverse(LinearAlgebra:-Transpose(ldata) . LinearAlgebra:-DiagonalMatrix(U) . ldata - (LinearAlgebra:-Transpose(ldata) . U) . LinearAlgebra:-Transpose((LinearAlgebra:-Transpose(ldata) . U)));
center := LinearAlgebra:-Transpose(ldata) . U;
evalues, evectors := LinearAlgebra:-Eigenvectors(A);
evectors := evectors(.., sort[index](1 /~ (sqrt~(Re~(evalues))), `>`, 'output' = 'permutation'));
evectors := evectors(.., sort[index](1 /~ (sqrt~(Re~(evalues))), `>`, ':-output' = ':-permutation'));
semiaxes := sort(1 /~ (sqrt~(Re~(evalues))), `>`);
alpha := arctan(Re(evectors[2,1]) / Re(evectors[1,1]));
if showplot = true then

x := t -> center[1] + semiaxes[1] * cos(t) * cos(alpha) - semiaxes[2] * sin(t) * sin(alpha);
y := t -> center[2] + semiaxes[1] * cos(t) * sin(alpha) + semiaxes[2] * sin(t) * cos(alpha);
if filled then
return plots:-display(subs(CURVES=POLYGONS, plot([x(t), y(t), t = 0..2*Pi], 'transparency' = 0.95, _rest)));
else
return plot([x(t), y(t), t = 0..2*Pi], _rest);

end if;
else
if l_output = [':-data'] then
return A, center, evectors, evalues;
elif has( l_output, ':-plot' ) then
x := t -> center[1] + semiaxes[1] * cos(t) * cos(alpha) - semiaxes[2] * sin(t) * sin(alpha);
y := t -> center[2] + semiaxes[1] * cos(t) * sin(alpha) + semiaxes[2] * sin(t) * cos(alpha);
if filled then
p1 := plots:-display(subs(CURVES=POLYGONS, plot([x(t), y(t), t = 0..2*Pi], ':-transparency' = 0.95, _rest)));
else
p1 := plot([x(t), y(t), t = 0..2*Pi], _rest);
end if;
return p1, `if`( has(l_output, ':-data'), op([A, center, evectors, evalues]), NULL );
end if;

end proc:

0 comments on commit 266e1b4

Please sign in to comment.