← Back to branch summary

~gabriel1984sibiu/octave/octave

« back to all changes in this revision

Viewing changes to doc/interpreter/splineimages.m

  • Committer: Grevutiu Gabriel
  • Date: 2014-01-02 13:05:54 UTC
  • Revision ID: gabriel1984sibiu@gmail.com-20140102130554-3r7ivdjln1ni6kcg
New version (3.8.0) from upstream.

Show diffs side-by-side

added added

removed removed

Lines of Context:
doc/interpreter/splineimages.m
 
1
## Copyright (C) 2012-2013 Ben Abbott, Jonas Lundgren
 
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
function splineimages (nm, typ)
 
20
  graphics_toolkit ("gnuplot");
 
21
  set_print_size ();
 
22
  hide_output ();
 
23
  if (strcmp (typ, "png"))
 
24
    set (0, "defaulttextfontname", "*");
 
25
  endif
 
26
  if (strcmp (typ, "eps"))
 
27
    d_typ = "-depsc2";
 
28
  else
 
29
    d_typ = ["-d" typ];
 
30
  endif
 
31
 
 
32
  if (strcmp (typ, "txt"))
 
33
    image_as_txt (nm);
 
34
  elseif (strcmp (nm, "splinefit1")) ## Breaks and Pieces
 
35
    x = 2 * pi * rand (1, 200);
 
36
    y = sin (x) + sin (2 * x) + 0.2 * randn (size (x));
 
37
    ## Uniform breaks
 
38
    breaks = linspace (0, 2 * pi, 41); ## 41 breaks, 40 pieces
 
39
    pp1 = splinefit (x, y, breaks);
 
40
    ## Breaks interpolated from data
 
41
    pp2 = splinefit (x, y, 10);  ## 11 breaks, 10 pieces
 
42
    ## Plot
 
43
    xx = linspace (0, 2 * pi, 400);
 
44
    y1 = ppval (pp1, xx);
 
45
    y2 = ppval (pp2, xx);
 
46
    plot (x, y, ".", xx, [y1; y2]);
 
47
    axis tight;
 
48
    ylim ([-2.5 2.5]);
 
49
    legend ("data", "41 breaks, 40 pieces", "11 breaks, 10 pieces");
 
50
    print ([nm "." typ], d_typ);
 
51
  elseif (strcmp (nm, "splinefit2")) ## Spline orders
 
52
    ## Data (200 points)
 
53
    x = 2 * pi * rand (1, 200);
 
54
    y = sin (x) + sin (2 * x) + 0.1 * randn (size (x));
 
55
    ## Splines
 
56
    pp1 = splinefit (x, y, 8, "order", 0);  ## Piecewise constant
 
57
    pp2 = splinefit (x, y, 8, "order", 1);  ## Piecewise linear
 
58
    pp3 = splinefit (x, y, 8, "order", 2);  ## Piecewise quadratic
 
59
    pp4 = splinefit (x, y, 8, "order", 3);  ## Piecewise cubic
 
60
    pp5 = splinefit (x, y, 8, "order", 4);  ## Etc.
 
61
    ## Plot
 
62
    xx = linspace (0, 2 * pi, 400);
 
63
    y1 = ppval (pp1, xx);
 
64
    y2 = ppval (pp2, xx);
 
65
    y3 = ppval (pp3, xx);
 
66
    y4 = ppval (pp4, xx);
 
67
    y5 = ppval (pp5, xx);
 
68
    plot (x, y, ".", xx, [y1; y2; y3; y4; y5]);
 
69
    axis tight;
 
70
    ylim ([-2.5 2.5]);
 
71
    legend ({"data", "order 0", "order 1", "order 2", "order 3", "order 4"});
 
72
    print ([nm, "." typ], d_typ);
 
73
  elseif (strcmp (nm, "splinefit3"))
 
74
    ## Data (100 points)
 
75
    x = 2 * pi * [0, (rand (1, 98)), 1];
 
76
    y = sin (x) - cos (2 * x) + 0.2 * randn (size (x));
 
77
    ## No constraints
 
78
    pp1 = splinefit (x, y, 10, "order", 5);
 
79
    ## Periodic boundaries
 
80
    pp2 = splinefit (x, y, 10, "order", 5, "periodic", true);
 
81
    ## Plot
 
82
    xx = linspace (0, 2 * pi, 400);
 
83
    y1 = ppval (pp1, xx);
 
84
    y2 = ppval (pp2, xx);
 
85
    plot (x, y, ".", xx, [y1; y2]);
 
86
    axis tight;
 
87
    ylim ([-2 3]);
 
88
    legend ({"data", "no constraints", "periodic"});
 
89
    print ([nm "." typ], d_typ);
 
90
  elseif (strcmp (nm, "splinefit4"))
 
91
    ## Data (200 points)
 
92
    x = 2 * pi * rand (1, 200);
 
93
    y = sin (2 * x) + 0.1 * randn (size (x));
 
94
    ## Breaks
 
95
    breaks = linspace (0, 2 * pi, 10);
 
96
    ## Clamped endpoints, y = y" = 0
 
97
    xc = [0, 0, 2*pi, 2*pi];
 
98
    cc = [(eye (2)), (eye (2))];
 
99
    con = struct ("xc", xc, "cc", cc);
 
100
    pp1 = splinefit (x, y, breaks, "constraints", con);
 
101
    ## Hinged periodic endpoints, y = 0
 
102
    con = struct ("xc", 0);
 
103
    pp2 = splinefit (x, y, breaks, "constraints", con, "periodic", true);
 
104
    ## Plot
 
105
    xx = linspace (0, 2 * pi, 400);
 
106
    y1 = ppval (pp1, xx);
 
107
    y2 = ppval (pp2, xx);
 
108
    plot (x, y, ".", xx, [y1; y2]);
 
109
    axis tight;
 
110
    ylim ([-1.5 1.5]);
 
111
    legend({"data", "clamped", "hinged periodic"});
 
112
    print ([nm "." typ], d_typ);
 
113
  elseif (strcmp (nm, "splinefit5"))
 
114
    ## Truncated data
 
115
    x = [0,  1,  2,  4,  8, 16, 24, 40, 56, 72, 80] / 80;
 
116
    y = [0, 28, 39, 53, 70, 86, 90, 79, 55, 22,  2] / 1000;
 
117
    xy = [x; y];
 
118
    ## Curve length parameter
 
119
    ds = sqrt (diff (x).^2 + diff (y).^2);
 
120
    s = [0, cumsum(ds)];
 
121
    ## Constraints at s = 0: (x,y) = (0,0), (dx/ds,dy/ds) = (0,1)
 
122
    con = struct ("xc", [0 0], "yc", [0 0; 0 1], "cc", eye (2));
 
123
    ## Fit a spline with 4 pieces
 
124
    pp = splinefit (s, xy, 4, "constraints", con);
 
125
    ## Plot
 
126
    ss = linspace (0, s(end), 400);
 
127
    xyfit = ppval (pp, ss);
 
128
    xyb = ppval (pp, pp.breaks);
 
129
    plot (x, y, ".", xyfit(1,:), xyfit(2,:), "r", xyb(1,:), xyb(2,:), "ro");
 
130
    legend ({"data", "spline", "breaks"});
 
131
    axis tight;
 
132
    ylim ([0 0.1]);
 
133
    print ([nm "." typ], d_typ);
 
134
  elseif (strcmp (nm, "splinefit6"))
 
135
    ## Data
 
136
    x = linspace (0, 2*pi, 200);
 
137
    y = sin (x) + sin (2 * x) + 0.05 * randn (size (x));
 
138
    ## Add outliers
 
139
    x = [x, linspace(0,2*pi,60)];
 
140
    y = [y, -ones(1,60)];
 
141
    ## Fit splines with hinged conditions
 
142
    con = struct ("xc", [0, 2*pi]);
 
143
    pp1 = splinefit (x, y, 8, "constraints", con, "beta", 0.25); ## Robust fitting
 
144
    pp2 = splinefit (x, y, 8, "constraints", con, "beta", 0.75); ## Robust fitting
 
145
    pp3 = splinefit (x, y, 8, "constraints", con); ## No robust fitting
 
146
    ## Plot
 
147
    xx = linspace (0, 2*pi, 400);
 
148
    y1 = ppval (pp1, xx);
 
149
    y2 = ppval (pp2, xx);
 
150
    y3 = ppval (pp3, xx);
 
151
    plot (x, y, ".", xx, [y1; y2; y3]);
 
152
    legend ({"data with outliers","robust, beta = 0.25", ...
 
153
             "robust, beta = 0.75", "no robust fitting"});
 
154
    axis tight;
 
155
    ylim ([-2 2]);
 
156
    print ([nm "." typ], d_typ);
 
157
  endif
 
158
  hide_output ();  
 
159
endfunction
 
160
 
 
161
function set_print_size ()
 
162
  image_size = [5.0, 3.5]; # in inches, 16:9 format
 
163
  border = 0;              # For postscript use 50/72
 
164
  set (0, "defaultfigurepapertype", "<custom>");
 
165
  set (0, "defaultfigurepaperorientation", "landscape");
 
166
  set (0, "defaultfigurepapersize", image_size + 2*border);
 
167
  set (0, "defaultfigurepaperposition", [border, border, image_size]);
 
168
endfunction
 
169
 
 
170
## Use this function before plotting commands and after every call to
 
171
## print since print() resets output to stdout (unfortunately, gnpulot
 
172
## can't pop output as it can the terminal type).
 
173
function hide_output ()
 
174
  f = figure (1);
 
175
  set (f, "visible", "off");
 
176
endfunction
 
177
 
 
178
## generate something for the texinfo @image command to process
 
179
function image_as_txt(nm)
 
180
  fid = fopen (sprintf ("%s.txt", nm), "wt");
 
181
  fputs (fid, "\n");
 
182
  fputs (fid, "+---------------------------------+\n");
 
183
  fputs (fid, "| Image unavailable in text mode. |\n");
 
184
  fputs (fid, "+---------------------------------+\n");
 
185
  fclose (fid);
 
186
endfunction
 
187
 
 
188
 
 
189
%!demo
 
190
%! for s = 1:6
 
191
%!   splineimages (sprintf ("splinefit##d", s), "pdf")
 
192
%! endfor
 
193