Rework of volume shading
[blender.git] / source / blender / render / intern / source / volumetric.c
1 /**
2  *
3  * ***** BEGIN GPL LICENSE BLOCK *****
4  *
5  * This program is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU General Public License
7  * as published by the Free Software Foundation; either version 2
8  * of the License, or (at your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program; if not, write to the Free Software Foundation,
17  * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
18  *
19  * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
20  * All rights reserved.
21  *
22  * The Original Code is: all of this file.
23  *
24  * Contributor(s): Matt Ebb, Raul Fernandez Hernandez (Farsthary)
25  *
26  * ***** END GPL LICENSE BLOCK *****
27  */
28
29 #include <math.h>
30 #include <stdlib.h>
31 #include <string.h>
32 #include <float.h>
33
34 #include "MEM_guardedalloc.h"
35
36 #include "BLI_blenlib.h"
37 #include "BLI_arithb.h"
38 #include "BLI_rand.h"
39 #include "BLI_voxel.h"
40
41 #include "RE_shader_ext.h"
42 #include "RE_raytrace.h"
43
44 #include "DNA_material_types.h"
45 #include "DNA_group_types.h"
46 #include "DNA_lamp_types.h"
47 #include "DNA_meta_types.h"
48
49 #include "BKE_global.h"
50
51 #include "render_types.h"
52 #include "pixelshading.h"
53 #include "shading.h"
54 #include "texture.h"
55 #include "volumetric.h"
56 #include "volume_precache.h"
57
58 #if defined( _MSC_VER ) && !defined( __cplusplus )
59 # define inline __inline
60 #endif // defined( _MSC_VER ) && !defined( __cplusplus )
61
62 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
63 /* defined in pipeline.c, is hardcopy of active dynamic allocated Render */
64 /* only to be used here in this file, it's for speed */
65 extern struct Render R;
66 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
67
68 /* luminance rec. 709 */
69 inline float luminance(float* col)
70 {
71         return (0.212671f*col[0] + 0.71516f*col[1] + 0.072169f*col[2]);
72 }
73
74 /* tracing */
75
76 static int vol_get_bounds(ShadeInput *shi, float *co, float *vec, float *hitco, Isect *isect, int intersect_type)
77 {
78         float maxsize = RE_ray_tree_max_size(R.raytree);
79         
80         /* XXX TODO - get raytrace max distance from object instance's bounding box */
81         /* need to account for scaling only, but keep coords in camera space...
82          * below code is WIP and doesn't work!
83         VecSubf(bb_dim, shi->obi->obr->boundbox[1], shi->obi->obr->boundbox[2]);
84         Mat3MulVecfl(shi->obi->nmat, bb_dim);
85         maxsize = VecLength(bb_dim);
86         */
87         
88         VECCOPY(isect->start, co);
89         isect->end[0] = co[0] + vec[0] * maxsize;
90         isect->end[1] = co[1] + vec[1] * maxsize;
91         isect->end[2] = co[2] + vec[2] * maxsize;
92         
93         isect->mode= RE_RAY_MIRROR;
94         isect->oborig= RAY_OBJECT_SET(&R, shi->obi);
95         isect->face_last= NULL;
96         isect->ob_last= 0;
97         isect->lay= -1;
98         
99         if (intersect_type == VOL_BOUNDS_DEPTH) isect->faceorig= (RayFace*)shi->vlr;
100         else if (intersect_type == VOL_BOUNDS_SS) isect->faceorig= NULL;
101         
102         if(RE_ray_tree_intersect(R.raytree, isect))
103         {
104                 hitco[0] = isect->start[0] + isect->labda*isect->vec[0];
105                 hitco[1] = isect->start[1] + isect->labda*isect->vec[1];
106                 hitco[2] = isect->start[2] + isect->labda*isect->vec[2];
107                 return 1;
108         } else {
109                 return 0;
110         }
111 }
112
113 static void shade_intersection(ShadeInput *shi, float *col, Isect *is)
114 {
115         ShadeInput shi_new;
116         ShadeResult shr_new;
117         
118         memset(&shi_new, 0, sizeof(ShadeInput)); 
119         
120         shi_new.mask= shi->mask;
121         shi_new.osatex= shi->osatex;
122         shi_new.thread= shi->thread;
123         shi_new.depth = shi->depth + 1;
124         shi_new.volume_depth= shi->volume_depth + 1;
125         shi_new.xs= shi->xs;
126         shi_new.ys= shi->ys;
127         shi_new.lay= shi->lay;
128         shi_new.passflag= SCE_PASS_COMBINED; /* result of tracing needs no pass info */
129         shi_new.combinedflag= 0xFFFFFF;          /* ray trace does all options */
130         shi_new.light_override= shi->light_override;
131         shi_new.mat_override= shi->mat_override;
132         
133         VECCOPY(shi_new.camera_co, is->start);
134         
135         memset(&shr_new, 0, sizeof(ShadeResult));
136         
137         /* hardcoded limit of 100 for now - prevents problems in weird geometry */
138         if (shi->volume_depth < 100) {
139                 shade_ray(is, &shi_new, &shr_new);
140         }
141         
142         VecCopyf(col, shr_new.combined);
143         col[3] = shr_new.alpha;
144 }
145
146 static void vol_trace_behind(ShadeInput *shi, VlakRen *vlr, float *co, float *col)
147 {
148         Isect isect;
149         float maxsize = RE_ray_tree_max_size(R.raytree);
150         
151         VECCOPY(isect.start, co);
152         isect.end[0] = isect.start[0] + shi->view[0] * maxsize;
153         isect.end[1] = isect.start[1] + shi->view[1] * maxsize;
154         isect.end[2] = isect.start[2] + shi->view[2] * maxsize;
155         
156         isect.faceorig= (RayFace *)vlr;
157         
158         isect.mode= RE_RAY_MIRROR;
159         isect.oborig= RAY_OBJECT_SET(&R, shi->obi);
160         isect.face_last= NULL;
161         isect.ob_last= 0;
162         isect.lay= -1;
163         
164         /* check to see if there's anything behind the volume, otherwise shade the sky */
165         if(RE_ray_tree_intersect(R.raytree, &isect)) {
166                 shade_intersection(shi, col, &isect);
167         } else {
168                 shadeSkyView(col, co, shi->view, NULL, shi->thread);
169                 shadeSunView(col, shi->view);
170         }
171 }
172
173
174 /* trilinear interpolation */
175 static void vol_get_precached_scattering(ShadeInput *shi, float *scatter_col, float *co)
176 {
177         VolumePrecache *vp = shi->obi->volume_precache;
178         float bbmin[3], bbmax[3], dim[3];
179         float sample_co[3];
180         
181         if (!vp) return;
182         
183         /* convert input coords to 0.0, 1.0 */
184         VECCOPY(bbmin, shi->obi->obr->boundbox[0]);
185         VECCOPY(bbmax, shi->obi->obr->boundbox[1]);
186         VecSubf(dim, bbmax, bbmin);
187
188         sample_co[0] = ((co[0] - bbmin[0]) / dim[0]);
189         sample_co[1] = ((co[1] - bbmin[1]) / dim[1]);
190         sample_co[2] = ((co[2] - bbmin[2]) / dim[2]);
191
192         scatter_col[0] = voxel_sample_triquadratic(vp->data_r, vp->res, sample_co);
193         scatter_col[1] = voxel_sample_triquadratic(vp->data_g, vp->res, sample_co);
194         scatter_col[2] = voxel_sample_triquadratic(vp->data_b, vp->res, sample_co);
195 }
196
197 /* Meta object density, brute force for now 
198  * (might be good enough anyway, don't need huge number of metaobs to model volumetric objects */
199 static float metadensity(Object* ob, float* co)
200 {
201         float mat[4][4], imat[4][4], dens = 0.f;
202         MetaBall* mb = (MetaBall*)ob->data;
203         MetaElem* ml;
204         
205         /* transform co to meta-element */
206         float tco[3] = {co[0], co[1], co[2]};
207         Mat4MulMat4(mat, ob->obmat, R.viewmat);
208         Mat4Invert(imat, mat);
209         Mat4MulVecfl(imat, tco);
210         
211         for (ml = mb->elems.first; ml; ml=ml->next) {
212                 float bmat[3][3], dist2;
213                 
214                 /* element rotation transform */
215                 float tp[3] = {ml->x - tco[0], ml->y - tco[1], ml->z - tco[2]};
216                 QuatToMat3(ml->quat, bmat);
217                 Mat3Transp(bmat);       // rot.only, so inverse == transpose
218                 Mat3MulVecfl(bmat, tp);
219                 
220                 /* MB_BALL default */
221                 switch (ml->type) {
222                         case MB_ELIPSOID:
223                                 tp[0] /= ml->expx, tp[1] /= ml->expy, tp[2] /= ml->expz;
224                                 break;
225                         case MB_CUBE:
226                                 tp[2] = (tp[2] > ml->expz) ? (tp[2] - ml->expz) : ((tp[2] < -ml->expz) ? (tp[2] + ml->expz) : 0.f);
227                                 // no break, xy as plane
228                         case MB_PLANE:
229                                 tp[1] = (tp[1] > ml->expy) ? (tp[1] - ml->expy) : ((tp[1] < -ml->expy) ? (tp[1] + ml->expy) : 0.f);
230                                 // no break, x as tube
231                         case MB_TUBE:
232                                 tp[0] = (tp[0] > ml->expx) ? (tp[0] - ml->expx) : ((tp[0] < -ml->expx) ? (tp[0] + ml->expx) : 0.f);
233                 }
234                 
235                 /* ml->rad2 is not set */
236                 dist2 = 1.f - ((tp[0]*tp[0] + tp[1]*tp[1] + tp[2]*tp[2]) / (ml->rad*ml->rad));
237                 if (dist2 > 0.f)
238                         dens += (ml->flag & MB_NEGATIVE) ? -ml->s*dist2*dist2*dist2 : ml->s*dist2*dist2*dist2;
239         }
240         
241         dens -= mb->thresh;
242         return (dens < 0.f) ? 0.f : dens;
243 }
244
245 float vol_get_density(struct ShadeInput *shi, float *co)
246 {
247         float density = shi->mat->vol.density;
248         float density_scale = shi->mat->vol.density_scale;
249                 
250         if (shi->mat->mapto_textured & MAP_DENSITY)
251                 do_volume_tex(shi, co, MAP_DENSITY, NULL, &density);
252         
253         // if meta-object, modulate by metadensity without increasing it
254         if (shi->obi->obr->ob->type == OB_MBALL) {
255                 const float md = metadensity(shi->obi->obr->ob, co);
256                 if (md < 1.f) density *= md;
257         }
258         
259         return density * density_scale;
260 }
261
262
263 /* Color of light that gets scattered out by the volume */
264 /* Uses same physically based scattering parameter as in transmission calculations, 
265  * along with artificial reflection scale/reflection color tint */
266 void vol_get_reflection_color(ShadeInput *shi, float *ref_col, float *co)
267 {
268         float scatter = shi->mat->vol.scattering;
269         float reflection= shi->mat->vol.reflection;
270         VECCOPY(ref_col, shi->mat->vol.reflection_col);
271         
272         if (shi->mat->mapto_textured & (MAP_SCATTERING+MAP_REFLECTION_COL))
273                 do_volume_tex(shi, co, MAP_SCATTERING+MAP_REFLECTION_COL, ref_col, &scatter);
274         
275         /* only one single float parameter at a time... :s */
276         if (shi->mat->mapto_textured & (MAP_REFLECTION))
277                 do_volume_tex(shi, co, MAP_REFLECTION, NULL, &reflection);
278         
279         ref_col[0] = reflection * ref_col[0] * scatter;
280         ref_col[1] = reflection * ref_col[1] * scatter;
281         ref_col[2] = reflection * ref_col[2] * scatter;
282 }
283
284 /* compute emission component, amount of radiance to add per segment
285  * can be textured with 'emit' */
286 void vol_get_emission(ShadeInput *shi, float *emission_col, float *co)
287 {
288         float emission = shi->mat->vol.emission;
289         VECCOPY(emission_col, shi->mat->vol.emission_col);
290         
291         if (shi->mat->mapto_textured & (MAP_EMISSION+MAP_EMISSION_COL))
292                 do_volume_tex(shi, co, MAP_EMISSION+MAP_EMISSION_COL, emission_col, &emission);
293         
294         emission_col[0] = emission_col[0] * emission;
295         emission_col[1] = emission_col[1] * emission;
296         emission_col[2] = emission_col[2] * emission;
297 }
298
299
300 /* A combination of scattering and absorption -> known as sigma T.
301  * This can possibly use a specific scattering colour, 
302  * and absorption multiplier factor too, but these parameters are left out for simplicity.
303  * It's easy enough to get a good wide range of results with just these two parameters. */
304 void vol_get_sigma_t(ShadeInput *shi, float *sigma_t, float *co)
305 {
306         /* technically absorption, but named transmission color 
307          * since it describes the effect of the coloring *after* absorption */
308         float transmission_col[3] = {shi->mat->vol.transmission_col[0], shi->mat->vol.transmission_col[1], shi->mat->vol.transmission_col[2]};
309         float scattering = shi->mat->vol.scattering;
310         
311         if (shi->mat->mapto_textured & (MAP_SCATTERING+MAP_TRANSMISSION_COL))
312                 do_volume_tex(shi, co, MAP_SCATTERING+MAP_TRANSMISSION_COL, transmission_col, &scattering);
313         
314         sigma_t[0] = (1.0f - transmission_col[0]) + scattering;
315         sigma_t[1] = (1.0f - transmission_col[1]) + scattering;
316         sigma_t[2] = (1.0f - transmission_col[2]) + scattering;
317 }
318
319 /* phase function - determines in which directions the light 
320  * is scattered in the volume relative to incoming direction 
321  * and view direction */
322 float vol_get_phasefunc(ShadeInput *shi, float g, float *w, float *wp)
323 {
324         const float normalize = 0.25f; // = 1.f/4.f = M_PI/(4.f*M_PI)
325         
326         /* normalization constant is 1/4 rather than 1/4pi, since
327          * Blender's shading system doesn't normalise for
328          * energy conservation - eg. multiplying by pdf ( 1/pi for a lambert brdf ).
329          * This means that lambert surfaces in Blender are pi times brighter than they 'should be'
330          * and therefore, with correct energy conservation, volumes will darker than other solid objects,
331          * for the same lighting intensity.
332          * To correct this, scale up the phase function values by pi
333          * until Blender's shading system supports this better. --matt
334          */
335         
336         if (g == 0.f) { /* isotropic */
337                 return normalize * 1.f;
338         } else {                /* schlick */
339                 const float k = 1.55f * g - .55f * g * g * g;
340                 const float kcostheta = k * Inpf(w, wp);
341                 return normalize * (1.f - k*k) / ((1.f - kcostheta) * (1.f - kcostheta));
342         }
343         
344         /*
345          * not used, but here for reference:
346         switch (phasefunc_type) {
347                 case MA_VOL_PH_MIEHAZY:
348                         return normalize * (0.5f + 4.5f * powf(0.5 * (1.f + costheta), 8.f));
349                 case MA_VOL_PH_MIEMURKY:
350                         return normalize * (0.5f + 16.5f * powf(0.5 * (1.f + costheta), 32.f));
351                 case MA_VOL_PH_RAYLEIGH:
352                         return normalize * 3.f/4.f * (1 + costheta * costheta);
353                 case MA_VOL_PH_HG:
354                         return normalize * (1.f - g*g) / powf(1.f + g*g - 2.f * g * costheta, 1.5f));
355                 case MA_VOL_PH_SCHLICK:
356                 {
357                         const float k = 1.55f * g - .55f * g * g * g;
358                         const float kcostheta = k * costheta;
359                         return normalize * (1.f - k*k) / ((1.f - kcostheta) * (1.f - kcostheta));
360                 }
361                 case MA_VOL_PH_ISOTROPIC:
362                 default:
363                         return normalize * 1.f;
364         }
365         */
366 }
367
368 /* Compute transmittance = e^(-attenuation) */
369 void vol_get_transmittance_seg(ShadeInput *shi, float *tr, float stepsize, float *co, float density)
370 {
371         /* input density = density at co */
372         float tau[3] = {0.f, 0.f, 0.f};
373         const float stepd = density * stepsize;
374         float sigma_t[3];
375         
376         vol_get_sigma_t(shi, sigma_t, co);
377         
378         /* homogenous volume within the sampled distance */
379         tau[0] += stepd * sigma_t[0];
380         tau[1] += stepd * sigma_t[1];
381         tau[2] += stepd * sigma_t[2];
382         
383         tr[0] *= exp(-tau[0]);
384         tr[1] *= exp(-tau[1]);
385         tr[2] *= exp(-tau[2]);
386 }
387
388 /* Compute transmittance = e^(-attenuation) */
389 static void vol_get_transmittance(ShadeInput *shi, float *tr, float *co, float *endco)
390 {
391         float p[3] = {co[0], co[1], co[2]};
392         float step_vec[3] = {endco[0] - co[0], endco[1] - co[1], endco[2] - co[2]};
393         float tau[3] = {0.f, 0.f, 0.f};
394
395         float t0 = 0.f;
396         float t1 = Normalize(step_vec);
397         float pt0 = t0;
398         
399         t0 += shi->mat->vol.stepsize * ((shi->mat->vol.stepsize_type == MA_VOL_STEP_CONSTANT) ? 0.5f : BLI_thread_frand(shi->thread));
400         p[0] += t0 * step_vec[0];
401         p[1] += t0 * step_vec[1];
402         p[2] += t0 * step_vec[2];
403         VecMulf(step_vec, shi->mat->vol.stepsize);
404
405         for (; t0 < t1; pt0 = t0, t0 += shi->mat->vol.stepsize) {
406                 const float d = vol_get_density(shi, p);
407                 const float stepd = (t0 - pt0) * d;
408                 float sigma_t[3];
409                 
410                 vol_get_sigma_t(shi, sigma_t, co);
411                 
412                 tau[0] += stepd * sigma_t[0];
413                 tau[1] += stepd * sigma_t[1];
414                 tau[2] += stepd * sigma_t[2];
415                 
416                 VecAddf(p, p, step_vec);
417         }
418         
419         /* return transmittance */
420         tr[0] = expf(-tau[0]);
421         tr[1] = expf(-tau[1]);
422         tr[2] = expf(-tau[2]);
423 }
424
425 void vol_shade_one_lamp(struct ShadeInput *shi, float *co, LampRen *lar, float *lacol)
426 {
427         float visifac, lv[3], lampdist;
428         float tr[3]={1.0,1.0,1.0};
429         float hitco[3], *atten_co;
430         float p, ref_col[3];
431         
432         if (lar->mode & LA_LAYER) if((lar->lay & shi->obi->lay)==0) return;
433         if ((lar->lay & shi->lay)==0) return;
434         if (lar->energy == 0.0) return;
435         
436         if ((visifac= lamp_get_visibility(lar, co, lv, &lampdist)) == 0.f) return;
437         
438         VecCopyf(lacol, &lar->r);
439         
440         if(lar->mode & LA_TEXTURE) {
441                 shi->osatex= 0;
442                 do_lamp_tex(lar, lv, shi, lacol, LA_TEXTURE);
443         }
444
445         VecMulf(lacol, visifac);
446
447         if (ELEM(lar->type, LA_SUN, LA_HEMI))
448                 VECCOPY(lv, lar->vec);
449         VecMulf(lv, -1.0f);
450         
451         if (shi->mat->vol.shade_type != MA_VOL_SHADE_NONE) {
452                 Isect is;
453                 
454                 /* find minimum of volume bounds, or lamp coord */
455                 if (vol_get_bounds(shi, co, lv, hitco, &is, VOL_BOUNDS_SS)) {
456                         float dist = VecLenf(co, hitco);
457                         VlakRen *vlr = (VlakRen *)is.face;
458                         
459                         /* simple internal shadowing */
460                         if (vlr->mat->material_type == MA_TYPE_SURFACE) {
461                                 lacol[0] = lacol[1] = lacol[2] = 0.0f;
462                                 return;
463                         }
464
465                         if (ELEM(lar->type, LA_SUN, LA_HEMI))
466                                 /* infinite lights, can never be inside volume */
467                                 atten_co = hitco;
468                         else if ( lampdist < dist ) {
469                                 atten_co = lar->co;
470                         } else
471                                 atten_co = hitco;
472                         
473                         vol_get_transmittance(shi, tr, co, atten_co);
474                         
475                         VecMulVecf(lacol, lacol, tr);
476                 }
477                 else {
478                         /* Point is on the outside edge of the volume,
479                          * therefore no attenuation, full transmission.
480                          * Radiance from lamp remains unchanged */
481                 }
482         }
483         
484         if (luminance(lacol) < 0.001f) return;
485         
486         p = vol_get_phasefunc(shi, shi->mat->vol.asymmetry, shi->view, lv);
487         
488         /* physically based scattering with non-physically based RGB gain */
489         vol_get_reflection_color(shi, ref_col, co);
490         
491         lacol[0] *= p * ref_col[0];
492         lacol[1] *= p * ref_col[1];
493         lacol[2] *= p * ref_col[2];
494 }
495
496 /* single scattering only for now */
497 void vol_get_scattering(ShadeInput *shi, float *scatter_col, float *co)
498 {
499         ListBase *lights;
500         GroupObject *go;
501         LampRen *lar;
502         
503         scatter_col[0] = scatter_col[1] = scatter_col[2] = 0.f;
504         
505         lights= get_lights(shi);
506         for(go=lights->first; go; go= go->next)
507         {
508                 float lacol[3] = {0.f, 0.f, 0.f};
509                 lar= go->lampren;
510                 
511                 if (lar) {
512                         vol_shade_one_lamp(shi, co, lar, lacol);
513                         VecAddf(scatter_col, scatter_col, lacol);
514                 }
515         }
516 }
517
518         
519 /*
520 The main volumetric integrator, using an emission/absorption/scattering model.
521
522 Incoming radiance = 
523
524 outgoing radiance from behind surface * beam transmittance/attenuation
525 + added radiance from all points along the ray due to participating media
526         --> radiance for each segment = 
527                 (radiance added by scattering + radiance added by emission) * beam transmittance/attenuation
528 */
529
530 /* For ease of use, I've also introduced a 'reflection' and 'reflection color' parameter, which isn't 
531  * physically correct. This works as an RGB tint/gain on out-scattered light, but doesn't affect the light 
532  * that is transmitted through the volume. While having wavelength dependent absorption/scattering is more correct,
533  * it also makes it harder to control the overall look of the volume since colouring the outscattered light results
534  * in the inverse colour being transmitted through the rest of the volume.
535  */
536 static void volumeintegrate(struct ShadeInput *shi, float *col, float *co, float *endco)
537 {
538         float radiance[3] = {0.f, 0.f, 0.f};
539         float tr[3] = {1.f, 1.f, 1.f};
540         float p[3] = {co[0], co[1], co[2]};
541         float step_vec[3] = {endco[0] - co[0], endco[1] - co[1], endco[2] - co[2]};
542         const float stepsize = shi->mat->vol.stepsize;
543         
544         float t0 = 0.f;
545         float pt0 = t0;
546         float t1 = Normalize(step_vec); /* returns vector length */
547         
548         t0 += stepsize * ((shi->mat->vol.stepsize_type == MA_VOL_STEP_CONSTANT) ? 0.5f : BLI_thread_frand(shi->thread));
549         p[0] += t0 * step_vec[0];
550         p[1] += t0 * step_vec[1];
551         p[2] += t0 * step_vec[2];
552         VecMulf(step_vec, stepsize);
553         
554         for (; t0 < t1; pt0 = t0, t0 += stepsize) {
555                 const float density = vol_get_density(shi, p);
556                 
557                 if (density > 0.01f) {
558                         float scatter_col[3], emit_col[3];
559                         const float stepd = (t0 - pt0) * density;
560                         
561                         /* transmittance component (alpha) */
562                         vol_get_transmittance_seg(shi, tr, stepsize, co, density);
563                         
564                         if (luminance(tr) < shi->mat->vol.depth_cutoff) break;
565                         
566                         vol_get_emission(shi, emit_col, p);
567                         
568                         if (shi->obi->volume_precache) {
569                                 float p2[3];
570                                 
571                                 p2[0] = p[0] + (step_vec[0] * 0.5);
572                                 p2[1] = p[1] + (step_vec[1] * 0.5);
573                                 p2[2] = p[2] + (step_vec[2] * 0.5);
574                                 
575                                 vol_get_precached_scattering(shi, scatter_col, p2);
576                         } else
577                                 vol_get_scattering(shi, scatter_col, p);
578                         
579                         radiance[0] += stepd * tr[0] * (emit_col[0] + scatter_col[0]);
580                         radiance[1] += stepd * tr[1] * (emit_col[1] + scatter_col[1]);
581                         radiance[2] += stepd * tr[2] * (emit_col[2] + scatter_col[2]);
582                 }
583                 VecAddf(p, p, step_vec);
584         }
585         
586         /* multiply original color (from behind volume) with transmittance over entire distance */
587         VecMulVecf(col, tr, col);
588         VecAddf(col, col, radiance);
589         
590         /* alpha <-- transmission luminance */
591         col[3] = 1.0f - luminance(tr);
592 }
593
594 /* the main entry point for volume shading */
595 static void volume_trace(struct ShadeInput *shi, struct ShadeResult *shr, int inside_volume)
596 {
597         float hitco[3], col[4] = {0.f,0.f,0.f,0.f};
598         float *startco, *endco;
599         int trace_behind = 1;
600         const int ztransp= ((shi->depth==0) && (shi->mat->mode & MA_TRANSP) && (shi->mat->mode & MA_ZTRANSP));
601         Isect is;
602
603         /* check for shading an internal face a volume object directly */
604         if (inside_volume == VOL_SHADE_INSIDE)
605                 trace_behind = 0;
606         else if (inside_volume == VOL_SHADE_OUTSIDE) {
607                 if (shi->flippednor)
608                         inside_volume = VOL_SHADE_INSIDE;
609         }
610         
611         if (ztransp && inside_volume == VOL_SHADE_INSIDE) {
612                 MatInside *mi;
613                 int render_this=0;
614                 
615                 /* don't render the backfaces of ztransp volume materials.
616                  
617                  * volume shading renders the internal volume from between the
618                  * ' view intersection of the solid volume to the
619                  * intersection on the other side, as part of the shading of
620                  * the front face.
621                  
622                  * Because ztransp renders both front and back faces independently
623                  * this will double up, so here we prevent rendering the backface as well, 
624                  * which would otherwise render the volume in between the camera and the backface
625                  * --matt */
626                 
627                 for (mi=R.render_volumes_inside.first; mi; mi=mi->next) {
628                         /* weak... */
629                         if (mi->ma == shi->mat) render_this=1;
630                 }
631                 if (!render_this) return;
632         }
633         
634
635         if (inside_volume == VOL_SHADE_INSIDE)
636         {
637                 startco = shi->camera_co;
638                 endco = shi->co;
639                 
640                 if (trace_behind) {
641                         if (!ztransp)
642                                 /* trace behind the volume object */
643                                 vol_trace_behind(shi, shi->vlr, endco, col);
644                 } else {
645                         /* we're tracing through the volume between the camera 
646                          * and a solid surface, so use that pre-shaded radiance */
647                         QUATCOPY(col, shr->combined);
648                 }
649                 
650                 /* shade volume from 'camera' to 1st hit point */
651                 volumeintegrate(shi, col, startco, endco);
652         }
653         /* trace to find a backface, the other side bounds of the volume */
654         /* (ray intersect ignores front faces here) */
655         else if (vol_get_bounds(shi, shi->co, shi->view, hitco, &is, VOL_BOUNDS_DEPTH))
656         {
657                 VlakRen *vlr = (VlakRen *)is.face;
658                 
659                 startco = shi->co;
660                 endco = hitco;
661                 
662                 if (!ztransp) {
663                         /* if it's another face in the same material */
664                         if (vlr->mat == shi->mat) {
665                                 /* trace behind the 2nd (raytrace) hit point */
666                                 vol_trace_behind(shi, (VlakRen *)is.face, endco, col);
667                         } else {
668                                 shade_intersection(shi, col, &is);
669                         }
670                 }
671                 
672                 /* shade volume from 1st hit point to 2nd hit point */
673                 volumeintegrate(shi, col, startco, endco);
674         }
675         
676         if (ztransp)
677                 col[3] = col[3]>1.f?1.f:col[3];
678         else
679                 col[3] = 1.f;
680         
681         VecCopyf(shr->combined, col);
682         shr->alpha = col[3];
683         
684         VECCOPY(shr->diff, shr->combined);
685 }
686
687 /* Traces a shadow through the object, 
688  * pretty much gets the transmission over a ray path */
689 void shade_volume_shadow(struct ShadeInput *shi, struct ShadeResult *shr, struct Isect *last_is)
690 {
691         float hitco[3];
692         float tr[3] = {1.0,1.0,1.0};
693         Isect is;
694         float *startco, *endco;
695         float density=0.f;
696
697         memset(shr, 0, sizeof(ShadeResult));
698         
699         /* if 1st hit normal is facing away from the camera, 
700          * then we're inside the volume already. */
701         if (shi->flippednor) {
702                 startco = last_is->start;
703                 endco = shi->co;
704         }
705         /* trace to find a backface, the other side bounds of the volume */
706         /* (ray intersect ignores front faces here) */
707         else if (vol_get_bounds(shi, shi->co, shi->view, hitco, &is, VOL_BOUNDS_DEPTH)) {
708                 startco = shi->co;
709                 endco = hitco;
710         }
711         else {
712                 shr->combined[0] = shr->combined[1] = shr->combined[2] = 0.f;
713                 shr->alpha = shr->combined[3] = 1.f;
714                 return;
715         }
716         
717         density = vol_get_density(shi, startco);
718         vol_get_transmittance(shi, tr, startco, endco);
719         
720         VecCopyf(shr->combined, tr);
721         shr->combined[3] = 1.0f - luminance(tr);
722 }
723
724
725 /* delivers a fully filled in ShadeResult, for all passes */
726 void shade_volume_outside(ShadeInput *shi, ShadeResult *shr)
727 {
728         memset(shr, 0, sizeof(ShadeResult));
729         volume_trace(shi, shr, VOL_SHADE_OUTSIDE);
730 }
731
732
733 void shade_volume_inside(ShadeInput *shi, ShadeResult *shr)
734 {
735         MatInside *m;
736         Material *mat_backup;
737         ObjectInstanceRen *obi_backup;
738         float prev_alpha = shr->alpha;
739         
740         //if (BLI_countlist(&R.render_volumes_inside) == 0) return;
741         
742         /* XXX: extend to multiple volumes perhaps later */
743         mat_backup = shi->mat;
744         obi_backup = shi->obi;
745         
746         m = R.render_volumes_inside.first;
747         shi->mat = m->ma;
748         shi->obi = m->obi;
749         shi->obr = m->obi->obr;
750         
751         volume_trace(shi, shr, VOL_SHADE_INSIDE);
752         shr->alpha += prev_alpha;
753         CLAMP(shr->alpha, 0.f, 1.f);
754         
755         shi->mat = mat_backup;
756         shi->obi = obi_backup;
757         shi->obr = obi_backup->obr;
758 }