1## Copyright (C) 2018 Martin Janda <janda.martin1@gmail.com>
2##
3## This program is free software: you can redistribute it and/or modify it
4## under the terms of the GNU General Public License as published by
5## the Free Software Foundation, either version 3 of the License, or
6## (at your option) any later version.
7##
8## This program is distributed in the hope that it will be useful, but
9## WITHOUT ANY WARRANTY; without even the implied warranty of
10## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
11## GNU General Public License for more details.
12##
13## You should have received a copy of the GNU General Public License
14## along with this program.  If not, see
15## <https://www.gnu.org/licenses/>.
16
17## -*- texinfo -*-
18## @deftypefn {} {[@var{i}, @var{j}, @var{k}] =} worldToSubscript (@var{r}, @var{xWorld}, @var{yWorld}, @var{zWorld})
19## Convert world coordinates to row, column and plane subscripts.
20##
21## Converts world coordinates to row, column and plane subscripts of an image
22## associated with the spatial referencing object @var{r}.  A point located at
23## (@var{xWorld}(i), @var{yWorld}(i), @var{zWorld}(i)) world coordinates maps
24## to row, column and plane subscripts @var{i}(i), @var{j}(i), @var{k}(i)
25## respectively.  Note the reversed order of the first two dimensions.  If the
26## point falls outside the bounds of the image, all of its subscripts are NaN.
27##
28## @seealso{imref2d, imref3d, worldToIntrinsic}
29## @end deftypefn
30
31function [i, j, k] = worldToSubscript (r, xWorld, yWorld, zWorld)
32  if (nargin != 4)
33    print_usage ();
34  endif
35
36  validateattributes (xWorld, {"numeric"}, ...
37  {"real"}, "imref3d", "xWorld");
38  validateattributes (yWorld, {"numeric"}, ...
39  {"real"}, "imref3d", "yWorld");
40  validateattributes (zWorld, {"numeric"}, ...
41  {"real"}, "imref3d", "zWorld");
42
43  if (! all (size (xWorld) == size (yWorld)) ...
44    || ! all (size (xWorld) == size (zWorld)))
45    error ("Octave:invalid-input-arg", ...
46    "xWorld, yWorld and zWorld must be of the same size");
47  endif
48
49  [xIntrinsic, yIntrinsic, zIntrinsic] ...
50  = worldToIntrinsic (r, xWorld, yWorld, zWorld);
51  xIntrinsicLimits = r.XIntrinsicLimits;
52  yIntrinsicLimits = r.YIntrinsicLimits;
53  zIntrinsicLimits = r.ZIntrinsicLimits;
54  inImage = contains (r, xWorld, yWorld, zWorld);
55  xIntrinsic(! inImage) = NaN;
56  yIntrinsic(! inImage) = NaN;
57  zIntrinsic(! inImage) = NaN;
58  i = round (yIntrinsic);
59  j = round (xIntrinsic);
60  k = round (zIntrinsic);
61endfunction
62
63%!error id=Octave:invalid-fun-call worldToSubscript (imref3d)
64%!error id=Octave:invalid-fun-call worldToSubscript (imref3d, 1)
65%!error id=Octave:invalid-fun-call worldToSubscript (imref3d, 1, 2)
66%!error id=Octave:invalid-fun-call worldToSubscript (imref3d, 1, 2, 3, 4)
67%!error id=Octave:expected-real worldToSubscript (imref3d, 1j, 2, 3)
68%!error id=Octave:expected-real worldToSubscript (imref3d, 1, 2j, 3)
69%!error id=Octave:expected-real worldToSubscript (imref3d, 1, 2, 3j)
70%!error id=Octave:invalid-input-arg worldToSubscript (imref3d, [1, 2], 3, 4)
71%!error id=Octave:invalid-input-arg worldToSubscript (imref3d, 1, [2, 3], 4)
72%!error id=Octave:invalid-input-arg worldToSubscript (imref3d, 1, 2, [3, 4])
73
74%!test
75%! r = imref3d ([128, 128, 27], 2, 2, 4);
76%! xW = [108, 108, 113.2, 2];
77%! yW = [92, 92, 92, -1];
78%! zW = [52, 55, 52, 0.33];
79%! [rS, cS, pS] = worldToSubscript (r, xW, yW, zW);
80%! assert (rS, [46, 46, 46, NaN])
81%! assert (cS, [54, 54, 57, NaN])
82%! assert (pS, [13, 14, 13, NaN])