svn merge -r 23207:23528 https://svn.blender.org/svnroot/bf-blender/trunk/blender
[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         /* XXX TODO - get raytrace max distance from object instance's bounding box */
79         /* need to account for scaling only, but keep coords in camera space...
80          * below code is WIP and doesn't work!
81         VecSubf(bb_dim, shi->obi->obr->boundbox[1], shi->obi->obr->boundbox[2]);
82         Mat3MulVecfl(shi->obi->nmat, bb_dim);
83         maxsize = VecLength(bb_dim);
84         */
85         
86         VECCOPY(isect->start, co);
87         VECCOPY(isect->vec, vec );
88         isect->labda = FLT_MAX;
89         /*
90         isect->end[0] = co[0] + vec[0] * maxsize;
91         isect->end[1] = co[1] + vec[1] * maxsize;
92         isect->end[2] = co[2] + vec[2] * maxsize;
93         */
94         
95         isect->mode= RE_RAY_MIRROR;
96         isect->orig.ob = (void*)shi->obi;
97         isect->last_hit = NULL;
98         isect->lay= -1;
99         
100         if (intersect_type == VOL_BOUNDS_DEPTH) isect->orig.face = (void*)shi->vlr;
101         else if (intersect_type == VOL_BOUNDS_SS) isect->orig.face= NULL;
102         
103         if(RE_rayobject_raycast(R.raytree, isect))
104         {
105                 hitco[0] = isect->start[0] + isect->labda*isect->vec[0];
106                 hitco[1] = isect->start[1] + isect->labda*isect->vec[1];
107                 hitco[2] = isect->start[2] + isect->labda*isect->vec[2];
108                 return 1;
109         } else {
110                 return 0;
111         }
112 }
113
114 static void shade_intersection(ShadeInput *shi, float *col, Isect *is)
115 {
116         ShadeInput shi_new;
117         ShadeResult shr_new;
118         
119         memset(&shi_new, 0, sizeof(ShadeInput)); 
120         
121         shi_new.mask= shi->mask;
122         shi_new.osatex= shi->osatex;
123         shi_new.thread= shi->thread;
124         shi_new.depth = shi->depth + 1;
125         shi_new.volume_depth= shi->volume_depth + 1;
126         shi_new.xs= shi->xs;
127         shi_new.ys= shi->ys;
128         shi_new.lay= shi->lay;
129         shi_new.passflag= SCE_PASS_COMBINED; /* result of tracing needs no pass info */
130         shi_new.combinedflag= 0xFFFFFF;          /* ray trace does all options */
131         shi_new.light_override= shi->light_override;
132         shi_new.mat_override= shi->mat_override;
133         
134         VECCOPY(shi_new.camera_co, is->start);
135         
136         memset(&shr_new, 0, sizeof(ShadeResult));
137         
138         /* hardcoded limit of 100 for now - prevents problems in weird geometry */
139         if (shi->volume_depth < 100) {
140                 shade_ray(is, &shi_new, &shr_new);
141         }
142         
143         VecCopyf(col, shr_new.combined);
144         col[3] = shr_new.alpha;
145 }
146
147 static void vol_trace_behind(ShadeInput *shi, VlakRen *vlr, float *co, float *col)
148 {
149         Isect isect;
150         
151         VECCOPY(isect.start, co);
152         VECCOPY(isect.vec, shi->view);
153         isect.labda = FLT_MAX;
154         
155         isect.mode= RE_RAY_MIRROR;
156         isect.orig.ob = (void*) shi->obi;
157         isect.orig.face = (void*)vlr;
158         isect.last_hit = NULL;
159         isect.lay= -1;
160         
161         /* check to see if there's anything behind the volume, otherwise shade the sky */
162         if(RE_rayobject_raycast(R.raytree, &isect)) {
163                 shade_intersection(shi, col, &isect);
164         } else {
165                 shadeSkyView(col, co, shi->view, NULL, shi->thread);
166                 shadeSunView(col, shi->view);
167         }
168 }
169
170 /* input shader data */
171
172 float vol_get_stepsize(struct ShadeInput *shi, int context)
173 {
174         if (shi->mat->vol.stepsize_type == MA_VOL_STEP_RANDOMIZED) {
175                 /* range between 0.75 and 1.25 */
176                 const float rnd = 0.5f * BLI_thread_frand(shi->thread) + 0.75f;
177         
178                 if (context == STEPSIZE_VIEW)
179                         return shi->mat->vol.stepsize * rnd;
180                 else if (context == STEPSIZE_SHADE)
181                         return shi->mat->vol.shade_stepsize * rnd;
182         }
183         else {  // MA_VOL_STEP_CONSTANT
184                 
185                 if (context == STEPSIZE_VIEW)
186                         return shi->mat->vol.stepsize;
187                 else if (context == STEPSIZE_SHADE)
188                         return shi->mat->vol.shade_stepsize;
189         }
190         
191         return shi->mat->vol.stepsize;
192 }
193
194 /* trilinear interpolation */
195 static void vol_get_precached_scattering(ShadeInput *shi, float *scatter_col, float *co)
196 {
197         VolumePrecache *vp = shi->obi->volume_precache;
198         float bbmin[3], bbmax[3], dim[3];
199         float sample_co[3];
200         
201         if (!vp) return;
202         
203         /* convert input coords to 0.0, 1.0 */
204         VECCOPY(bbmin, shi->obi->obr->boundbox[0]);
205         VECCOPY(bbmax, shi->obi->obr->boundbox[1]);
206         VecSubf(dim, bbmax, bbmin);
207
208         sample_co[0] = ((co[0] - bbmin[0]) / dim[0]);
209         sample_co[1] = ((co[1] - bbmin[1]) / dim[1]);
210         sample_co[2] = ((co[2] - bbmin[2]) / dim[2]);
211
212         scatter_col[0] = voxel_sample_trilinear(vp->data_r, vp->res, sample_co);
213         scatter_col[1] = voxel_sample_trilinear(vp->data_g, vp->res, sample_co);
214         scatter_col[2] = voxel_sample_trilinear(vp->data_b, vp->res, sample_co);
215 }
216
217 /* Meta object density, brute force for now 
218  * (might be good enough anyway, don't need huge number of metaobs to model volumetric objects */
219 static float metadensity(Object* ob, float* co)
220 {
221         float mat[4][4], imat[4][4], dens = 0.f;
222         MetaBall* mb = (MetaBall*)ob->data;
223         MetaElem* ml;
224         
225         /* transform co to meta-element */
226         float tco[3] = {co[0], co[1], co[2]};
227         Mat4MulMat4(mat, ob->obmat, R.viewmat);
228         Mat4Invert(imat, mat);
229         Mat4MulVecfl(imat, tco);
230         
231         for (ml = mb->elems.first; ml; ml=ml->next) {
232                 float bmat[3][3], dist2;
233                 
234                 /* element rotation transform */
235                 float tp[3] = {ml->x - tco[0], ml->y - tco[1], ml->z - tco[2]};
236                 QuatToMat3(ml->quat, bmat);
237                 Mat3Transp(bmat);       // rot.only, so inverse == transpose
238                 Mat3MulVecfl(bmat, tp);
239                 
240                 /* MB_BALL default */
241                 switch (ml->type) {
242                         case MB_ELIPSOID:
243                                 tp[0] /= ml->expx, tp[1] /= ml->expy, tp[2] /= ml->expz;
244                                 break;
245                         case MB_CUBE:
246                                 tp[2] = (tp[2] > ml->expz) ? (tp[2] - ml->expz) : ((tp[2] < -ml->expz) ? (tp[2] + ml->expz) : 0.f);
247                                 // no break, xy as plane
248                         case MB_PLANE:
249                                 tp[1] = (tp[1] > ml->expy) ? (tp[1] - ml->expy) : ((tp[1] < -ml->expy) ? (tp[1] + ml->expy) : 0.f);
250                                 // no break, x as tube
251                         case MB_TUBE:
252                                 tp[0] = (tp[0] > ml->expx) ? (tp[0] - ml->expx) : ((tp[0] < -ml->expx) ? (tp[0] + ml->expx) : 0.f);
253                 }
254                 
255                 /* ml->rad2 is not set */
256                 dist2 = 1.f - ((tp[0]*tp[0] + tp[1]*tp[1] + tp[2]*tp[2]) / (ml->rad*ml->rad));
257                 if (dist2 > 0.f)
258                         dens += (ml->flag & MB_NEGATIVE) ? -ml->s*dist2*dist2*dist2 : ml->s*dist2*dist2*dist2;
259         }
260         
261         dens -= mb->thresh;
262         return (dens < 0.f) ? 0.f : dens;
263 }
264
265 float vol_get_density(struct ShadeInput *shi, float *co)
266 {
267         float density = shi->mat->vol.density;
268         float density_scale = shi->mat->vol.density_scale;
269                 
270         do_volume_tex(shi, co, MAP_DENSITY, NULL, &density);
271         
272         // if meta-object, modulate by metadensity without increasing it
273         if (shi->obi->obr->ob->type == OB_MBALL) {
274                 const float md = metadensity(shi->obi->obr->ob, co);
275                 if (md < 1.f) density *= md;
276         }
277         
278         return density * density_scale;
279 }
280
281 /* scattering multiplier, values above 1.0 are non-physical, 
282  * but can be useful to tweak lighting */
283 float vol_get_scattering_fac(ShadeInput *shi, float *co)
284 {
285         float scatter = shi->mat->vol.scattering;
286         float col[3] = {0.0, 0.0, 0.0};
287         
288         do_volume_tex(shi, co, MAP_SCATTERING, col, &scatter);
289         
290         return scatter;
291 }
292
293 /* compute emission component, amount of radiance to add per segment
294  * can be textured with 'emit' */
295 void vol_get_emission(ShadeInput *shi, float *emission_col, float *co, float density)
296 {
297         float emission = shi->mat->vol.emission;
298         VECCOPY(emission_col, shi->mat->vol.emission_col);
299         
300         do_volume_tex(shi, co, MAP_EMISSION+MAP_EMISSION_COL, emission_col, &emission);
301         
302         emission_col[0] = emission_col[0] * emission * density;
303         emission_col[1] = emission_col[1] * emission * density;
304         emission_col[2] = emission_col[2] * emission * density;
305 }
306
307 void vol_get_absorption(ShadeInput *shi, float *absorb_col, float *co)
308 {
309         float absorption = shi->mat->vol.absorption;
310         VECCOPY(absorb_col, shi->mat->vol.absorption_col);
311         
312         do_volume_tex(shi, co, MAP_ABSORPTION+MAP_ABSORPTION_COL, absorb_col, &absorption);
313         
314         absorb_col[0] = (1.0f - absorb_col[0]) * absorption;
315         absorb_col[1] = (1.0f - absorb_col[1]) * absorption;
316         absorb_col[2] = (1.0f - absorb_col[2]) * absorption;
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, short phasefunc_type, float g, float *w, float *wp)
323 {
324         const float costheta = Inpf(w, wp);
325         const float scale = M_PI;
326         
327         /*
328          * Scale constant is required, since Blender's shading system doesn't normalise for
329          * energy conservation - eg. scaling by 1/pi for a lambert shader.
330          * This makes volumes darker than other solid objects, for the same lighting intensity.
331          * To correct this, scale up the phase function values
332          * until Blender's shading system supports this better. --matt
333          */
334         
335         switch (phasefunc_type) {
336                 case MA_VOL_PH_MIEHAZY:
337                         return scale * (0.5f + 4.5f * powf(0.5 * (1.f + costheta), 8.f)) / (4.f*M_PI);
338                 case MA_VOL_PH_MIEMURKY:
339                         return scale * (0.5f + 16.5f * powf(0.5 * (1.f + costheta), 32.f)) / (4.f*M_PI);
340                 case MA_VOL_PH_RAYLEIGH:
341                         return scale * 3.f/(16.f*M_PI) * (1 + costheta * costheta);
342                 case MA_VOL_PH_HG:
343                         return scale * (1.f / (4.f * M_PI) * (1.f - g*g) / powf(1.f + g*g - 2.f * g * costheta, 1.5f));
344                 case MA_VOL_PH_SCHLICK:
345                 {
346                         const float k = 1.55f * g - .55f * g * g * g;
347                         const float kcostheta = k * costheta;
348                         return scale * (1.f / (4.f * M_PI) * (1.f - k*k) / ((1.f - kcostheta) * (1.f - kcostheta)));
349                 }
350                 case MA_VOL_PH_ISOTROPIC:
351                 default:
352                         return scale * (1.f / (4.f * M_PI));
353         }
354 }
355
356 /* Compute transmittance = e^(-attenuation) */
357 void vol_get_transmittance_seg(ShadeInput *shi, float *tr, float stepsize, float *co, float density)
358 {
359         /* input density = density at co */
360         float tau[3] = {0.f, 0.f, 0.f};
361         float absorb[3];
362         const float scatter_dens = vol_get_scattering_fac(shi, co) * density * stepsize;
363
364         vol_get_absorption(shi, absorb, co);
365         
366         /* homogenous volume within the sampled distance */
367         tau[0] += scatter_dens * absorb[0];
368         tau[1] += scatter_dens * absorb[1];
369         tau[2] += scatter_dens * absorb[2];
370         
371         tr[0] *= exp(-tau[0]);
372         tr[1] *= exp(-tau[1]);
373         tr[2] *= exp(-tau[2]);
374 }
375
376 /* Compute transmittance = e^(-attenuation) */
377 static void vol_get_transmittance(ShadeInput *shi, float *tr, float *co, float *endco)
378 {
379         float p[3] = {co[0], co[1], co[2]};
380         float step_vec[3] = {endco[0] - co[0], endco[1] - co[1], endco[2] - co[2]};
381         //const float ambtau = -logf(shi->mat->vol.depth_cutoff);       // never zero
382         float tau[3] = {0.f, 0.f, 0.f};
383
384         float t0 = 0.f;
385         float t1 = Normalize(step_vec);
386         float pt0 = t0;
387         
388         t0 += shi->mat->vol.shade_stepsize * ((shi->mat->vol.stepsize_type == MA_VOL_STEP_CONSTANT) ? 0.5f : BLI_thread_frand(shi->thread));
389         p[0] += t0 * step_vec[0];
390         p[1] += t0 * step_vec[1];
391         p[2] += t0 * step_vec[2];
392         VecMulf(step_vec, shi->mat->vol.shade_stepsize);
393
394         for (; t0 < t1; pt0 = t0, t0 += shi->mat->vol.shade_stepsize) {
395                 float absorb[3];
396                 const float d = vol_get_density(shi, p);
397                 const float stepd = (t0 - pt0) * d;
398                 const float scatter_dens = vol_get_scattering_fac(shi, p) * stepd;
399                 vol_get_absorption(shi, absorb, p);
400                 
401                 tau[0] += scatter_dens * absorb[0];
402                 tau[1] += scatter_dens * absorb[1];
403                 tau[2] += scatter_dens * absorb[2];
404                 
405                 //if (luminance(tau) >= ambtau) break;
406                 VecAddf(p, p, step_vec);
407         }
408         
409         /* return transmittance */
410         tr[0] = expf(-tau[0]);
411         tr[1] = expf(-tau[1]);
412         tr[2] = expf(-tau[2]);
413 }
414
415 void vol_shade_one_lamp(struct ShadeInput *shi, float *co, LampRen *lar, float *lacol)
416 {
417         float visifac, lv[3], lampdist;
418         float tr[3]={1.0,1.0,1.0};
419         float hitco[3], *atten_co;
420         float p;
421         float scatter_fac;
422         
423         if (lar->mode & LA_LAYER) if((lar->lay & shi->obi->lay)==0) return;
424         if ((lar->lay & shi->lay)==0) return;
425         if (lar->energy == 0.0) return;
426         
427         if ((visifac= lamp_get_visibility(lar, co, lv, &lampdist)) == 0.f) return;
428         
429         VecCopyf(lacol, &lar->r);
430         
431         if(lar->mode & LA_TEXTURE) {
432                 shi->osatex= 0;
433                 do_lamp_tex(lar, lv, shi, lacol, LA_TEXTURE);
434         }
435
436         VecMulf(lacol, visifac);
437
438         if (ELEM(lar->type, LA_SUN, LA_HEMI))
439                 VECCOPY(lv, lar->vec);
440         VecMulf(lv, -1.0f);
441         
442         if (shi->mat->vol.shade_type != MA_VOL_SHADE_NONE) {
443                 Isect is;
444                 
445                 /* find minimum of volume bounds, or lamp coord */
446                 if (vol_get_bounds(shi, co, lv, hitco, &is, VOL_BOUNDS_SS)) {
447                         float dist = VecLenf(co, hitco);
448                         VlakRen *vlr = (VlakRen *)is.hit.face;
449                         
450                         /* simple internal shadowing */
451                         if (vlr->mat->material_type == MA_TYPE_SURFACE) {
452                                 lacol[0] = lacol[1] = lacol[2] = 0.0f;
453                                 return;
454                         }
455
456                         if (ELEM(lar->type, LA_SUN, LA_HEMI))
457                                 /* infinite lights, can never be inside volume */
458                                 atten_co = hitco;
459                         else if ( lampdist < dist ) {
460                                 atten_co = lar->co;
461                         } else
462                                 atten_co = hitco;
463                         
464                         vol_get_transmittance(shi, tr, co, atten_co);
465                         
466                         VecMulVecf(lacol, lacol, tr);
467                 }
468                 else {
469                         /* Point is on the outside edge of the volume,
470                          * therefore no attenuation, full transmission.
471                          * Radiance from lamp remains unchanged */
472                 }
473         }
474         
475         p = vol_get_phasefunc(shi, shi->mat->vol.phasefunc_type, shi->mat->vol.phasefunc_g, shi->view, lv);
476         VecMulf(lacol, p);
477         
478         scatter_fac = vol_get_scattering_fac(shi, co);
479         VecMulf(lacol, scatter_fac);
480 }
481
482 /* single scattering only for now */
483 void vol_get_scattering(ShadeInput *shi, float *scatter_col, float *co, float stepsize, float density)
484 {
485         ListBase *lights;
486         GroupObject *go;
487         LampRen *lar;
488         
489         scatter_col[0] = scatter_col[1] = scatter_col[2] = 0.f;
490         
491         lights= get_lights(shi);
492         for(go=lights->first; go; go= go->next)
493         {
494                 float lacol[3] = {0.f, 0.f, 0.f};
495                 lar= go->lampren;
496                 
497                 if (lar) {
498                         vol_shade_one_lamp(shi, co, lar, lacol);
499                         VecAddf(scatter_col, scatter_col, lacol);
500                 }
501         }
502 }
503
504         
505 /*
506 The main volumetric integrator, using an emission/absorption/scattering model.
507
508 Incoming radiance = 
509
510 outgoing radiance from behind surface * beam transmittance/attenuation
511 + added radiance from all points along the ray due to participating media
512         --> radiance for each segment = 
513                 (radiance added by scattering + radiance added by emission) * beam transmittance/attenuation
514 */
515 static void volumeintegrate(struct ShadeInput *shi, float *col, float *co, float *endco)
516 {
517         float tr[3] = {1.0f, 1.0f, 1.0f};
518         float radiance[3] = {0.f, 0.f, 0.f}, d_radiance[3] = {0.f, 0.f, 0.f};
519         float stepsize = vol_get_stepsize(shi, STEPSIZE_VIEW);
520         int nsteps, s;
521         float emit_col[3], scatter_col[3] = {0.0, 0.0, 0.0};
522         float stepvec[3], step_sta[3], step_end[3], step_mid[3];
523         float density;
524         const float depth_cutoff = shi->mat->vol.depth_cutoff;
525
526         /* ray marching */
527         nsteps = (int)((VecLenf(co, endco) / stepsize) + 0.5);
528         
529         VecSubf(stepvec, endco, co);
530         VecMulf(stepvec, 1.0f / nsteps);
531         VecCopyf(step_sta, co);
532         VecAddf(step_end, step_sta, stepvec);
533         
534         /* get radiance from all points along the ray due to participating media */
535         for (s = 0; s < nsteps; s++) {
536
537                 density = vol_get_density(shi, step_sta);
538                 
539                 /* there's only any use in shading here if there's actually some density to shade! */
540                 if (density > 0.01f) {
541                 
542                         /* transmittance component (alpha) */
543                         vol_get_transmittance_seg(shi, tr, stepsize, co, density);
544
545                         step_mid[0] = step_sta[0] + (stepvec[0] * 0.5);
546                         step_mid[1] = step_sta[1] + (stepvec[1] * 0.5);
547                         step_mid[2] = step_sta[2] + (stepvec[2] * 0.5);
548                 
549                         /* incoming light via emission or scattering (additive) */
550                         vol_get_emission(shi, emit_col, step_mid, density);
551                         
552                         if (shi->obi->volume_precache)
553                                 vol_get_precached_scattering(shi, scatter_col, step_mid);
554                         else
555                                 vol_get_scattering(shi, scatter_col, step_mid, stepsize, density);
556                         
557                         VecMulf(scatter_col, density);
558                         VecAddf(d_radiance, emit_col, scatter_col);
559                         
560                         /*   Lv += Tr * (Lve() + Ld) */
561                         VecMulVecf(d_radiance, tr, d_radiance);
562                         VecMulf(d_radiance, stepsize);
563                         
564                         VecAddf(radiance, radiance, d_radiance);        
565                 }
566
567                 VecCopyf(step_sta, step_end);
568                 VecAddf(step_end, step_end, stepvec);
569                 
570                 /* luminance rec. 709 */
571                 if ((0.2126*tr[0] + 0.7152*tr[1] + 0.0722*tr[2]) < depth_cutoff) break; 
572         }
573         
574         /* multiply original color (behind volume) with beam transmittance over entire distance */
575         VecMulVecf(col, tr, col);       
576         VecAddf(col, col, radiance);
577         
578         /* alpha <-- transmission luminance */
579         col[3] = 1.0f -(0.2126*tr[0] + 0.7152*tr[1] + 0.0722*tr[2]);
580 }
581
582 /* the main entry point for volume shading */
583 static void volume_trace(struct ShadeInput *shi, struct ShadeResult *shr, int inside_volume)
584 {
585         float hitco[3], col[4] = {0.f,0.f,0.f,0.f};
586         float *startco, *endco;
587         int trace_behind = 1;
588         const int ztransp= ((shi->depth==0) && (shi->mat->mode & MA_TRANSP) && (shi->mat->mode & MA_ZTRANSP));
589         Isect is;
590
591         /* check for shading an internal face a volume object directly */
592         if (inside_volume == VOL_SHADE_INSIDE)
593                 trace_behind = 0;
594         else if (inside_volume == VOL_SHADE_OUTSIDE) {
595                 if (shi->flippednor)
596                         inside_volume = VOL_SHADE_INSIDE;
597         }
598         
599         if (ztransp && inside_volume == VOL_SHADE_INSIDE) {
600                 MatInside *mi;
601                 int render_this=0;
602                 
603                 /* don't render the backfaces of ztransp volume materials.
604                  
605                  * volume shading renders the internal volume from between the
606                  * near view intersection of the solid volume to the
607                  * intersection on the other side, as part of the shading of
608                  * the front face.
609                  
610                  * Because ztransp renders both front and back faces independently
611                  * this will double up, so here we prevent rendering the backface as well, 
612                  * which would otherwise render the volume in between the camera and the backface
613                  * --matt */
614                 
615                 for (mi=R.render_volumes_inside.first; mi; mi=mi->next) {
616                         /* weak... */
617                         if (mi->ma == shi->mat) render_this=1;
618                 }
619                 if (!render_this) return;
620         }
621         
622
623         if (inside_volume == VOL_SHADE_INSIDE)
624         {
625                 startco = shi->camera_co;
626                 endco = shi->co;
627                 
628                 if (trace_behind) {
629                         if (!ztransp)
630                                 /* trace behind the volume object */
631                                 vol_trace_behind(shi, shi->vlr, endco, col);
632                 } else {
633                         /* we're tracing through the volume between the camera 
634                          * and a solid surface, so use that pre-shaded radiance */
635                         QUATCOPY(col, shr->combined);
636                 }
637                 
638                 /* shade volume from 'camera' to 1st hit point */
639                 volumeintegrate(shi, col, startco, endco);
640         }
641         /* trace to find a backface, the other side bounds of the volume */
642         /* (ray intersect ignores front faces here) */
643         else if (vol_get_bounds(shi, shi->co, shi->view, hitco, &is, VOL_BOUNDS_DEPTH))
644         {
645                 VlakRen *vlr = (VlakRen *)is.hit.face;
646                 
647                 startco = shi->co;
648                 endco = hitco;
649                 
650                 if (!ztransp) {
651                         /* if it's another face in the same material */
652                         if (vlr->mat == shi->mat) {
653                                 /* trace behind the 2nd (raytrace) hit point */
654                                 vol_trace_behind(shi, (VlakRen *)is.hit.face, endco, col);
655                         } else {
656                                 shade_intersection(shi, col, &is);
657                         }
658                 }
659                 
660                 /* shade volume from 1st hit point to 2nd hit point */
661                 volumeintegrate(shi, col, startco, endco);
662         }
663         
664         if (ztransp)
665                 col[3] = col[3]>1.f?1.f:col[3];
666         else
667                 col[3] = 1.f;
668         
669         VecCopyf(shr->combined, col);
670         shr->alpha = col[3];
671         
672         VECCOPY(shr->diff, shr->combined);
673 }
674
675 /* Traces a shadow through the object, 
676  * pretty much gets the transmission over a ray path */
677 void shade_volume_shadow(struct ShadeInput *shi, struct ShadeResult *shr, struct Isect *last_is)
678 {
679         float hitco[3];
680         float tr[3] = {1.0,1.0,1.0};
681         Isect is;
682         float *startco, *endco;
683         float density=0.f;
684
685         memset(shr, 0, sizeof(ShadeResult));
686         
687         /* if 1st hit normal is facing away from the camera, 
688          * then we're inside the volume already. */
689         if (shi->flippednor) {
690                 startco = last_is->start;
691                 endco = shi->co;
692         }
693         /* trace to find a backface, the other side bounds of the volume */
694         /* (ray intersect ignores front faces here) */
695         else if (vol_get_bounds(shi, shi->co, shi->view, hitco, &is, VOL_BOUNDS_DEPTH)) {
696                 startco = shi->co;
697                 endco = hitco;
698         }
699         else {
700                 shr->combined[0] = shr->combined[1] = shr->combined[2] = 0.f;
701                 shr->alpha = shr->combined[3] = 1.f;
702                 return;
703         }
704         
705         density = vol_get_density(shi, startco);
706         vol_get_transmittance(shi, tr, startco, endco);
707         
708         VecCopyf(shr->combined, tr);
709         shr->combined[3] = 1.0f -(0.2126*tr[0] + 0.7152*tr[1] + 0.0722*tr[2]);
710         shr->alpha = shr->combined[3];
711 }
712
713
714 /* delivers a fully filled in ShadeResult, for all passes */
715 void shade_volume_outside(ShadeInput *shi, ShadeResult *shr)
716 {
717         memset(shr, 0, sizeof(ShadeResult));
718         volume_trace(shi, shr, VOL_SHADE_OUTSIDE);
719 }
720
721
722 void shade_volume_inside(ShadeInput *shi, ShadeResult *shr)
723 {
724         MatInside *m;
725         Material *mat_backup;
726         ObjectInstanceRen *obi_backup;
727         float prev_alpha = shr->alpha;
728         
729         //if (BLI_countlist(&R.render_volumes_inside) == 0) return;
730         
731         /* XXX: extend to multiple volumes perhaps later */
732         mat_backup = shi->mat;
733         obi_backup = shi->obi;
734         
735         m = R.render_volumes_inside.first;
736         shi->mat = m->ma;
737         shi->obi = m->obi;
738         shi->obr = m->obi->obr;
739         
740         volume_trace(shi, shr, VOL_SHADE_INSIDE);
741         shr->alpha += prev_alpha;
742         CLAMP(shr->alpha, 0.f, 1.f);
743         
744         shi->mat = mat_backup;
745         shi->obi = obi_backup;
746         shi->obr = obi_backup->obr;
747 }