1 ## Copyright (C) 2010 VZLU Prague, a.s., Czech Republic
3 ## This program is free software; you can redistribute it and/or modify it under
4 ## the terms of the GNU General Public License as published by the Free Software
5 ## Foundation; either version 3 of the License, or (at your option) any later
8 ## This program is distributed in the hope that it will be useful, but WITHOUT
9 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
10 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
13 ## You should have received a copy of the GNU General Public License along with
14 ## this program; if not, see <http://www.gnu.org/licenses/>.
17 ## @deftypefn{Function File} {@var{y} =} ndcovlt (@var{x}, @var{t1}, @var{t2}, @dots{})
18 ## Computes an n-dimensional covariant linear transform of an n-d tensor, given a
19 ## transformation matrix for each dimension. The number of columns of each transformation
20 ## matrix must match the corresponding extent of @var{x}, and the number of rows determines
21 ## the corresponding extent of @var{y}. For example:
24 ## size (@var{x}, 2) == columns (@var{t2})
25 ## size (@var{y}, 2) == rows (@var{t2})
28 ## The element @code{@var{y}(i1, i2, @dots{})} is defined as a sum of
31 ## @var{x}(j1, j2, @dots{}) * @var{t1}(i1, j1) * @var{t2}(i2, j2) * @dots{}
34 ## over all j1, j2, @dots{}. For two dimensions, this reduces to
36 ## @var{y} = @var{t1} * @var{x} * @var{t2}.'
39 ## [] passed as a transformation matrix is converted to identity matrix for
40 ## the corresponding dimension.
44 ## Author: Jaroslav Hajek <highegg@gmail.com>
46 function y = ndcovlt (x, varargin)
47 nd = max (ndims (x), nargin - 1);
48 varargin = resize (varargin, 1, nd);
53 if (isnumeric (ti) && ndims (ti) == 2)
56 varargin{i} = eye (size (x, i));
57 elseif (c != size (x, i))
58 error ("ndcovt: dimension mismatch for x-th transformation matrix");
61 error ("ndcovt: transformation matrices must be numeric 2d matrices");
66 szy = cellfun (@rows, varargin);
72 ## First transformation.
73 y = ldtrans (x, varargin{1});
75 ## Always shift one dimension.
77 y = ldtrans (permute (y, ldp), varargin{i});
80 ## Permute to normal order now to save one permutation.
82 y = ipermute (y, [nd-1:nd, 1:nd-2]);
85 ## Now multiply from the right.
90 y = reshape (y, [], size (y, nd));
91 y = reshape (y * m.', szy);
95 function y = ldtrans (x, m)
98 y = reshape (m * x(:,:), sz);