Instructions for using the program.
The code uses only low level input/output routines so it can
run in several environments. But a unix like operating system
with built-in input/output redirects will make it easier to
use the program.
Save the source code in a file. A good file name would be
"find_center.c". The program needs to be linked with the
standard math library during compilation. Under a unix like
operating system the command is:
% cc find_center.c -lm -o find_center
While the program can read and write to a terminal window it is
convenient to place the data to be read into a file and use
input redirect. For example suppose the data is in a file
called "perturbC". The command to run the program with this
data would be:
% find_center < perturbC
The program will then proceed to look for the minimal variation
center. The input itself is simply a list of x and
y coordinates for the data points separated only by
white space characters. A sample data file looks like this:
42.23944217 -0.4917291616
-20.87269253 18.59472133
2.253700853 -20.23047981
8.656859039 11.37945103
-10.44990849 -1.977633283
6.291805073 -5.339319904
Each line has the x and y coordinates
for a data point. The program reads input until it reaches
the end of the file, it receives a cntrl-D character, or it
has read in 500 x and y coordinates.
Once it is done reading in data it proceeds to compute the
minimal variation center.
The ouput is displayed on the terminal window. The program
outputs the starting estimate for the minimal varaition center
(which is the center of mass of the data set). Once every
million iterations it outputs the most updated estimate for
the center. When it is finished it displays the final
estimate for the center. A typical run would look something
like this:
% find_center < perturbC
Starting estimate for center: (4.686534, 0.322502)
Estimate for center: (3.116342, -0.066155) after 1000000 iterations
Estimate for center: (2.174011, -0.328605) after 2000000 iterations
Estimate for center: (1.558815, -0.466284) after 3000000 iterations
Estimate for center: (1.131538, -0.526998) after 4000000 iterations
Estimate for center: (0.822636, -0.547612) after 5000000 iterations
Estimate for center: (0.593575, -0.549673) after 6000000 iterations
Estimate for center: (0.420937, -0.544380) after 7000000 iterations
...................
...................
...................
Estimate for center: (-0.161164, -0.504152) after 39000000 iterations
Estimate for center: (-0.161211, -0.504149) after 40000000 iterations
Estimate for center: (-0.161248, -0.504147) after 41000000 iterations
Estimate for center: (-0.161276, -0.504146) after 42000000 iterations
Estimate for center: (-0.161299, -0.504144) after 43000000 iterations
Estimate for center: (-0.161316, -0.504144) after 44000000 iterations
Estimate for center: (-0.161330, -0.504143) after 45000000 iterations
Final estimate for center: (-0.161339, -0.504142)
%
The variable "epsilon" affects the distance by which the
estimate for the center is incremented. The larger epsilon is
the further the estimate will be incremented each iteration.
A larger value of epsilon can allow the program to find the
center more quickly. However an excessively large value for
epsilon can cause the algorithm to fall into a specious
periodic orbit in which case the program will not exit on its
own. The variable "tolerance" is the minimal change in both
the x and y coordinates for the
estimated position of the center. So long as the change in
at least one of the coordinates is above the tolerance level
the program will continue to try improve the estimate. The
smaller tolerance is the more precise the estimate for the
center will be when the program is done. Of course a small
value for tolerance will cause the program to run for a longer
amount of time. But very precise values can generally be
obtained in a matter of minutes.
An explanation of how the algorithm implemented by this program
works is in:
Finding the Center of a Phyllotactic Pattern (Scott Hotton,
J. Theor. Biology , Vol. 225,
Issue 1 (2003) pp. 15-32) (pdf 341 K)
Data Files
Files for the data used in the paper are provided below. There
are ready for use by the program.
Very similar to the previous code (see instructions above). The
input is now a list of u , v , and w
coordinates for the data points separated only by white
space characters. The output is the set of four numbers
(φ0, ψ0,x0,y0)
that specify which line is the minimal variation axis of the
the data points.
(φ0 and ψ0 are expressed in degrees.)
A sample data set is provided. The data points lie on the
surface of a cone whose axis is the w-axis of the
coordinate system. The coordintes of the points are given by
the formulas
uj = (200-t/4)cos(δ t)
vj = (200-t/4)sin(δ t)
wj = t
where where t takes on integer values 0 to 99 and
δ=137.5 degrees. These coordinates are then rotated by
the matrices
( 1
0
0
)
( sin(φ0)
0
cos(φ0)
)
( 0
sin(ψ0)
cos(ψ0)
)
( 0
1
0
)
( 0
-cos(ψ0)
sin(ψ0)
)
( -cos(φ0)
0
sin(φ0)
)
Where φ0=π/2-0.2=78.540844 degrees and
ψ0=&pi/2+0.3 = 107.188734 degrees.
This makes the axis of the cone point in the direction
(cos(φ0), sin(φ0)
cos(ψ0), sin(φ0),
sin(ψ0)). Next all the points are translated by
the vector (-0.9,0.7,0). The reader can verify that the axis of
the cone now passes through the point 0.8409623589
ξ0-0.6687355420η0.
This data set can be downloaded by clicking here:
