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., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, 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, Ra˙l Fern·ndez Hern·ndez (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_math.h"
  38 #include "BLI_threads.h"
  39 #include "BLI_voxel.h"
  40
  41 #include "PIL_time.h"
  42
  43 #include "RE_shader_ext.h"
  44 #include "RE_raytrace.h"
  45
  46 #include "DNA_material_types.h"
  47
  48 #include "render_types.h"
  49 #include "rendercore.h"
  50 #include "renderdatabase.h"
  51 #include "volumetric.h"
  52 #include "volume_precache.h"
  53
  54 #if defined( _MSC_VER ) && !defined( __cplusplus )
  55 # define inline __inline
  56 #endif // defined( _MSC_VER ) && !defined( __cplusplus )
  57
  58 #include "BKE_global.h"
  59
  60 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
  61 /* defined in pipeline.c, is hardcopy of active dynamic allocated Render */
  62 /* only to be used here in this file, it's for speed */
  63 extern struct Render R;
  64 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
  65
  66 /* *** utility code to set up an individual raytree for objectinstance, for checking inside/outside *** */
  67
  68 /* Recursive test for intersections, from a point inside the mesh, to outside
  69  * Number of intersections (depth) determine if a point is inside or outside the mesh */
  70 int intersect_outside_volume(RayObject *tree, Isect *isect, float *offset, int limit, int depth)
  71 {
  72         if (limit == 0) return depth;
  73
  74         if (RE_rayobject_raycast(tree, isect)) {
  75
  76                 isect->start[0] = isect->start[0] + isect->labda*isect->vec[0];
  77                 isect->start[1] = isect->start[1] + isect->labda*isect->vec[1];
  78                 isect->start[2] = isect->start[2] + isect->labda*isect->vec[2];
  79
  80                 isect->labda = FLT_MAX;
  81                 isect->skip = RE_SKIP_VLR_NEIGHBOUR;
  82                 isect->orig.face= isect->hit.face;
  83                 isect->orig.ob= isect->hit.ob;
  84
  85                 return intersect_outside_volume(tree, isect, offset, limit-1, depth+1);
  86         } else {
  87                 return depth;
  88         }
  89 }
  90
  91 /* Uses ray tracing to check if a point is inside or outside an ObjectInstanceRen */
  92 int point_inside_obi(RayObject *tree, ObjectInstanceRen *obi, float *co)
  93 {
  94         Isect isect;
  95         float vec[3] = {0.0f,0.0f,1.0f};
  96         int final_depth=0, depth=0, limit=20;
  97
  98         /* set up the isect */
  99         memset(&isect, 0, sizeof(isect));
 100         VECCOPY(isect.start, co);
 101         VECCOPY(isect.vec, vec);
 102         isect.mode= RE_RAY_MIRROR;
 103         isect.last_hit= NULL;
 104         isect.lay= -1;
 105
 106         isect.labda = FLT_MAX;
 107         isect.orig.face= NULL;
 108         isect.orig.ob = NULL;
 109
 110         final_depth = intersect_outside_volume(tree, &isect, vec, limit, depth);
 111
 112         /* even number of intersections: point is outside
 113          * odd number: point is inside */
 114         if (final_depth % 2 == 0) return 0;
 115         else return 1;
 116 }
 117
 118 /* find the bounding box of an objectinstance in global space */
 119 void global_bounds_obi(Render *re, ObjectInstanceRen *obi, float *bbmin, float *bbmax)
 120 {
 121         ObjectRen *obr = obi->obr;
 122         VolumePrecache *vp = obi->volume_precache;
 123         VertRen *ver= NULL;
 124         float co[3];
 125         int a;
 126
 127         if (vp->bbmin != NULL && vp->bbmax != NULL) {
 128                 copy_v3_v3(bbmin, vp->bbmin);
 129                 copy_v3_v3(bbmax, vp->bbmax);
 130                 return;
 131         }
 132
 133         vp->bbmin = MEM_callocN(sizeof(float)*3, "volume precache min boundbox corner");
 134         vp->bbmax = MEM_callocN(sizeof(float)*3, "volume precache max boundbox corner");
 135
 136         INIT_MINMAX(bbmin, bbmax);
 137
 138         for(a=0; a<obr->totvert; a++) {
 139                 if((a & 255)==0) ver= obr->vertnodes[a>>8].vert;
 140                 else ver++;
 141
 142                 copy_v3_v3(co, ver->co);
 143
 144                 /* transformed object instance in camera space */
 145                 if(obi->flag & R_TRANSFORMED)
 146                         mul_m4_v3(obi->mat, co);
 147
 148                 /* convert to global space */
 149                 mul_m4_v3(re->viewinv, co);
 150
 151                 DO_MINMAX(co, vp->bbmin, vp->bbmax);
 152         }
 153
 154         copy_v3_v3(bbmin, vp->bbmin);
 155         copy_v3_v3(bbmax, vp->bbmax);
 156
 157 }
 158
 159 /* *** light cache filtering *** */
 160
 161 static float get_avg_surrounds(float *cache, int *res, int xx, int yy, int zz)
 162 {
 163         int x, y, z, x_, y_, z_;
 164         int added=0;
 165         float tot=0.0f;
 166
 167         for (z=-1; z <= 1; z++) {
 168                 z_ = zz+z;
 169                 if (z_ >= 0 && z_ <= res[2]-1) {
 170
 171                         for (y=-1; y <= 1; y++) {
 172                                 y_ = yy+y;
 173                                 if (y_ >= 0 && y_ <= res[1]-1) {
 174
 175                                         for (x=-1; x <= 1; x++) {
 176                                                 x_ = xx+x;
 177                                                 if (x_ >= 0 && x_ <= res[0]-1) {
 178                                                         const int i= V_I(x_, y_, z_, res);
 179
 180                                                         if (cache[i] > 0.0f) {
 181                                                                 tot += cache[i];
 182                                                                 added++;
 183                                                         }
 184
 185                                                 }
 186                                         }
 187                                 }
 188                         }
 189                 }
 190         }
 191
 192         if (added > 0) tot /= added;
 193
 194         return tot;
 195 }
 196
 197 /* function to filter the edges of the light cache, where there was no volume originally.
 198  * For each voxel which was originally external to the mesh, it finds the average values of
 199  * the surrounding internal voxels and sets the original external voxel to that average amount.
 200  * Works almost a bit like a 'dilate' filter */
 201 static void lightcache_filter(VolumePrecache *vp)
 202 {
 203         int x, y, z;
 204
 205         for (z=0; z < vp->res[2]; z++) {
 206                 for (y=0; y < vp->res[1]; y++) {
 207                         for (x=0; x < vp->res[0]; x++) {
 208                                 /* trigger for outside mesh */
 209                                 const int i= V_I(x, y, z, vp->res);
 210
 211                                 if (vp->data_r[i] < -0.f)
 212                                         vp->data_r[i] = get_avg_surrounds(vp->data_r, vp->res, x, y, z);
 213                                 if (vp->data_g[i] < -0.f)
 214                                         vp->data_g[i] = get_avg_surrounds(vp->data_g, vp->res, x, y, z);
 215                                 if (vp->data_b[i] < -0.f)
 216                                         vp->data_b[i] = get_avg_surrounds(vp->data_b, vp->res, x, y, z);
 217                         }
 218                 }
 219         }
 220 }
 221
 222 #if 0
 223 static void lightcache_filter2(VolumePrecache *vp)
 224 {
 225         int x, y, z;
 226         float *new_r, *new_g, *new_b;
 227         int field_size = vp->res[0]*vp->res[1]*vp->res[2]*sizeof(float);
 228
 229         new_r = MEM_mallocN(field_size, "temp buffer for light cache filter r channel");
 230         new_g = MEM_mallocN(field_size, "temp buffer for light cache filter g channel");
 231         new_b = MEM_mallocN(field_size, "temp buffer for light cache filter b channel");
 232
 233         memcpy(new_r, vp->data_r, field_size);
 234         memcpy(new_g, vp->data_g, field_size);
 235         memcpy(new_b, vp->data_b, field_size);
 236
 237         for (z=0; z < vp->res[2]; z++) {
 238                 for (y=0; y < vp->res[1]; y++) {
 239                         for (x=0; x < vp->res[0]; x++) {
 240                                 /* trigger for outside mesh */
 241                                 const int i= V_I(x, y, z, vp->res);
 242                                 if (vp->data_r[i] < -0.f)
 243                                         new_r[i] = get_avg_surrounds(vp->data_r, vp->res, x, y, z);
 244                                 if (vp->data_g[i] < -0.f)
 245                                         new_g[i] = get_avg_surrounds(vp->data_g, vp->res, x, y, z);
 246                                 if (vp->data_b[i] < -0.f)
 247                                         new_b[i] = get_avg_surrounds(vp->data_b, vp->res, x, y, z);
 248                         }
 249                 }
 250         }
 251
 252         SWAP(float *, vp->data_r, new_r);
 253         SWAP(float *, vp->data_g, new_g);
 254         SWAP(float *, vp->data_b, new_b);
 255
 256         if (new_r) { MEM_freeN(new_r); new_r=NULL; }
 257         if (new_g) { MEM_freeN(new_g); new_g=NULL; }
 258         if (new_b) { MEM_freeN(new_b); new_b=NULL; }
 259 }
 260 #endif
 261
 262 static inline int ms_I(int x, int y, int z, int *n) //has a pad of 1 voxel surrounding the core for boundary simulation
 263 {
 264         /* different ordering to light cache */
 265         return x*(n[1]+2)*(n[2]+2) + y*(n[2]+2) + z;
 266 }
 267
 268 static inline int v_I_pad(int x, int y, int z, int *n) //has a pad of 1 voxel surrounding the core for boundary simulation
 269 {
 270         /* same ordering to light cache, with padding */
 271         return z*(n[1]+2)*(n[0]+2) + y*(n[0]+2) + x;
 272 }
 273
 274 static inline int lc_to_ms_I(int x, int y, int z, int *n)
 275 {
 276         /* converting light cache index to multiple scattering index */
 277         return (x-1)*(n[1]*n[2]) + (y-1)*(n[2]) + z-1;
 278 }
 279
 280 /* *** multiple scattering approximation *** */
 281
 282 /* get the total amount of light energy in the light cache. used to normalise after multiple scattering */
 283 static float total_ss_energy(VolumePrecache *vp)
 284 {
 285         int x, y, z;
 286         int *res = vp->res;
 287         float energy=0.f;
 288
 289         for (z=0; z < res[2]; z++) {
 290                 for (y=0; y < res[1]; y++) {
 291                         for (x=0; x < res[0]; x++) {
 292                                 const int i=V_I(x, y, z, res);
 293
 294                                 if (vp->data_r[i] > 0.f) energy += vp->data_r[i];
 295                                 if (vp->data_g[i] > 0.f) energy += vp->data_g[i];
 296                                 if (vp->data_b[i] > 0.f) energy += vp->data_b[i];
 297                         }
 298                 }
 299         }
 300
 301         return energy;
 302 }
 303
 304 static float total_ms_energy(float *sr, float *sg, float *sb, int *res)
 305 {
 306         int x, y, z;
 307         float energy=0.f;
 308
 309         for (z=1;z<=res[2];z++) {
 310                 for (y=1;y<=res[1];y++) {
 311                         for (x=1;x<=res[0];x++) {
 312                                 const int i = ms_I(x,y,z,res);
 313
 314                                 if (sr[i] > 0.f) energy += sr[i];
 315                                 if (sg[i] > 0.f) energy += sg[i];
 316                                 if (sb[i] > 0.f) energy += sb[i];
 317                         }
 318                 }
 319         }
 320
 321         return energy;
 322 }
 323
 324 static void ms_diffuse(float *x0, float *x, float diff, int *n) //n is the unpadded resolution
 325 {
 326         int i, j, k, l;
 327         const float dt = VOL_MS_TIMESTEP;
 328         const float a = dt*diff*n[0]*n[1]*n[2];
 329
 330         for (l=0; l<20; l++)
 331         {
 332                 for (k=1; k<=n[2]; k++)
 333                 {
 334                         for (j=1; j<=n[1]; j++)
 335                         {
 336                                 for (i=1; i<=n[0]; i++)
 337                                 {
 338                                    x[v_I_pad(i,j,k,n)] = (x0[v_I_pad(i,j,k,n)]) + a*(   x0[v_I_pad(i-1,j,k,n)]+ x0[v_I_pad(i+1,j,k,n)]+ x0[v_I_pad(i,j-1,k,n)]+
 339                                                                                                                                                 x0[v_I_pad(i,j+1,k,n)]+ x0[v_I_pad(i,j,k-1,n)]+x0[v_I_pad(i,j,k+1,n)]
 340                                                                                                                                                 ) / (1+6*a);
 341                                 }
 342                         }
 343                 }
 344         }
 345 }
 346
 347 void multiple_scattering_diffusion(Render *re, VolumePrecache *vp, Material *ma)
 348 {
 349         const float diff = ma->vol.ms_diff * 0.001f;    /* compensate for scaling for a nicer UI range */
 350         const int simframes = (int)(ma->vol.ms_spread * (float)MAX3(vp->res[0], vp->res[1], vp->res[2]));
 351         const int shade_type = ma->vol.shade_type;
 352         float fac = ma->vol.ms_intensity;
 353
 354         int x, y, z, m;
 355         int *n = vp->res;
 356         const int size = (n[0]+2)*(n[1]+2)*(n[2]+2);
 357         double time, lasttime= PIL_check_seconds_timer();
 358         float total;
 359         float c=1.0f;
 360         float origf;    /* factor for blending in original light cache */
 361         float energy_ss, energy_ms;
 362
 363         float *sr0=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 364         float *sr=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 365         float *sg0=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 366         float *sg=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 367         float *sb0=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 368         float *sb=(float *)MEM_callocN(size*sizeof(float), "temporary multiple scattering buffer");
 369
 370         total = (float)(n[0]*n[1]*n[2]*simframes);
 371
 372         energy_ss = total_ss_energy(vp);
 373
 374         /* Scattering as diffusion pass */
 375         for (m=0; m<simframes; m++)
 376         {
 377                 /* add sources */
 378                 for (z=1; z<=n[2]; z++)
 379                 {
 380                         for (y=1; y<=n[1]; y++)
 381                         {
 382                                 for (x=1; x<=n[0]; x++)
 383                                 {
 384                                         const int i = lc_to_ms_I(x, y ,z, n);   //lc index
 385                                         const int j = ms_I(x, y, z, n);                 //ms index
 386
 387                                         time= PIL_check_seconds_timer();
 388                                         c++;
 389                                         if (vp->data_r[i] > 0.0f)
 390                                                 sr[j] += vp->data_r[i];
 391                                         if (vp->data_g[i] > 0.0f)
 392                                                 sg[j] += vp->data_g[i];
 393                                         if (vp->data_b[i] > 0.0f)
 394                                                 sb[j] += vp->data_b[i];
 395
 396                                         /* Displays progress every second */
 397                                         if(time-lasttime>1.0f) {
 398                                                 char str[64];
 399                                                 sprintf(str, "Simulating multiple scattering: %d%%", (int)(100.0f * (c / total)));
 400                                                 re->i.infostr= str;
 401                                                 re->stats_draw(re->sdh, &re->i);
 402                                                 re->i.infostr= NULL;
 403                                                 lasttime= time;
 404                                         }
 405                                 }
 406                         }
 407                 }
 408                 SWAP(float *,sr,sr0);
 409                 SWAP(float *,sg,sg0);
 410                 SWAP(float *,sb,sb0);
 411
 412                 /* main diffusion simulation */
 413                 ms_diffuse(sr0, sr, diff, n);
 414                 ms_diffuse(sg0, sg, diff, n);
 415                 ms_diffuse(sb0, sb, diff, n);
 416
 417                 if (re->test_break(re->tbh)) break;
 418         }
 419
 420         /* normalisation factor to conserve energy */
 421         energy_ms = total_ms_energy(sr, sg, sb, n);
 422         fac *= (energy_ss / energy_ms);
 423
 424         /* blend multiple scattering back in the light cache */
 425         if (shade_type == MA_VOL_SHADE_SHADEDPLUSMULTIPLE) {
 426                 /* conserve energy - half single, half multiple */
 427                 origf = 0.5f;
 428                 fac *= 0.5f;
 429         } else {
 430                 origf = 0.0f;
 431         }
 432
 433         for (z=1;z<=n[2];z++)
 434         {
 435                 for (y=1;y<=n[1];y++)
 436                 {
 437                         for (x=1;x<=n[0];x++)
 438                         {
 439                                 const int i = lc_to_ms_I(x, y ,z, n);   //lc index
 440                                 const int j = ms_I(x, y, z, n);                 //ms index
 441
 442                                 vp->data_r[i] = origf * vp->data_r[i] + fac * sr[j];
 443                                 vp->data_g[i] = origf * vp->data_g[i] + fac * sg[j];
 444                                 vp->data_b[i] = origf * vp->data_b[i] + fac * sb[j];
 445                         }
 446                 }
 447         }
 448
 449         MEM_freeN(sr0);
 450         MEM_freeN(sr);
 451         MEM_freeN(sg0);
 452         MEM_freeN(sg);
 453         MEM_freeN(sb0);
 454         MEM_freeN(sb);
 455 }
 456
 457
 458
 459 #if 0 // debug stuff
 460 static void *vol_precache_part_test(void *data)
 461 {
 462         VolPrecachePart *pa = data;
 463
 464         printf("part number: %d \n", pa->num);
 465         printf("done: %d \n", pa->done);
 466         printf("x min: %d   x max: %d \n", pa->minx, pa->maxx);
 467         printf("y min: %d   y max: %d \n", pa->miny, pa->maxy);
 468         printf("z min: %d   z max: %d \n", pa->minz, pa->maxz);
 469
 470         return NULL;
 471 }
 472 #endif
 473
 474 /* Iterate over the 3d voxel grid, and fill the voxels with scattering information
 475  *
 476  * It's stored in memory as 3 big float grids next to each other, one for each RGB channel.
 477  * I'm guessing the memory alignment may work out better this way for the purposes
 478  * of doing linear interpolation, but I haven't actually tested this theory! :)
 479  */
 480 static void *vol_precache_part(void *data)
 481 {
 482         VolPrecachePart *pa =  (VolPrecachePart *)data;
 483         ObjectInstanceRen *obi = pa->obi;
 484         RayObject *tree = pa->tree;
 485         ShadeInput *shi = pa->shi;
 486         float scatter_col[3] = {0.f, 0.f, 0.f};
 487         float co[3], cco[3];
 488         int x, y, z, i;
 489         const int res[3]= {pa->res[0], pa->res[1], pa->res[2]};
 490
 491         for (z= pa->minz; z < pa->maxz; z++) {
 492                 co[2] = pa->bbmin[2] + (pa->voxel[2] * (z + 0.5f));
 493
 494                 for (y= pa->miny; y < pa->maxy; y++) {
 495                         co[1] = pa->bbmin[1] + (pa->voxel[1] * (y + 0.5f));
 496
 497                         for (x=pa->minx; x < pa->maxx; x++) {
 498                                 co[0] = pa->bbmin[0] + (pa->voxel[0] * (x + 0.5f));
 499
 500                                 /* convert from world->camera space for shading */
 501                                 mul_v3_m4v3(cco, pa->viewmat, co);
 502
 503                                 i= V_I(x, y, z, res);
 504
 505                                 // don't bother if the point is not inside the volume mesh
 506                                 if (!point_inside_obi(tree, obi, cco)) {
 507                                         obi->volume_precache->data_r[i] = -1.0f;
 508                                         obi->volume_precache->data_g[i] = -1.0f;
 509                                         obi->volume_precache->data_b[i] = -1.0f;
 510                                         continue;
 511                                 }
 512
 513                                 /* this view coordinate is very wrong! */
 514                                 copy_v3_v3(shi->view, cco);
 515                                 normalize_v3(shi->view);
 516                                 vol_get_scattering(shi, scatter_col, cco);
 517
 518                                 obi->volume_precache->data_r[i] = scatter_col[0];
 519                                 obi->volume_precache->data_g[i] = scatter_col[1];
 520                                 obi->volume_precache->data_b[i] = scatter_col[2];
 521
 522                         }
 523                 }
 524         }
 525
 526         pa->done = 1;
 527
 528         return 0;
 529 }
 530
 531
 532 static void precache_setup_shadeinput(Render *re, ObjectInstanceRen *obi, Material *ma, ShadeInput *shi)
 533 {
 534         memset(shi, 0, sizeof(ShadeInput));
 535         shi->depth= 1;
 536         shi->mask= 1;
 537         shi->mat = ma;
 538         shi->vlr = NULL;
 539         memcpy(&shi->r, &shi->mat->r, 23*sizeof(float));        // note, keep this synced with render_types.h
 540         shi->har= shi->mat->har;
 541         shi->obi= obi;
 542         shi->obr= obi->obr;
 543         shi->lay = re->lay;
 544 }
 545
 546 static void precache_init_parts(Render *re, RayObject *tree, ShadeInput *shi, ObjectInstanceRen *obi, int totthread, int *parts)
 547 {
 548         VolumePrecache *vp = obi->volume_precache;
 549         int i=0, x, y, z;
 550         float voxel[3];
 551         int sizex, sizey, sizez;
 552         float bbmin[3], bbmax[3];
 553         int *res;
 554         int minx, maxx;
 555         int miny, maxy;
 556         int minz, maxz;
 557
 558         if (!vp) return;
 559
 560         BLI_freelistN(&re->volume_precache_parts);
 561
 562         /* currently we just subdivide the box, number of threads per side */
 563         parts[0] = parts[1] = parts[2] = totthread;
 564         res = vp->res;
 565
 566         /* using boundbox in worldspace */
 567         global_bounds_obi(re, obi, bbmin, bbmax);
 568         sub_v3_v3v3(voxel, bbmax, bbmin);
 569
 570         voxel[0] /= (float)res[0];
 571         voxel[1] /= (float)res[1];
 572         voxel[2] /= (float)res[2];
 573
 574         for (x=0; x < parts[0]; x++) {
 575                 sizex = ceil(res[0] / (float)parts[0]);
 576                 minx = x * sizex;
 577                 maxx = minx + sizex;
 578                 maxx = (maxx>res[0])?res[0]:maxx;
 579
 580                 for (y=0; y < parts[1]; y++) {
 581                         sizey = ceil(res[1] / (float)parts[1]);
 582                         miny = y * sizey;
 583                         maxy = miny + sizey;
 584                         maxy = (maxy>res[1])?res[1]:maxy;
 585
 586                         for (z=0; z < parts[2]; z++) {
 587                                 VolPrecachePart *pa= MEM_callocN(sizeof(VolPrecachePart), "new precache part");
 588
 589                                 sizez = ceil(res[2] / (float)parts[2]);
 590                                 minz = z * sizez;
 591                                 maxz = minz + sizez;
 592                                 maxz = (maxz>res[2])?res[2]:maxz;
 593
 594                                 pa->done = 0;
 595                                 pa->working = 0;
 596
 597                                 pa->num = i;
 598                                 pa->tree = tree;
 599                                 pa->shi = shi;
 600                                 pa->obi = obi;
 601                                 copy_m4_m4(pa->viewmat, re->viewmat);
 602
 603                                 copy_v3_v3(pa->bbmin, bbmin);
 604                                 copy_v3_v3(pa->voxel, voxel);
 605                                 VECCOPY(pa->res, res);
 606
 607                                 pa->minx = minx; pa->maxx = maxx;
 608                                 pa->miny = miny; pa->maxy = maxy;
 609                                 pa->minz = minz; pa->maxz = maxz;
 610
 611
 612                                 BLI_addtail(&re->volume_precache_parts, pa);
 613
 614                                 i++;
 615                         }
 616                 }
 617         }
 618 }
 619
 620 static VolPrecachePart *precache_get_new_part(Render *re)
 621 {
 622         VolPrecachePart *pa, *nextpa=NULL;
 623
 624         for (pa = re->volume_precache_parts.first; pa; pa=pa->next)
 625         {
 626                 if (pa->done==0 && pa->working==0) {
 627                         nextpa = pa;
 628                         break;
 629                 }
 630         }
 631
 632         return nextpa;
 633 }
 634
 635 /* calculate resolution from bounding box in world space */
 636 static int precache_resolution(Render *re, VolumePrecache *vp, ObjectInstanceRen *obi, int res)
 637 {
 638         float dim[3], div;
 639         float bbmin[3], bbmax[3];
 640
 641         /* bound box in global space */
 642         global_bounds_obi(re, obi, bbmin, bbmax);
 643         sub_v3_v3v3(dim, bbmax, bbmin);
 644
 645         div = MAX3(dim[0], dim[1], dim[2]);
 646         dim[0] /= div;
 647         dim[1] /= div;
 648         dim[2] /= div;
 649
 650         vp->res[0] = ceil(dim[0] * res);
 651         vp->res[1] = ceil(dim[1] * res);
 652         vp->res[2] = ceil(dim[2] * res);
 653
 654         if ((vp->res[0] < 1) || (vp->res[1] < 1) || (vp->res[2] < 1))
 655                 return 0;
 656
 657         return 1;
 658 }
 659
 660 /* Precache a volume into a 3D voxel grid.
 661  * The voxel grid is stored in the ObjectInstanceRen,
 662  * in camera space, aligned with the ObjectRen's bounding box.
 663  * Resolution is defined by the user.
 664  */
 665 void vol_precache_objectinstance_threads(Render *re, ObjectInstanceRen *obi, Material *ma)
 666 {
 667         VolumePrecache *vp;
 668         VolPrecachePart *nextpa, *pa;
 669         RayObject *tree;
 670         ShadeInput shi;
 671         ListBase threads;
 672         int parts[3] = {1, 1, 1}, totparts;
 673
 674         int caching=1, counter=0;
 675         int totthread = re->r.threads;
 676
 677         double time, lasttime= PIL_check_seconds_timer();
 678
 679         R = *re;
 680
 681         /* create a raytree with just the faces of the instanced ObjectRen,
 682          * used for checking if the cached point is inside or outside. */
 683         tree = makeraytree_object(&R, obi);
 684         if (!tree) return;
 685
 686         vp = MEM_callocN(sizeof(VolumePrecache), "volume light cache");
 687         obi->volume_precache = vp;
 688
 689         if (!precache_resolution(re, vp, obi, ma->vol.precache_resolution)) {
 690                 MEM_freeN(vp);
 691                 vp = NULL;
 692                 return;
 693         }
 694
 695         vp->data_r = MEM_callocN(sizeof(float)*vp->res[0]*vp->res[1]*vp->res[2], "volume light cache data red channel");
 696         vp->data_g = MEM_callocN(sizeof(float)*vp->res[0]*vp->res[1]*vp->res[2], "volume light cache data green channel");
 697         vp->data_b = MEM_callocN(sizeof(float)*vp->res[0]*vp->res[1]*vp->res[2], "volume light cache data blue channel");
 698         if (vp->data_r==0 || vp->data_g==0 || vp->data_b==0) {
 699                 MEM_freeN(vp);
 700                 vp = NULL;
 701                 return;
 702         }
 703         //obi->volume_precache = vp;
 704
 705         /* Need a shadeinput to calculate scattering */
 706         precache_setup_shadeinput(re, obi, ma, &shi);
 707
 708         precache_init_parts(re, tree, &shi, obi, totthread, parts);
 709         totparts = parts[0] * parts[1] * parts[2];
 710
 711         BLI_init_threads(&threads, vol_precache_part, totthread);
 712
 713         while(caching) {
 714
 715                 if(BLI_available_threads(&threads) && !(re->test_break(re->tbh))) {
 716                         nextpa = precache_get_new_part(re);
 717                         if (nextpa) {
 718                                 nextpa->working = 1;
 719                                 BLI_insert_thread(&threads, nextpa);
 720                         }
 721                 }
 722                 else PIL_sleep_ms(50);
 723
 724                 caching=0;
 725                 counter=0;
 726                 for(pa= re->volume_precache_parts.first; pa; pa= pa->next) {
 727
 728                         if(pa->done) {
 729                                 counter++;
 730                                 BLI_remove_thread(&threads, pa);
 731                         } else
 732                                 caching = 1;
 733                 }
 734
 735                 if (re->test_break(re->tbh) && BLI_available_threads(&threads)==totthread)
 736                         caching=0;
 737
 738                 time= PIL_check_seconds_timer();
 739                 if(time-lasttime>1.0f) {
 740                         char str[64];
 741                         sprintf(str, "Precaching volume: %d%%", (int)(100.0f * ((float)counter / (float)totparts)));
 742                         re->i.infostr= str;
 743                         re->stats_draw(re->sdh, &re->i);
 744                         re->i.infostr= NULL;
 745                         lasttime= time;
 746                 }
 747         }
 748
 749         BLI_end_threads(&threads);
 750         BLI_freelistN(&re->volume_precache_parts);
 751
 752         if(tree) {
 753                 //TODO: makeraytree_object creates a tree and saves it on OBI, if we free this tree we should also clear other pointers to it
 754                 //RE_rayobject_free(tree);
 755                 //tree= NULL;
 756         }
 757
 758         if (ELEM(ma->vol.shade_type, MA_VOL_SHADE_MULTIPLE, MA_VOL_SHADE_SHADEDPLUSMULTIPLE))
 759         {
 760                 /* this should be before the filtering */
 761                 multiple_scattering_diffusion(re, obi->volume_precache, ma);
 762         }
 763
 764         lightcache_filter(obi->volume_precache);
 765 }
 766
 767 static int using_lightcache(Material *ma)
 768 {
 769         return (((ma->vol.shadeflag & MA_VOL_PRECACHESHADING) && (ma->vol.shade_type == MA_VOL_SHADE_SHADED))
 770                 || (ELEM(ma->vol.shade_type, MA_VOL_SHADE_MULTIPLE, MA_VOL_SHADE_SHADEDPLUSMULTIPLE)));
 771 }
 772
 773 /* loop through all objects (and their associated materials)
 774  * marked for pre-caching in convertblender.c, and pre-cache them */
 775 void volume_precache(Render *re)
 776 {
 777         ObjectInstanceRen *obi;
 778         VolumeOb *vo;
 779
 780         for(vo= re->volumes.first; vo; vo= vo->next) {
 781                 if (using_lightcache(vo->ma)) {
 782                         for(obi= re->instancetable.first; obi; obi= obi->next) {
 783                                 if (obi->obr == vo->obr) {
 784                                         vol_precache_objectinstance_threads(re, obi, vo->ma);
 785                                 }
 786                         }
 787                 }
 788         }
 789
 790         re->i.infostr= NULL;
 791         re->stats_draw(re->sdh, &re->i);
 792 }
 793
 794 void free_volume_precache(Render *re)
 795 {
 796         ObjectInstanceRen *obi;
 797
 798         for(obi= re->instancetable.first; obi; obi= obi->next) {
 799                 if (obi->volume_precache != NULL) {
 800                         MEM_freeN(obi->volume_precache->data_r);
 801                         MEM_freeN(obi->volume_precache->data_g);
 802                         MEM_freeN(obi->volume_precache->data_b);
 803                         MEM_freeN(obi->volume_precache->bbmin);
 804                         MEM_freeN(obi->volume_precache->bbmax);
 805                         MEM_freeN(obi->volume_precache);
 806                         obi->volume_precache = NULL;
 807                 }
 808         }
 809
 810         BLI_freelistN(&re->volumes);
 811 }
 812
 813 int point_inside_volume_objectinstance(Render *re, ObjectInstanceRen *obi, float *co)
 814 {
 815         RayObject *tree;
 816         int inside=0;
 817
 818         tree = makeraytree_object(re, obi);
 819         if (!tree) return 0;
 820
 821         inside = point_inside_obi(tree, obi, co);
 822
 823         //TODO: makeraytree_object creates a tree and saves it on OBI, if we free this tree we should also clear other pointers to it
 824         //RE_rayobject_free(tree);
 825         //tree= NULL;
 826
 827         return inside;
 828 }
 829