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