7183
|
1 ## Copyright (C) 2007 Kai Habel, David Bateman |
|
2 ## |
|
3 ## This file is part of Octave. |
|
4 ## |
|
5 ## Octave is free software; you can redistribute it and/or modify it |
|
6 ## under the terms of the GNU General Public License as published by |
|
7 ## the Free Software Foundation; either version 3 of the License, or (at |
|
8 ## your option) any later version. |
|
9 ## |
|
10 ## Octave is distributed in the hope that it will be useful, but |
|
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
|
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
|
13 ## General Public License for more details. |
|
14 ## |
|
15 ## You should have received a copy of the GNU General Public License |
|
16 ## along with Octave; see the file COPYING. If not, see |
|
17 ## <http://www.gnu.org/licenses/>. |
|
18 |
|
19 ## -*- texinfo -*- |
|
20 ## @deftypefn {Function File} {} slice (@var{X}, @var{Y}, @var{Z}, @var{V}, @var{SX}, @var{SY}, @var{SZ}) |
|
21 ## @deftypefnx {Function File} {} slice (@var{X}, @var{Y}, @var{Z}, @var{V}, @var{XI}, @var{YI}, @var{ZI}) |
|
22 ## @deftypefnx {Function File} {} slice (@var{V}, @var{SX}, @var{SY}, @var{SZ}) |
|
23 ## @deftypefnx {Function File} {} slice (@var{V}, @var{XI}, @var{YI}, @var{ZI}) |
|
24 ## @deftypefnx {Function File} {@var{H} =} slice (...) |
|
25 ## @deftypefnx {Function File} {@var{H} =} slice (...,@var{METHOD}) |
|
26 ## Plots slice(s) of 3D data/scalar fields. Each element of then 3-dimensional |
|
27 ## array @var{v} represents a scalar value at a location given by the parameters |
|
28 ## @var{x}, @var{y}, and @var{z}. The parameters @var{x}, @var{x}, and |
|
29 ## @var{z} are either 3-dimensional arrays of the same size as the array |
|
30 ## @var{v} in the 'meshgrid' format or vectors. The parameters @var{xi}, etc |
|
31 ## respect a similar format to @var{x}, etc, and they represent the points |
|
32 ## at which the array @var{vi} is interpolated using interp3. The vectors |
|
33 ## @var{sx}, @var{sy}, and @var{sz} contain points of orthogonal slices of |
|
34 ## the respective axes. |
|
35 ## |
|
36 ## If @var{x}, @var{y}, @var{z} are omitted, they are assumed to be |
|
37 ## @code{x = 1 : size (@var{v}, 2)}, @code{y = 1 : size (@var{v}, 1)} and |
|
38 ## @code{z = 1 : size (@var{v}, 3)}. |
|
39 ## |
|
40 ## @var{Method} is one of: |
|
41 ## |
|
42 ## @table @asis |
|
43 ## @item 'nearest' |
|
44 ## Return the nearest neighbour. |
|
45 ## @item 'linear' |
|
46 ## Linear interpolation from nearest neighbours. |
|
47 ## @item 'cubic' |
|
48 ## Cubic interpolation from four nearest neighbours (not implemented yet). |
|
49 ## @item 'spline' |
|
50 ## Cubic spline interpolation--smooth first and second derivatives |
|
51 ## throughout the curve. |
|
52 ## @end table |
|
53 ## |
|
54 ## The default method is 'linear'. |
|
55 ## The optional return value @var{H} is a vector of handles to the surface graphic |
|
56 ## objects. |
|
57 ## |
|
58 ## Examples: |
|
59 ## @example |
|
60 ## [X,Y,Z] = meshgrid(linspace(-8,8,32)); |
|
61 ## V = sin (sqrt (X.^2 + Y.^2 + Z.^2)) ./ (sqrt (X.^2 + Y.^2 + Z.^2)) |
|
62 ## slice(X,Y,Z,V,[],0,[]) |
|
63 ## [XI,YI]=meshgrid(linspace(-7,7)); |
|
64 ## ZI=XI+YI; |
|
65 ## slice(X,Y,Z,V,XI,YI,ZI) |
|
66 ## @end example |
|
67 ## @seealso{interp3, surface, pcolor} |
|
68 ## @end deftypefn |
|
69 |
|
70 ## Author: Kai Habel <kai.habel at gmx.de> |
|
71 |
|
72 function h = slice(varargin) |
|
73 |
|
74 method = "linear"; |
|
75 extrapval = NA; |
|
76 nargs = nargin; |
|
77 |
|
78 if (ischar (varargin{end})) |
|
79 method = varargin{end}; |
|
80 nargs -= 1; |
|
81 endif |
|
82 |
|
83 if (nargs == 4) |
|
84 V = varargin{1}; |
|
85 if (ndims (V) != 3) |
|
86 error ("slice: expect 3-dimensional array of values"); |
|
87 endif |
|
88 [nx, ny, nz] = size(V); |
|
89 [X,Y,Z] = meshgrid(1:nx,1:ny,1:nz); |
|
90 sx = varargin{2}; |
|
91 sy = varargin{3}; |
|
92 sz = varargin{4}; |
|
93 elseif (nargs == 7) |
|
94 V = varargin{4}; |
|
95 if (ndims (V) != 3) |
|
96 error ("slice: expect 3-dimensional array of values"); |
|
97 endif |
|
98 X = varargin{1}; |
|
99 Y = varargin{2}; |
|
100 Z = varargin{3}; |
|
101 if (all([isvector(X) isvector(Y) isvector(Z)])) |
|
102 [X,Y,Z] = meshgrid(X,Y,Z); |
|
103 elseif ((ndims(X) == 3) && size_equal(X,Y) && size_equal(X,Z)) |
|
104 ##do nothing |
|
105 else |
|
106 error("slice: X,Y,Z size mismatch") |
|
107 endif |
|
108 sx = varargin{5}; |
|
109 sy = varargin{6}; |
|
110 sz = varargin{7}; |
|
111 else |
|
112 print_usage(); |
|
113 endif |
|
114 |
|
115 if (any([isvector(sx), isvector(sy), isvector(sz)])) |
|
116 have_sval = true(); |
|
117 elseif ((ndims(sx) == 2) && size_equal(sx,sy) && size_equal(sx,sz)) |
|
118 have_sval = false(); |
|
119 else |
|
120 error ("slice: dimensional mismatch for (XI,YI,ZI) or (sx,sy,sz)"); |
|
121 endif |
|
122 |
|
123 newplot (); |
|
124 ax = gca; |
|
125 sidx = 1; |
|
126 maxv = max(V(:)); |
|
127 minv = min(V(:)); |
|
128 set(ax, "CLim", [minv, maxv]); |
|
129 |
|
130 if (have_sval) |
|
131 ns = length(sx) + length(sy) + length(sz); |
|
132 hs = zeros(ns,1); |
|
133 [ny, nx, nz] = size(V); |
|
134 if (length(sz) > 0) |
|
135 for i=1:length(sz) |
|
136 [XI,YI,ZI] = meshgrid(squeeze(X(1,:,1)),squeeze(Y(:,1,1)),sz(i)); |
|
137 Vz = squeeze(interp3(X,Y,Z,V,XI,YI,ZI,method)); |
|
138 tmp(sidx++) = surface(XI,YI,sz(i)*ones(size(YI)),Vz); |
|
139 endfor |
|
140 endif |
|
141 |
|
142 if (length(sy) > 0) |
|
143 for i=length(sy):-1:1 |
|
144 [XI,YI,ZI] = meshgrid(squeeze(X(1,:,1)),sy(i),squeeze(Z(1,1,:))); |
|
145 Vy = squeeze(interp3(X,Y,Z,V,XI,YI,ZI,method)); |
|
146 tmp(sidx++) = surface(squeeze(XI),squeeze(sy(i)*ones(size(ZI))),squeeze(ZI),Vy); |
|
147 endfor |
|
148 endif |
|
149 |
|
150 if (length(sx) > 0) |
|
151 for i=length(sx):-1:1 |
|
152 [XI,YI,ZI] = meshgrid(sx(i),squeeze(Y(:,1,1)),squeeze(Z(1,1,:))); |
|
153 Vx = squeeze(interp3(X,Y,Z,V,XI,YI,ZI,method)); |
|
154 tmp(sidx++) = surface(squeeze(sx(i)*ones(size(ZI))),squeeze(YI),squeeze(ZI),Vx); |
|
155 endfor |
|
156 endif |
|
157 else |
|
158 VI = interp3(X,Y,Z,V,sx,sy,sz); |
|
159 tmp(sidx++) = surface(sx,sy,sz,VI); |
|
160 endif |
|
161 |
|
162 if (! ishold ()) |
|
163 set (ax, "view", [-37.5, 30.0]); |
|
164 endif |
|
165 |
|
166 if (nargout > 0) |
|
167 h = tmp; |
|
168 endif |
|
169 |
|
170 end |