I Hinnerup Net har vi, mest bare fordi vi kan, udviklet en WebGL shader der kører på computerens GPU (bedre kendt som grafikprocessor). Shaderen beregner konstant afstanden fra et vilkårligt punkt til en dynamisk skiftende kvadratisk Bézier-kurve.

Advarsel: Formler følger herefter i stride og voldsomme mængder. Hvis du vil direkte til demoen så hop til bunden af indlægget! 😉

Indenfor computer-graphic området er Bézier-kurver ikke til at komme udenom. Oprindeligt blev Bernsteins polynomie introduceret som en del af teoretisk matematik i 1800-tallet af Charles Hermite og Sergei Bernstein. Bézier-kurven, der er en udvidelse heraf, blev en realitet i 1900 tallet, da Pierre Bézier og Paul de Casteljau, fandt på at bruge del-summene i Bernsteins polynomie til at skalerere et array af vektorer (også kaldet kontrolpunkter).

Ligningen for en Bézier-kurve af n^te-grad er som følger:

P(t) er et punkt P, som ændrer sig som en funktion af t.
t er en floating point værdi, som bruges i intervallet 0.0 til 1.0.
n er en integer, som har værdien antal kontrol punkter minus 1.
i er en integer, som bevæger sig fra 0 til n.
P er et array af vektorer eller kontrol punkter, som indekseres af integeren i.
! er fakultet tegnet, som er produktet af en talrække af de positive hele tal fra 1 til og med tallet selv. Fakultet-funktionen n! kan dermed udtrykkes som:

En Bézier-kurve af n^te-grad har n+1 kontrolpunkter. I dette tilfælde ønsker vi at have 3 kontrolpunkter. Derfor skal ligningen blot løses for n = 2, hvilket giver følgende:

 
Herunder er koden til at beregne et punkt på Bézier-kurven.

vec2 getPositionOnBezierCurve(float t, vec2 p0, vec2 p1, vec2 p2) 
{
	float fOneMinusT = 1.0-t;
	vec2 pos = fOneMinusT*fOneMinusT*p0+2.0*fOneMinusT*t*p1+t*t*p2;
	return pos;
}

 Hvis det undlades at gange kontrolpunkterne på hver af del-summene i Bernsteins polynomial, kommer de 3 kurver af t til at se således ud:

Hvis de 3 kontrolpunkter bliver ganget på de 3 kurver af t, så bliver resultatet en kvadratisk Bézier-kurve.

At finde afstanden imellem et vilkårligt punkt Q og denne kvadratiske Bézier-kurve er overkommeligt hvis t er 0.0 eller 1.0, da P(t) jo her vil give henholdvis det første eller sidste iforvejen kendte kontrolpunkt. Men da shaderen skal virke med et vilkårligt punkt Q og en vilkårlig kvadratisk Bézier kurve skal der lidt mere formel-galop til. For at kunne beregne den korteste afstand imellem et vilkårligt punkt Q og et punkt på Bézier-kurven P(t) skal tangenten til P(t) anvendes. Tangenten beregnes ved at finde den afledte af funktionen P(t).

Ligningen for tangenten til P(t) kan omskrives til:

Dette kan reduceres yderligere til:

Den korteste afstand imellem et vilkårligt punkt Q og et punkt på Bézier-kurve P(t) vil være, der hvor skalarproduktet imellem vektoren Q-P(t) og tangenten til punktet P(t) er lig med 0.

 
. Skalarproduktet giver som bekendt cosinus til vinkelen imellem vektorene. Cosinus giver nul ved 90 grader, så med andre ord finder ligningen den værdi af t hvor vektoren Q-P(t) er ortogonal med tangenten til punktet P(t). Det er netop dette punkt, som vil have den korteste afstand til Bézier-kurven fra punktet Q.

Når ligningen bliver løst for t bliver det til en trediegradsligning:

 

Udtrykt i kode, ser overstående sådan ud:

vec2 dP0Q = p0-q;
vec2 A = p1-p0;
vec2 B = p0+p2-p1*2.0;
float a = dot(B,B);
float b = dot(A,B)*3.0;
float c = dot(A,A)*2.0+dot(dP0Q, B);
float d = dot(dP0Q,A);

Der findes forskellige måder til at løse en trediegradsligning på. Vi valgte at bruge Cardano’s metode.

Udtrykt i kode bliver Cardano’s metode til følgende:

#define EPSILON 0.000000001
#define PI 3.14159265358979
int findRoots(float a, float b, float c, float d, out float r[3]) {
	if (abs(a) > EPSILON) {
		float z = 1.0/a;
		float d3 = 1.0/3.0;
		float d27 = 1.0/27.0;
		a = b*z;
		b = c*z;
		c = d*z;
		float p = b-a*a*d3;
		float q = a*(2.0*a*a-9.0*b)*d27+c;
		float ppp = p*p*p;
		float D = q*q+4.0*ppp*d27;
		float delta = -a*d3;
		if (D > EPSILON) {
			z = sqrt(D);
			float u = (-q+z)*0.5;
			float v = (-q-z)*0.5;
			u = sign(u)*pow(abs(u),d3);
			v = sign(v)*pow(abs(v),d3);
			r[0] = u+v+delta;
			return 1;
		}
		else if (D < -EPSILON) {
			float u = sqrt(-p*d3)*2.0;
			float v = acos(-sqrt(-27.0/ppp)*q*0.5)*d3;
			r[0] = u*cos(v)+delta;
			r[1] = u*cos(v+2.0*PI*d3)+delta;
			r[2] = u*cos(v+4.0*PI*d3)+delta;
			return 3;
		}		
		else {
			q = sign(q)*pow(abs(q)*0.5,d3);
			r[0] = 2.0*-q+delta;
			r[1] = q+delta;
			return 2;
		}
	}
	else {
		if (abs(b) <= EPSILON && abs(c) > EPSILON) {
			r[0] = -d/c;
			return 1;
		}
		else {
			float D = c*c-4.0*b*d;
			float z = 1.0/(2.0*b);
			if (D > EPSILON) {
				D = sqrt(D);
				r[0] = (-c-D)*z;
				r[1] = (-c+D)*z;
				return 2;
			}
			else if (D > -EPSILON) {
				r[0] = -c*z;
				return 1;
			}
		}
	}
	return 0;
}

Løsningen af trediegrads ligningen vil give et sted imellem 0 og 3 løsninger. I tilfælde af mere end 1 løsning, så bruges kun den løsning hvor t er indenfor intervallet 0.0 til 1.0. Når den rigtige værdi af t således er fundet, skal denne sættes ind i den oprindelige ligning for en kvadratisk Bézier-kurve for at beregne punktet K, som vil være punktet på kurven der er tættest på Q. Til sidst beregnes selve afstanden g imellem K og Q med Pythagoras.

 

Koden til henholdvis at beregne K og sikre at t er i intervallet 0.0 til 1.0 og samt beregne afstanden imellem K og Q ser sådan ud:

float g = distance(q,getPositionOnBezierCurve(clamp(r[0],0.0,1.0),p0,p1,p2));
g = min(g, distance(q,getPositionOnBezierCurve(clamp(r[1],0.0,1.0),p0,p1,p2)));
g = min(g, distance(q,getPositionOnBezierCurve(clamp(r[2],0.0,1.0),p0,p1,p2)));

Afstandsværdien g bruges til at indeksere hele farvespektret og hermed oversat til en enkelt farve for hver pixel på skærmen. Koden til at gøre dette er:

float white = max(5.0-pow(d,2.0),0.0);
float red = max(-0.0007*pow(d-40.0,2.0)+1.0,0.0);
float yellow = max(-0.0007*pow(d-75.0,2.0)+1.0,0.0);
float green = max(-0.0007*pow(d-110.0,2.0)+1.0,0.0);
float cyan = max(-0.0007*pow(d-145.0,2.0)+1.0,0.0);
float blue = max(-0.0007*pow(d-180.0,2.0)+1.0,0.0);
vec4 c = vec4(white,white,white,0);
c += vec4(red,0,0,0);
c += vec4(yellow,yellow,0,0);
c += vec4(0,green,0,0);
c += vec4(0,cyan,cyan,0);
c += vec4(0,0,blue,0);

Hvis din browser understøtter WebGL kan du se en realtime version herunder.