polytope.glsl

precision highp float;

uniform vec2 resolution;

uniform mat4 viewMatrix;
uniform vec3 cameraPosition;

uniform mat4 cameraWorldMatrix;
uniform mat4 cameraProjectionMatrixInverse;


// Original DE from Knighty's Fragmentarium frag: http://www.fractalforums.com/fragmentarium/solids-many-many-solids/


#define PI 3.141592


// Return rotation matrix for rotating around vector v by angle
mat3  rotationMatrix3(vec3 v, float angle)
{
	float c = cos(radians(angle));
	float s = sin(radians(angle));
	
	return mat3(c + (1.0 - c) * v.x * v.x, (1.0 - c) * v.x * v.y - s * v.z, (1.0 - c) * v.x * v.z + s * v.y,
		(1.0 - c) * v.x * v.y + s * v.z, c + (1.0 - c) * v.y * v.y, (1.0 - c) * v.y * v.z - s * v.x,
		(1.0 - c) * v.x * v.z - s * v.y, (1.0 - c) * v.y * v.z + s * v.x, c + (1.0 - c) * v.z * v.z
		);
}

vec2 rotate(vec2 v, float a) {
	return vec2(cos(a)*v.x + sin(a)*v.y, -sin(a)*v.x + cos(a)*v.y);
}

// control-group: coordinate
uniform int Degree; // control[4, 3-5]
uniform float U; // control[0, 0-1]
uniform float V; // control[0, 0-1]
uniform float W; // control[0, 0-1]
uniform float T; // control[1, 0-1]

// control-group: rotation
uniform float RotationX; // control[1, 0-1]
uniform float RotationY; // control[0, 0-1]
uniform float RotationZ; // control[0, 0-1]
uniform float Angle; // control[0, 0-360]

// control-group: style
uniform float VRadius; // control[0.04, 0-0.2]
uniform float SRadius; // control[0.03, 0-0.2]
 
uniform bool displaySegments; // control[true]
uniform bool displayVertices; // control[true]




mat3 rot;
vec4 nc,nd,p;
void init() {

	float cospin=cos(PI/float(Degree)), isinpin=1./sin(PI/float(Degree));
	float scospin=sqrt(2./3.-cospin*cospin), issinpin=1./sqrt(3.-4.*cospin*cospin);

	nc=0.5*vec4(0,-1,sqrt(3.),0.);
	nd=vec4(-cospin,-0.5,-0.5/sqrt(3.),scospin);

	vec4 pabc,pbdc,pcda,pdba;
	pabc=vec4(0.,0.,0.,1.);
	pbdc=0.5*sqrt(3.)*vec4(scospin,0.,0.,cospin);
	pcda=isinpin*vec4(0.,0.5*sqrt(3.)*scospin,0.5*scospin,1./sqrt(3.));
	pdba=issinpin*vec4(0.,0.,2.*scospin,1./sqrt(3.));
	
	p=normalize(V*pabc+U*pbdc+W*pcda+T*pdba);

	rot = rotationMatrix3(normalize(vec3(RotationX,RotationY,RotationZ)), Angle);//in reality we need a 4D rotation
}

vec4 fold(vec4 pos) {
	for(int i=0;i<25;i++){
		vec4 tmp = pos;
		pos.xy=abs(pos.xy);
		float t=-2.*min(0.,dot(pos,nc)); pos+=t*nc;
		t=-2.*min(0.,dot(pos,nd)); pos+=t*nd;
		if (tmp == pos) { return pos; }
	}
	return pos;
}

float DD(float ca, float sa, float r){
	//magic formula to convert from spherical distance to planar distance.
	//involves transforming from 3-plane to 3-sphere, getting the distance
	//on the sphere then going back to the 3-plane.
	return r-(2.*r*ca-(1.-r*r)*sa)/((1.-r*r)*ca+2.*r*sa+1.+r*r);
}

float dist2Vertex(vec4 z, float r){
	float ca=dot(z,p), sa=0.5*length(p-z)*length(p+z);//sqrt(1.-ca*ca);//
	return DD(ca,sa,r)-VRadius;
}

float dist2Segment(vec4 z, vec4 n, float r){
	//pmin is the orthogonal projection of z onto the plane defined by p and n
	//then pmin is projected onto the unit sphere
	float zn=dot(z,n),zp=dot(z,p),np=dot(n,p);
	float alpha=zp-zn*np, beta=zn-zp*np;
	vec4 pmin=normalize(alpha*p+min(0.,beta)*n);
	//ca and sa are the cosine and sine of the angle between z and pmin. This is the spherical distance.
	float ca=dot(z,pmin), sa=0.5*length(pmin-z)*length(pmin+z);//sqrt(1.-ca*ca);//
	float factor = 1.0;// DD(ca,sa,r)/DD(ca+0.01,sa,r);
	return (DD(ca,sa,r)-SRadius*factor)*min(1.0/factor,1.0);
	
}
//it is possible to compute the distance to a face just as for segments: pmin will be the orthogonal projection
// of z onto the 3-plane defined by p and two n's (na and nb, na and nc, na and and, nb and nd... and so on).
//that involves solving a system of 3 linear equations.
//it's not implemented here because it is better with transparency

float dist2Segments(vec4 z, float r){
	float da=dist2Segment(z, vec4(1.,0.,0.,0.), r);
	float db=dist2Segment(z, vec4(0.,1.,0.,0.), r);
	float dc=dist2Segment(z, nc, r);
	float dd=dist2Segment(z, nd, r);
	
	return min(min(da,db),min(dc,dd));
}

float DE(vec3 pos) {
	
	float r=length(pos);
	vec4 z4=vec4(2.*pos,1.-r*r)*1./(1.+r*r);//Inverse stereographic projection of pos: z4 lies onto the unit 3-sphere centered at 0.
	z4.xyw=rot*z4.xyw;
	z4=fold(z4);//fold it
	float d=10000.;
	if(displayVertices ) d=min(d,dist2Vertex(z4,r));
	if(displaySegments) d=min(d,dist2Segments(z4, r));
	return d ;
}

vec4 baseColor(vec3 pos, vec3 normal){
	return vec4(1.0);
}