#include <maps.h>

/*----------------------------------------------------------------------
 *
 * threads puts a thread pattern on the surface
 * of a "cylindrical workpiece" with a cylindrical projection.
 *
 * freq - "frequency" of the threads - number per unit length of cylinder
 * depth - depth the threads are cut into the surface of the cylinder
 * phase - rotates the threads about the central axis as this goes from
 *	0 to 1, to connect threads from different cylinders
 * zmin,zmax - extent of threads. No threads are cut outside here.
 * dampzone - distance inward from zmin,zmax to do a smoothstep diminishing
 * 	of depth
 * maporigin - the origin of the cylindrical projection
 * xaxis, yaxis, zaxis - 3 points defining the coordinate system of
 *        the cylinder.  The central axis is along zaxis-maporigin.
 *	  The point on the equator where (uu,vv)==(0,0) is mapped is at xaxis.
 *--------------------------------------------------------------------*/
displacement
threads ( float freq = 5.0, float depth = 0.1, float phase = 0.0,
	float zmin = -1e20, float zmax = 1e20, float dampzone = 0,
	point maporigin=point ("shader", 0,0,0),
	point xaxis=point ("shader", 1,0,0), point yaxis=point ("shader", 0,1,0),
	point zaxis=point ("shader", 0,0,1))
{
    point	PP,O,X,Y,Z;
    normal      sN;
    float	uu, vv, T,
			scale;


    /* do projection in shader space */
	PP=transform("shader",P);
	O=transform("shader",maporigin);
	X=transform("shader",xaxis);
	Y=transform("shader",yaxis);
	Z=transform("shader",zaxis);
	cylinmap(PP,uu,vv,O,X,Y,Z);

    if ((vv >= zmin) && (vv <= zmax)) {
    /* calculate new surface point and normal */
	if ((vv >= zmin+dampzone) && (vv <= zmax-dampzone)) {
		scale = depth;
	} else {
		scale = depth * min(smoothstep(zmin,zmin+dampzone,vv),
			1-smoothstep(zmax-dampzone,zmax,vv));
	}
	T = freq*vv + uu + phase;
	T = T - floor(T);
	scale = (2 * abs(T-0.5) - 1) * scale;
    /* do scaling in shader space */
	sN = normalize(transform("shader",N));
	P = transform("shader","common",PP + sN*scale);
	N = calculatenormal(P);
    }
}
