Skip to content

Instantly share code, notes, and snippets.

@mwalczyk
Created May 9, 2020 22:23
Show Gist options
  • Save mwalczyk/448da6b8c828b7dce7f8bc52e88bf9f9 to your computer and use it in GitHub Desktop.
Save mwalczyk/448da6b8c828b7dce7f8bc52e88bf9f9 to your computer and use it in GitHub Desktop.
layout (triangles) in;
layout (triangle_strip, max_vertices = 512) out;
in VS_OUT
{
vec3 normal;
vec3 view;
} gs_in[];
out GS_OUT
{
vec3 color;
} gs_out;
struct Edge
{
int a;
int b;
vec3 tangent_a;
vec3 tangent_b;
vec3 silhouette_point;
vec3 n;
};
vec3 triangle_normal()
{
vec3 a = vec3(gl_in[0].gl_Position) - vec3(gl_in[1].gl_Position);
vec3 b = vec3(gl_in[2].gl_Position) - vec3(gl_in[1].gl_Position);
return normalize(cross(a, b));
}
vec3 triangle_center()
{
return (vec3(gl_in[0].gl_Position) +
vec3(gl_in[1].gl_Position) +
vec3(gl_in[2].gl_Position)) / 3.0;
}
vec3 point(int index)
{
return vec3(gl_in[index].gl_Position);
}
vec3 normal(int index)
{
return gs_in[index].normal;
}
vec3 view(int index)
{
return gs_in[index].view;
}
vec3 tangent_at_point(in vec3 p, in vec3 n)
{
return vec3(0.0);
}
vec3 hermite(in vec3 v1, in vec3 v2, in vec3 t1, in vec3 t2, float u)
{
// https://www.cubic.org/docs/hermite.htm
return (2.0 * v1 - 2.0 * v2 + t1 + t2) * pow(u, 3.0) -
(3.0 * v1 - 3.0 * v2 + 2.0 * t1 + t2) * pow(u, 2.0) +
t1 * u +
v1;
}
vec3 evaluate_silhouette_bridge(in Edge edge, float u)
{
vec3 v1 = point(edge.a);
vec3 v2 = point(edge.b);
vec3 t1 = edge.tangent_a;
vec3 t2 = edge.tangent_b;
return hermite(v1, v2, t1, t2, u);
}
vec3 evaluate_silhouette_segment(in Edge edge_a, in Edge edge_b, float u)
{
// Interpolate across the silhouette segment connecting the two
// edges using each edge's silhouette point and normal vector (that is,
// the normal vector at the location of the silhouette point along the
// corresponding silhouette bridge)
vec3 v1 = edge_a.silhouette_point;
vec3 v2 = edge_b.silhouette_point;
// The two normal vectors at each of the silhouette points
vec3 n1 = edge_a.n;
vec3 n2 = edge_b.n;
// Calculate the tangent vectors at each of the silhouette points
vec3 diff = v2 - v1;
vec3 t1 = normalize(diff - dot(diff, n1) * n1);
vec3 t2 = normalize(diff - dot(diff, n2) * n2);
float cos_theta_a = dot(normalize(diff), t1);
float cos_theta_b = dot(normalize(diff), t2);
float len_a = (2.0 * length(diff)) / (1.0 + cos_theta_a);
float len_b = (2.0 * length(diff)) / (1.0 + cos_theta_b);
t1 *= len_a;
t2 *= len_b;
return hermite(v1, v2, t1, t2, u);
}
void calculate_tangents(inout Edge edge)
{
// First, grab the two endpoints of the edge
vec3 a = point(edge.a);
vec3 b = point(edge.b);
// Calculate the vector pointing towards `b` from `a`
vec3 diff = b - a;
// Grab the vertex normals at each endpoint of this edge
vec3 n_a = normal(edge.a);
vec3 n_b = normal(edge.b);
// Calculate the projection of `b - a` onto the tangent plane
// of each normal vector and its associated vertex
vec3 tangent_a = normalize(diff - dot(diff, n_a) * n_a);
vec3 tangent_b = normalize(diff - dot(diff, n_b) * n_b);
// Calculate the angles between `b - a` and each tangent vector
float cos_theta_a = dot(normalize(diff), tangent_a);
float cos_theta_b = dot(normalize(diff), tangent_b);
// Calculate the required length of each tanget vector (so as to
// approximate an arc with the silhouette bridge)
float len_a = (2.0 * length(diff)) / (1.0 + cos_theta_a);
float len_b = (2.0 * length(diff)) / (1.0 + cos_theta_b);
edge.tangent_a = tangent_a * len_a;
edge.tangent_b = tangent_b * len_b;
}
void calculate_silhouette_point(inout Edge edge)
{
// Now, solve the quadratic equation to find `u0`:
//
// d_tilde = (1 - u) * d1 + u * d2
// n_tilde = (1 - u) * n1 + u * n2
//
// Find the value of `u` such that `dot(d_tilde, n_tilde) = 0`
vec3 n1 = normalize(normal(edge.a)); // Vertex normal at `a`
vec3 n2 = normalize(normal(edge.b)); // Vertex normal at `b`
vec3 d1 = normalize(view(edge.a)); // View direction vector at `a`
vec3 d2 = normalize(view(edge.b)); // View direction vector at `b`
float coeffA = dot(d1, n1);
float coeffB = dot(d1, n2);
float coeffC = dot(d2, n1);
float coeffD = dot(d2, n2);
float a = (coeffA - coeffB - coeffC + coeffD);
float b = (coeffB + coeffC - 2 * coeffA);
float c = coeffA;
const float root = sqrt(b*b - 4.0 * a * c);
const float denom = 2.0 * a;
// The interpolant `u` should be between 0..1
float u0 = (-b + root) / denom;
if (u0 < 0.0 || u0 > 1.0)
{
u0 = (-b - root) / denom;
}
vec3 n_tilde = (1.0 - u0) * n1 + u0 * n2;
n_tilde = normalize(n_tilde);
vec3 silhouette_point = evaluate_silhouette_bridge(edge, u0);
edge.silhouette_point = silhouette_point;
edge.n = n_tilde;
}
Edge[2] calculate_silhouette_bridges(in bvec3 visibility, int index)
{
// Use the XOR operation to determine whether each pair of
// vertices had different visibility statuses (in which case,
// the containing edge will be part of the silhouette bridge)
// ...
// X -> Y - edge 0, 1
// Y -> Z - edge 1, 2
// Z -> X - edge 2, 0
bvec3 edge_parity = visibility.xyz ^^ visibility.yzx;
// The vertex at `index` will always be the "odd-one-out", i.e. the triangle
// vertex with differing visibility status
Edge edges[2] =
{
Edge(index, (index + 1) % 3, vec3(0.0), vec3(0.0), vec3(0.0), vec3(0.0)),
Edge(index, (index + 2) % 3, vec3(0.0), vec3(0.0), vec3(0.0), vec3(0.0))
};
for (int i = 0; i < edges.length(); ++i)
{
calculate_tangents(edges[i]);
calculate_silhouette_point(edges[i]);
}
return edges;
}
void main()
{
// Steps:
//
// 1) Calculate triangle visibility: if one or more of the triangle's vertices
// differ in "sign", continue - otherwise, break
// 2) Find the two edges of the triangle whose vertices differ in "sign" - by
// definition, there will be exactly 2 such edges
// 3) For each edge, calculate the silhouette bridge:
// a. Find the tangent vectors `t1` and `t2` at either end of the edge
// b. Find the lengths of the two tangent vectors
// c. Use the tangent vectors to calculate points along the resulting
// silhouette bridge
// 4) Calculate the silhouette points (i.e. the point on each of the silhouette
// bridges whose surface normal vector is perpendicular to the viewing direction) -
// this is done by solving a simple quadratic equation
// 5) Calculate the silhouette segment, which is the Hermite interpolation between
// the two silhouette points
// 6) Determine the number of sample points to generate along the silhouette segment
// based on that segment's projected length
// 7) Remesh by adding new triangles
// If `dot(n, v)` is greater than (or equal to) 0, then the vertex is said to be "visible" -
// note that this is different from occlusion
bvec3 visibility =
{
dot(gs_in[0].normal, gs_in[0].view) >= 0.0,
dot(gs_in[1].normal, gs_in[1].view) >= 0.0,
dot(gs_in[2].normal, gs_in[2].view) >= 0.0
};
// Silhouette refinement and local remeshing only needs to be done
// for triangles whose vertices do not all have the same visibility
// status
if (all(not(visibility)) || all(visibility))
{
// Pass-through geometry
for(int i = 0; i < gl_in.length; ++i)
{
vec3 position = point(i);
gl_Position = TDWorldToProj(vec4(position, 1.0), 0);
// Simple RGB colors for now
gs_out.color = (position / 2.0) * 0.5 + 0.5;
EmitVertex();
}
EndPrimitive();
}
else
{
// First, find the index of the triangle vertex whose visibility
// status differs from the other two vertices that make up this
// triangle - it *should* only be one such vertex, and its index
// will be 0, 1, or 2
int index = 0;
bvec3 test = visibility;
test.x = !test.x;
if (all(not(test)) || all(test))
{
index = 0;
}
test = visibility;
test.y = !test.y;
if (all(not(test)) || all(test))
{
index = 1;
}
test = visibility;
test.z = !test.z;
if (all(not(test)) || all(test))
{
index = 2;
}
Edge edges[2] = calculate_silhouette_bridges(visibility, index);
const int number_of_samples = 8;
const float step_size = 1.0 / float(number_of_samples);
for (uint i = 0; i < number_of_samples; ++i)
{
// The first point will be the point at `index` (calculate above)
// whose status differs from the other 2 triangle vertices
vec3 position = point(index);
gl_Position = TDWorldToProj(vec4(position, 1.0), 0);
// Simple RGB colors for now
gs_out.color = (position / 2.0) * 0.5 + 0.5;
EmitVertex();
// The next 2 points will be neighbors along the silhouette segment
vec3 segment_sample;
float t;
t = step_size * (i + 0);
segment_sample = evaluate_silhouette_segment(edges[0], edges[1], t);
gl_Position = TDWorldToProj(vec4(segment_sample, 1.0), 0);
gs_out.color = (segment_sample / 2.0) * 0.5 + 0.5;
EmitVertex();
t = step_size * (i + 1);
segment_sample = evaluate_silhouette_segment(edges[0], edges[1], t);
gl_Position = TDWorldToProj(vec4(segment_sample, 1.0), 0);
gs_out.color = (segment_sample / 2.0) * 0.5 + 0.5;
EmitVertex();
// Emit the triangle and reset
EndPrimitive();
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment