1##  Drawing "spirograph" curves - epitrochoids, cycolids, roulettes
2##
3##  Copyright (C) 2007 Arjen Markus
4##  Copyright (C) 2008 Andrew Ross
5##
6##  This file is part of PLplot.
7##
8##  PLplot is free software; you can redistribute it and/or modify
9##  it under the terms of the GNU Library General Public License as published
10##  by the Free Software Foundation; either version 2 of the License, or
11##  (at your option) any later version.
12##
13##  PLplot is distributed in the hope that it will be useful,
14##  but WITHOUT ANY WARRANTY; without even the implied warranty of
15##  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16##  GNU Library General Public License for more details.
17##
18##  You should have received a copy of the GNU Library General Public License
19##  along with PLplot; if not, write to the Free Software
20##  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
21##
22
23
24
251;
26
27function ix27c
28
29##--------------------------------------------------------------------------*\
30## Generates two kinds of plots:
31##   - construction of a cycloid (animated)
32##   - series of epitrochoids and hypotrochoids
33##--------------------------------------------------------------------------*/
34
35  ## R, r, p, N
36  ## R and r should be integers to give correct termination of the
37  ## angle loop using gcd.
38  ## N.B. N is just a place holder since it is no longer used
39  ## (because we now have proper termination of the angle loop).
40  params = [
41	    21.0,  7.0,  7.0,  3.0;  ## Deltoid
42	    21.0,  7.0, 10.0,  3.0;
43	    21.0, -7.0, 10.0,  3.0;
44	    20.0,  3.0,  7.0, 20.0;
45	    20.0,  3.0, 10.0, 20.0;
46	    20.0, -3.0, 10.0, 20.0;
47	    20.0, 13.0,  7.0, 20.0;
48	    20.0, 13.0, 20.0, 20.0;
49	    20.0,-13.0, 20.0, 20.0];
50
51  ## Parse and process command line arguments
52
53  ## (void) plparseopts(&argc, argv, PL_PARSE_FULL);
54
55  ## Initialize plplot
56  plinit();
57
58  ## Illustrate the construction of a cycloid
59
60  cycloid();
61
62  ## Loop over the various curves
63  ## First an overview, then all curves one by one
64
65  plssub(3, 3); ## Three by three window
66
67  fill = 0;
68  for i = 1:9
69    pladv(0);
70    plvpor( 0.0, 1.0, 0.0, 1.0 );
71    spiro( params(i,:), fill );
72  endfor
73
74  pladv(0);
75  plssub(1, 1); ## One window per curve
76
77  for i=1:9
78    pladv(0);
79    plvpor( 0.0, 1.0, 0.0, 1.0 );
80    spiro( params(i,:), fill );
81  endfor
82
83  ## Fill the curves.
84  fill = 1;
85  pladv( 0 );
86  plssub( 1, 1 ); ## One window per curve
87  for i=1:9
88    pladv(0);
89    plvpor( 0.0, 1.0, 0.0, 1.0 );
90    spiro( params(i,:), fill );
91  endfor
92
93  ## Finally, an example to test out plarc capabilities
94  arcs();
95
96
97  ## Don't forget to call plend() to finish off!
98
99  plend1();
100end
101
102##------------------------------------------------------------------------
103## Calculate greatest common divisor following pseudo-code for the
104## Euclidian algorithm at http://en.wikipedia.org/wiki/Euclidean_algorithm
105
106function [value] = gcd (a,  b)
107    a = floor(abs(a));
108    b = floor(abs(b));
109    while(b!=0)
110        t = b;
111        b = mod(a,b);
112        a = t;
113    endwhile
114    value = a;
115end
116
117## ===============================================================
118
119function cycloid()
120    ## TODO
121endfunction
122
123## ===============================================================
124
125function spiro(params, fill)
126
127  NPNT=2000;
128
129  ## Fill the coordinates
130
131  ## Proper termination of the angle loop very near the beginning
132  ## point, see
133  ## http://mathforum.org/mathimages/index.php/Hypotrochoid.
134  windings = floor(abs(params(2))/gcd(params(1), params(2)));
135  steps    = floor(NPNT/windings);
136  dphi     = 2.0*pi/steps;
137
138  i = (0:windings*steps)';
139  phi = i*dphi;
140  phiw = (params(1)-params(2))/params(2)*phi;
141  xcoord = (params(1)-params(2))*cos(phi) + params(3)*cos(phiw);
142  ycoord = (params(1)-params(2))*sin(phi) - params(3)*sin(phiw);
143
144  xmin = min(xcoord);
145  xmax = max(xcoord);
146  ymin = min(ycoord);
147  ymax = max(ycoord);
148
149  xrange_adjust = 0.15 * (xmax - xmin);
150  xmin -= xrange_adjust;
151  xmax += xrange_adjust;
152  yrange_adjust = 0.15 * (ymax - ymin);
153  ymin -= yrange_adjust;
154  ymax += yrange_adjust;
155
156  plwind( xmin, xmax, ymin, ymax );
157
158  plcol0(1);
159  if ( fill )
160    plfill( xcoord, ycoord );
161  else
162    plline( xcoord, ycoord );
163  endif
164
165endfunction
166
167function arcs
168  NSEG = 8;
169
170  theta = 0.0;
171  dtheta = 360.0 / NSEG;
172  plenv( -10.0, 10.0, -10.0, 10.0, 1, 0 );
173
174  ## Plot segments of circle in different colors
175  for i = 0:NSEG-1
176    plcol0( mod(i,2) + 1 );
177    plarc(0.0, 0.0, 8.0, 8.0, theta, theta + dtheta, 0.0, 0);
178    theta = theta + dtheta;
179  endfor
180
181  ## Draw several filled ellipses inside the circle at different
182  ## angles.
183  a = 3.0;
184  b = a * tan( (dtheta/180.0*pi)/2.0 );
185  theta = dtheta/2.0;
186  for i = 0:NSEG-1
187    plcol0( 2 - mod(i,2) );
188    plarc( a*cos(theta/180.0*pi), a*sin(theta/180.0*pi), a, b, 0.0, 360.0, theta, 1);
189    theta = theta + dtheta;
190  endfor
191
192endfunction
193
194ix27c
195