function ans = moebius(w,k,c,n);
% moebius(w,k,c,n) - Called by moeb.m
%                    Creates the moebius strip.
%   w is the strip width, try w = 1/3, then w = 1.
%   k is the number of half twists, usually integer.  Default is 1.
%   c is the curvature, vary from 0. to 1.  Default is 1.
%   n determines number of grid points, try n = 30

%   Copyright (c) 1984-98 by The MathWorks, Inc. All Rights Reserved
%   $Revision: 1.3 $
%   $Date: 1997/12/03 17:29:58 $

if nargin < 4, n = 64; end
if nargin < 3, c = 1; end
if nargin < 2, k = 1; end
if nargin < 1, w = 1/3; end

if c==0, c=eps; end
m = round(w*n/pi);
if m==0, m=1; end
[s,t] = meshgrid(-1:1/m:1, -1:1/n:1);
r = (1-w) + w*s.*sin(k/2*pi*t);
x = r.*sin(c*pi*t)/c;
y = r.*(1 - (1-cos(c*pi*t))/c);
z = w*s.*cos(k/2*pi*t);
q =[min(min(x)) max(max(x))
    min(min(y)) max(max(y))
    min(min(z)) max(max(z))];
Q = max(max(abs(q)));
if Q > 1, x = x/Q; y = y/Q; z = z/Q; end
surf(x,y,z,t)
colormap(hsv)
axis equal
axis off
