ACG: Ray Marching

This is a full example of doing ray marching from scratch in GLSL. Everything is done in a single fragment shader using all the advanced tricks and techniques we learned in class.

/* Constants defined at the top. */
const float NEAR_CLIPPING_PLANE  = 0.1;
const float FAR_CLIPPING_PLANE   = 10.0;
const int MARCHING_STEPS         = 100;
const float EPSILON              = 0.0001;
const float DISTANCE_BIAS        = 1.0;
const vec3 LIGHT_COL             = vec3( 2.5 );

float sdSphere( vec3 p, float r ) {
   return length(p) - r;
}
float sdBox( vec3 p, in vec3 b ) {
   vec3 q = abs(p) - b;
   return length(max(q,0.0)) + min(max(q.x,max(q.y,q.z)),0.0);
}
float sdTorus( vec3 p, vec2 t ) {
  vec2 q = vec2(length(p.xy)-t.x,p.z);
  return length(q)-t.y;
}
float sdPlane( vec3 p, vec4 n ) {
   return dot( p, n.xyz ) + n.w;
}
vec2 sdUnion( vec2 a, vec2 b ) {
   return (a.x < b.x) ? a : b;
}
float sdSmoothUnion( float a, float b, float k ) {
   float h = max( k-abs(a-b), 0.0 )/k;
   return min( a, b ) - h*h*k*(1.0/4.0);
}

/* The scene is defined here, and is a function of time (u_time uniform) */
vec2 scene( vec3 pos )
{
   float ds = sdSphere( vec3(pos) + vec3(0.0,-0.05,0.0), 0.25 );
   float dx = sdTorus( vec3(pos) + vec3(0.0,-0.05,0.0), vec2(0.2, 0.05));
   float db = sdBox(pos - vec3(0.0, -0.25,0.0), vec3(0.2,0.05,0.2));
   float dp = sdPlane( vec3(pos), vec4( 0.0, 1.0, 0.0, 0.25 ) );

   vec2 dtop = mix( vec2(ds,0.0), vec2(dx,1.0), sin(u_time)*0.5+0.5 );
   float dbase = sdSmoothUnion(db,dp,0.1);
   return sdUnion( vec2(dbase,2.0), dtop );
}

/* Computes the normals using a tetrahedron for efficiency. */
vec3 normal( vec3 pos )
{
   const float e = 0.0001;
   const vec2 k = vec2(1.0,-1.0);
   return normalize(
      k.xyy*scene( pos + k.xyy*e ).x +
      k.yyx*scene( pos + k.yyx*e ).x +
      k.yxy*scene( pos + k.yxy*e ).x +
      k.xxx*scene( pos + k.xxx*e ).x
   );
}

/* The heart of our raymarching approach. */
vec2 raymarch( vec3 ro, vec3 rd )
{
   vec2 t = vec2( NEAR_CLIPPING_PLANE, 0.0 );
   for (int i=0; i < MARCHING_STEPS; i++) {
      vec3 pos = ro + t.x*rd;

      vec2 h = scene(pos);
      if (h.x < EPSILON)
         break;
      t.x += h.x * DISTANCE_BIAS;
      t.y = h.y;
      if (t.x > FAR_CLIPPING_PLANE)
         break;
   }
   if (t.x > FAR_CLIPPING_PLANE)
      t = vec2( -1.0, 0.0 );
   return t;
}

/* We are approximating soft shadows for "free" with raymarching. */
float softshadow( vec3 ro, vec3 rd, float k )
{
   float res = 1.0;
   float ph = 1e20;
   float t = NEAR_CLIPPING_PLANE;
   for (int i=0; i < MARCHING_STEPS; i++) {
      float h = scene( ro+t*rd ).x;
      if (h < EPSILON)
         return 0.0;
      float y = h*h/(2.0*ph);
      float d = sqrt(h*h-y*y);
      res = min( res, k*d/max(0.0,t-y) );
      ph = h;
      t += h * DISTANCE_BIAS;
      if (t > FAR_CLIPPING_PLANE)
         break;
   }
   return res;
}

/* This is an analytically anti-aliased checkerboard pattern. */
float checkerboard(vec2 p)
{
   vec2 w = fwidth(p) + 0.001;
   vec2 i = 2.0*(abs(fract((p-0.5*w)*0.5)-0.5)-abs(fract((p+0.5*w)*0.5)-0.5))/w;
   return 0.5 - 0.5*i.x*i.y;
}

void main ()
{
   /* Light is a point light with no falloff that moves around based on time. */
   vec3 LIGHT_POS = 3.0 * vec3( sin(u_time*0.7), 1.0, cos(u_time*0.7) );

   /* Very simple orthogonal camera model. Rays go in the negative z-axis. */
   vec2 uv = (2.0*gl_FragCoord.xy -  u_resolution.xy) / u_resolution.y;
   vec3 ro = vec3(0.0,0.0,1.0);
   vec3 rd = normalize( vec3(uv, -1.5) );

   /* Raymarch the scene here. */
   vec2 t = raymarch( ro, rd );

   vec3 col; /* We determine the colour from the material. */
   if (t.x > 0.0) {
      /* Determine the colour mased on material (stored in t.y) */
      vec3 pos = ro + t.x * rd;
      vec3 nor = normal( pos );
      vec3 light_dir = normalize(LIGHT_POS - pos);
      vec3 mat;
      /* For the floor we use the checkerboard pattern. */
      if (t.y == 2.0)
         mat = vec3(checkerboard( pos.xz*20.0 ) * 0.2 + 0.7) * 0.2;
      /* Otherwise we interpolate between two different colours. */
      else
         mat = mix( vec3(0.2,0.0,0.0), vec3(0.0,0.2,0.0), t.y );

      /* Shadows computed with raymarching. */
      float sha = softshadow( pos + nor*0.001, light_dir, 8.0 );
      /* Lighting is simple based on a point light. */
      float sun = max( dot( nor, light_dir ), 0.0 );
      vec3 lin  = sun*LIGHT_COL * pow(vec3(sha),vec3(1.5));

      /* Ilumination and material are multiplied together. */
      col = lin * mat;
   }
   /* Ray did not hit anything so we draw "the sky" */
   else
      col = vec3(0.4,0.5,0.75) - 0.4*rd.y;

   /* Gamma correction. */
   col = pow( col, vec3(1.0/2.2) );
   colour_out = vec4(col,1.0);
}
149
 
1
/* Constants defined at the top. */
2
const float NEAR_CLIPPING_PLANE  = 0.1;
3
const float FAR_CLIPPING_PLANE   = 10.0;
4
const int MARCHING_STEPS         = 100;
5
const float EPSILON              = 0.0001;
6
const float DISTANCE_BIAS        = 1.0;
7
const vec3 LIGHT_COL             = vec3( 2.5 );
8
​
9
float sdSphere( vec3 p, float r ) {
10
   return length(p) - r;
11
}
12
float sdBox( vec3 p, in vec3 b ) {
13
   vec3 q = abs(p) - b;
14
   return length(max(q,0.0)) + min(max(q.x,max(q.y,q.z)),0.0);
15
}
16
float sdTorus( vec3 p, vec2 t ) {
17
  vec2 q = vec2(length(p.xy)-t.x,p.z);
18
  return length(q)-t.y;
19
}
20
float sdPlane( vec3 p, vec4 n ) {
21
   return dot( p, n.xyz ) + n.w;
22
}
23
vec2 sdUnion( vec2 a, vec2 b ) {
24
   return (a.x < b.x) ? a : b;
25
}
26
float sdSmoothUnion( float a, float b, float k ) {
27
   float h = max( k-abs(a-b), 0.0 )/k;
28
   return min( a, b ) - h*h*k*(1.0/4.0);
29
}
30
​
31
/* The scene is defined here, and is a function of time (u_time uniform) */
32
vec2 scene( vec3 pos )
33
{
34
   float ds = sdSphere( vec3(pos) + vec3(0.0,-0.05,0.0), 0.25 );
35
   float dx = sdTorus( vec3(pos) + vec3(0.0,-0.05,0.0), vec2(0.2, 0.05));
36
   float db = sdBox(pos - vec3(0.0, -0.25,0.0), vec3(0.2,0.05,0.2));
37
   float dp = sdPlane( vec3(pos), vec4( 0.0, 1.0, 0.0, 0.25 ) );
38
​
39
   vec2 dtop = mix( vec2(ds,0.0), vec2(dx,1.0), sin(u_time)*0.5+0.5 );
40
   float dbase = sdSmoothUnion(db,dp,0.1);
41
   return sdUnion( vec2(dbase,2.0), dtop );
42
}
43
​
44
/* Computes the normals using a tetrahedron for efficiency. */
45
vec3 normal( vec3 pos )
46
{
47
   const float e = 0.0001;
48
   const vec2 k = vec2(1.0,-1.0);
49
   return normalize(
50
      k.xyy*scene( pos + k.xyy*e ).x +
51
      k.yyx*scene( pos + k.yyx*e ).x +
52
      k.yxy*scene( pos + k.yxy*e ).x +
53
      k.xxx*scene( pos + k.xxx*e ).x
54
   );
55
}
56
​
57
/* The heart of our raymarching approach. */
58
vec2 raymarch( vec3 ro, vec3 rd )
59
{
60
   vec2 t = vec2( NEAR_CLIPPING_PLANE, 0.0 );
61
   for (int i=0; i < MARCHING_STEPS; i++) {
62
      vec3 pos = ro + t.x*rd;
63
​
64
      vec2 h = scene(pos);
65
      if (h.x < EPSILON)
66
         break;
67
      t.x += h.x * DISTANCE_BIAS;
68
      t.y = h.y;
69
      if (t.x > FAR_CLIPPING_PLANE)
70
         break;
71
   }
72
   if (t.x > FAR_CLIPPING_PLANE)
73
      t = vec2( -1.0, 0.0 );
74
   return t;
75
}
76
​
77
/* We are approximating soft shadows for "free" with raymarching. */
78
float softshadow( vec3 ro, vec3 rd, float k )
79
{
80
   float res = 1.0;
81
   float ph = 1e20;
82
   float t = NEAR_CLIPPING_PLANE;
83
   for (int i=0; i < MARCHING_STEPS; i++) {
84
      float h = scene( ro+t*rd ).x;
85
      if (h < EPSILON)
86
         return 0.0;
87
      float y = h*h/(2.0*ph);
88
      float d = sqrt(h*h-y*y);
89
      res = min( res, k*d/max(0.0,t-y) );
90
      ph = h;
91
      t += h * DISTANCE_BIAS;
92
      if (t > FAR_CLIPPING_PLANE)
93
         break;
94
   }
95
   return res;
96
}
97
​
98
/* This is an analytically anti-aliased checkerboard pattern. */
99
float checkerboard(vec2 p)
100
{
101
   vec2 w = fwidth(p) + 0.001;
102
   vec2 i = 2.0*(abs(fract((p-0.5*w)*0.5)-0.5)-abs(fract((p+0.5*w)*0.5)-0.5))/w;
103
   return 0.5 - 0.5*i.x*i.y;
104
}
105
​
106
void main ()
107
{
108
   /* Light is a point light with no falloff that moves around based on time. */
109
   vec3 LIGHT_POS = 3.0 * vec3( sin(u_time*0.7), 1.0, cos(u_time*0.7) );
110
​
111
   /* Very simple orthogonal camera model. Rays go in the negative z-axis. */
112
   vec2 uv = (2.0*gl_FragCoord.xy -  u_resolution.xy) / u_resolution.y;
113
   vec3 ro = vec3(0.0,0.0,1.0);
114
   vec3 rd = normalize( vec3(uv, -1.5) );
115
​
116
   /* Raymarch the scene here. */
117
   vec2 t = raymarch( ro, rd );
118
​
119
   vec3 col; /* We determine the colour from the material. */
120
   if (t.x > 0.0) {
121
      /* Determine the colour mased on material (stored in t.y) */
122
      vec3 pos = ro + t.x * rd;
123
      vec3 nor = normal( pos );
124
      vec3 light_dir = normalize(LIGHT_POS - pos);
125
      vec3 mat;
126
      /* For the floor we use the checkerboard pattern. */
127
      if (t.y == 2.0)
128
         mat = vec3(checkerboard( pos.xz*20.0 ) * 0.2 + 0.7) * 0.2;
129
      /* Otherwise we interpolate between two different colours. */
130
      else
131
         mat = mix( vec3(0.2,0.0,0.0), vec3(0.0,0.2,0.0), t.y );
132
​
133
      /* Shadows computed with raymarching. */
134
      float sha = softshadow( pos + nor*0.001, light_dir, 8.0 );
135
      /* Lighting is simple based on a point light. */
136
      float sun = max( dot( nor, light_dir ), 0.0 );
137
      vec3 lin  = sun*LIGHT_COL * pow(vec3(sha),vec3(1.5));
138
​
139
      /* Ilumination and material are multiplied together. */
140
      col = lin * mat;
141
   }
142
   /* Ray did not hit anything so we draw "the sky" */
143
   else
144
      col = vec3(0.4,0.5,0.75) - 0.4*rd.y;
145
​
146
   /* Gamma correction. */
147
   col = pow( col, vec3(1.0/2.2) );
148
   colour_out = vec4(col,1.0);
149
}