% This matlab script generates an Open Inventor file
% with a bunch of small ellipsoids aligned along a spiral path.
% By changing the definitions of x, y, z and dx, dy, dz,
% any other path can be taken instead.

% Copyright (c)  Alexei Aganin
%       aaganin@math.fsu.edu
%       02/21/2003


clear all;

% number of beads per one loop of the spiral
n = 46;

% radius of the beads
radius = 0.14;

% elongation of the ellipsoids
ellps = 1.8;

% range of the parameter (angle)
t0 = 0;
tend = pi*7.3;

% discrete parameter
t = t0 : (tend-t0)/n : tend;
t = t(:);

% coordinates of the points on the spiral 
% (i.e centers of the beads)
x = cos(t);
y = sin(t);
z = t/7;

% tangent directions at these points
dx = -sin(t);
dy = cos(t);
dz = 1/7 * ones( size(t) );

% normalize the tangent vectors
for n = 1:length(t),
   dl = sqrt( dx(n)*dx(n) + dy(n)*dy(n) + dz(n)*dz(n) );
   dx(n) = dx(n) / dl;
   dy(n) = dy(n) / dl;
   dz(n) = dz(n) / dl;
end

% initially an ellipsoid is aligned along the y-axis
% let's call v0 the vector along the y-axis
v0 = [ 0 1 0 ]';

% open a file for writing
fid = fopen( 'beads.iv', 'w' );

% print some Open Inventor stuff
fprintf( fid, '#Inventor V2.1 ascii\n\n' );
fprintf( fid, 'Separator {\n\n' );

% compute and print the parameters for each bead
for n = 1:length(t),
   % for each bead we want to find a rotation that 
   % aligns v0 with the tangent at the center of the bead;
   % rotation is specified by the axis and the angle
   
    % compute angle
   a = acos( [dx(n) dy(n) dz(n)] * v0 );
   
   % compute axis
   r = cross( v0, [dx(n) dy(n) dz(n)] );
   
   % print the Open Inventor stuff
   
   fprintf( fid, 'Separator {\n' );
   
   % print the coordinates of the center [ x(n) y(n) z(n) ]
   fprintf( fid, ' Translation {\n' );
   fprintf( fid, '      translation %12.8f  %12.8f  %12.8f\n', [x(n) y(n) z(n)] );
   fprintf( fid, '  }\n' );
   
   % print the rotation parameters [ r a ] (four numbers)
   fprintf( fid, ' Rotation {\n' );
   fprintf( fid, '      rotation %12.8f  %12.8f  %12.8f  %12.8f\n', [r a] );
   fprintf( fid, '  }\n' );
   
   % scaling to make an ellipsoid from a sphere
   fprintf( fid, '  Scale {\n' );
   fprintf( fid, '      scaleFactor 1 %12.8f 1\n', ellps );
   fprintf( fid, '  }\n' );
   
   fprintf( fid, '  Material {\n' );
   % RGB color of the beads
   rr = 0.2; gg = 0.4; bb = 0.4;
   fprintf( fid, '      diffuseColor    %8.2f  %8.2f  %8.2f\n', [rr gg bb] );     
   fprintf( fid, '      specularColor   0.1 0.1 0.1\n' );
   fprintf( fid, '      shininess   1.0\n' );
   fprintf( fid, '  }\n' );
   
   fprintf( fid, '  Sphere {\n' );
   fprintf( fid, '      radius %12.8f\n', radius ); % radius of the beads
   fprintf( fid, '  }\n' );
   fprintf( fid, '}\n\n' );
            
 end

% print some Open Inventor stuff
fprintf( fid, '}\n\n' );

% close the file
fclose( fid );