wip commit to work on at home, nothing to see
[blender.git] / source / blender / render / intern / source / volume_precache.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.
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
39 #include "PIL_time.h"
40
41 #include "RE_shader_ext.h"
42 #include "RE_raytrace.h"
43
44 #include "DNA_material_types.h"
45
46 #include "render_types.h"
47 #include "renderdatabase.h"
48 #include "volumetric.h"
49
50 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
51 /* defined in pipeline.c, is hardcopy of active dynamic allocated Render */
52 /* only to be used here in this file, it's for speed */
53 extern struct Render R;
54 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
55
56
57 /* Recursive test for intersections, from a point inside the mesh, to outside
58  * Number of intersections (depth) determine if a point is inside or outside the mesh */
59 int intersect_outside_volume(RayTree *tree, Isect *isect, float *offset, int limit, int depth)
60 {
61         if (limit == 0) return depth;
62         
63         if (RE_ray_tree_intersect(tree, isect)) {
64                 float hitco[3];
65                 
66                 hitco[0] = isect->start[0] + isect->labda*isect->vec[0];
67                 hitco[1] = isect->start[1] + isect->labda*isect->vec[1];
68                 hitco[2] = isect->start[2] + isect->labda*isect->vec[2];
69                 VecAddf(isect->start, hitco, offset);
70
71                 return intersect_outside_volume(tree, isect, offset, limit-1, depth+1);
72         } else {
73                 return depth;
74         }
75 }
76
77 /* Uses ray tracing to check if a point is inside or outside an ObjectInstanceRen */
78 int point_inside_obi(RayTree *tree, ObjectInstanceRen *obi, float *co)
79 {
80         float maxsize = RE_ray_tree_max_size(tree);
81         Isect isect;
82         float vec[3] = {0.0f,0.0f,1.0f};
83         int final_depth=0, depth=0, limit=20;
84         
85         /* set up the isect */
86         memset(&isect, 0, sizeof(isect));
87         VECCOPY(isect.start, co);
88         isect.end[0] = co[0] + vec[0] * maxsize;
89         isect.end[1] = co[1] + vec[1] * maxsize;
90         isect.end[2] = co[2] + vec[2] * maxsize;
91         
92         /* and give it a little offset to prevent self-intersections */
93         VecMulf(vec, 1e-5);
94         VecAddf(isect.start, isect.start, vec);
95         
96         isect.mode= RE_RAY_MIRROR;
97         isect.face_last= NULL;
98         isect.lay= -1;
99         
100         final_depth = intersect_outside_volume(tree, &isect, vec, limit, depth);
101         
102         /* even number of intersections: point is outside
103          * odd number: point is inside */
104         if (final_depth % 2 == 0) return 0;
105         else return 1;
106 }
107
108 static int inside_check_func(Isect *is, int ob, RayFace *face)
109 {
110         return 1;
111 }
112 static void vlr_face_coords(RayFace *face, float **v1, float **v2, float **v3, float **v4)
113 {
114         VlakRen *vlr= (VlakRen*)face;
115
116         *v1 = (vlr->v1)? vlr->v1->co: NULL;
117         *v2 = (vlr->v2)? vlr->v2->co: NULL;
118         *v3 = (vlr->v3)? vlr->v3->co: NULL;
119         *v4 = (vlr->v4)? vlr->v4->co: NULL;
120 }
121
122 RayTree *create_raytree_obi(ObjectInstanceRen *obi, float *bbmin, float *bbmax)
123 {
124         int v;
125         VlakRen *vlr= NULL;
126         
127         /* create empty raytree */
128         RayTree *tree = RE_ray_tree_create(64, obi->obr->totvlak, bbmin, bbmax,
129                 vlr_face_coords, inside_check_func, NULL, NULL);
130         
131         /* fill it with faces */
132         for(v=0; v<obi->obr->totvlak; v++) {
133                 if((v & 255)==0)
134                         vlr= obi->obr->vlaknodes[v>>8].vlak;
135                 else
136                         vlr++;
137         
138                 RE_ray_tree_add_face(tree, 0, vlr);
139         }
140         
141         RE_ray_tree_done(tree);
142         
143         return tree;
144 }
145
146 static float get_avg_surrounds(float *cache, int res, int res_2, int res_3, int rgb, int xx, int yy, int zz)
147 {
148         int x, y, z, x_, y_, z_;
149         int added=0;
150         float tot=0.0f;
151         int i;
152         
153         for (x=-1; x <= 1; x++) {
154                 x_ = xx+x;
155                 if (x_ >= 0 && x_ <= res-1) {
156                 
157                         for (y=-1; y <= 1; y++) {
158                                 y_ = yy+y;
159                                 if (y_ >= 0 && y_ <= res-1) {
160                                 
161                                         for (z=-1; z <= 1; z++) {
162                                                 z_ = zz+z;
163                                                 if (z_ >= 0 && z_ <= res-1) {
164                                                 
165                                                         i = rgb*res_3 + x_*res_2 + y_*res + z_;
166                                                         if (cache[i] > 0.0f) {
167                                                                 tot += cache[i];
168                                                                 added++;
169                                                         }
170                                                         
171                                                 }
172                                         }
173                                 }
174                         }
175                 }
176         }
177         
178         tot /= added;
179         
180         return ((added>0)?tot:0.0f);
181 }
182
183 /* function to filter the edges of the light cache, where there was no volume originally.
184  * For each voxel which was originally external to the mesh, it finds the average values of
185  * the surrounding internal voxels and sets the original external voxel to that average amount.
186  * Works almost a bit like a 'dilate' filter */
187 static void lightcache_filter(float *cache, int res)
188 {
189         int x, y, z, rgb;
190         int res_2, res_3;
191         int i;
192         
193         res_2 = res*res;
194         res_3 = res*res*res;
195
196         for (x=0; x < res; x++) {
197                 for (y=0; y < res; y++) {
198                         for (z=0; z < res; z++) {
199                                 for (rgb=0; rgb < 3; rgb++) {
200                                         i = rgb*res_3 + x*res_2 + y*res + z;
201
202                                         /* trigger for outside mesh */
203                                         if (cache[i] < 0.5f) cache[i] = get_avg_surrounds(cache, res, res_2, res_3, rgb, x, y, z);
204                                 }
205                         }
206                 }
207         }
208 }
209
210 /* Precache a volume into a 3D voxel grid.
211  * The voxel grid is stored in the ObjectInstanceRen, 
212  * in camera space, aligned with the ObjectRen's bounding box.
213  * Resolution is defined by the user.
214  */
215 void vol_precache_objectinstance(Render *re, ObjectInstanceRen *obi, Material *ma, float *bbmin, float *bbmax)
216 {
217         int x, y, z;
218
219         float co[3], voxel[3], scatter_col[3];
220         ShadeInput shi;
221         float view[3] = {0.0,0.0,-1.0};
222         float density;
223         float stepsize;
224         
225         float resf, res_3f;
226         int res_2, res_3;
227         
228         float i = 1.0f;
229         double time, lasttime= PIL_check_seconds_timer();
230         const int res = ma->vol_precache_resolution;
231         RayTree *tree;
232         
233         R = *re;
234
235         /* create a raytree with just the faces of the instanced ObjectRen, 
236          * used for checking if the cached point is inside or outside. */
237         tree = create_raytree_obi(obi, bbmin, bbmax);
238         if (!tree) return;
239
240         /* Need a shadeinput to calculate scattering */
241         memset(&shi, 0, sizeof(ShadeInput)); 
242         shi.depth= 1;
243         shi.mask= 1;
244         shi.mat = ma;
245         shi.vlr = NULL;
246         memcpy(&shi.r, &shi.mat->r, 23*sizeof(float));  // note, keep this synced with render_types.h
247         shi.har= shi.mat->har;
248         shi.obi= obi;
249         shi.obr= obi->obr;
250         shi.lay = re->scene->lay;
251         VECCOPY(shi.view, view);
252         
253         stepsize = vol_get_stepsize(&shi, STEPSIZE_VIEW);
254
255         resf = (float)res;
256         res_2 = res*res;
257         res_3 = res*res*res;
258         res_3f = (float)res_3;
259         
260         VecSubf(voxel, bbmax, bbmin);
261         if ((voxel[0] < FLT_EPSILON) || (voxel[1] < FLT_EPSILON) || (voxel[2] < FLT_EPSILON))
262                 return;
263         VecMulf(voxel, 1.0f/res);
264         
265         obi->volume_precache = MEM_callocN(sizeof(float)*res_3*3, "volume light cache");
266         
267         /* Iterate over the 3d voxel grid, and fill the voxels with scattering information
268          *
269          * It's stored in memory as 3 big float grids next to each other, one for each RGB channel.
270          * I'm guessing the memory alignment may work out better this way for the purposes
271          * of doing linear interpolation, but I haven't actually tested this theory! :)
272          */
273         for (x=0; x < res; x++) {
274                 co[0] = bbmin[0] + (voxel[0] * x);
275                 
276                 for (y=0; y < res; y++) {
277                         co[1] = bbmin[1] + (voxel[1] * y);
278                         
279                         for (z=0; z < res; z++) {
280                                 co[2] = bbmin[2] + (voxel[2] * z);
281                         
282                                 time= PIL_check_seconds_timer();
283                                 i++;
284                         
285                                 /* display progress every second */
286                                 if(re->test_break()) {
287                                         if(tree) {
288                                                 RE_ray_tree_free(tree);
289                                                 tree= NULL;
290                                         }
291                                         return;
292                                 }
293                                 if(time-lasttime>1.0f) {
294                                         char str[64];
295                                         sprintf(str, "Precaching volume: %d%%", (int)(100.0f * (i / res_3f)));
296                                         re->i.infostr= str;
297                                         re->stats_draw(&re->i);
298                                         re->i.infostr= NULL;
299                                         lasttime= time;
300                                 }
301                                 
302                                 /* don't bother if the point is not inside the volume mesh */
303                                 if (!point_inside_obi(tree, obi, co)) {
304                                         obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = -1.0f;
305                                         obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = -1.0f;
306                                         obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = -1.0f;
307                                         continue;
308                                 }
309                                 density = vol_get_density(&shi, co);
310                                 vol_get_scattering(&shi, scatter_col, co, stepsize, density);
311                         
312                                 obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = scatter_col[0];
313                                 obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = scatter_col[1];
314                                 obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = scatter_col[2];
315                         }
316                 }
317         }
318
319         if(tree) {
320                 RE_ray_tree_free(tree);
321                 tree= NULL;
322         }
323         
324         lightcache_filter(obi->volume_precache, res);
325
326 }
327
328 #if 0
329 typedef struct VolPrecachePart {
330         struct VolPrecachePart *next, *prev;
331         int num;
332         int minx, maxx;
333         int miny, maxy;
334         int minz, maxz;
335         int res;
336         float bbmin[3], voxel[3];
337         struct RayTree *tree;
338         struct ShadeInput *shi;
339         struct ObjectInstanceRen *obi;
340         int done;
341 } VolPrecachePart;
342
343 static void *vol_precache_part_test(void *data)
344 {
345         VolPrecachePart *vpt =  (VolPrecachePart *)data;
346
347         printf("part number: %d \n", vpt->num);
348
349         return 0;
350 }
351
352 /* Iterate over the 3d voxel grid, and fill the voxels with scattering information
353  *
354  * It's stored in memory as 3 big float grids next to each other, one for each RGB channel.
355  * I'm guessing the memory alignment may work out better this way for the purposes
356  * of doing linear interpolation, but I haven't actually tested this theory! :)
357  */
358 static void *vol_precache_part(void *data)
359 {
360         VolPrecachePart *vpt =  (VolPrecachePart *)data;
361         ObjectInstanceRen *obi = vpt->obi;
362         RayTree *tree = vpt->tree;
363         ShadeInput *shi = vpt->shi;
364         float scatter_col[3] = {0.f, 0.f, 0.f};
365         float co[3];
366         int x, y, z;
367         const int res=vpt->res, res_2=vpt->res*vpt->res, res_3=vpt->res*vpt->res*vpt->res;
368         const float stepsize = vol_get_stepsize(shi, STEPSIZE_VIEW);
369         
370         res = vpt->res;
371         res_2 = res*res;
372         res_3 = res*res*res;
373         
374         for (x= vpt->minx; x < vpt->maxx; x++) {
375                 co[0] = vpt->bbmin[0] + (vpt->voxel[0] * x);
376                 
377                 for (y= vpt->miny; y < vpt->maxy; y++) {
378                         co[1] = vpt->bbmin[1] + (vpt->voxel[1] * y);
379                         
380                         for (z=vpt->minz; z < vpt->maxz; z++) {
381                                 co[2] = vpt->bbmin[2] + (vpt->voxel[2] * z);
382                         
383                                 // don't bother if the point is not inside the volume mesh
384                                 if (!point_inside_obi(tree, obi, co)) {
385                                         obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = -1.0f;
386                                         obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = -1.0f;
387                                         obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = -1.0f;
388                                         continue;
389                                 }
390                                 density = vol_get_density(shi, co);
391                                 vol_get_scattering(shi, scatter_col, co, stepsize, density);
392                         
393                                 obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = scatter_col[0];
394                                 obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = scatter_col[1];
395                                 obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = scatter_col[2];
396                         }
397                 }
398         }
399         
400         return 0;
401 }
402
403 void precache_setup_shadeinput(Render *re, ObjectInstanceRen *obi, Material *ma, ShadeInput *shi)
404 {
405         float view[3] = {0.0,0.0,-1.0};
406         
407         memset(&shi, 0, sizeof(ShadeInput)); 
408         shi->depth= 1;
409         shi->mask= 1;
410         shi->mat = ma;
411         shi->vlr = NULL;
412         memcpy(&shi->r, &shi->mat->r, 23*sizeof(float));        // note, keep this synced with render_types.h
413         shi->har= shi->mat->har;
414         shi->obi= obi;
415         shi->obr= obi->obr;
416         shi->lay = re->scene->lay;
417         VECCOPY(shi->view, view);
418 }
419
420 void vol_precache_objectinstance_threads(Render *re, ObjectInstanceRen *obi, Material *ma, float *bbmin, float *bbmax)
421 {
422         int x, y, z;
423
424         float co[3], voxel[3], scatter_col[3];
425         ShadeInput shi;
426         
427         float density;
428         float stepsize;
429         
430         float resf, res_3f;
431         int res_2, res_3;
432         
433         int edgeparts=2;
434         int totparts;
435         ListBase threads;
436         int cont= 1;
437         int xparts, yparts, zparts;
438         float part[3];
439         int totthread = re->r.threads;
440         ListBase precache_parts;
441         VolPrecachePart *nextpa;
442         int j;
443         
444         float i = 1.0f;
445         double time, lasttime= PIL_check_seconds_timer();
446         const int res = ma->vol_precache_resolution;
447         RayTree *tree;
448         
449         R = *re;
450
451         /* create a raytree with just the faces of the instanced ObjectRen, 
452          * used for checking if the cached point is inside or outside. */
453         tree = create_raytree_obi(obi, bbmin, bbmax);
454         if (!tree) return;
455
456         /* Need a shadeinput to calculate scattering */
457         precache_setup_shadeinput(re, obi, ma, &shi);
458
459         VecSubf(voxel, bbmax, bbmin);
460         if ((voxel[0] < FLT_EPSILON) || (voxel[1] < FLT_EPSILON) || (voxel[2] < FLT_EPSILON))
461                 return;
462         VecMulf(voxel, 1.0f/res);
463         
464         part[0] = ceil(res/(float)xparts); 
465         part[1] = ceil(res/(float)yparts);
466         part[2] = ceil(res/(float)zparts);
467         
468         obi->volume_precache = MEM_callocN(sizeof(float)*res_3*3, "volume light cache");
469         
470         totparts = edgeparts*edgeparts*edgeparts;
471         precache_parts= MEM_callocN(sizeof(VolPrecachePart)*totparts, "VolPrecachePart");
472         memset(precache_parts, 0, sizeof(VolPrecachePart)*totparts);
473
474         precache_init_parts(precache_parts);
475
476         for(j=0; j < totparts; j++) {
477                 VolPrecachePart *pa= MEM_callocN(sizeof(VolPrecachePart), "new precache part");
478         
479                 pa->done = 0;
480                 pa->num = j;
481                 
482                 pa->res = res;
483                 VECCOPY(pa->bbmin, bbmin);
484                 VECCOPY(precache_parts[j].voxel, voxel);
485                 precache_parts[j].tree = tree;
486                 precache_parts[j].shi = shi;
487                 precache_parts[j].obi = obi;
488         }
489         
490         BLI_init_threads(&threads, vol_precache_part, totthread);
491         
492         nextpa = precache_get_new_part(precache_threads);
493         
494         while(cont) {
495
496                 if(BLI_available_threads(&threads) && !(re->test_break())) {
497                         
498                         precache_get_new_part(
499                         // get new job (data pointer)
500                         for(j=0; j < totparts; j++) {
501                                 if (!precache_threads[j].done) {
502                                         // tag job 'processed
503                                         precache_threads[j].done = 1;
504                                 }
505                         }
506                 
507                         BLI_insert_thread(&threads, precache_get_new_part(precache_threads));
508                 }
509                 else PIL_sleep_ms(50);
510
511                 // find if a job is ready, this the do_something_func() should write in job somewhere
512                 cont= 0;
513                 for(go over all jobs)
514                         if(job is ready) {
515                                 if(job was not removed) {
516                                         BLI_remove_thread(&lb, job);
517                                 }
518                         }
519                         else cont= 1;
520                 }
521                 // conditions to exit loop
522                 if(if escape loop event) {
523                 if(BLI_available_threadslots(&lb)==maxthreads)
524                         break;
525                 }
526         }
527
528         BLI_end_threads(&threads);
529         
530         //
531         
532         /* Iterate over the 3d voxel grid, and fill the voxels with scattering information
533          *
534          * It's stored in memory as 3 big float grids next to each other, one for each RGB channel.
535          * I'm guessing the memory alignment may work out better this way for the purposes
536          * of doing linear interpolation, but I haven't actually tested this theory! :)
537          */
538          /*
539         for (x=0; x < res; x++) {
540                 co[0] = bbmin[0] + (voxel[0] * x);
541                 
542                 for (y=0; y < res; y++) {
543                         co[1] = bbmin[1] + (voxel[1] * y);
544                         
545                         for (z=0; z < res; z++) {
546                                 co[2] = bbmin[2] + (voxel[2] * z);
547                         
548                                 time= PIL_check_seconds_timer();
549                                 i++;
550                         
551                                 // display progress every second
552                                 if(re->test_break()) {
553                                         if(tree) {
554                                                 RE_ray_tree_free(tree);
555                                                 tree= NULL;
556                                         }
557                                         return;
558                                 }
559                                 if(time-lasttime>1.0f) {
560                                         char str[64];
561                                         sprintf(str, "Precaching volume: %d%%", (int)(100.0f * (i / res_3f)));
562                                         re->i.infostr= str;
563                                         re->stats_draw(&re->i);
564                                         re->i.infostr= NULL;
565                                         lasttime= time;
566                                 }
567                                 
568                                 // don't bother if the point is not inside the volume mesh
569                                 
570                                 if (!point_inside_obi(tree, obi, co)) {
571                                         obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = -1.0f;
572                                         obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = -1.0f;
573                                         obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = -1.0f;
574                                         continue;
575                                 }
576                                 density = vol_get_density(&shi, co);
577                                 vol_get_scattering(&shi, scatter_col, co, stepsize, density);
578                         
579                                 obi->volume_precache[0*res_3 + x*res_2 + y*res + z] = scatter_col[0];
580                                 obi->volume_precache[1*res_3 + x*res_2 + y*res + z] = scatter_col[1];
581                                 obi->volume_precache[2*res_3 + x*res_2 + y*res + z] = scatter_col[2];
582                                 
583                         }
584                 }
585         }
586         */
587
588         if(tree) {
589                 RE_ray_tree_free(tree);
590                 tree= NULL;
591         }
592         
593         lightcache_filter(obi->volume_precache, res);
594
595 }
596 #endif
597
598 /* loop through all objects (and their associated materials)
599  * marked for pre-caching in convertblender.c, and pre-cache them */
600 void volume_precache(Render *re)
601 {
602         ObjectInstanceRen *obi;
603         VolPrecache *vp;
604
605         for(vp= re->vol_precache_obs.first; vp; vp= vp->next) {
606                 for(obi= re->instancetable.first; obi; obi= obi->next) {
607                         if (obi->obr == vp->obr) 
608                                 vol_precache_objectinstance(re, obi, vp->ma, obi->obr->boundbox[0], obi->obr->boundbox[1]);
609                 }
610         }
611         
612         re->i.infostr= NULL;
613         re->stats_draw(&re->i);
614 }
615
616 void free_volume_precache(Render *re)
617 {
618         ObjectInstanceRen *obi;
619         
620         for(obi= re->instancetable.first; obi; obi= obi->next) {
621                 if (obi->volume_precache)
622                         MEM_freeN(obi->volume_precache);
623         }
624         
625         BLI_freelistN(&re->vol_precache_obs);
626 }