Workspace: Move engines to workspace and Properties Editor cleanup
[blender.git] / source / blender / render / intern / source / pointdensity.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  * Contributors: Matt Ebb
22  *
23  * ***** END GPL LICENSE BLOCK *****
24  */
25
26 /** \file blender/render/intern/source/pointdensity.c
27  *  \ingroup render
28  */
29
30
31 #include <math.h>
32 #include <stdlib.h>
33 #include <stdio.h>
34
35 #include "MEM_guardedalloc.h"
36
37 #include "BLI_math.h"
38 #include "BLI_blenlib.h"
39 #include "BLI_noise.h"
40 #include "BLI_kdopbvh.h"
41 #include "BLI_utildefines.h"
42 #include "BLI_task.h"
43
44 #include "BLT_translation.h"
45
46 #include "DNA_meshdata_types.h"
47 #include "DNA_object_types.h"
48 #include "DNA_particle_types.h"
49 #include "DNA_texture_types.h"
50
51 #include "BKE_deform.h"
52 #include "BKE_DerivedMesh.h"
53 #include "BKE_lattice.h"
54 #include "BKE_main.h"
55 #include "BKE_object.h"
56 #include "BKE_particle.h"
57 #include "BKE_scene.h"
58 #include "BKE_texture.h"
59 #include "BKE_colortools.h"
60
61 #include "DEG_depsgraph.h"
62
63 #include "render_types.h"
64 #include "texture.h"
65 #include "pointdensity.h"
66
67 #include "RE_render_ext.h"
68
69 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
70 /* defined in pipeline.c, is hardcopy of active dynamic allocated Render */
71 /* only to be used here in this file, it's for speed */
72 extern struct Render R;
73 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
74
75 static ThreadMutex sample_mutex = PTHREAD_MUTEX_INITIALIZER;
76
77 static int point_data_used(PointDensity *pd)
78 {
79         int pd_bitflag = 0;
80
81         if (pd->source == TEX_PD_PSYS) {
82                 if ((pd->noise_influence == TEX_PD_NOISE_VEL) ||
83                     (pd->falloff_type == TEX_PD_FALLOFF_PARTICLE_VEL) ||
84                     (pd->color_source == TEX_PD_COLOR_PARTVEL) ||
85                     (pd->color_source == TEX_PD_COLOR_PARTSPEED))
86                 {
87                         pd_bitflag |= POINT_DATA_VEL;
88                 }
89                 if ((pd->noise_influence == TEX_PD_NOISE_AGE) ||
90                     (pd->color_source == TEX_PD_COLOR_PARTAGE) ||
91                     (pd->falloff_type == TEX_PD_FALLOFF_PARTICLE_AGE))
92                 {
93                         pd_bitflag |= POINT_DATA_LIFE;
94                 }
95         }
96         else if (pd->source == TEX_PD_OBJECT) {
97                 if (ELEM(pd->ob_color_source, TEX_PD_COLOR_VERTCOL, TEX_PD_COLOR_VERTWEIGHT, TEX_PD_COLOR_VERTNOR)) {
98                         pd_bitflag |= POINT_DATA_COLOR;
99                 }
100         }
101
102         return pd_bitflag;
103 }
104
105 static void point_data_pointers(PointDensity *pd,
106                                 float **r_data_velocity, float **r_data_life, float **r_data_color)
107 {
108         const int data_used = point_data_used(pd);
109         const int totpoint = pd->totpoints;
110         float *data = pd->point_data;
111         int offset = 0;
112         
113         if (data_used & POINT_DATA_VEL) {
114                 if (r_data_velocity)
115                         *r_data_velocity = data + offset;
116                 offset += 3 * totpoint;
117         }
118         else {
119                 if (r_data_velocity)
120                         *r_data_velocity = NULL;
121         }
122         
123         if (data_used & POINT_DATA_LIFE) {
124                 if (r_data_life)
125                         *r_data_life = data + offset;
126                 offset += totpoint;
127         }
128         else {
129                 if (r_data_life)
130                         *r_data_life = NULL;
131         }
132         
133         if (data_used & POINT_DATA_COLOR) {
134                 if (r_data_color)
135                         *r_data_color = data + offset;
136                 offset += 3 * totpoint;
137         }
138         else {
139                 if (r_data_color)
140                         *r_data_color = NULL;
141         }
142 }
143
144 /* additional data stored alongside the point density BVH,
145  * accessible by point index number to retrieve other information
146  * such as particle velocity or lifetime */
147 static void alloc_point_data(PointDensity *pd)
148 {
149         const int totpoints = pd->totpoints;
150         int data_used = point_data_used(pd);
151         int data_size = 0;
152
153         if (data_used & POINT_DATA_VEL) {
154                 /* store 3 channels of velocity data */
155                 data_size += 3;
156         }
157         if (data_used & POINT_DATA_LIFE) {
158                 /* store 1 channel of lifetime data */
159                 data_size += 1;
160         }
161         if (data_used & POINT_DATA_COLOR) {
162                 /* store 3 channels of RGB data */
163                 data_size += 3;
164         }
165
166         if (data_size) {
167                 pd->point_data = MEM_callocN(sizeof(float) * data_size * totpoints,
168                                              "particle point data");
169         }
170 }
171
172 static void pointdensity_cache_psys(EvaluationContext *eval_ctx, Scene *scene,
173                                     PointDensity *pd,
174                                     Object *ob,
175                                     ParticleSystem *psys,
176                                     float viewmat[4][4],
177                                     float winmat[4][4],
178                                     int winx, int winy,
179                                     const bool use_render_params)
180 {
181         DerivedMesh *dm;
182         ParticleKey state;
183         ParticleCacheKey *cache;
184         ParticleSimulationData sim = {NULL};
185         ParticleData *pa = NULL;
186         float cfra = BKE_scene_frame_get(scene);
187         int i /*, childexists*/ /* UNUSED */;
188         int total_particles;
189         int data_used;
190         float *data_vel, *data_life;
191         float partco[3];
192
193         /* init everything */
194         if (!psys || !ob || !pd) {
195                 return;
196         }
197
198         data_used = point_data_used(pd);
199
200         /* Just to create a valid rendering context for particles */
201         if (use_render_params) {
202                 psys_render_set(ob, psys, viewmat, winmat, winx, winy, 0);
203         }
204
205         if (use_render_params) {
206                 dm = mesh_create_derived_render(eval_ctx, scene,
207                                                 ob,
208                                                 CD_MASK_BAREMESH | CD_MASK_MTFACE | CD_MASK_MCOL);
209         }
210         else {
211                 dm = mesh_get_derived_final(eval_ctx, scene,
212                                             ob,
213                                             CD_MASK_BAREMESH | CD_MASK_MTFACE | CD_MASK_MCOL);
214         }
215
216         if (!psys_check_enabled(ob, psys, use_render_params)) {
217                 psys_render_restore(ob, psys);
218                 return;
219         }
220
221         sim.eval_ctx = eval_ctx;
222         sim.scene = scene;
223         sim.ob = ob;
224         sim.psys = psys;
225         sim.psmd = psys_get_modifier(ob, psys);
226
227         /* in case ob->imat isn't up-to-date */
228         invert_m4_m4(ob->imat, ob->obmat);
229
230         total_particles = psys->totpart + psys->totchild;
231         psys->lattice_deform_data = psys_create_lattice_deform_data(&sim);
232
233         pd->point_tree = BLI_bvhtree_new(total_particles, 0.0, 4, 6);
234         pd->totpoints = total_particles;
235         alloc_point_data(pd);
236         point_data_pointers(pd, &data_vel, &data_life, NULL);
237
238 #if 0 /* UNUSED */
239         if (psys->totchild > 0 && !(psys->part->draw & PART_DRAW_PARENT))
240                 childexists = 1;
241 #endif
242
243         for (i = 0, pa = psys->particles; i < total_particles; i++, pa++) {
244
245                 if (psys->part->type == PART_HAIR) {
246                         /* hair particles */
247                         if (i < psys->totpart && psys->pathcache)
248                                 cache = psys->pathcache[i];
249                         else if (i >= psys->totpart && psys->childcache)
250                                 cache = psys->childcache[i - psys->totpart];
251                         else
252                                 continue;
253
254                         cache += cache->segments; /* use endpoint */
255
256                         copy_v3_v3(state.co, cache->co);
257                         zero_v3(state.vel);
258                         state.time = 0.0f;
259                 }
260                 else {
261                         /* emitter particles */
262                         state.time = cfra;
263
264                         if (!psys_get_particle_state(&sim, i, &state, 0))
265                                 continue;
266
267                         if (data_used & POINT_DATA_LIFE) {
268                                 if (i < psys->totpart) {
269                                         state.time = (cfra - pa->time) / pa->lifetime;
270                                 }
271                                 else {
272                                         ChildParticle *cpa = (psys->child + i) - psys->totpart;
273                                         float pa_birthtime, pa_dietime;
274
275                                         state.time = psys_get_child_time(psys, cpa, cfra, &pa_birthtime, &pa_dietime);
276                                 }
277                         }
278                 }
279
280                 copy_v3_v3(partco, state.co);
281
282                 if (pd->psys_cache_space == TEX_PD_OBJECTSPACE)
283                         mul_m4_v3(ob->imat, partco);
284                 else if (pd->psys_cache_space == TEX_PD_OBJECTLOC) {
285                         sub_v3_v3(partco, ob->loc);
286                 }
287                 else {
288                         /* TEX_PD_WORLDSPACE */
289                 }
290
291                 BLI_bvhtree_insert(pd->point_tree, i, partco, 1);
292
293                 if (data_vel) {
294                         data_vel[i*3 + 0] = state.vel[0];
295                         data_vel[i*3 + 1] = state.vel[1];
296                         data_vel[i*3 + 2] = state.vel[2];
297                 }
298                 if (data_life) {
299                         data_life[i] = state.time;
300                 }
301         }
302
303         BLI_bvhtree_balance(pd->point_tree);
304         dm->release(dm);
305
306         if (psys->lattice_deform_data) {
307                 end_latt_deform(psys->lattice_deform_data);
308                 psys->lattice_deform_data = NULL;
309         }
310
311         if (use_render_params) {
312                 psys_render_restore(ob, psys);
313         }
314 }
315
316
317 static void pointdensity_cache_vertex_color(PointDensity *pd, Object *UNUSED(ob), DerivedMesh *dm, float *data_color)
318 {
319         const MLoop *mloop = dm->getLoopArray(dm);
320         const int totloop = dm->getNumLoops(dm);
321         const MLoopCol *mcol;
322         char layername[MAX_CUSTOMDATA_LAYER_NAME];
323         int i;
324         
325         BLI_assert(data_color);
326         
327         if (!CustomData_has_layer(&dm->loopData, CD_MLOOPCOL))
328                 return;
329         CustomData_validate_layer_name(&dm->loopData, CD_MLOOPCOL, pd->vertex_attribute_name, layername);
330         mcol = CustomData_get_layer_named(&dm->loopData, CD_MLOOPCOL, layername);
331         if (!mcol)
332                 return;
333         
334         /* Stores the number of MLoops using the same vertex, so we can normalize colors. */
335         int *mcorners = MEM_callocN(sizeof(int) * pd->totpoints, "point density corner count");
336         
337         for (i = 0; i < totloop; i++) {
338                 int v = mloop[i].v;
339
340                 if (mcorners[v] == 0) {
341                         rgb_uchar_to_float(&data_color[v * 3], &mcol[i].r);
342                 }
343                 else {
344                         float col[3];
345                         rgb_uchar_to_float(col, &mcol[i].r);
346                         add_v3_v3(&data_color[v * 3], col);
347                 }
348
349                 ++mcorners[v];
350         }
351         
352         /* Normalize colors by averaging over mcorners.
353          * All the corners share the same vertex, ie. occupy the same point in space.
354          */
355         for (i = 0; i < pd->totpoints; i++) {
356                 if (mcorners[i] > 0)
357                         mul_v3_fl(&data_color[i*3], 1.0f / mcorners[i]);
358         }
359         
360         MEM_freeN(mcorners);
361 }
362
363 static void pointdensity_cache_vertex_weight(PointDensity *pd, Object *ob, DerivedMesh *dm, float *data_color)
364 {
365         const int totvert = dm->getNumVerts(dm);
366         const MDeformVert *mdef, *dv;
367         int mdef_index;
368         int i;
369         
370         BLI_assert(data_color);
371         
372         mdef = CustomData_get_layer(&dm->vertData, CD_MDEFORMVERT);
373         if (!mdef)
374                 return;
375         mdef_index = defgroup_name_index(ob, pd->vertex_attribute_name);
376         if (mdef_index < 0)
377                 mdef_index = ob->actdef - 1;
378         if (mdef_index < 0)
379                 return;
380         
381         for (i = 0, dv = mdef; i < totvert; ++i, ++dv, data_color += 3) {
382                 MDeformWeight *dw;
383                 int j;
384                 
385                 for (j = 0, dw = dv->dw; j < dv->totweight; ++j, ++dw) {
386                         if (dw->def_nr == mdef_index) {
387                                 copy_v3_fl(data_color, dw->weight);
388                                 break;
389                         }
390                 }
391         }
392 }
393
394 static void pointdensity_cache_vertex_normal(PointDensity *pd, Object *UNUSED(ob), DerivedMesh *dm, float *data_color)
395 {
396         MVert *mvert = dm->getVertArray(dm), *mv;
397         int i;
398         
399         BLI_assert(data_color);
400         
401         for (i = 0, mv = mvert; i < pd->totpoints; i++, mv++, data_color += 3) {
402                 normal_short_to_float_v3(data_color, mv->no);
403         }
404 }
405
406 static void pointdensity_cache_object(EvaluationContext *eval_ctx, Scene *scene,
407                                       PointDensity *pd,
408                                       Object *ob,
409                                       const bool use_render_params)
410 {
411         float *data_color;
412         int i;
413         DerivedMesh *dm;
414         CustomDataMask mask = CD_MASK_BAREMESH | CD_MASK_MTFACE | CD_MASK_MCOL;
415         MVert *mvert = NULL, *mv;
416
417         switch (pd->ob_color_source) {
418                 case TEX_PD_COLOR_VERTCOL:
419                         mask |= CD_MASK_MLOOPCOL;
420                         break;
421                 case TEX_PD_COLOR_VERTWEIGHT:
422                         mask |= CD_MASK_MDEFORMVERT;
423                         break;
424         }
425
426         if (use_render_params) {
427                 dm = mesh_create_derived_render(eval_ctx, scene, ob, mask);
428         }
429         else {
430                 dm = mesh_get_derived_final(eval_ctx, scene, ob, mask);
431         }
432
433         mvert = dm->getVertArray(dm);   /* local object space */
434         pd->totpoints = dm->getNumVerts(dm);
435         if (pd->totpoints == 0) {
436                 return;
437         }
438
439         pd->point_tree = BLI_bvhtree_new(pd->totpoints, 0.0, 4, 6);
440         alloc_point_data(pd);
441         point_data_pointers(pd, NULL, NULL, &data_color);
442
443         for (i = 0, mv = mvert; i < pd->totpoints; i++, mv++) {
444                 float co[3];
445
446                 copy_v3_v3(co, mv->co);
447
448                 switch (pd->ob_cache_space) {
449                         case TEX_PD_OBJECTSPACE:
450                                 break;
451                         case TEX_PD_OBJECTLOC:
452                                 mul_m4_v3(ob->obmat, co);
453                                 sub_v3_v3(co, ob->loc);
454                                 break;
455                         case TEX_PD_WORLDSPACE:
456                         default:
457                                 mul_m4_v3(ob->obmat, co);
458                                 break;
459                 }
460
461                 BLI_bvhtree_insert(pd->point_tree, i, co, 1);
462         }
463         
464         switch (pd->ob_color_source) {
465                 case TEX_PD_COLOR_VERTCOL:
466                         pointdensity_cache_vertex_color(pd, ob, dm, data_color);
467                         break;
468                 case TEX_PD_COLOR_VERTWEIGHT:
469                         pointdensity_cache_vertex_weight(pd, ob, dm, data_color);
470                         break;
471                 case TEX_PD_COLOR_VERTNOR:
472                         pointdensity_cache_vertex_normal(pd, ob, dm, data_color);
473                         break;
474         }
475
476         BLI_bvhtree_balance(pd->point_tree);
477         dm->release(dm);
478
479 }
480
481 static void cache_pointdensity_ex(EvaluationContext *eval_ctx, Scene *scene,
482                                   PointDensity *pd,
483                                   float viewmat[4][4],
484                                   float winmat[4][4],
485                                   int winx, int winy,
486                                   const bool use_render_params)
487 {
488         if (pd == NULL) {
489                 return;
490         }
491
492         if (pd->point_tree) {
493                 BLI_bvhtree_free(pd->point_tree);
494                 pd->point_tree = NULL;
495         }
496
497         if (pd->source == TEX_PD_PSYS) {
498                 Object *ob = pd->object;
499                 ParticleSystem *psys;
500
501                 if (!ob || !pd->psys) {
502                         return;
503                 }
504
505                 psys = BLI_findlink(&ob->particlesystem, pd->psys - 1);
506                 if (!psys) {
507                         return;
508                 }
509
510                 pointdensity_cache_psys(eval_ctx,
511                                         scene,
512                                         pd,
513                                         ob,
514                                         psys,
515                                         viewmat, winmat,
516                                         winx, winy,
517                                         use_render_params);
518         }
519         else if (pd->source == TEX_PD_OBJECT) {
520                 Object *ob = pd->object;
521                 if (ob && ob->type == OB_MESH)
522                         pointdensity_cache_object(eval_ctx, scene, pd, ob, use_render_params);
523         }
524 }
525
526 void cache_pointdensity(Render *re, PointDensity *pd)
527 {
528         cache_pointdensity_ex(re->eval_ctx,
529                               re->scene,
530                               pd,
531                               re->viewmat, re->winmat,
532                               re->winx, re->winy,
533                               true);
534 }
535
536 void free_pointdensity(PointDensity *pd)
537 {
538         if (pd == NULL) {
539                 return;
540         }
541
542         if (pd->point_tree) {
543                 BLI_bvhtree_free(pd->point_tree);
544                 pd->point_tree = NULL;
545         }
546
547         if (pd->point_data) {
548                 MEM_freeN(pd->point_data);
549                 pd->point_data = NULL;
550         }
551         pd->totpoints = 0;
552 }
553
554 void make_pointdensities(Render *re)
555 {
556         Tex *tex;
557
558         if (re->scene->r.scemode & R_BUTS_PREVIEW) {
559                 return;
560         }
561
562         re->i.infostr = IFACE_("Caching Point Densities");
563         re->stats_draw(re->sdh, &re->i);
564
565         for (tex = re->main->tex.first; tex != NULL; tex = tex->id.next) {
566                 if (tex->id.us && tex->type == TEX_POINTDENSITY) {
567                         cache_pointdensity(re, tex->pd);
568                 }
569         }
570
571         re->i.infostr = NULL;
572         re->stats_draw(re->sdh, &re->i);
573 }
574
575 void free_pointdensities(Render *re)
576 {
577         Tex *tex;
578
579         if (re->scene->r.scemode & R_BUTS_PREVIEW)
580                 return;
581
582         for (tex = re->main->tex.first; tex != NULL; tex = tex->id.next) {
583                 if (tex->id.us && tex->type == TEX_POINTDENSITY) {
584                         free_pointdensity(tex->pd);
585                 }
586         }
587 }
588
589 typedef struct PointDensityRangeData {
590         float *density;
591         float squared_radius;
592         float *point_data_life;
593         float *point_data_velocity;
594         float *point_data_color;
595         float *vec;
596         float *col;
597         float softness;
598         short falloff_type;
599         short noise_influence;
600         float *age;
601         struct CurveMapping *density_curve;
602         float velscale;
603 } PointDensityRangeData;
604
605 static float density_falloff(PointDensityRangeData *pdr, int index, float squared_dist)
606 {
607         const float dist = (pdr->squared_radius - squared_dist) / pdr->squared_radius * 0.5f;
608         float density = 0.0f;
609         
610         switch (pdr->falloff_type) {
611                 case TEX_PD_FALLOFF_STD:
612                         density = dist;
613                         break;
614                 case TEX_PD_FALLOFF_SMOOTH:
615                         density = 3.0f * dist * dist - 2.0f * dist * dist * dist;
616                         break;
617                 case TEX_PD_FALLOFF_SOFT:
618                         density = pow(dist, pdr->softness);
619                         break;
620                 case TEX_PD_FALLOFF_CONSTANT:
621                         density = pdr->squared_radius;
622                         break;
623                 case TEX_PD_FALLOFF_ROOT:
624                         density = sqrtf(dist);
625                         break;
626                 case TEX_PD_FALLOFF_PARTICLE_AGE:
627                         if (pdr->point_data_life)
628                                 density = dist * MIN2(pdr->point_data_life[index], 1.0f);
629                         else
630                                 density = dist;
631                         break;
632                 case TEX_PD_FALLOFF_PARTICLE_VEL:
633                         if (pdr->point_data_velocity)
634                                 density = dist * len_v3(&pdr->point_data_velocity[index * 3]) * pdr->velscale;
635                         else
636                                 density = dist;
637                         break;
638         }
639         
640         if (pdr->density_curve && dist != 0.0f) {
641                 curvemapping_initialize(pdr->density_curve);
642                 density = curvemapping_evaluateF(pdr->density_curve, 0, density / dist) * dist;
643         }
644         
645         return density;
646 }
647
648 static void accum_density(void *userdata, int index, const float co[3], float squared_dist)
649 {
650         PointDensityRangeData *pdr = (PointDensityRangeData *)userdata;
651         float density = 0.0f;
652
653         UNUSED_VARS(co);
654
655         if (pdr->point_data_velocity) {
656                 pdr->vec[0] += pdr->point_data_velocity[index * 3 + 0]; // * density;
657                 pdr->vec[1] += pdr->point_data_velocity[index * 3 + 1]; // * density;
658                 pdr->vec[2] += pdr->point_data_velocity[index * 3 + 2]; // * density;
659         }
660         if (pdr->point_data_life) {
661                 *pdr->age += pdr->point_data_life[index]; // * density;
662         }
663         if (pdr->point_data_color) {
664                 add_v3_v3(pdr->col, &pdr->point_data_color[index * 3]); // * density;
665         }
666
667         density = density_falloff(pdr, index, squared_dist);
668
669         *pdr->density += density;
670 }
671
672
673 static void init_pointdensityrangedata(PointDensity *pd, PointDensityRangeData *pdr,
674         float *density, float *vec, float *age, float *col, struct CurveMapping *density_curve, float velscale)
675 {
676         pdr->squared_radius = pd->radius * pd->radius;
677         pdr->density = density;
678         pdr->falloff_type = pd->falloff_type;
679         pdr->vec = vec;
680         pdr->age = age;
681         pdr->col = col;
682         pdr->softness = pd->falloff_softness;
683         pdr->noise_influence = pd->noise_influence;
684         point_data_pointers(pd, &pdr->point_data_velocity, &pdr->point_data_life, &pdr->point_data_color);
685         pdr->density_curve = density_curve;
686         pdr->velscale = velscale;
687 }
688
689
690 static int pointdensity(PointDensity *pd,
691                         const float texvec[3],
692                         TexResult *texres,
693                         float r_vec[3],
694                         float *r_age,
695                         float r_col[3])
696 {
697         int retval = TEX_INT;
698         PointDensityRangeData pdr;
699         float density = 0.0f, age = 0.0f, time = 0.0f;
700         float vec[3] = {0.0f, 0.0f, 0.0f}, col[3] = {0.0f, 0.0f, 0.0f}, co[3];
701         float turb, noise_fac;
702         int num = 0;
703
704         texres->tin = 0.0f;
705
706         if ((!pd) || (!pd->point_tree))
707                 return 0;
708
709         init_pointdensityrangedata(pd, &pdr, &density, vec, &age, col,
710                 (pd->flag & TEX_PD_FALLOFF_CURVE ? pd->falloff_curve : NULL),
711                 pd->falloff_speed_scale * 0.001f);
712         noise_fac = pd->noise_fac * 0.5f;       /* better default */
713
714         copy_v3_v3(co, texvec);
715
716         if (point_data_used(pd)) {
717                 /* does a BVH lookup to find accumulated density and additional point data *
718                  * stores particle velocity vector in 'vec', and particle lifetime in 'time' */
719                 num = BLI_bvhtree_range_query(pd->point_tree, co, pd->radius, accum_density, &pdr);
720                 if (num > 0) {
721                         age /= num;
722                         mul_v3_fl(vec, 1.0f / num);
723                         mul_v3_fl(col, 1.0f / num);
724                 }
725
726                 /* reset */
727                 density = 0.0f;
728                 zero_v3(vec);
729                 zero_v3(col);
730         }
731
732         if (pd->flag & TEX_PD_TURBULENCE) {
733
734                 if (pd->noise_influence == TEX_PD_NOISE_AGE) {
735                         turb = BLI_gTurbulence(pd->noise_size, texvec[0] + age, texvec[1] + age, texvec[2] + age,
736                                                pd->noise_depth, 0, pd->noise_basis);
737                 }
738                 else if (pd->noise_influence == TEX_PD_NOISE_TIME) {
739                         time = R.r.cfra / (float)R.r.efra;
740                         turb = BLI_gTurbulence(pd->noise_size, texvec[0] + time, texvec[1] + time, texvec[2] + time,
741                                                pd->noise_depth, 0, pd->noise_basis);
742                         //turb = BLI_turbulence(pd->noise_size, texvec[0]+time, texvec[1]+time, texvec[2]+time, pd->noise_depth);
743                 }
744                 else {
745                         turb = BLI_gTurbulence(pd->noise_size, texvec[0] + vec[0], texvec[1] + vec[1], texvec[2] + vec[2],
746                                                pd->noise_depth, 0, pd->noise_basis);
747                 }
748
749                 turb -= 0.5f;   /* re-center 0.0-1.0 range around 0 to prevent offsetting result */
750
751                 /* now we have an offset coordinate to use for the density lookup */
752                 co[0] = texvec[0] + noise_fac * turb;
753                 co[1] = texvec[1] + noise_fac * turb;
754                 co[2] = texvec[2] + noise_fac * turb;
755         }
756
757         /* BVH query with the potentially perturbed coordinates */
758         num = BLI_bvhtree_range_query(pd->point_tree, co, pd->radius, accum_density, &pdr);
759         if (num > 0) {
760                 age /= num;
761                 mul_v3_fl(vec, 1.0f / num);
762                 mul_v3_fl(col, 1.0f / num);
763         }
764
765         texres->tin = density;
766         if (r_age != NULL) {
767                 *r_age = age;
768         }
769         if (r_vec != NULL) {
770                 copy_v3_v3(r_vec, vec);
771         }
772         if (r_col != NULL) {
773                 copy_v3_v3(r_col, col);
774         }
775
776         return retval;
777 }
778
779 static int pointdensity_color(PointDensity *pd, TexResult *texres, float age, const float vec[3], const float col[3])
780 {
781         int retval = TEX_RGB;
782
783         if (pd->source == TEX_PD_PSYS) {
784                 float rgba[4];
785                 
786                 switch (pd->color_source) {
787                         case TEX_PD_COLOR_PARTAGE:
788                                 if (pd->coba) {
789                                         if (do_colorband(pd->coba, age, rgba)) {
790                                                 texres->talpha = true;
791                                                 copy_v3_v3(&texres->tr, rgba);
792                                                 texres->tin *= rgba[3];
793                                                 texres->ta = texres->tin;
794                                         }
795                                 }
796                                 break;
797                         case TEX_PD_COLOR_PARTSPEED:
798                         {
799                                 float speed = len_v3(vec) * pd->speed_scale;
800                                 
801                                 if (pd->coba) {
802                                         if (do_colorband(pd->coba, speed, rgba)) {
803                                                 texres->talpha = true;
804                                                 copy_v3_v3(&texres->tr, rgba);
805                                                 texres->tin *= rgba[3];
806                                                 texres->ta = texres->tin;
807                                         }
808                                 }
809                                 break;
810                         }
811                         case TEX_PD_COLOR_PARTVEL:
812                                 texres->talpha = true;
813                                 mul_v3_v3fl(&texres->tr, vec, pd->speed_scale);
814                                 texres->ta = texres->tin;
815                                 break;
816                         case TEX_PD_COLOR_CONSTANT:
817                         default:
818                                 texres->tr = texres->tg = texres->tb = texres->ta = 1.0f;
819                                 retval = TEX_INT;
820                                 break;
821                 }
822         }
823         else {
824                 float rgba[4];
825                 
826                 switch (pd->ob_color_source) {
827                         case TEX_PD_COLOR_VERTCOL:
828                                 texres->talpha = true;
829                                 copy_v3_v3(&texres->tr, col);
830                                 texres->ta = texres->tin;
831                                 break;
832                         case TEX_PD_COLOR_VERTWEIGHT:
833                                 texres->talpha = true;
834                                 if (pd->coba && do_colorband(pd->coba, col[0], rgba)) {
835                                         copy_v3_v3(&texres->tr, rgba);
836                                         texres->tin *= rgba[3];
837                                 }
838                                 else {
839                                         copy_v3_v3(&texres->tr, col);
840                                 }
841                                 texres->ta = texres->tin;
842                                 break;
843                         case TEX_PD_COLOR_VERTNOR:
844                                 texres->talpha = true;
845                                 copy_v3_v3(&texres->tr, col);
846                                 texres->ta = texres->tin;
847                                 break;
848                         case TEX_PD_COLOR_CONSTANT:
849                         default:
850                                 texres->tr = texres->tg = texres->tb = texres->ta = 1.0f;
851                                 retval = TEX_INT;
852                                 break;
853                 }
854         }
855         
856         return retval;
857 }
858
859 int pointdensitytex(Tex *tex, const float texvec[3], TexResult *texres)
860 {
861         PointDensity *pd = tex->pd;
862         float age = 0.0f;
863         float vec[3] = {0.0f, 0.0f, 0.0f};
864         float col[3] = {0.0f, 0.0f, 0.0f};
865         int retval = pointdensity(pd, texvec, texres, vec, &age, col);
866
867         retval |= pointdensity_color(pd, texres, age, vec, col);
868         BRICONTRGB;
869         
870         return retval;
871
872 #if 0
873         if (texres->nor!=NULL) {
874                 texres->nor[0] = texres->nor[1] = texres->nor[2] = 0.0f;
875         }
876 #endif
877 }
878
879 static void sample_dummy_point_density(int resolution, float *values)
880 {
881         memset(values, 0, sizeof(float) * 4 * resolution * resolution * resolution);
882 }
883
884 static void particle_system_minmax(EvaluationContext *eval_ctx,
885                                    Scene *scene,
886                                    Object *object,
887                                    ParticleSystem *psys,
888                                    float radius,
889                                    const bool use_render_params,
890                                    float min[3], float max[3])
891 {
892         const float size[3] = {radius, radius, radius};
893         const float cfra = BKE_scene_frame_get(scene);
894         ParticleSettings *part = psys->part;
895         ParticleSimulationData sim = {NULL};
896         ParticleData *pa = NULL;
897         int i;
898         int total_particles;
899         float mat[4][4], imat[4][4];
900
901         INIT_MINMAX(min, max);
902         if (part->type == PART_HAIR) {
903                 /* TOOD(sergey): Not supported currently. */
904                 return;
905         }
906
907         unit_m4(mat);
908         if (use_render_params) {
909                 psys_render_set(object, psys, mat, mat, 1, 1, 0);
910         }
911
912         sim.eval_ctx = eval_ctx;
913         sim.scene = scene;
914         sim.ob = object;
915         sim.psys = psys;
916         sim.psmd = psys_get_modifier(object, psys);
917
918         invert_m4_m4(imat, object->obmat);
919         total_particles = psys->totpart + psys->totchild;
920         psys->lattice_deform_data = psys_create_lattice_deform_data(&sim);
921
922         for (i = 0, pa = psys->particles; i < total_particles; i++, pa++) {
923                 float co_object[3], co_min[3], co_max[3];
924                 ParticleKey state;
925                 state.time = cfra;
926                 if (!psys_get_particle_state(&sim, i, &state, 0)) {
927                         continue;
928                 }
929                 mul_v3_m4v3(co_object, imat, state.co);
930                 sub_v3_v3v3(co_min, co_object, size);
931                 add_v3_v3v3(co_max, co_object, size);
932                 minmax_v3v3_v3(min, max, co_min);
933                 minmax_v3v3_v3(min, max, co_max);
934         }
935
936         if (psys->lattice_deform_data) {
937                 end_latt_deform(psys->lattice_deform_data);
938                 psys->lattice_deform_data = NULL;
939         }
940
941         if (use_render_params) {
942                 psys_render_restore(object, psys);
943         }
944 }
945
946 void RE_point_density_cache(
947         Scene *scene,
948         SceneLayer *scene_layer,
949         PointDensity *pd,
950         const bool use_render_params)
951 {
952         EvaluationContext eval_ctx = {0};
953         float mat[4][4];
954
955         DEG_evaluation_context_init(&eval_ctx, use_render_params ? DAG_EVAL_RENDER
956                                                                  : DAG_EVAL_VIEWPORT);
957
958         eval_ctx.scene_layer = scene_layer;
959
960         /* Same matricies/resolution as dupli_render_particle_set(). */
961         unit_m4(mat);
962         BLI_mutex_lock(&sample_mutex);
963         cache_pointdensity_ex(&eval_ctx, scene, pd, mat, mat, 1, 1, use_render_params);
964         BLI_mutex_unlock(&sample_mutex);
965 }
966
967 void RE_point_density_minmax(
968         struct Scene *scene,
969         SceneLayer *sl,
970         struct PointDensity *pd,
971         const bool use_render_params,
972         float r_min[3], float r_max[3])
973 {
974         Object *object = pd->object;
975         if (object == NULL) {
976                 zero_v3(r_min);
977                 zero_v3(r_max);
978                 return;
979         }
980         if (pd->source == TEX_PD_PSYS) {
981                 ParticleSystem *psys;
982                 EvaluationContext eval_ctx = {0};
983
984                 if (pd->psys == 0) {
985                         zero_v3(r_min);
986                         zero_v3(r_max);
987                         return;
988                 }
989                 psys = BLI_findlink(&object->particlesystem, pd->psys - 1);
990                 if (psys == NULL) {
991                         zero_v3(r_min);
992                         zero_v3(r_max);
993                         return;
994                 }
995
996                 DEG_evaluation_context_init(&eval_ctx, use_render_params ? DAG_EVAL_RENDER :
997                                                                            DAG_EVAL_VIEWPORT);
998
999                 eval_ctx.ctime = (float)scene->r.cfra + scene->r.subframe;
1000                 eval_ctx.scene_layer = sl;
1001
1002                 particle_system_minmax(&eval_ctx,
1003                                        scene,
1004                                        object,
1005                                        psys,
1006                                        pd->radius,
1007                                        use_render_params,
1008                                        r_min, r_max);
1009         }
1010         else {
1011                 float radius[3] = {pd->radius, pd->radius, pd->radius};
1012                 BoundBox *bb = BKE_object_boundbox_get(object);
1013
1014                 if (bb != NULL) {
1015                         BLI_assert((bb->flag & BOUNDBOX_DIRTY) == 0);
1016                         copy_v3_v3(r_min, bb->vec[0]);
1017                         copy_v3_v3(r_max, bb->vec[6]);
1018                         /* Adjust texture space to include density points on the boundaries. */
1019                         sub_v3_v3(r_min, radius);
1020                         add_v3_v3(r_max, radius);
1021                 }
1022                 else {
1023                         zero_v3(r_min);
1024                         zero_v3(r_max);
1025                 }
1026         }
1027 }
1028
1029 typedef struct SampleCallbackData {
1030         PointDensity *pd;
1031         int resolution;
1032         float *min, *dim;
1033         float *values;
1034 } SampleCallbackData;
1035
1036 static void point_density_sample_func(void *data_v, const int iter)
1037 {
1038         SampleCallbackData *data = (SampleCallbackData *)data_v;
1039
1040         const int resolution = data->resolution;
1041         const int resolution2 = resolution * resolution;
1042         const float *min = data->min, *dim = data->dim;
1043         PointDensity *pd = data->pd;
1044         float *values = data->values;
1045
1046         size_t z = (size_t)iter;
1047         for (size_t y = 0; y < resolution; ++y) {
1048                 for (size_t x = 0; x < resolution; ++x) {
1049                         size_t index = z * resolution2 + y * resolution + x;
1050                         float texvec[3];
1051                         float age, vec[3], col[3];
1052                         TexResult texres;
1053
1054                         copy_v3_v3(texvec, min);
1055                         texvec[0] += dim[0] * (float)x / (float)resolution;
1056                         texvec[1] += dim[1] * (float)y / (float)resolution;
1057                         texvec[2] += dim[2] * (float)z / (float)resolution;
1058
1059                         pointdensity(pd, texvec, &texres, vec, &age, col);
1060                         pointdensity_color(pd, &texres, age, vec, col);
1061
1062                         copy_v3_v3(&values[index*4 + 0], &texres.tr);
1063                         values[index*4 + 3] = texres.tin;
1064                 }
1065         }
1066 }
1067
1068 /* NOTE 1: Requires RE_point_density_cache() to be called first.
1069  * NOTE 2: Frees point density structure after sampling.
1070  */
1071 void RE_point_density_sample(
1072         Scene *scene,
1073         SceneLayer *sl,
1074         PointDensity *pd,
1075         const int resolution,
1076         const bool use_render_params,
1077         float *values)
1078 {
1079         Object *object = pd->object;
1080         float min[3], max[3], dim[3];
1081
1082         /* TODO(sergey): Implement some sort of assert() that point density
1083          * was cached already.
1084          */
1085
1086         if (object == NULL) {
1087                 sample_dummy_point_density(resolution, values);
1088                 return;
1089         }
1090
1091         BLI_mutex_lock(&sample_mutex);
1092         RE_point_density_minmax(scene,
1093                                 sl,
1094                                 pd,
1095                                 use_render_params,
1096                                 min,
1097                                 max);
1098         BLI_mutex_unlock(&sample_mutex);
1099         sub_v3_v3v3(dim, max, min);
1100         if (dim[0] <= 0.0f || dim[1] <= 0.0f || dim[2] <= 0.0f) {
1101                 sample_dummy_point_density(resolution, values);
1102                 return;
1103         }
1104
1105         SampleCallbackData data;
1106         data.pd = pd;
1107         data.resolution = resolution;
1108         data.min = min;
1109         data.dim = dim;
1110         data.values = values;
1111         BLI_task_parallel_range(0,
1112                                 resolution,
1113                                 &data,
1114                                 point_density_sample_func,
1115                                 resolution > 32);
1116
1117         free_pointdensity(pd);
1118 }
1119
1120 void RE_point_density_free(struct PointDensity *pd)
1121 {
1122         free_pointdensity(pd);
1123 }