1########################################################################
2##
3## Copyright (C) 1993-2021 The Octave Project Developers
4##
5## See the file COPYRIGHT.md in the top-level directory of this
6## distribution or <https://octave.org/copyright/>.
7##
8## This file is part of Octave.
9##
10## Octave is free software: you can redistribute it and/or modify it
11## under the terms of the GNU General Public License as published by
12## the Free Software Foundation, either version 3 of the License, or
13## (at your option) any later version.
14##
15## Octave is distributed in the hope that it will be useful, but
16## WITHOUT ANY WARRANTY; without even the implied warranty of
17## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18## GNU General Public License for more details.
19##
20## You should have received a copy of the GNU General Public License
21## along with Octave; see the file COPYING.  If not, see
22## <https://www.gnu.org/licenses/>.
23##
24########################################################################
25
26## -*- texinfo -*-
27## @deftypefn  {} {} hankel (@var{c})
28## @deftypefnx {} {} hankel (@var{c}, @var{r})
29## Return the Hankel matrix constructed from the first column @var{c}, and
30## (optionally) the last row @var{r}.
31##
32## If the last element of @var{c} is not the same as the first element of
33## @var{r}, the last element of @var{c} is used.  If the second argument is
34## omitted, it is assumed to be a vector of zeros with the same size as
35## @var{c}.
36##
37## A Hankel matrix formed from an m-vector @var{c}, and an n-vector @var{r},
38## has the elements
39## @tex
40## $$
41## H(i, j) = \cases{c_{i+j-1},&$i+j-1\le m$;\cr r_{i+j-m},&otherwise.\cr}
42## $$
43## @end tex
44## @ifnottex
45##
46## @example
47## @group
48## H(i,j) = c(i+j-1),  i+j-1 <= m;
49## H(i,j) = r(i+j-m),  otherwise
50## @end group
51## @end example
52##
53## @end ifnottex
54## @seealso{hadamard, toeplitz}
55## @end deftypefn
56
57function retval = hankel (c, r)
58
59  if (nargin < 1 || nargin > 2)
60    print_usage ();
61  endif
62
63  if (nargin == 1)
64
65    if (! isvector (c))
66      error ("hankel: C must be a vector");
67    endif
68
69    nr = length (c);
70    nc = nr;
71    data = [c(:) ; zeros(nr, 1)];
72
73  else
74
75    if (! (isvector (c) && isvector (r)))
76      error ("hankel: C and R must be vectors");
77    elseif (r(1) != c(end))
78      warning ("hankel: column wins anti-diagonal conflict");
79    endif
80
81    nr = length (c);
82    nc = length (r);
83    data = [c(:) ; r(2:end)(:)];
84
85  endif
86
87  slices = cellslices (data, 1:nc, nr:1:nc+nr-1);
88  retval = horzcat (slices{:});
89
90endfunction
91
92
93%!assert (hankel (1), [1])
94%!assert (hankel ([1, 2]), [1, 2; 2, 0])
95%!assert (hankel ([1, 2], [2; -1; -3]), [1, 2, -1; 2, -1, -3])
96%!assert (hankel (1:3), [1,2,3;2,3,0;3,0,0])
97%!assert (hankel (1:3,3:6), [1,2,3,4;2,3,4,5;3,4,5,6])
98%!assert (hankel (1:3,3:4), [1,2;2,3;3,4])
99%!warning <column wins anti-diagonal conflict>
100%!  assert (hankel (1:3,4:6), [1,2,3;2,3,5;3,5,6]);
101
102%!error hankel ()
103%!error hankel (1, 2, 3)
104%!error <C must be a vector> hankel ([1, 2; 3, 4])
105%!error <C and R must be vectors> hankel (1:4, [1, 2; 3, 4])
106