annotate scripts/general/cumtrapz.m @ 11523:fd0a3ac60b0e

update copyright notices
author John W. Eaton <jwe@octave.org>
date Fri, 14 Jan 2011 05:47:45 -0500
parents 35adf2a71f3f
children c792872f8942
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
11523
fd0a3ac60b0e update copyright notices
John W. Eaton <jwe@octave.org>
parents: 10690
diff changeset
1 ## Copyright (C) 2000-2011 Kai Habel
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
2 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
3 ## This file is part of Octave.
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
4 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
5 ## Octave is free software; you can redistribute it and/or modify it
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
6 ## under the terms of the GNU General Public License as published by
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
7 ## the Free Software Foundation; either version 3 of the License, or (at
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
8 ## your option) any later version.
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
9 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
10 ## Octave is distributed in the hope that it will be useful, but
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
13 ## General Public License for more details.
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
14 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
15 ## You should have received a copy of the GNU General Public License
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
16 ## along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 7001
diff changeset
17 ## <http://www.gnu.org/licenses/>.
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
18
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
19 ## -*- texinfo -*-
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
20 ## @deftypefn {Function File} {@var{z} =} cumtrapz (@var{y})
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
21 ## @deftypefnx {Function File} {@var{z} =} cumtrapz (@var{x}, @var{y})
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
22 ## @deftypefnx {Function File} {@var{z} =} cumtrapz (@dots{}, @var{dim})
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
23 ##
7001
8b0cfeb06365 [project @ 2007-10-10 18:02:59 by jwe]
jwe
parents: 6741
diff changeset
24 ## Cumulative numerical integration using trapezoidal method.
8b0cfeb06365 [project @ 2007-10-10 18:02:59 by jwe]
jwe
parents: 6741
diff changeset
25 ## @code{cumtrapz (@var{y})} computes the cumulative integral of the
9051
1bf0ce0930be Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents: 8920
diff changeset
26 ## @var{y} along the first non-singleton dimension. If the argument
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
27 ## @var{x} is omitted an equally spaced vector is assumed. @code{cumtrapz
7001
8b0cfeb06365 [project @ 2007-10-10 18:02:59 by jwe]
jwe
parents: 6741
diff changeset
28 ## (@var{x}, @var{y})} evaluates the cumulative integral with respect
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
29 ## to @var{x}.
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
30 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
31 ## @seealso{trapz,cumsum}
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
32 ## @end deftypefn
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
33
10549
95c3e38098bf Untabify .m scripts
Rik <code@nomad.inbox5.com>
parents: 9051
diff changeset
34 ## Author: Kai Habel <kai.habel@gmx.de>
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
35 ##
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
36 ## also: June 2000 Paul Kienzle (fixes,suggestions)
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
37 ## 2006-05-12 David Bateman - Modified for NDArrays
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
38
10549
95c3e38098bf Untabify .m scripts
Rik <code@nomad.inbox5.com>
parents: 9051
diff changeset
39 function z = cumtrapz (x, y, dim)
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
40
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
41 if (nargin < 1) || (nargin > 3)
6046
34f96dd5441b [project @ 2006-10-10 16:10:25 by jwe]
jwe
parents: 5820
diff changeset
42 print_usage ();
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
43 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
44
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
45 nd = ndims (x);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
46 sz = size (x);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
47
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
48 have_x = false;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
49 have_dim = false;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
50 if (nargin == 3)
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
51 have_x = true;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
52 have_dim = true;
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
53 elseif (nargin == 2)
6157
045038e0108a [project @ 2006-11-13 22:22:53 by jwe]
jwe
parents: 6046
diff changeset
54 if (! size_equal (x, y) && isscalar (y))
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
55 dim = y;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
56 have_dim = true;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
57 else
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
58 have_x = true;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
59 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
60 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
61
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
62 if (! have_dim)
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
63 ## Find the first non-singleton dimension.
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
64 dim = find (sz > 1, 1);
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
65 if (isempty (dim))
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
66 dim = 1;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
67 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
68 else
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
69 dim = floor (dim);
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
70 if (! (isscalar (dim) && 1 <= dim && dim <= nd))
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
71 error ("cumtrapz: invalid dimension DIM");
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
72 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
73 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
74
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
75 n = sz(dim);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
76 idx1 = cell ();
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
77 for i = 1:nd
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
78 idx1{i} = 1:sz(i);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
79 endfor
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
80 idx2 = idx1;
10549
95c3e38098bf Untabify .m scripts
Rik <code@nomad.inbox5.com>
parents: 9051
diff changeset
81 idx1{dim} = 2 : n;
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
82 idx2{dim} = 1 : (n - 1);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
83
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
84 if (! have_x)
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
85 z = 0.5 * cumsum (x(idx1{:}) + x(idx2{:}), dim);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
86 else
6157
045038e0108a [project @ 2006-11-13 22:22:53 by jwe]
jwe
parents: 6046
diff changeset
87 if (! size_equal (x, y))
10690
35adf2a71f3f Use common code block to find first non-singleton dimension.
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
88 error ("cumtrapz: X and Y must have the same shape");
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
89 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
90 z = 0.5 * cumsum ((x(idx1{:}) - x(idx2{:})) .*
10549
95c3e38098bf Untabify .m scripts
Rik <code@nomad.inbox5.com>
parents: 9051
diff changeset
91 (y(idx1{:}) + y(idx2{:})), dim);
5820
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
92 endif
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
93
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
94 sz(dim) = 1;
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
95 z = cat (dim, zeros (sz), z);
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
96
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
97 endfunction
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
98
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
99 %!shared x1,x2,y
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
100 %! x1 = [0,0,0;2,2,2];
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
101 %! x2 = [0,2,4;0,2,4];
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
102 %! y = [1,2,3;4,5,6];
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
103 %!assert (cumtrapz(y),[0,0,0;2.5,3.5,4.5])
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
104 %!assert (cumtrapz(x1,y),[0,0,0;5,7,9])
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
105 %!assert (cumtrapz(y,1),[0,0,0;2.5,3.5,4.5])
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
106 %!assert (cumtrapz(x1,y,1),[0,0,0;5,7,9])
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
107 %!assert (cumtrapz(y,2),[0,1.5,4;0,4.5,10])
27c966e4b2dc [project @ 2006-05-17 21:00:54 by jwe]
jwe
parents:
diff changeset
108 %!assert (cumtrapz(x2,y,2),[0,3,8;0,9,20])