1d86ed7fbStbbdev /*
2*b15aabb3Stbbdev Copyright (c) 2005-2021 Intel Corporation
3d86ed7fbStbbdev
4d86ed7fbStbbdev Licensed under the Apache License, Version 2.0 (the "License");
5d86ed7fbStbbdev you may not use this file except in compliance with the License.
6d86ed7fbStbbdev You may obtain a copy of the License at
7d86ed7fbStbbdev
8d86ed7fbStbbdev http://www.apache.org/licenses/LICENSE-2.0
9d86ed7fbStbbdev
10d86ed7fbStbbdev Unless required by applicable law or agreed to in writing, software
11d86ed7fbStbbdev distributed under the License is distributed on an "AS IS" BASIS,
12d86ed7fbStbbdev WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13d86ed7fbStbbdev See the License for the specific language governing permissions and
14d86ed7fbStbbdev limitations under the License.
15d86ed7fbStbbdev */
16d86ed7fbStbbdev
17d86ed7fbStbbdev /*
18d86ed7fbStbbdev The original source for this example is
19d86ed7fbStbbdev Copyright (c) 1994-2008 John E. Stone
20d86ed7fbStbbdev All rights reserved.
21d86ed7fbStbbdev
22d86ed7fbStbbdev Redistribution and use in source and binary forms, with or without
23d86ed7fbStbbdev modification, are permitted provided that the following conditions
24d86ed7fbStbbdev are met:
25d86ed7fbStbbdev 1. Redistributions of source code must retain the above copyright
26d86ed7fbStbbdev notice, this list of conditions and the following disclaimer.
27d86ed7fbStbbdev 2. Redistributions in binary form must reproduce the above copyright
28d86ed7fbStbbdev notice, this list of conditions and the following disclaimer in the
29d86ed7fbStbbdev documentation and/or other materials provided with the distribution.
30d86ed7fbStbbdev 3. The name of the author may not be used to endorse or promote products
31d86ed7fbStbbdev derived from this software without specific prior written permission.
32d86ed7fbStbbdev
33d86ed7fbStbbdev THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS
34d86ed7fbStbbdev OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
35d86ed7fbStbbdev WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
36d86ed7fbStbbdev ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY
37d86ed7fbStbbdev DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
38d86ed7fbStbbdev DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
39d86ed7fbStbbdev OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
40d86ed7fbStbbdev HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
41d86ed7fbStbbdev LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
42d86ed7fbStbbdev OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
43d86ed7fbStbbdev SUCH DAMAGE.
44d86ed7fbStbbdev */
45d86ed7fbStbbdev
46d86ed7fbStbbdev /*
47d86ed7fbStbbdev * extvol.cpp - Volume rendering helper routines etc.
48d86ed7fbStbbdev */
49d86ed7fbStbbdev
50d86ed7fbStbbdev #include <cstdio>
51d86ed7fbStbbdev
52d86ed7fbStbbdev #include "machine.hpp"
53d86ed7fbStbbdev #include "types.hpp"
54d86ed7fbStbbdev #include "macros.hpp"
55d86ed7fbStbbdev #include "vector.hpp"
56d86ed7fbStbbdev #include "util.hpp"
57d86ed7fbStbbdev #include "box.hpp"
58d86ed7fbStbbdev #include "extvol.hpp"
59d86ed7fbStbbdev #include "trace.hpp"
60d86ed7fbStbbdev #include "sphere.hpp"
61d86ed7fbStbbdev #include "light.hpp"
62d86ed7fbStbbdev #include "shade.hpp"
63d86ed7fbStbbdev #include "global.hpp"
64d86ed7fbStbbdev
extvol_bbox(void * obj,vector * min,vector * max)65d86ed7fbStbbdev int extvol_bbox(void *obj, vector *min, vector *max) {
66d86ed7fbStbbdev box *b = (box *)obj;
67d86ed7fbStbbdev
68d86ed7fbStbbdev *min = b->min;
69d86ed7fbStbbdev *max = b->max;
70d86ed7fbStbbdev
71d86ed7fbStbbdev return 1;
72d86ed7fbStbbdev }
73d86ed7fbStbbdev
74d86ed7fbStbbdev static object_methods extvol_methods = { (void (*)(void *, void *))(box_intersect),
75d86ed7fbStbbdev (void (*)(void *, void *, void *, void *))(box_normal),
76d86ed7fbStbbdev extvol_bbox,
77d86ed7fbStbbdev free };
78d86ed7fbStbbdev
newextvol(void * voidtex,vector min,vector max,int samples,flt (* evaluator)(flt,flt,flt))79d86ed7fbStbbdev extvol *newextvol(void *voidtex,
80d86ed7fbStbbdev vector min,
81d86ed7fbStbbdev vector max,
82d86ed7fbStbbdev int samples,
83d86ed7fbStbbdev flt (*evaluator)(flt, flt, flt)) {
84d86ed7fbStbbdev extvol *xvol;
85d86ed7fbStbbdev texture *tex;
86d86ed7fbStbbdev
87d86ed7fbStbbdev tex = (texture *)voidtex;
88d86ed7fbStbbdev
89d86ed7fbStbbdev xvol = (extvol *)rt_getmem(sizeof(extvol));
90d86ed7fbStbbdev memset(xvol, 0, sizeof(extvol));
91d86ed7fbStbbdev
92d86ed7fbStbbdev xvol->methods = &extvol_methods;
93d86ed7fbStbbdev
94d86ed7fbStbbdev xvol->min = min;
95d86ed7fbStbbdev xvol->max = max;
96d86ed7fbStbbdev xvol->evaluator = evaluator;
97d86ed7fbStbbdev xvol->ambient = tex->ambient;
98d86ed7fbStbbdev xvol->diffuse = tex->diffuse;
99d86ed7fbStbbdev xvol->opacity = tex->opacity;
100d86ed7fbStbbdev xvol->samples = samples;
101d86ed7fbStbbdev
102d86ed7fbStbbdev xvol->tex = (texture *)rt_getmem(sizeof(texture));
103d86ed7fbStbbdev memset(xvol->tex, 0, sizeof(texture));
104d86ed7fbStbbdev
105d86ed7fbStbbdev xvol->tex->ctr.x = 0.0;
106d86ed7fbStbbdev xvol->tex->ctr.y = 0.0;
107d86ed7fbStbbdev xvol->tex->ctr.z = 0.0;
108d86ed7fbStbbdev xvol->tex->rot = xvol->tex->ctr;
109d86ed7fbStbbdev xvol->tex->scale = xvol->tex->ctr;
110d86ed7fbStbbdev xvol->tex->uaxs = xvol->tex->ctr;
111d86ed7fbStbbdev xvol->tex->vaxs = xvol->tex->ctr;
112d86ed7fbStbbdev xvol->tex->islight = 0;
113d86ed7fbStbbdev xvol->tex->shadowcast = 0;
114d86ed7fbStbbdev
115d86ed7fbStbbdev xvol->tex->col = tex->col;
116d86ed7fbStbbdev xvol->tex->ambient = 1.0;
117d86ed7fbStbbdev xvol->tex->diffuse = 0.0;
118d86ed7fbStbbdev xvol->tex->specular = 0.0;
119d86ed7fbStbbdev xvol->tex->opacity = 1.0;
120d86ed7fbStbbdev xvol->tex->img = nullptr;
121d86ed7fbStbbdev xvol->tex->texfunc = (color(*)(void *, void *, void *))(ext_volume_texture);
122d86ed7fbStbbdev xvol->tex->obj = (void *)xvol; /* XXX hack! */
123d86ed7fbStbbdev
124d86ed7fbStbbdev return xvol;
125d86ed7fbStbbdev }
126d86ed7fbStbbdev
ExtVoxelColor(flt scalar)127d86ed7fbStbbdev color ExtVoxelColor(flt scalar) {
128d86ed7fbStbbdev color col;
129d86ed7fbStbbdev
130d86ed7fbStbbdev if (scalar > 1.0)
131d86ed7fbStbbdev scalar = 1.0;
132d86ed7fbStbbdev
133d86ed7fbStbbdev if (scalar < 0.0)
134d86ed7fbStbbdev scalar = 0.0;
135d86ed7fbStbbdev
136d86ed7fbStbbdev if (scalar < 0.5) {
137d86ed7fbStbbdev col.g = 0.0;
138d86ed7fbStbbdev }
139d86ed7fbStbbdev else {
140d86ed7fbStbbdev col.g = (scalar - 0.5) * 2.0;
141d86ed7fbStbbdev }
142d86ed7fbStbbdev
143d86ed7fbStbbdev col.r = scalar;
144d86ed7fbStbbdev col.b = 1.0 - (scalar / 2.0);
145d86ed7fbStbbdev
146d86ed7fbStbbdev return col;
147d86ed7fbStbbdev }
148d86ed7fbStbbdev
ext_volume_texture(vector * hit,texture * tex,ray * ry)149d86ed7fbStbbdev color ext_volume_texture(vector *hit, texture *tex, ray *ry) {
150d86ed7fbStbbdev color col, col2;
151d86ed7fbStbbdev box *bx;
152d86ed7fbStbbdev extvol *xvol;
153d86ed7fbStbbdev flt a, tx1, tx2, ty1, ty2, tz1, tz2;
154d86ed7fbStbbdev flt tnear, tfar;
155d86ed7fbStbbdev flt t, tdist, dt, ddt, sum, tt;
156d86ed7fbStbbdev vector pnt, bln;
157d86ed7fbStbbdev flt scalar, transval;
158d86ed7fbStbbdev int i;
159d86ed7fbStbbdev point_light *li;
160d86ed7fbStbbdev color diffint;
161d86ed7fbStbbdev vector N, L;
162d86ed7fbStbbdev flt inten;
163d86ed7fbStbbdev
164d86ed7fbStbbdev col.r = 0.0;
165d86ed7fbStbbdev col.g = 0.0;
166d86ed7fbStbbdev col.b = 0.0;
167d86ed7fbStbbdev
168d86ed7fbStbbdev bx = (box *)tex->obj;
169d86ed7fbStbbdev xvol = (extvol *)tex->obj;
170d86ed7fbStbbdev
171d86ed7fbStbbdev tnear = -FHUGE;
172d86ed7fbStbbdev tfar = FHUGE;
173d86ed7fbStbbdev
174d86ed7fbStbbdev if (ry->d.x == 0.0) {
175d86ed7fbStbbdev if ((ry->o.x < bx->min.x) || (ry->o.x > bx->max.x))
176d86ed7fbStbbdev return col;
177d86ed7fbStbbdev }
178d86ed7fbStbbdev else {
179d86ed7fbStbbdev tx1 = (bx->min.x - ry->o.x) / ry->d.x;
180d86ed7fbStbbdev tx2 = (bx->max.x - ry->o.x) / ry->d.x;
181d86ed7fbStbbdev if (tx1 > tx2) {
182d86ed7fbStbbdev a = tx1;
183d86ed7fbStbbdev tx1 = tx2;
184d86ed7fbStbbdev tx2 = a;
185d86ed7fbStbbdev }
186d86ed7fbStbbdev if (tx1 > tnear)
187d86ed7fbStbbdev tnear = tx1;
188d86ed7fbStbbdev if (tx2 < tfar)
189d86ed7fbStbbdev tfar = tx2;
190d86ed7fbStbbdev }
191d86ed7fbStbbdev if (tnear > tfar)
192d86ed7fbStbbdev return col;
193d86ed7fbStbbdev if (tfar < 0.0)
194d86ed7fbStbbdev return col;
195d86ed7fbStbbdev
196d86ed7fbStbbdev if (ry->d.y == 0.0) {
197d86ed7fbStbbdev if ((ry->o.y < bx->min.y) || (ry->o.y > bx->max.y))
198d86ed7fbStbbdev return col;
199d86ed7fbStbbdev }
200d86ed7fbStbbdev else {
201d86ed7fbStbbdev ty1 = (bx->min.y - ry->o.y) / ry->d.y;
202d86ed7fbStbbdev ty2 = (bx->max.y - ry->o.y) / ry->d.y;
203d86ed7fbStbbdev if (ty1 > ty2) {
204d86ed7fbStbbdev a = ty1;
205d86ed7fbStbbdev ty1 = ty2;
206d86ed7fbStbbdev ty2 = a;
207d86ed7fbStbbdev }
208d86ed7fbStbbdev if (ty1 > tnear)
209d86ed7fbStbbdev tnear = ty1;
210d86ed7fbStbbdev if (ty2 < tfar)
211d86ed7fbStbbdev tfar = ty2;
212d86ed7fbStbbdev }
213d86ed7fbStbbdev if (tnear > tfar)
214d86ed7fbStbbdev return col;
215d86ed7fbStbbdev if (tfar < 0.0)
216d86ed7fbStbbdev return col;
217d86ed7fbStbbdev
218d86ed7fbStbbdev if (ry->d.z == 0.0) {
219d86ed7fbStbbdev if ((ry->o.z < bx->min.z) || (ry->o.z > bx->max.z))
220d86ed7fbStbbdev return col;
221d86ed7fbStbbdev }
222d86ed7fbStbbdev else {
223d86ed7fbStbbdev tz1 = (bx->min.z - ry->o.z) / ry->d.z;
224d86ed7fbStbbdev tz2 = (bx->max.z - ry->o.z) / ry->d.z;
225d86ed7fbStbbdev if (tz1 > tz2) {
226d86ed7fbStbbdev a = tz1;
227d86ed7fbStbbdev tz1 = tz2;
228d86ed7fbStbbdev tz2 = a;
229d86ed7fbStbbdev }
230d86ed7fbStbbdev if (tz1 > tnear)
231d86ed7fbStbbdev tnear = tz1;
232d86ed7fbStbbdev if (tz2 < tfar)
233d86ed7fbStbbdev tfar = tz2;
234d86ed7fbStbbdev }
235d86ed7fbStbbdev if (tnear > tfar)
236d86ed7fbStbbdev return col;
237d86ed7fbStbbdev if (tfar < 0.0)
238d86ed7fbStbbdev return col;
239d86ed7fbStbbdev
240d86ed7fbStbbdev if (tnear < 0.0)
241d86ed7fbStbbdev tnear = 0.0;
242d86ed7fbStbbdev
243d86ed7fbStbbdev tdist = xvol->samples;
244d86ed7fbStbbdev
245d86ed7fbStbbdev tt = (xvol->opacity / tdist);
246d86ed7fbStbbdev
247d86ed7fbStbbdev bln.x = fabs(bx->min.x - bx->max.x);
248d86ed7fbStbbdev bln.y = fabs(bx->min.y - bx->max.y);
249d86ed7fbStbbdev bln.z = fabs(bx->min.z - bx->max.z);
250d86ed7fbStbbdev
251d86ed7fbStbbdev dt = 1.0 / tdist;
252d86ed7fbStbbdev sum = 0.0;
253d86ed7fbStbbdev
254d86ed7fbStbbdev /* Accumulate color as the ray passes through the voxels */
255d86ed7fbStbbdev for (t = tnear; t <= tfar; t += dt) {
256d86ed7fbStbbdev if (sum < 1.0) {
257d86ed7fbStbbdev pnt.x = ((ry->o.x + (ry->d.x * t)) - bx->min.x) / bln.x;
258d86ed7fbStbbdev pnt.y = ((ry->o.y + (ry->d.y * t)) - bx->min.y) / bln.y;
259d86ed7fbStbbdev pnt.z = ((ry->o.z + (ry->d.z * t)) - bx->min.z) / bln.z;
260d86ed7fbStbbdev
261d86ed7fbStbbdev /* call external evaluator assume 0.0 -> 1.0 range.. */
262d86ed7fbStbbdev scalar = xvol->evaluator(pnt.x, pnt.y, pnt.z);
263d86ed7fbStbbdev
264d86ed7fbStbbdev transval = tt * scalar;
265d86ed7fbStbbdev sum += transval;
266d86ed7fbStbbdev
267d86ed7fbStbbdev col2 = ExtVoxelColor(scalar);
268d86ed7fbStbbdev
269d86ed7fbStbbdev col.r += transval * col2.r * xvol->ambient;
270d86ed7fbStbbdev col.g += transval * col2.g * xvol->ambient;
271d86ed7fbStbbdev col.b += transval * col2.b * xvol->ambient;
272d86ed7fbStbbdev
273d86ed7fbStbbdev ddt = dt;
274d86ed7fbStbbdev
275d86ed7fbStbbdev /* Add in diffuse shaded light sources (no shadows) */
276d86ed7fbStbbdev if (xvol->diffuse > 0.0) {
277d86ed7fbStbbdev /* Calculate the Volume gradient at the voxel */
278d86ed7fbStbbdev N.x = (xvol->evaluator(pnt.x - ddt, pnt.y, pnt.z) -
279d86ed7fbStbbdev xvol->evaluator(pnt.x + ddt, pnt.y, pnt.z)) *
280d86ed7fbStbbdev 8.0 * tt;
281d86ed7fbStbbdev
282d86ed7fbStbbdev N.y = (xvol->evaluator(pnt.x, pnt.y - ddt, pnt.z) -
283d86ed7fbStbbdev xvol->evaluator(pnt.x, pnt.y + ddt, pnt.z)) *
284d86ed7fbStbbdev 8.0 * tt;
285d86ed7fbStbbdev
286d86ed7fbStbbdev N.z = (xvol->evaluator(pnt.x, pnt.y, pnt.z - ddt) -
287d86ed7fbStbbdev xvol->evaluator(pnt.x, pnt.y, pnt.z + ddt)) *
288d86ed7fbStbbdev 8.0 * tt;
289d86ed7fbStbbdev
290d86ed7fbStbbdev /* only light surfaces with enough of a normal.. */
291d86ed7fbStbbdev if ((N.x * N.x + N.y * N.y + N.z * N.z) > 0.0) {
292d86ed7fbStbbdev diffint.r = 0.0;
293d86ed7fbStbbdev diffint.g = 0.0;
294d86ed7fbStbbdev diffint.b = 0.0;
295d86ed7fbStbbdev
296d86ed7fbStbbdev /* add the contribution of each of the lights.. */
297d86ed7fbStbbdev for (i = 0; i < numlights; i++) {
298d86ed7fbStbbdev li = lightlist[i];
299d86ed7fbStbbdev VSUB(li->ctr, (*hit), L)
300d86ed7fbStbbdev VNorm(&L);
301d86ed7fbStbbdev VDOT(inten, N, L)
302d86ed7fbStbbdev
303d86ed7fbStbbdev /* only add light if its from the front of the surface */
304d86ed7fbStbbdev /* could add back-lighting if we wanted to later.. */
305d86ed7fbStbbdev if (inten > 0.0) {
306d86ed7fbStbbdev diffint.r += inten * li->tex->col.r;
307d86ed7fbStbbdev diffint.g += inten * li->tex->col.g;
308d86ed7fbStbbdev diffint.b += inten * li->tex->col.b;
309d86ed7fbStbbdev }
310d86ed7fbStbbdev }
311d86ed7fbStbbdev col.r += col2.r * diffint.r * xvol->diffuse;
312d86ed7fbStbbdev col.g += col2.g * diffint.g * xvol->diffuse;
313d86ed7fbStbbdev col.b += col2.b * diffint.b * xvol->diffuse;
314d86ed7fbStbbdev }
315d86ed7fbStbbdev }
316d86ed7fbStbbdev }
317d86ed7fbStbbdev else {
318d86ed7fbStbbdev sum = 1.0;
319d86ed7fbStbbdev }
320d86ed7fbStbbdev }
321d86ed7fbStbbdev
322d86ed7fbStbbdev /* Add in transmitted ray from outside environment */
323d86ed7fbStbbdev if (sum < 1.0) { /* spawn transmission rays / refraction */
324d86ed7fbStbbdev color transcol;
325d86ed7fbStbbdev
326d86ed7fbStbbdev transcol = shade_transmission(ry, hit, 1.0 - sum);
327d86ed7fbStbbdev
328d86ed7fbStbbdev col.r += transcol.r; /* add the transmitted ray */
329d86ed7fbStbbdev col.g += transcol.g; /* to the diffuse and */
330d86ed7fbStbbdev col.b += transcol.b; /* transmission total.. */
331d86ed7fbStbbdev }
332d86ed7fbStbbdev
333d86ed7fbStbbdev return col;
334d86ed7fbStbbdev }
335