capuaf/convertv.m at master · celsoa/capuaf · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
function [vout,T] = convertv(i1,i2,v)
%CONVERTV convert vectors among different bases
%
% INPUT
%   i1      index of input moment tensor basis (see convert_MT.m)
%   i2      index of output moment tensor basis (see convert_MT.m)
%   V       3 x n set of vectors (or 3 x 3 x n set of bases)
%
% OUTPUT
%   vout    3 x n set of vectors in basis of i2 (or 3 x 3 x n set of bases in i2)
%   T       transformation matrix to change basis of v from i1 to i2:
%              vout = T*v
%
% See also convert_MT.m
%
% calls convert_getbasis.m
%
% Carl Tape, 10/2013
%

% get transformation matrix
T = convert_getbasis(i1,i2);

[a,b,c] = size(v);
if c==1     % v is 3 x n matrix
    if a~=3
        v = v'; n = b; disp('convertv.m: taking transpose of input v');
    else
        n = a;
        disp(sprintf('convertv.m: transforming a set of %i vectors by T',n));
    end
    vout = T*v;
else        % v is 3 x 3 x n
    if and(a==3,b==3)
        disp(sprintf('convertv.m: transforming a set of %i 3 x 3 bases by T',c));
        vout = NaN(a,b,c);
        for ii=1:c
           vout(:,:,ii) = T*v(:,:,ii);
        end
    else
        whos v
        disp('check dimensions of input v');
    end
end

%==========================================================================
% EXAMPLES

if 0==1
    % input standard basis for up-south-east
    % in this case we convert from up-south-east to north-east-down
    % (1,0,0) points up    for v and is (0,0,-1) for vout
    % (0,1,0) points south for v and is (-1,0,0) for vout
    % (0,0,1) points east  for v and is (0,1,0) for vout
    i1 = 1; i2 = 2;
    v = eye(3,3)
    [vout,T] = convertv(i1,i2,v)

    % what about the other way?
    % in this case we convert from north-east-down to up-south-east
    % (1,0,0) points north for v and is (0,-1,0) for vout
    % (0,1,0) points east  for v and is (0,0,1) for vout
    % (0,0,1) points down  for v and is (-1,0,0) for vout
    i1 = 2; i2 = 1;
    v = eye(3,3)
    [vout,T] = convertv(i1,i2,v)

    % check that converting the basis vectors is the same as converting
    % the moment tensor
    % Mout = T Min T' = T U D U' T' = (T U) D (T U)' = Uout D Uout'
    i1 = 1; i2 = 2;
    M = rand(6,1);
    Mmat = Mvec2Mmat(M,1);
    [U,D] = eig(Mmat);
    % convert the basis vectors
    [Uout,T] = convertv(i1,i2,U)
    Mcheck = Uout*D*Uout'
    % convert the moment tensor
    Mout = convert_MT(i1,i2,M);
    Mvec2Mmat(Mout,1)
end


%==========================================================================