rayfract

view sdr/julia.p.glsl @ 1:03022062c464

lalal
author John Tsiombikas <nuclear@siggraph.org>
date Tue, 26 Oct 2010 08:45:37 +0300
parents 09bb67c000bc
children dfe7c98cf373
line source
1 /* vim: set ft=glsl:ts=4:sw=4 */
2 uniform vec4 seed;
3 uniform sampler2D ray_tex;
4 uniform float err_thres;
5 uniform int iter;
6 uniform float reflectivity;
7 uniform vec3 diffuse_color;
9 #define quat(s, x, y, z) vec4(x, y, z, s)
10 #define quat_identity() vec4(0.0, 0.0, 0.0, 1.0)
12 #define vec2quat(v) (v).wxyz
14 struct Ray {
15 vec3 origin;
16 vec3 dir;
17 };
19 struct Julia {
20 bool inside;
21 vec4 q;
22 vec4 qprime;
23 };
25 struct Material {
26 vec3 kd, ks;
27 float spow;
28 float kr;
29 };
31 struct ISect {
32 bool hit;
33 float t;
34 vec3 pos;
35 vec3 normal;
36 Material mat;
37 };
39 ISect find_intersection(Ray ray);
40 vec3 shade(Ray ray, ISect isect);
41 float amboc(ISect isect);
42 vec3 sky(Ray ray);
43 Julia julia(vec4 q, vec4 c);
44 float julia_dist(vec4 z);
45 vec3 julia_grad(vec4 z);
46 vec4 quat_mul(vec4 q1, vec4 q2);
47 vec4 quat_sq(vec4 q);
48 float quat_length_sq(vec4 q);
49 ISect ray_julia(Ray ray);
50 ISect ray_sphere(Ray ray, float rad);
51 ISect ray_floor(Ray ray);
52 Ray get_primary_ray();
54 vec3 steps_color(int steps);
56 void main()
57 {
58 Ray ray = get_primary_ray();
60 float energy = 1.0;
61 vec3 color = vec3(0.0, 0.0, 0.0);
63 while(energy > 0.001) {
64 ISect res = find_intersection(ray);
66 if(res.hit) {
67 color += shade(ray, res) * energy;
68 energy *= res.mat.kr;
70 ray.origin = res.pos;
71 ray.dir = reflect(ray.dir, res.normal);
72 } else {
73 color += sky(ray) * energy;
74 break;
75 }
76 }
78 gl_FragColor = vec4(color, 1.0);
79 }
82 ISect find_intersection(Ray ray)
83 {
84 ISect res;
85 res.hit = false;
87 ISect bhit = ray_sphere(ray, 2.0);
88 if(bhit.hit) {
89 ray.origin = bhit.pos;
90 res = ray_julia(ray);
91 }
93 if(!res.hit) {
94 res = ray_floor(ray);
95 }
96 return res;
97 }
99 vec3 shade(Ray ray, ISect isect)
100 {
101 vec3 ldir = normalize(vec3(10.0, 10.0, -10.0) - isect.pos);
102 vec3 vdir = -ray.dir;
103 vec3 hdir = normalize(ldir + vdir);
105 float ndotl = dot(ldir, isect.normal);
106 float ndoth = dot(hdir, isect.normal);
108 vec3 dcol = isect.mat.kd;// * abs(ndotl);
109 vec3 scol = isect.mat.ks * pow(abs(ndoth), isect.mat.spow);
111 return vec3(0.05, 0.05, 0.05) + dcol + scol;
112 }
114 #define AO_STEP 0.04
115 #define AO_MAGIC 8.0
116 float amboc(ISect isect)
117 {
118 float sum = 0.0;
120 for(float fi=0.0; fi<5.0; fi+=1.0) {
121 float sample_dist = fi * AO_STEP;
122 vec3 pt = isect.pos + isect.normal * sample_dist;
123 float jdist = julia_dist(quat(pt.x, pt.y, pt.z, 0.0));
125 sum += 1.0 / pow(2.0, fi) * (sample_dist - jdist);
126 }
128 float res = 1.0 - AO_MAGIC * sum;
129 return clamp(res, 0.0, 1.0);
130 }
132 vec3 sky(Ray ray)
133 {
134 vec3 col1 = vec3(0.75, 0.78, 0.8);
135 vec3 col2 = vec3(0.56, 0.7, 1.0);
137 float t = max(ray.dir.y, -0.5);
138 return mix(col1, col2, t);
139 }
141 Julia julia(vec4 q, vec4 c)
142 {
143 Julia res;
144 res.inside = true;
146 res.q = q;
147 res.qprime = quat_identity();
149 for(int i=0; i<iter; i++) {
150 res.qprime = 2.0 * quat_mul(res.q, res.qprime);
151 res.q = quat_sq(res.q) + c;
153 if(dot(res.q, res.q) > 8.0) {
154 res.inside = false;
155 break;
156 }
157 }
158 return res;
159 }
161 float julia_dist(vec4 z)
162 {
163 Julia jres = julia(z, seed);
165 float lenq = length(jres.q);
166 float lenqprime = length(jres.qprime);
168 return 0.5 * lenq * log(lenq) / lenqprime;
169 }
171 #define OFFS 1e-4
172 vec3 julia_grad(vec4 z)
173 {
174 vec3 grad;
175 grad.x = julia_dist(z + quat(OFFS, 0.0, 0.0, 0.0)) - julia_dist(z - quat(OFFS, 0.0, 0.0, 0.0));
176 grad.y = julia_dist(z + quat(0.0, OFFS, 0.0, 0.0)) - julia_dist(z - quat(0.0, OFFS, 0.0, 0.0));
177 grad.z = julia_dist(z + quat(0.0, 0.0, OFFS, 0.0)) - julia_dist(z - quat(0.0, 0.0, OFFS, 0.0));
178 return grad;
179 }
181 vec4 quat_mul(vec4 q1, vec4 q2)
182 {
183 vec4 res;
184 res.w = q1.w * q2.w - dot(q1.xyz, q2.xyz);
185 res.xyz = q1.w * q2.xyz + q2.w * q1.xyz + cross(q1.xyz, q2.xyz);
186 return res;
187 }
189 vec4 quat_sq(vec4 q)
190 {
191 vec4 res;
192 res.w = q.w * q.w - dot(q.xyz, q.xyz);
193 res.xyz = 2.0 * q.w * q.xyz;
194 return res;
195 }
197 #define MIN_STEP 0.001
198 ISect ray_julia(Ray inray)
199 {
200 float dist_acc = 0.0;
201 Ray ray = inray;
202 ISect res;
204 int i = 0;
205 for(float fi=0.0; ; fi+=0.1) {
206 i++;
207 vec4 q = quat(ray.origin.x, ray.origin.y, ray.origin.z, 0.0);
209 float dist = max(julia_dist(q), MIN_STEP);
211 ray.origin += ray.dir * dist;
212 dist_acc += dist;
214 if(dist < err_thres) {
215 res.hit = true;
216 res.t = dist_acc;
217 res.pos = ray.origin;
218 res.normal = normalize(julia_grad(quat(res.pos.x, res.pos.y, res.pos.z, 0.0)));
219 res.mat.kr = reflectivity;
220 res.mat.kd = diffuse_color * amboc(res) * (1.0 - res.mat.kr);
221 res.mat.ks = vec3(0.4, 0.4, 0.4);//vec3(res.mat.kr, res.mat.kr, res.mat.kr);
222 res.mat.spow = 50.0;
223 break;
224 }
226 if(dot(ray.origin, ray.origin) > 100.0) {
227 res.hit = false;
228 break;
229 }
230 }
232 return res;
233 }
235 ISect ray_sphere(Ray ray, float rad)
236 {
237 ISect res;
238 res.hit = false;
240 float a = dot(ray.dir, ray.dir);
241 float b = 2.0 * dot(ray.dir, ray.origin);
242 float c = dot(ray.origin, ray.origin) - rad * rad;
244 float d = b * b - 4.0 * a * c;
245 if(d < 0.0) return res;
247 float sqrt_d = sqrt(d);
248 float t1 = (-b + sqrt_d) / (2.0 * a);
249 float t2 = (-b - sqrt_d) / (2.0 * a);
251 if((t1 >= 0.0 || t2 >= 0.0)) {
252 if(t1 < 0.0) t1 = t2;
253 if(t2 < 0.0) t2 = t1;
255 res.hit = true;
256 res.t = min(t1, t2);
257 res.pos = ray.origin + ray.dir * res.t;
258 //res.mat.kd = vec3(1.0, 0.3, 0.2);
259 //res.normal = res.pos / rad;
260 }
262 return res;
263 }
265 #define FLOOR_HEIGHT (-2.0)
267 ISect ray_floor(Ray ray)
268 {
269 ISect res;
270 res.hit = false;
272 if(ray.origin.y < FLOOR_HEIGHT || ray.dir.y >= 0.0) {
273 return res;
274 }
276 res.normal = vec3(0.0, 1.0, 0.0);
277 float ndotdir = dot(res.normal, ray.dir);
279 float t = (FLOOR_HEIGHT - ray.origin.y) / ndotdir;
280 res.pos = ray.origin + ray.dir * t;
282 if(abs(res.pos.x) > 8.0 || abs(res.pos.z) > 8.0) {
283 res.hit = false;
284 } else {
285 res.hit = true;
287 float chess = mod(floor(res.pos.x) + floor(res.pos.z), 2.0);
288 res.mat.kd = mix(vec3(0.498, 0.165, 0.149), vec3(0.776, 0.851, 0.847), chess);
289 res.mat.ks = vec3(0.0, 0.0, 0.0);
290 res.mat.spow = 1.0;
291 res.mat.kr = 0.0;
292 }
293 return res;
294 }
296 Ray get_primary_ray()
297 {
298 Ray ray;
299 vec2 tc = gl_TexCoord[0].xy;
300 ray.dir = gl_NormalMatrix * normalize(texture2D(ray_tex, tc).xyz);
301 ray.origin = (gl_ModelViewMatrix * vec4(0.0, 0.0, 0.0, 1.0)).xyz;
302 return ray;
303 }
306 vec3 steps_color(int steps)
307 {
308 if(steps <= 1) {
309 return vec3(0.0, 0.5, 0.0);
310 } else if(steps == 2) {
311 return vec3(0.0, 1.0, 0.0);
312 } else if(steps == 3) {
313 return vec3(0.0, 0.0, 0.5);
314 } else if(steps == 4) {
315 return vec3(0.0, 0.0, 1.0);
316 } else if(steps == 5) {
317 return vec3(0.0, 0.5, 0.5);
318 } else if(steps == 6) {
319 return vec3(0.0, 1.0, 1.0);
320 } else if(steps == 7) {
321 return vec3(0.5, 0.0, 0.5);
322 } else if(steps == 8) {
323 return vec3(1.0, 0.0, 1.0);
324 } else if(steps == 9) {
325 return vec3(0.5, 0.0, 0.0);
326 }
327 return vec3(0.5 + float(steps - 9) / 10.0, 0.0, 0.0);
328 }