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