Moved particle code for distributions into own file to make hacking
authorLukas Tönne <lukas.toenne@gmail.com>
Tue, 28 Oct 2014 15:29:33 +0000 (16:29 +0100)
committerLukas Tönne <lukas.toenne@gmail.com>
Tue, 20 Jan 2015 08:30:04 +0000 (09:30 +0100)
easier.

This code is badly broken and needs to be replaced, but at least having
a workable code structure might help with quick hacks to fix the worst
cases.

source/blender/blenkernel/BKE_particle.h
source/blender/blenkernel/CMakeLists.txt
source/blender/blenkernel/intern/particle_distribute.c [new file with mode: 0644]
source/blender/blenkernel/intern/particle_system.c

index 15d3c2f5ee8f5924bc2cc82d44ac55c1d51fa2a9..31c0125584a5127acc5804078991570697f60f7d 100644 (file)
@@ -279,6 +279,9 @@ BLI_INLINE void psys_frand_vec(ParticleSystem *psys, unsigned int seed, float ve
 int count_particles(struct ParticleSystem *psys);
 int count_particles_mod(struct ParticleSystem *psys, int totgr, int cur);
 
+int psys_get_child_number(struct Scene *scene, struct ParticleSystem *psys);
+int psys_get_tot_child(struct Scene *scene, struct ParticleSystem *psys);
+
 struct ParticleSystem *psys_get_current(struct Object *ob);
 /* for rna */
 short psys_get_current_num(struct Object *ob);
@@ -404,6 +407,7 @@ void psys_particle_on_dm(struct DerivedMesh *dm, int from, int index, int index_
                          float orco[3], float ornor[3]);
 
 /* particle_system.c */
+void distribute_particles(struct ParticleSimulationData *sim, int from);
 void initialize_particle(struct ParticleSimulationData *sim, struct ParticleData *pa);
 void psys_calc_dmcache(struct Object *ob, struct DerivedMesh *dm, struct ParticleSystem *psys);
 int psys_particle_dm_face_lookup(struct Object *ob, struct DerivedMesh *dm, int index, const float fw[4], struct LinkNode *node);
index 28b6aaa3b7c8abbaf9722880139592c792db8e1a..d94710a7373c7ab23739a9d19e096557592272c8 100644 (file)
@@ -137,6 +137,7 @@ set(SRC
        intern/packedFile.c
        intern/paint.c
        intern/particle.c
+       intern/particle_distribute.c
        intern/particle_system.c
        intern/pbvh.c
        intern/pbvh_bmesh.c
diff --git a/source/blender/blenkernel/intern/particle_distribute.c b/source/blender/blenkernel/intern/particle_distribute.c
new file mode 100644 (file)
index 0000000..6e16435
--- /dev/null
@@ -0,0 +1,1104 @@
+/*
+ * ***** BEGIN GPL LICENSE BLOCK *****
+ *
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ *
+ * The Original Code is Copyright (C) 2007 by Janne Karhu.
+ * All rights reserved.
+ *
+ * The Original Code is: all of this file.
+ *
+ * Contributor(s): Raul Fernandez Hernandez (Farsthary),
+ *                 Stephen Swhitehorn,
+ *                 Lukas Toenne
+ *
+ * ***** END GPL LICENSE BLOCK *****
+ */
+
+/** \file blender/blenkernel/intern/particle_distribute.c
+ *  \ingroup bke
+ */
+
+#include <string.h>
+
+#include "MEM_guardedalloc.h"
+
+#include "BLI_utildefines.h"
+#include "BLI_jitter.h"
+#include "BLI_kdtree.h"
+#include "BLI_math.h"
+#include "BLI_rand.h"
+#include "BLI_sort.h"
+#include "BLI_task.h"
+
+#include "DNA_mesh_types.h"
+#include "DNA_meshdata_types.h"
+#include "DNA_modifier_types.h"
+#include "DNA_particle_types.h"
+#include "DNA_scene_types.h"
+
+#include "BKE_cdderivedmesh.h"
+#include "BKE_DerivedMesh.h"
+#include "BKE_global.h"
+#include "BKE_mesh.h"
+#include "BKE_object.h"
+#include "BKE_particle.h"
+
+static void alloc_child_particles(ParticleSystem *psys, int tot)
+{
+       if (psys->child) {
+               /* only re-allocate if we have to */
+               if (psys->part->childtype && psys->totchild == tot) {
+                       memset(psys->child, 0, tot*sizeof(ChildParticle));
+                       return;
+               }
+
+               MEM_freeN(psys->child);
+               psys->child=NULL;
+               psys->totchild=0;
+       }
+
+       if (psys->part->childtype) {
+               psys->totchild= tot;
+               if (psys->totchild)
+                       psys->child= MEM_callocN(psys->totchild*sizeof(ChildParticle), "child_particles");
+       }
+}
+
+static void distribute_simple_children(Scene *scene, Object *ob, DerivedMesh *finaldm, ParticleSystem *psys)
+{
+       ChildParticle *cpa = NULL;
+       int i, p;
+       int child_nbr= psys_get_child_number(scene, psys);
+       int totpart= psys_get_tot_child(scene, psys);
+
+       alloc_child_particles(psys, totpart);
+
+       cpa = psys->child;
+       for (i=0; i<child_nbr; i++) {
+               for (p=0; p<psys->totpart; p++,cpa++) {
+                       float length=2.0;
+                       cpa->parent=p;
+                                       
+                       /* create even spherical distribution inside unit sphere */
+                       while (length>=1.0f) {
+                               cpa->fuv[0]=2.0f*BLI_frand()-1.0f;
+                               cpa->fuv[1]=2.0f*BLI_frand()-1.0f;
+                               cpa->fuv[2]=2.0f*BLI_frand()-1.0f;
+                               length=len_v3(cpa->fuv);
+                       }
+
+                       cpa->num=-1;
+               }
+       }
+       /* dmcache must be updated for parent particles if children from faces is used */
+       psys_calc_dmcache(ob, finaldm, psys);
+}
+static void distribute_grid(DerivedMesh *dm, ParticleSystem *psys)
+{
+       ParticleData *pa=NULL;
+       float min[3], max[3], delta[3], d;
+       MVert *mv, *mvert = dm->getVertDataArray(dm,0);
+       int totvert=dm->getNumVerts(dm), from=psys->part->from;
+       int i, j, k, p, res=psys->part->grid_res, size[3], axis;
+
+       /* find bounding box of dm */
+       if (totvert > 0) {
+               mv=mvert;
+               copy_v3_v3(min, mv->co);
+               copy_v3_v3(max, mv->co);
+               mv++;
+               for (i = 1; i < totvert; i++, mv++) {
+                       minmax_v3v3_v3(min, max, mv->co);
+               }
+       }
+       else {
+               zero_v3(min);
+               zero_v3(max);
+       }
+
+       sub_v3_v3v3(delta, max, min);
+
+       /* determine major axis */
+       axis = axis_dominant_v3_single(delta);
+        
+       d = delta[axis]/(float)res;
+
+       size[axis] = res;
+       size[(axis+1)%3] = (int)ceil(delta[(axis+1)%3]/d);
+       size[(axis+2)%3] = (int)ceil(delta[(axis+2)%3]/d);
+
+       /* float errors grrr.. */
+       size[(axis+1)%3] = MIN2(size[(axis+1)%3],res);
+       size[(axis+2)%3] = MIN2(size[(axis+2)%3],res);
+
+       size[0] = MAX2(size[0], 1);
+       size[1] = MAX2(size[1], 1);
+       size[2] = MAX2(size[2], 1);
+
+       /* no full offset for flat/thin objects */
+       min[0]+= d < delta[0] ? d/2.f : delta[0]/2.f;
+       min[1]+= d < delta[1] ? d/2.f : delta[1]/2.f;
+       min[2]+= d < delta[2] ? d/2.f : delta[2]/2.f;
+
+       for (i=0,p=0,pa=psys->particles; i<res; i++) {
+               for (j=0; j<res; j++) {
+                       for (k=0; k<res; k++,p++,pa++) {
+                               pa->fuv[0] = min[0] + (float)i*d;
+                               pa->fuv[1] = min[1] + (float)j*d;
+                               pa->fuv[2] = min[2] + (float)k*d;
+                               pa->flag |= PARS_UNEXIST;
+                               pa->hair_index = 0; /* abused in volume calculation */
+                       }
+               }
+       }
+
+       /* enable particles near verts/edges/faces/inside surface */
+       if (from==PART_FROM_VERT) {
+               float vec[3];
+
+               pa=psys->particles;
+
+               min[0] -= d/2.0f;
+               min[1] -= d/2.0f;
+               min[2] -= d/2.0f;
+
+               for (i=0,mv=mvert; i<totvert; i++,mv++) {
+                       sub_v3_v3v3(vec,mv->co,min);
+                       vec[0]/=delta[0];
+                       vec[1]/=delta[1];
+                       vec[2]/=delta[2];
+                       pa[((int)(vec[0] * (size[0] - 1))  * res +
+                           (int)(vec[1] * (size[1] - 1))) * res +
+                           (int)(vec[2] * (size[2] - 1))].flag &= ~PARS_UNEXIST;
+               }
+       }
+       else if (ELEM(from,PART_FROM_FACE,PART_FROM_VOLUME)) {
+               float co1[3], co2[3];
+
+               MFace *mface= NULL, *mface_array;
+               float v1[3], v2[3], v3[3], v4[4], lambda;
+               int a, a1, a2, a0mul, a1mul, a2mul, totface;
+               int amax= from==PART_FROM_FACE ? 3 : 1;
+
+               totface=dm->getNumTessFaces(dm);
+               mface=mface_array=dm->getTessFaceDataArray(dm,CD_MFACE);
+               
+               for (a=0; a<amax; a++) {
+                       if (a==0) { a0mul=res*res; a1mul=res; a2mul=1; }
+                       else if (a==1) { a0mul=res; a1mul=1; a2mul=res*res; }
+                       else { a0mul=1; a1mul=res*res; a2mul=res; }
+
+                       for (a1=0; a1<size[(a+1)%3]; a1++) {
+                               for (a2=0; a2<size[(a+2)%3]; a2++) {
+                                       mface= mface_array;
+
+                                       pa = psys->particles + a1*a1mul + a2*a2mul;
+                                       copy_v3_v3(co1, pa->fuv);
+                                       co1[a] -= d < delta[a] ? d/2.f : delta[a]/2.f;
+                                       copy_v3_v3(co2, co1);
+                                       co2[a] += delta[a] + 0.001f*d;
+                                       co1[a] -= 0.001f*d;
+                                       
+                                       /* lets intersect the faces */
+                                       for (i=0; i<totface; i++,mface++) {
+                                               copy_v3_v3(v1, mvert[mface->v1].co);
+                                               copy_v3_v3(v2, mvert[mface->v2].co);
+                                               copy_v3_v3(v3, mvert[mface->v3].co);
+
+                                               if (isect_axial_line_tri_v3(a, co1, co2, v2, v3, v1, &lambda)) {
+                                                       if (from==PART_FROM_FACE)
+                                                               (pa+(int)(lambda*size[a])*a0mul)->flag &= ~PARS_UNEXIST;
+                                                       else /* store number of intersections */
+                                                               (pa+(int)(lambda*size[a])*a0mul)->hair_index++;
+                                               }
+                                               else if (mface->v4) {
+                                                       copy_v3_v3(v4, mvert[mface->v4].co);
+
+                                                       if (isect_axial_line_tri_v3(a, co1, co2, v4, v1, v3, &lambda)) {
+                                                               if (from==PART_FROM_FACE)
+                                                                       (pa+(int)(lambda*size[a])*a0mul)->flag &= ~PARS_UNEXIST;
+                                                               else
+                                                                       (pa+(int)(lambda*size[a])*a0mul)->hair_index++;
+                                                       }
+                                               }
+                                       }
+
+                                       if (from==PART_FROM_VOLUME) {
+                                               int in=pa->hair_index%2;
+                                               if (in) pa->hair_index++;
+                                               for (i=0; i<size[0]; i++) {
+                                                       if (in || (pa+i*a0mul)->hair_index%2)
+                                                               (pa+i*a0mul)->flag &= ~PARS_UNEXIST;
+                                                       /* odd intersections == in->out / out->in */
+                                                       /* even intersections -> in stays same */
+                                                       in=(in + (pa+i*a0mul)->hair_index) % 2;
+                                               }
+                                       }
+                               }
+                       }
+               }
+       }
+
+       if (psys->part->flag & PART_GRID_HEXAGONAL) {
+               for (i=0,p=0,pa=psys->particles; i<res; i++) {
+                       for (j=0; j<res; j++) {
+                               for (k=0; k<res; k++,p++,pa++) {
+                                       if (j%2)
+                                               pa->fuv[0] += d/2.f;
+
+                                       if (k%2) {
+                                               pa->fuv[0] += d/2.f;
+                                               pa->fuv[1] += d/2.f;
+                                       }
+                               }
+                       }
+               }
+       }
+
+       if (psys->part->flag & PART_GRID_INVERT) {
+               for (i=0; i<size[0]; i++) {
+                       for (j=0; j<size[1]; j++) {
+                               pa=psys->particles + res*(i*res + j);
+                               for (k=0; k<size[2]; k++, pa++) {
+                                       pa->flag ^= PARS_UNEXIST;
+                               }
+                       }
+               }
+       }
+
+       if (psys->part->grid_rand > 0.f) {
+               float rfac = d * psys->part->grid_rand;
+               for (p=0,pa=psys->particles; p<psys->totpart; p++,pa++) {
+                       if (pa->flag & PARS_UNEXIST)
+                               continue;
+
+                       pa->fuv[0] += rfac * (psys_frand(psys, p + 31) - 0.5f);
+                       pa->fuv[1] += rfac * (psys_frand(psys, p + 32) - 0.5f);
+                       pa->fuv[2] += rfac * (psys_frand(psys, p + 33) - 0.5f);
+               }
+       }
+}
+
+/* modified copy from rayshade.c */
+static void hammersley_create(float *out, int n, int seed, float amount)
+{
+       RNG *rng;
+       double p, t, offs[2];
+       int k, kk;
+
+       rng = BLI_rng_new(31415926 + n + seed);
+       offs[0] = BLI_rng_get_double(rng) + (double)amount;
+       offs[1] = BLI_rng_get_double(rng) + (double)amount;
+       BLI_rng_free(rng);
+
+       for (k = 0; k < n; k++) {
+               t = 0;
+               for (p = 0.5, kk = k; kk; p *= 0.5, kk >>= 1)
+                       if (kk & 1) /* kk mod 2 = 1 */
+                               t += p;
+
+               out[2*k + 0] = fmod((double)k/(double)n + offs[0], 1.0);
+               out[2*k + 1] = fmod(t + offs[1], 1.0);
+       }
+}
+
+/* almost exact copy of BLI_jitter_init */
+static void init_mv_jit(float *jit, int num, int seed2, float amount)
+{
+       RNG *rng;
+       float *jit2, x, rad1, rad2, rad3;
+       int i, num2;
+
+       if (num==0) return;
+
+       rad1= (float)(1.0f/sqrtf((float)num));
+       rad2= (float)(1.0f/((float)num));
+       rad3= (float)sqrtf((float)num)/((float)num);
+
+       rng = BLI_rng_new(31415926 + num + seed2);
+       x= 0;
+               num2 = 2 * num;
+       for (i=0; i<num2; i+=2) {
+       
+               jit[i] = x + amount*rad1*(0.5f - BLI_rng_get_float(rng));
+               jit[i+1] = i/(2.0f*num) + amount*rad1*(0.5f - BLI_rng_get_float(rng));
+               
+               jit[i]-= (float)floor(jit[i]);
+               jit[i+1]-= (float)floor(jit[i+1]);
+               
+               x+= rad3;
+               x -= (float)floor(x);
+       }
+
+       jit2= MEM_mallocN(12 + 2*sizeof(float)*num, "initjit");
+
+       for (i=0 ; i<4 ; i++) {
+               BLI_jitterate1((float (*)[2])jit, (float (*)[2])jit2, num, rad1);
+               BLI_jitterate1((float (*)[2])jit, (float (*)[2])jit2, num, rad1);
+               BLI_jitterate2((float (*)[2])jit, (float (*)[2])jit2, num, rad2);
+       }
+       MEM_freeN(jit2);
+       BLI_rng_free(rng);
+}
+
+static void psys_uv_to_w(float u, float v, int quad, float *w)
+{
+       float vert[4][3], co[3];
+
+       if (!quad) {
+               if (u+v > 1.0f)
+                       v= 1.0f-v;
+               else
+                       u= 1.0f-u;
+       }
+
+       vert[0][0] = 0.0f; vert[0][1] = 0.0f; vert[0][2] = 0.0f;
+       vert[1][0] = 1.0f; vert[1][1] = 0.0f; vert[1][2] = 0.0f;
+       vert[2][0] = 1.0f; vert[2][1] = 1.0f; vert[2][2] = 0.0f;
+
+       co[0] = u;
+       co[1] = v;
+       co[2] = 0.0f;
+
+       if (quad) {
+               vert[3][0] = 0.0f; vert[3][1] = 1.0f; vert[3][2] = 0.0f;
+               interp_weights_poly_v3( w,vert, 4, co);
+       }
+       else {
+               interp_weights_poly_v3( w,vert, 3, co);
+               w[3] = 0.0f;
+       }
+}
+
+/* Find the index in "sum" array before "value" is crossed. */
+static int distribute_binary_search(float *sum, int n, float value)
+{
+       int mid, low=0, high=n;
+
+       if (value == 0.f)
+               return 0;
+
+       while (low <= high) {
+               mid= (low + high)/2;
+               
+               if (sum[mid] < value && value <= sum[mid+1])
+                       return mid;
+               
+               if (sum[mid] >= value)
+                       high= mid - 1;
+               else if (sum[mid] < value)
+                       low= mid + 1;
+               else
+                       return mid;
+       }
+
+       return low;
+}
+
+/* the max number if calls to rng_* funcs within psys_thread_distribute_particle
+ * be sure to keep up to date if this changes */
+#define PSYS_RND_DIST_SKIP 2
+
+/* note: this function must be thread safe, for from == PART_FROM_CHILD */
+#define ONLY_WORKING_WITH_PA_VERTS 0
+static void distribute_threads_exec(ParticleTask *thread, ParticleData *pa, ChildParticle *cpa, int p)
+{
+       ParticleThreadContext *ctx= thread->ctx;
+       Object *ob= ctx->sim.ob;
+       DerivedMesh *dm= ctx->dm;
+       float *v1, *v2, *v3, *v4, nor[3], orco1[3], co1[3], co2[3], nor1[3];
+       float cur_d, min_d, randu, randv;
+       int from= ctx->from;
+       int cfrom= ctx->cfrom;
+       int distr= ctx->distr;
+       int i, intersect, tot;
+       int rng_skip_tot= PSYS_RND_DIST_SKIP; /* count how many rng_* calls wont need skipping */
+
+       if (from == PART_FROM_VERT) {
+               /* TODO_PARTICLE - use original index */
+               pa->num= ctx->index[p];
+               pa->fuv[0] = 1.0f;
+               pa->fuv[1] = pa->fuv[2] = pa->fuv[3] = 0.0;
+
+#if ONLY_WORKING_WITH_PA_VERTS
+               if (ctx->tree) {
+                       KDTreeNearest ptn[3];
+                       int w, maxw;
+
+                       psys_particle_on_dm(ctx->dm,from,pa->num,pa->num_dmcache,pa->fuv,pa->foffset,co1,0,0,0,orco1,0);
+                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco1, 1, 1);
+                       maxw = BLI_kdtree_find_nearest_n(ctx->tree,orco1,ptn,3);
+
+                       for (w=0; w<maxw; w++) {
+                               pa->verts[w]=ptn->num;
+                       }
+               }
+#endif
+       }
+       else if (from == PART_FROM_FACE || from == PART_FROM_VOLUME) {
+               MFace *mface;
+
+               pa->num = i = ctx->index[p];
+               mface = dm->getTessFaceData(dm,i,CD_MFACE);
+               
+               switch (distr) {
+               case PART_DISTR_JIT:
+                       if (ctx->jitlevel == 1) {
+                               if (mface->v4)
+                                       psys_uv_to_w(0.5f, 0.5f, mface->v4, pa->fuv);
+                               else
+                                       psys_uv_to_w(1.0f / 3.0f, 1.0f / 3.0f, mface->v4, pa->fuv);
+                       }
+                       else {
+                               ctx->jitoff[i] = fmod(ctx->jitoff[i],(float)ctx->jitlevel);
+                               if (!isnan(ctx->jitoff[i])) {
+                                       psys_uv_to_w(ctx->jit[2*(int)ctx->jitoff[i]], ctx->jit[2*(int)ctx->jitoff[i]+1], mface->v4, pa->fuv);
+                                       ctx->jitoff[i]++;
+                               }
+                       }
+                       break;
+               case PART_DISTR_RAND:
+                       randu= BLI_rng_get_float(thread->rng);
+                       randv= BLI_rng_get_float(thread->rng);
+                       rng_skip_tot -= 2;
+
+                       psys_uv_to_w(randu, randv, mface->v4, pa->fuv);
+                       break;
+               }
+               pa->foffset= 0.0f;
+               
+               /* experimental */
+               if (from==PART_FROM_VOLUME) {
+                       MVert *mvert=dm->getVertDataArray(dm,CD_MVERT);
+
+                       tot=dm->getNumTessFaces(dm);
+
+                       psys_interpolate_face(mvert,mface,0,0,pa->fuv,co1,nor,0,0,0,0);
+
+                       normalize_v3(nor);
+                       mul_v3_fl(nor,-100.0);
+
+                       add_v3_v3v3(co2,co1,nor);
+
+                       min_d=2.0;
+                       intersect=0;
+
+                       for (i=0,mface=dm->getTessFaceDataArray(dm,CD_MFACE); i<tot; i++,mface++) {
+                               if (i==pa->num) continue;
+
+                               v1=mvert[mface->v1].co;
+                               v2=mvert[mface->v2].co;
+                               v3=mvert[mface->v3].co;
+
+                               if (isect_line_tri_v3(co1, co2, v2, v3, v1, &cur_d, 0)) {
+                                       if (cur_d<min_d) {
+                                               min_d=cur_d;
+                                               pa->foffset=cur_d*50.0f; /* to the middle of volume */
+                                               intersect=1;
+                                       }
+                               }
+                               if (mface->v4) {
+                                       v4=mvert[mface->v4].co;
+
+                                       if (isect_line_tri_v3(co1, co2, v4, v1, v3, &cur_d, 0)) {
+                                               if (cur_d<min_d) {
+                                                       min_d=cur_d;
+                                                       pa->foffset=cur_d*50.0f; /* to the middle of volume */
+                                                       intersect=1;
+                                               }
+                                       }
+                               }
+                       }
+                       if (intersect==0)
+                               pa->foffset=0.0;
+                       else {
+                               switch (distr) {
+                                       case PART_DISTR_JIT:
+                                               pa->foffset *= ctx->jit[p % (2 * ctx->jitlevel)];
+                                               break;
+                                       case PART_DISTR_RAND:
+                                               pa->foffset *= BLI_frand();
+                                               break;
+                               }
+                       }
+               }
+       }
+       else if (from == PART_FROM_CHILD) {
+               MFace *mf;
+
+               if (ctx->index[p] < 0) {
+                       cpa->num=0;
+                       cpa->fuv[0]=cpa->fuv[1]=cpa->fuv[2]=cpa->fuv[3]=0.0f;
+                       cpa->pa[0]=cpa->pa[1]=cpa->pa[2]=cpa->pa[3]=0;
+                       return;
+               }
+
+               mf= dm->getTessFaceData(dm, ctx->index[p], CD_MFACE);
+
+               randu= BLI_rng_get_float(thread->rng);
+               randv= BLI_rng_get_float(thread->rng);
+               rng_skip_tot -= 2;
+
+               psys_uv_to_w(randu, randv, mf->v4, cpa->fuv);
+
+               cpa->num = ctx->index[p];
+
+               if (ctx->tree) {
+                       KDTreeNearest ptn[10];
+                       int w,maxw;//, do_seams;
+                       float maxd /*, mind,dd */, totw= 0.0f;
+                       int parent[10];
+                       float pweight[10];
+
+                       psys_particle_on_dm(dm,cfrom,cpa->num,DMCACHE_ISCHILD,cpa->fuv,cpa->foffset,co1,nor1,NULL,NULL,orco1,NULL);
+                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco1, 1, 1);
+                       maxw = BLI_kdtree_find_nearest_n(ctx->tree,orco1,ptn,3);
+
+                       maxd=ptn[maxw-1].dist;
+                       /* mind=ptn[0].dist; */ /* UNUSED */
+                       
+                       /* the weights here could be done better */
+                       for (w=0; w<maxw; w++) {
+                               parent[w]=ptn[w].index;
+                               pweight[w]=(float)pow(2.0,(double)(-6.0f*ptn[w].dist/maxd));
+                       }
+                       for (;w<10; w++) {
+                               parent[w]=-1;
+                               pweight[w]=0.0f;
+                       }
+
+                       for (w=0,i=0; w<maxw && i<4; w++) {
+                               if (parent[w]>=0) {
+                                       cpa->pa[i]=parent[w];
+                                       cpa->w[i]=pweight[w];
+                                       totw+=pweight[w];
+                                       i++;
+                               }
+                       }
+                       for (;i<4; i++) {
+                               cpa->pa[i]=-1;
+                               cpa->w[i]=0.0f;
+                       }
+
+                       if (totw>0.0f) for (w=0; w<4; w++)
+                               cpa->w[w]/=totw;
+
+                       cpa->parent=cpa->pa[0];
+               }
+       }
+
+       if (rng_skip_tot > 0) /* should never be below zero */
+               BLI_rng_skip(thread->rng, rng_skip_tot);
+}
+
+static void exec_distribute_parent(TaskPool *UNUSED(pool), void *taskdata, int UNUSED(threadid))
+{
+       ParticleTask *task = taskdata;
+       ParticleSystem *psys= task->ctx->sim.psys;
+       ParticleData *pa;
+       int p;
+       
+       pa= psys->particles + task->begin;
+       for (p = task->begin; p < task->end; ++p, ++pa)
+               distribute_threads_exec(task, pa, NULL, p);
+}
+
+static void exec_distribute_child(TaskPool *UNUSED(pool), void *taskdata, int UNUSED(threadid))
+{
+       ParticleTask *task = taskdata;
+       ParticleSystem *psys = task->ctx->sim.psys;
+       ChildParticle *cpa;
+       int p;
+       
+       /* RNG skipping at the beginning */
+       cpa = psys->child;
+       for (p = 0; p < task->begin; ++p, ++cpa) {
+               if (task->ctx->skip) /* simplification skip */
+                       BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP * task->ctx->skip[p]);
+               
+               BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP);
+       }
+               
+       for (; p < task->end; ++p, ++cpa) {
+               if (task->ctx->skip) /* simplification skip */
+                       BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP * task->ctx->skip[p]);
+               
+               distribute_threads_exec(task, NULL, cpa, p);
+       }
+}
+
+static int distribute_compare_orig_index(const void *p1, const void *p2, void *user_data)
+{
+       int *orig_index = (int *) user_data;
+       int index1 = orig_index[*(const int *)p1];
+       int index2 = orig_index[*(const int *)p2];
+
+       if (index1 < index2)
+               return -1;
+       else if (index1 == index2) {
+               /* this pointer comparison appears to make qsort stable for glibc,
+                * and apparently on solaris too, makes the renders reproducible */
+               if (p1 < p2)
+                       return -1;
+               else if (p1 == p2)
+                       return 0;
+               else
+                       return 1;
+       }
+       else
+               return 1;
+}
+
+static void distribute_invalid(Scene *scene, ParticleSystem *psys, int from)
+{
+       if (from == PART_FROM_CHILD) {
+               ChildParticle *cpa;
+               int p, totchild = psys_get_tot_child(scene, psys);
+
+               if (psys->child && totchild) {
+                       for (p=0,cpa=psys->child; p<totchild; p++,cpa++) {
+                               cpa->fuv[0]=cpa->fuv[1]=cpa->fuv[2]=cpa->fuv[3] = 0.0;
+                               cpa->foffset= 0.0f;
+                               cpa->parent=0;
+                               cpa->pa[0]=cpa->pa[1]=cpa->pa[2]=cpa->pa[3]=0;
+                               cpa->num= -1;
+                       }
+               }
+       }
+       else {
+               PARTICLE_P;
+               LOOP_PARTICLES {
+                       pa->fuv[0] = pa->fuv[1] = pa->fuv[2] = pa->fuv[3] = 0.0;
+                       pa->foffset= 0.0f;
+                       pa->num= -1;
+               }
+       }
+}
+
+/* Creates a distribution of coordinates on a DerivedMesh      */
+/* This is to denote functionality that does not yet work with mesh - only derived mesh */
+static int psys_thread_context_init_distribute(ParticleThreadContext *ctx, ParticleSimulationData *sim, int from)
+{
+       Scene *scene = sim->scene;
+       DerivedMesh *finaldm = sim->psmd->dm;
+       Object *ob = sim->ob;
+       ParticleSystem *psys= sim->psys;
+       ParticleData *pa=0, *tpars= 0;
+       ParticleSettings *part;
+       ParticleSeam *seams= 0;
+       KDTree *tree=0;
+       DerivedMesh *dm= NULL;
+       float *jit= NULL;
+       int i, p=0;
+       int cfrom=0;
+       int totelem=0, totpart, *particle_element=0, children=0, totseam=0;
+       int jitlevel= 1, distr;
+       float *element_weight=NULL,*element_sum=NULL,*jitter_offset=NULL, *vweight=NULL;
+       float cur, maxweight=0.0, tweight, totweight, inv_totweight, co[3], nor[3], orco[3];
+       
+       if (ELEM(NULL, ob, psys, psys->part))
+               return 0;
+       
+       part=psys->part;
+       totpart=psys->totpart;
+       if (totpart==0)
+               return 0;
+       
+       if (!finaldm->deformedOnly && !finaldm->getTessFaceDataArray(finaldm, CD_ORIGINDEX)) {
+               printf("Can't create particles with the current modifier stack, disable destructive modifiers\n");
+// XXX         error("Can't paint with the current modifier stack, disable destructive modifiers");
+               return 0;
+       }
+       
+       psys_thread_context_init(ctx, sim);
+       
+       /* First handle special cases */
+       if (from == PART_FROM_CHILD) {
+               /* Simple children */
+               if (part->childtype != PART_CHILD_FACES) {
+                       BLI_srandom(31415926 + psys->seed + psys->child_seed);
+                       distribute_simple_children(scene, ob, finaldm, psys);
+                       return 0;
+               }
+       }
+       else {
+               /* Grid distribution */
+               if (part->distr==PART_DISTR_GRID && from != PART_FROM_VERT) {
+                       BLI_srandom(31415926 + psys->seed);
+                       dm= CDDM_from_mesh((Mesh*)ob->data);
+                       DM_ensure_tessface(dm);
+                       distribute_grid(dm,psys);
+                       dm->release(dm);
+                       return 0;
+               }
+       }
+       
+       /* Create trees and original coordinates if needed */
+       if (from == PART_FROM_CHILD) {
+               distr=PART_DISTR_RAND;
+               BLI_srandom(31415926 + psys->seed + psys->child_seed);
+               dm= finaldm;
+
+               /* BMESH ONLY */
+               DM_ensure_tessface(dm);
+
+               children=1;
+
+               tree=BLI_kdtree_new(totpart);
+
+               for (p=0,pa=psys->particles; p<totpart; p++,pa++) {
+                       psys_particle_on_dm(dm,part->from,pa->num,pa->num_dmcache,pa->fuv,pa->foffset,co,nor,0,0,orco,NULL);
+                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco, 1, 1);
+                       BLI_kdtree_insert(tree, p, orco);
+               }
+
+               BLI_kdtree_balance(tree);
+
+               totpart = psys_get_tot_child(scene, psys);
+               cfrom = from = PART_FROM_FACE;
+       }
+       else {
+               distr = part->distr;
+               BLI_srandom(31415926 + psys->seed);
+               
+               if (psys->part->use_modifier_stack)
+                       dm = finaldm;
+               else
+                       dm= CDDM_from_mesh((Mesh*)ob->data);
+
+               /* BMESH ONLY, for verts we don't care about tessfaces */
+               if (from != PART_FROM_VERT) {
+                       DM_ensure_tessface(dm);
+               }
+
+               /* we need orco for consistent distributions */
+               if (!CustomData_has_layer(&dm->vertData, CD_ORCO))
+                       DM_add_vert_layer(dm, CD_ORCO, CD_ASSIGN, BKE_mesh_orco_verts_get(ob));
+
+               if (from == PART_FROM_VERT) {
+                       MVert *mv= dm->getVertDataArray(dm, CD_MVERT);
+                       float (*orcodata)[3] = dm->getVertDataArray(dm, CD_ORCO);
+                       int totvert = dm->getNumVerts(dm);
+
+                       tree=BLI_kdtree_new(totvert);
+
+                       for (p=0; p<totvert; p++) {
+                               if (orcodata) {
+                                       copy_v3_v3(co,orcodata[p]);
+                                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co, 1, 1);
+                               }
+                               else
+                                       copy_v3_v3(co,mv[p].co);
+                               BLI_kdtree_insert(tree, p, co);
+                       }
+
+                       BLI_kdtree_balance(tree);
+               }
+       }
+
+       /* Get total number of emission elements and allocate needed arrays */
+       totelem = (from == PART_FROM_VERT) ? dm->getNumVerts(dm) : dm->getNumTessFaces(dm);
+
+       if (totelem == 0) {
+               distribute_invalid(scene, psys, children ? PART_FROM_CHILD : 0);
+
+               if (G.debug & G_DEBUG)
+                       fprintf(stderr,"Particle distribution error: Nothing to emit from!\n");
+
+               if (dm != finaldm) dm->release(dm);
+
+               BLI_kdtree_free(tree);
+
+               return 0;
+       }
+
+       element_weight  = MEM_callocN(sizeof(float)*totelem, "particle_distribution_weights");
+       particle_element= MEM_callocN(sizeof(int)*totpart, "particle_distribution_indexes");
+       element_sum             = MEM_callocN(sizeof(float)*(totelem+1), "particle_distribution_sum");
+       jitter_offset   = MEM_callocN(sizeof(float)*totelem, "particle_distribution_jitoff");
+
+       /* Calculate weights from face areas */
+       if ((part->flag&PART_EDISTR || children) && from != PART_FROM_VERT) {
+               MVert *v1, *v2, *v3, *v4;
+               float totarea=0.f, co1[3], co2[3], co3[3], co4[3];
+               float (*orcodata)[3];
+               
+               orcodata= dm->getVertDataArray(dm, CD_ORCO);
+
+               for (i=0; i<totelem; i++) {
+                       MFace *mf=dm->getTessFaceData(dm,i,CD_MFACE);
+
+                       if (orcodata) {
+                               copy_v3_v3(co1, orcodata[mf->v1]);
+                               copy_v3_v3(co2, orcodata[mf->v2]);
+                               copy_v3_v3(co3, orcodata[mf->v3]);
+                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co1, 1, 1);
+                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co2, 1, 1);
+                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co3, 1, 1);
+                               if (mf->v4) {
+                                       copy_v3_v3(co4, orcodata[mf->v4]);
+                                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co4, 1, 1);
+                               }
+                       }
+                       else {
+                               v1= (MVert*)dm->getVertData(dm,mf->v1,CD_MVERT);
+                               v2= (MVert*)dm->getVertData(dm,mf->v2,CD_MVERT);
+                               v3= (MVert*)dm->getVertData(dm,mf->v3,CD_MVERT);
+                               copy_v3_v3(co1, v1->co);
+                               copy_v3_v3(co2, v2->co);
+                               copy_v3_v3(co3, v3->co);
+                               if (mf->v4) {
+                                       v4= (MVert*)dm->getVertData(dm,mf->v4,CD_MVERT);
+                                       copy_v3_v3(co4, v4->co);
+                               }
+                       }
+
+                       cur = mf->v4 ? area_quad_v3(co1, co2, co3, co4) : area_tri_v3(co1, co2, co3);
+                       
+                       if (cur > maxweight)
+                               maxweight = cur;
+
+                       element_weight[i] = cur;
+                       totarea += cur;
+               }
+
+               for (i=0; i<totelem; i++)
+                       element_weight[i] /= totarea;
+
+               maxweight /= totarea;
+       }
+       else {
+               float min=1.0f/(float)(MIN2(totelem,totpart));
+               for (i=0; i<totelem; i++)
+                       element_weight[i]=min;
+               maxweight=min;
+       }
+
+       /* Calculate weights from vgroup */
+       vweight = psys_cache_vgroup(dm,psys,PSYS_VG_DENSITY);
+
+       if (vweight) {
+               if (from==PART_FROM_VERT) {
+                       for (i=0;i<totelem; i++)
+                               element_weight[i]*=vweight[i];
+               }
+               else { /* PART_FROM_FACE / PART_FROM_VOLUME */
+                       for (i=0;i<totelem; i++) {
+                               MFace *mf=dm->getTessFaceData(dm,i,CD_MFACE);
+                               tweight = vweight[mf->v1] + vweight[mf->v2] + vweight[mf->v3];
+                               
+                               if (mf->v4) {
+                                       tweight += vweight[mf->v4];
+                                       tweight /= 4.0f;
+                               }
+                               else {
+                                       tweight /= 3.0f;
+                               }
+
+                               element_weight[i]*=tweight;
+                       }
+               }
+               MEM_freeN(vweight);
+       }
+
+       /* Calculate total weight of all elements */
+       totweight= 0.0f;
+       for (i=0;i<totelem; i++)
+               totweight += element_weight[i];
+
+       inv_totweight = (totweight > 0.f ? 1.f/totweight : 0.f);
+
+       /* Calculate cumulative weights */
+       element_sum[0] = 0.0f;
+       for (i=0; i<totelem; i++)
+               element_sum[i+1] = element_sum[i] + element_weight[i] * inv_totweight;
+       
+       /* Finally assign elements to particles */
+       if ((part->flag&PART_TRAND) || (part->simplify_flag&PART_SIMPLIFY_ENABLE)) {
+               float pos;
+
+               for (p=0; p<totpart; p++) {
+                       /* In theory element_sum[totelem] should be 1.0, but due to float errors this is not necessarily always true, so scale pos accordingly. */
+                       pos= BLI_frand() * element_sum[totelem];
+                       particle_element[p] = distribute_binary_search(element_sum, totelem, pos);
+                       particle_element[p] = MIN2(totelem-1, particle_element[p]);
+                       jitter_offset[particle_element[p]] = pos;
+               }
+       }
+       else {
+               double step, pos;
+               
+               step= (totpart < 2) ? 0.5 : 1.0/(double)totpart;
+               pos= 1e-6; /* tiny offset to avoid zero weight face */
+               i= 0;
+
+               for (p=0; p<totpart; p++, pos+=step) {
+                       while ((i < totelem) && (pos > (double)element_sum[i + 1]))
+                               i++;
+
+                       particle_element[p] = MIN2(totelem-1, i);
+
+                       /* avoid zero weight face */
+                       if (p == totpart-1 && element_weight[particle_element[p]] == 0.0f)
+                               particle_element[p] = particle_element[p-1];
+
+                       jitter_offset[particle_element[p]] = pos;
+               }
+       }
+
+       MEM_freeN(element_sum);
+
+       /* For hair, sort by origindex (allows optimization's in rendering), */
+       /* however with virtual parents the children need to be in random order. */
+       if (part->type == PART_HAIR && !(part->childtype==PART_CHILD_FACES && part->parents!=0.0f)) {
+               int *orig_index = NULL;
+
+               if (from == PART_FROM_VERT) {
+                       if (dm->numVertData)
+                               orig_index = dm->getVertDataArray(dm, CD_ORIGINDEX);
+               }
+               else {
+                       if (dm->numTessFaceData)
+                               orig_index = dm->getTessFaceDataArray(dm, CD_ORIGINDEX);
+               }
+
+               if (orig_index) {
+                       BLI_qsort_r(particle_element, totpart, sizeof(int), distribute_compare_orig_index, orig_index);
+               }
+       }
+
+       /* Create jittering if needed */
+       if (distr==PART_DISTR_JIT && ELEM(from,PART_FROM_FACE,PART_FROM_VOLUME)) {
+               jitlevel= part->userjit;
+               
+               if (jitlevel == 0) {
+                       jitlevel= totpart/totelem;
+                       if (part->flag & PART_EDISTR) jitlevel*= 2;     /* looks better in general, not very scietific */
+                       if (jitlevel<3) jitlevel= 3;
+               }
+               
+               jit= MEM_callocN((2+ jitlevel*2)*sizeof(float), "jit");
+
+               /* for small amounts of particles we use regular jitter since it looks
+                * a bit better, for larger amounts we switch to hammersley sequence 
+                * because it is much faster */
+               if (jitlevel < 25)
+                       init_mv_jit(jit, jitlevel, psys->seed, part->jitfac);
+               else
+                       hammersley_create(jit, jitlevel+1, psys->seed, part->jitfac);
+               BLI_array_randomize(jit, 2*sizeof(float), jitlevel, psys->seed); /* for custom jit or even distribution */
+       }
+
+       /* Setup things for threaded distribution */
+       ctx->tree= tree;
+       ctx->seams= seams;
+       ctx->totseam= totseam;
+       ctx->sim.psys= psys;
+       ctx->index= particle_element;
+       ctx->jit= jit;
+       ctx->jitlevel= jitlevel;
+       ctx->jitoff= jitter_offset;
+       ctx->weight= element_weight;
+       ctx->maxweight= maxweight;
+       ctx->from= (children) ? PART_FROM_CHILD : from;
+       ctx->cfrom= cfrom;
+       ctx->distr= distr;
+       ctx->dm= dm;
+       ctx->tpars= tpars;
+
+       if (children) {
+               totpart= psys_render_simplify_distribution(ctx, totpart);
+               alloc_child_particles(psys, totpart);
+       }
+
+       return 1;
+}
+
+static void psys_task_init_distribute(ParticleTask *task, ParticleSimulationData *sim)
+{
+       /* init random number generator */
+       int seed = 31415926 + sim->psys->seed;
+       
+       task->rng = BLI_rng_new(seed);
+}
+
+static void distribute_particles_on_dm(ParticleSimulationData *sim, int from)
+{
+       TaskScheduler *task_scheduler;
+       TaskPool *task_pool;
+       ParticleThreadContext ctx;
+       ParticleTask *tasks;
+       DerivedMesh *finaldm = sim->psmd->dm;
+       int i, totpart, numtasks;
+       
+       /* create a task pool for distribution tasks */
+       if (!psys_thread_context_init_distribute(&ctx, sim, from))
+               return;
+       
+       task_scheduler = BLI_task_scheduler_get();
+       task_pool = BLI_task_pool_create(task_scheduler, &ctx);
+       
+       totpart = (from == PART_FROM_CHILD ? sim->psys->totchild : sim->psys->totpart);
+       psys_tasks_create(&ctx, totpart, &tasks, &numtasks);
+       for (i = 0; i < numtasks; ++i) {
+               ParticleTask *task = &tasks[i];
+               
+               psys_task_init_distribute(task, sim);
+               if (from == PART_FROM_CHILD)
+                       BLI_task_pool_push(task_pool, exec_distribute_child, task, false, TASK_PRIORITY_LOW);
+               else
+                       BLI_task_pool_push(task_pool, exec_distribute_parent, task, false, TASK_PRIORITY_LOW);
+       }
+       BLI_task_pool_work_and_wait(task_pool);
+       
+       BLI_task_pool_free(task_pool);
+       
+       psys_calc_dmcache(sim->ob, finaldm, sim->psys);
+       
+       if (ctx.dm != finaldm)
+               ctx.dm->release(ctx.dm);
+       
+       psys_tasks_free(tasks, numtasks);
+}
+
+/* ready for future use, to emit particles without geometry */
+static void distribute_particles_on_shape(ParticleSimulationData *sim, int UNUSED(from))
+{
+       distribute_invalid(sim->scene, sim->psys, 0);
+
+       fprintf(stderr,"Shape emission not yet possible!\n");
+}
+
+void distribute_particles(ParticleSimulationData *sim, int from)
+{
+       PARTICLE_PSMD;
+       int distr_error=0;
+
+       if (psmd) {
+               if (psmd->dm)
+                       distribute_particles_on_dm(sim, from);
+               else
+                       distr_error=1;
+       }
+       else
+               distribute_particles_on_shape(sim, from);
+
+       if (distr_error) {
+               distribute_invalid(sim->scene, sim->psys, from);
+
+               fprintf(stderr,"Particle distribution error!\n");
+       }
+}
index 877c2ccdec9875e9363b9dd562d2568ff5cf4ecb..6956766f5eac628d30f737ad433f52be11783a2e 100644 (file)
@@ -285,7 +285,7 @@ static void realloc_particles(ParticleSimulationData *sim, int new_totpart)
        }
 }
 
-static int get_psys_child_number(struct Scene *scene, ParticleSystem *psys)
+int psys_get_child_number(Scene *scene, ParticleSystem *psys)
 {
        int nbr;
 
@@ -300,30 +300,9 @@ static int get_psys_child_number(struct Scene *scene, ParticleSystem *psys)
        return get_render_child_particle_number(&scene->r, nbr);
 }
 
-static int get_psys_tot_child(struct Scene *scene, ParticleSystem *psys)
+int psys_get_tot_child(Scene *scene, ParticleSystem *psys)
 {
-       return psys->totpart*get_psys_child_number(scene, psys);
-}
-
-static void alloc_child_particles(ParticleSystem *psys, int tot)
-{
-       if (psys->child) {
-               /* only re-allocate if we have to */
-               if (psys->part->childtype && psys->totchild == tot) {
-                       memset(psys->child, 0, tot*sizeof(ChildParticle));
-                       return;
-               }
-
-               MEM_freeN(psys->child);
-               psys->child=NULL;
-               psys->totchild=0;
-       }
-
-       if (psys->part->childtype) {
-               psys->totchild= tot;
-               if (psys->totchild)
-                       psys->child= MEM_callocN(psys->totchild*sizeof(ChildParticle), "child_particles");
-       }
+       return psys->totpart*psys_get_child_number(scene, psys);
 }
 
 /************************************************/
@@ -453,1032 +432,6 @@ void psys_calc_dmcache(Object *ob, DerivedMesh *dm, ParticleSystem *psys)
        }
 }
 
-static void distribute_simple_children(Scene *scene, Object *ob, DerivedMesh *finaldm, ParticleSystem *psys)
-{
-       ChildParticle *cpa = NULL;
-       int i, p;
-       int child_nbr= get_psys_child_number(scene, psys);
-       int totpart= get_psys_tot_child(scene, psys);
-
-       alloc_child_particles(psys, totpart);
-
-       cpa = psys->child;
-       for (i=0; i<child_nbr; i++) {
-               for (p=0; p<psys->totpart; p++,cpa++) {
-                       float length=2.0;
-                       cpa->parent=p;
-                                       
-                       /* create even spherical distribution inside unit sphere */
-                       while (length>=1.0f) {
-                               cpa->fuv[0]=2.0f*BLI_frand()-1.0f;
-                               cpa->fuv[1]=2.0f*BLI_frand()-1.0f;
-                               cpa->fuv[2]=2.0f*BLI_frand()-1.0f;
-                               length=len_v3(cpa->fuv);
-                       }
-
-                       cpa->num=-1;
-               }
-       }
-       /* dmcache must be updated for parent particles if children from faces is used */
-       psys_calc_dmcache(ob, finaldm, psys);
-}
-static void distribute_grid(DerivedMesh *dm, ParticleSystem *psys)
-{
-       ParticleData *pa=NULL;
-       float min[3], max[3], delta[3], d;
-       MVert *mv, *mvert = dm->getVertDataArray(dm,0);
-       int totvert=dm->getNumVerts(dm), from=psys->part->from;
-       int i, j, k, p, res=psys->part->grid_res, size[3], axis;
-
-       /* find bounding box of dm */
-       if (totvert > 0) {
-               mv=mvert;
-               copy_v3_v3(min, mv->co);
-               copy_v3_v3(max, mv->co);
-               mv++;
-               for (i = 1; i < totvert; i++, mv++) {
-                       minmax_v3v3_v3(min, max, mv->co);
-               }
-       }
-       else {
-               zero_v3(min);
-               zero_v3(max);
-       }
-
-       sub_v3_v3v3(delta, max, min);
-
-       /* determine major axis */
-       axis = axis_dominant_v3_single(delta);
-        
-       d = delta[axis]/(float)res;
-
-       size[axis] = res;
-       size[(axis+1)%3] = (int)ceil(delta[(axis+1)%3]/d);
-       size[(axis+2)%3] = (int)ceil(delta[(axis+2)%3]/d);
-
-       /* float errors grrr.. */
-       size[(axis+1)%3] = MIN2(size[(axis+1)%3],res);
-       size[(axis+2)%3] = MIN2(size[(axis+2)%3],res);
-
-       size[0] = MAX2(size[0], 1);
-       size[1] = MAX2(size[1], 1);
-       size[2] = MAX2(size[2], 1);
-
-       /* no full offset for flat/thin objects */
-       min[0]+= d < delta[0] ? d/2.f : delta[0]/2.f;
-       min[1]+= d < delta[1] ? d/2.f : delta[1]/2.f;
-       min[2]+= d < delta[2] ? d/2.f : delta[2]/2.f;
-
-       for (i=0,p=0,pa=psys->particles; i<res; i++) {
-               for (j=0; j<res; j++) {
-                       for (k=0; k<res; k++,p++,pa++) {
-                               pa->fuv[0] = min[0] + (float)i*d;
-                               pa->fuv[1] = min[1] + (float)j*d;
-                               pa->fuv[2] = min[2] + (float)k*d;
-                               pa->flag |= PARS_UNEXIST;
-                               pa->hair_index = 0; /* abused in volume calculation */
-                       }
-               }
-       }
-
-       /* enable particles near verts/edges/faces/inside surface */
-       if (from==PART_FROM_VERT) {
-               float vec[3];
-
-               pa=psys->particles;
-
-               min[0] -= d/2.0f;
-               min[1] -= d/2.0f;
-               min[2] -= d/2.0f;
-
-               for (i=0,mv=mvert; i<totvert; i++,mv++) {
-                       sub_v3_v3v3(vec,mv->co,min);
-                       vec[0]/=delta[0];
-                       vec[1]/=delta[1];
-                       vec[2]/=delta[2];
-                       pa[((int)(vec[0] * (size[0] - 1))  * res +
-                           (int)(vec[1] * (size[1] - 1))) * res +
-                           (int)(vec[2] * (size[2] - 1))].flag &= ~PARS_UNEXIST;
-               }
-       }
-       else if (ELEM(from,PART_FROM_FACE,PART_FROM_VOLUME)) {
-               float co1[3], co2[3];
-
-               MFace *mface= NULL, *mface_array;
-               float v1[3], v2[3], v3[3], v4[4], lambda;
-               int a, a1, a2, a0mul, a1mul, a2mul, totface;
-               int amax= from==PART_FROM_FACE ? 3 : 1;
-
-               totface=dm->getNumTessFaces(dm);
-               mface=mface_array=dm->getTessFaceDataArray(dm,CD_MFACE);
-               
-               for (a=0; a<amax; a++) {
-                       if (a==0) { a0mul=res*res; a1mul=res; a2mul=1; }
-                       else if (a==1) { a0mul=res; a1mul=1; a2mul=res*res; }
-                       else { a0mul=1; a1mul=res*res; a2mul=res; }
-
-                       for (a1=0; a1<size[(a+1)%3]; a1++) {
-                               for (a2=0; a2<size[(a+2)%3]; a2++) {
-                                       mface= mface_array;
-
-                                       pa = psys->particles + a1*a1mul + a2*a2mul;
-                                       copy_v3_v3(co1, pa->fuv);
-                                       co1[a] -= d < delta[a] ? d/2.f : delta[a]/2.f;
-                                       copy_v3_v3(co2, co1);
-                                       co2[a] += delta[a] + 0.001f*d;
-                                       co1[a] -= 0.001f*d;
-                                       
-                                       /* lets intersect the faces */
-                                       for (i=0; i<totface; i++,mface++) {
-                                               copy_v3_v3(v1, mvert[mface->v1].co);
-                                               copy_v3_v3(v2, mvert[mface->v2].co);
-                                               copy_v3_v3(v3, mvert[mface->v3].co);
-
-                                               if (isect_axial_line_tri_v3(a, co1, co2, v2, v3, v1, &lambda)) {
-                                                       if (from==PART_FROM_FACE)
-                                                               (pa+(int)(lambda*size[a])*a0mul)->flag &= ~PARS_UNEXIST;
-                                                       else /* store number of intersections */
-                                                               (pa+(int)(lambda*size[a])*a0mul)->hair_index++;
-                                               }
-                                               else if (mface->v4) {
-                                                       copy_v3_v3(v4, mvert[mface->v4].co);
-
-                                                       if (isect_axial_line_tri_v3(a, co1, co2, v4, v1, v3, &lambda)) {
-                                                               if (from==PART_FROM_FACE)
-                                                                       (pa+(int)(lambda*size[a])*a0mul)->flag &= ~PARS_UNEXIST;
-                                                               else
-                                                                       (pa+(int)(lambda*size[a])*a0mul)->hair_index++;
-                                                       }
-                                               }
-                                       }
-
-                                       if (from==PART_FROM_VOLUME) {
-                                               int in=pa->hair_index%2;
-                                               if (in) pa->hair_index++;
-                                               for (i=0; i<size[0]; i++) {
-                                                       if (in || (pa+i*a0mul)->hair_index%2)
-                                                               (pa+i*a0mul)->flag &= ~PARS_UNEXIST;
-                                                       /* odd intersections == in->out / out->in */
-                                                       /* even intersections -> in stays same */
-                                                       in=(in + (pa+i*a0mul)->hair_index) % 2;
-                                               }
-                                       }
-                               }
-                       }
-               }
-       }
-
-       if (psys->part->flag & PART_GRID_HEXAGONAL) {
-               for (i=0,p=0,pa=psys->particles; i<res; i++) {
-                       for (j=0; j<res; j++) {
-                               for (k=0; k<res; k++,p++,pa++) {
-                                       if (j%2)
-                                               pa->fuv[0] += d/2.f;
-
-                                       if (k%2) {
-                                               pa->fuv[0] += d/2.f;
-                                               pa->fuv[1] += d/2.f;
-                                       }
-                               }
-                       }
-               }
-       }
-
-       if (psys->part->flag & PART_GRID_INVERT) {
-               for (i=0; i<size[0]; i++) {
-                       for (j=0; j<size[1]; j++) {
-                               pa=psys->particles + res*(i*res + j);
-                               for (k=0; k<size[2]; k++, pa++) {
-                                       pa->flag ^= PARS_UNEXIST;
-                               }
-                       }
-               }
-       }
-
-       if (psys->part->grid_rand > 0.f) {
-               float rfac = d * psys->part->grid_rand;
-               for (p=0,pa=psys->particles; p<psys->totpart; p++,pa++) {
-                       if (pa->flag & PARS_UNEXIST)
-                               continue;
-
-                       pa->fuv[0] += rfac * (psys_frand(psys, p + 31) - 0.5f);
-                       pa->fuv[1] += rfac * (psys_frand(psys, p + 32) - 0.5f);
-                       pa->fuv[2] += rfac * (psys_frand(psys, p + 33) - 0.5f);
-               }
-       }
-}
-
-/* modified copy from rayshade.c */
-static void hammersley_create(float *out, int n, int seed, float amount)
-{
-       RNG *rng;
-       double p, t, offs[2];
-       int k, kk;
-
-       rng = BLI_rng_new(31415926 + n + seed);
-       offs[0] = BLI_rng_get_double(rng) + (double)amount;
-       offs[1] = BLI_rng_get_double(rng) + (double)amount;
-       BLI_rng_free(rng);
-
-       for (k = 0; k < n; k++) {
-               t = 0;
-               for (p = 0.5, kk = k; kk; p *= 0.5, kk >>= 1)
-                       if (kk & 1) /* kk mod 2 = 1 */
-                               t += p;
-
-               out[2*k + 0] = fmod((double)k/(double)n + offs[0], 1.0);
-               out[2*k + 1] = fmod(t + offs[1], 1.0);
-       }
-}
-
-/* almost exact copy of BLI_jitter_init */
-static void init_mv_jit(float *jit, int num, int seed2, float amount)
-{
-       RNG *rng;
-       float *jit2, x, rad1, rad2, rad3;
-       int i, num2;
-
-       if (num==0) return;
-
-       rad1= (float)(1.0f/sqrtf((float)num));
-       rad2= (float)(1.0f/((float)num));
-       rad3= (float)sqrtf((float)num)/((float)num);
-
-       rng = BLI_rng_new(31415926 + num + seed2);
-       x= 0;
-               num2 = 2 * num;
-       for (i=0; i<num2; i+=2) {
-       
-               jit[i] = x + amount*rad1*(0.5f - BLI_rng_get_float(rng));
-               jit[i+1] = i/(2.0f*num) + amount*rad1*(0.5f - BLI_rng_get_float(rng));
-               
-               jit[i]-= (float)floor(jit[i]);
-               jit[i+1]-= (float)floor(jit[i+1]);
-               
-               x+= rad3;
-               x -= (float)floor(x);
-       }
-
-       jit2= MEM_mallocN(12 + 2*sizeof(float)*num, "initjit");
-
-       for (i=0 ; i<4 ; i++) {
-               BLI_jitterate1((float (*)[2])jit, (float (*)[2])jit2, num, rad1);
-               BLI_jitterate1((float (*)[2])jit, (float (*)[2])jit2, num, rad1);
-               BLI_jitterate2((float (*)[2])jit, (float (*)[2])jit2, num, rad2);
-       }
-       MEM_freeN(jit2);
-       BLI_rng_free(rng);
-}
-
-static void psys_uv_to_w(float u, float v, int quad, float *w)
-{
-       float vert[4][3], co[3];
-
-       if (!quad) {
-               if (u+v > 1.0f)
-                       v= 1.0f-v;
-               else
-                       u= 1.0f-u;
-       }
-
-       vert[0][0] = 0.0f; vert[0][1] = 0.0f; vert[0][2] = 0.0f;
-       vert[1][0] = 1.0f; vert[1][1] = 0.0f; vert[1][2] = 0.0f;
-       vert[2][0] = 1.0f; vert[2][1] = 1.0f; vert[2][2] = 0.0f;
-
-       co[0] = u;
-       co[1] = v;
-       co[2] = 0.0f;
-
-       if (quad) {
-               vert[3][0] = 0.0f; vert[3][1] = 1.0f; vert[3][2] = 0.0f;
-               interp_weights_poly_v3( w,vert, 4, co);
-       }
-       else {
-               interp_weights_poly_v3( w,vert, 3, co);
-               w[3] = 0.0f;
-       }
-}
-
-/* Find the index in "sum" array before "value" is crossed. */
-static int distribute_binary_search(float *sum, int n, float value)
-{
-       int mid, low=0, high=n;
-
-       if (value == 0.f)
-               return 0;
-
-       while (low <= high) {
-               mid= (low + high)/2;
-               
-               if (sum[mid] < value && value <= sum[mid+1])
-                       return mid;
-               
-               if (sum[mid] >= value)
-                       high= mid - 1;
-               else if (sum[mid] < value)
-                       low= mid + 1;
-               else
-                       return mid;
-       }
-
-       return low;
-}
-
-/* the max number if calls to rng_* funcs within psys_thread_distribute_particle
- * be sure to keep up to date if this changes */
-#define PSYS_RND_DIST_SKIP 2
-
-/* note: this function must be thread safe, for from == PART_FROM_CHILD */
-#define ONLY_WORKING_WITH_PA_VERTS 0
-static void distribute_threads_exec(ParticleTask *thread, ParticleData *pa, ChildParticle *cpa, int p)
-{
-       ParticleThreadContext *ctx= thread->ctx;
-       Object *ob= ctx->sim.ob;
-       DerivedMesh *dm= ctx->dm;
-       float *v1, *v2, *v3, *v4, nor[3], orco1[3], co1[3], co2[3], nor1[3];
-       float cur_d, min_d, randu, randv;
-       int from= ctx->from;
-       int cfrom= ctx->cfrom;
-       int distr= ctx->distr;
-       int i, intersect, tot;
-       int rng_skip_tot= PSYS_RND_DIST_SKIP; /* count how many rng_* calls wont need skipping */
-
-       if (from == PART_FROM_VERT) {
-               /* TODO_PARTICLE - use original index */
-               pa->num= ctx->index[p];
-               pa->fuv[0] = 1.0f;
-               pa->fuv[1] = pa->fuv[2] = pa->fuv[3] = 0.0;
-
-#if ONLY_WORKING_WITH_PA_VERTS
-               if (ctx->tree) {
-                       KDTreeNearest ptn[3];
-                       int w, maxw;
-
-                       psys_particle_on_dm(ctx->dm,from,pa->num,pa->num_dmcache,pa->fuv,pa->foffset,co1,0,0,0,orco1,0);
-                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco1, 1, 1);
-                       maxw = BLI_kdtree_find_nearest_n(ctx->tree,orco1,ptn,3);
-
-                       for (w=0; w<maxw; w++) {
-                               pa->verts[w]=ptn->num;
-                       }
-               }
-#endif
-       }
-       else if (from == PART_FROM_FACE || from == PART_FROM_VOLUME) {
-               MFace *mface;
-
-               pa->num = i = ctx->index[p];
-               mface = dm->getTessFaceData(dm,i,CD_MFACE);
-               
-               switch (distr) {
-               case PART_DISTR_JIT:
-                       if (ctx->jitlevel == 1) {
-                               if (mface->v4)
-                                       psys_uv_to_w(0.5f, 0.5f, mface->v4, pa->fuv);
-                               else
-                                       psys_uv_to_w(1.0f / 3.0f, 1.0f / 3.0f, mface->v4, pa->fuv);
-                       }
-                       else {
-                               ctx->jitoff[i] = fmod(ctx->jitoff[i],(float)ctx->jitlevel);
-                               if (!isnan(ctx->jitoff[i])) {
-                                       psys_uv_to_w(ctx->jit[2*(int)ctx->jitoff[i]], ctx->jit[2*(int)ctx->jitoff[i]+1], mface->v4, pa->fuv);
-                                       ctx->jitoff[i]++;
-                               }
-                       }
-                       break;
-               case PART_DISTR_RAND:
-                       randu= BLI_rng_get_float(thread->rng);
-                       randv= BLI_rng_get_float(thread->rng);
-                       rng_skip_tot -= 2;
-
-                       psys_uv_to_w(randu, randv, mface->v4, pa->fuv);
-                       break;
-               }
-               pa->foffset= 0.0f;
-               
-               /* experimental */
-               if (from==PART_FROM_VOLUME) {
-                       MVert *mvert=dm->getVertDataArray(dm,CD_MVERT);
-
-                       tot=dm->getNumTessFaces(dm);
-
-                       psys_interpolate_face(mvert,mface,0,0,pa->fuv,co1,nor,0,0,0,0);
-
-                       normalize_v3(nor);
-                       mul_v3_fl(nor,-100.0);
-
-                       add_v3_v3v3(co2,co1,nor);
-
-                       min_d=2.0;
-                       intersect=0;
-
-                       for (i=0,mface=dm->getTessFaceDataArray(dm,CD_MFACE); i<tot; i++,mface++) {
-                               if (i==pa->num) continue;
-
-                               v1=mvert[mface->v1].co;
-                               v2=mvert[mface->v2].co;
-                               v3=mvert[mface->v3].co;
-
-                               if (isect_line_tri_v3(co1, co2, v2, v3, v1, &cur_d, 0)) {
-                                       if (cur_d<min_d) {
-                                               min_d=cur_d;
-                                               pa->foffset=cur_d*50.0f; /* to the middle of volume */
-                                               intersect=1;
-                                       }
-                               }
-                               if (mface->v4) {
-                                       v4=mvert[mface->v4].co;
-
-                                       if (isect_line_tri_v3(co1, co2, v4, v1, v3, &cur_d, 0)) {
-                                               if (cur_d<min_d) {
-                                                       min_d=cur_d;
-                                                       pa->foffset=cur_d*50.0f; /* to the middle of volume */
-                                                       intersect=1;
-                                               }
-                                       }
-                               }
-                       }
-                       if (intersect==0)
-                               pa->foffset=0.0;
-                       else {
-                               switch (distr) {
-                                       case PART_DISTR_JIT:
-                                               pa->foffset *= ctx->jit[p % (2 * ctx->jitlevel)];
-                                               break;
-                                       case PART_DISTR_RAND:
-                                               pa->foffset *= BLI_frand();
-                                               break;
-                               }
-                       }
-               }
-       }
-       else if (from == PART_FROM_CHILD) {
-               MFace *mf;
-
-               if (ctx->index[p] < 0) {
-                       cpa->num=0;
-                       cpa->fuv[0]=cpa->fuv[1]=cpa->fuv[2]=cpa->fuv[3]=0.0f;
-                       cpa->pa[0]=cpa->pa[1]=cpa->pa[2]=cpa->pa[3]=0;
-                       return;
-               }
-
-               mf= dm->getTessFaceData(dm, ctx->index[p], CD_MFACE);
-
-               randu= BLI_rng_get_float(thread->rng);
-               randv= BLI_rng_get_float(thread->rng);
-               rng_skip_tot -= 2;
-
-               psys_uv_to_w(randu, randv, mf->v4, cpa->fuv);
-
-               cpa->num = ctx->index[p];
-
-               if (ctx->tree) {
-                       KDTreeNearest ptn[10];
-                       int w,maxw;//, do_seams;
-                       float maxd /*, mind,dd */, totw= 0.0f;
-                       int parent[10];
-                       float pweight[10];
-
-                       psys_particle_on_dm(dm,cfrom,cpa->num,DMCACHE_ISCHILD,cpa->fuv,cpa->foffset,co1,nor1,NULL,NULL,orco1,NULL);
-                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco1, 1, 1);
-                       maxw = BLI_kdtree_find_nearest_n(ctx->tree,orco1,ptn,3);
-
-                       maxd=ptn[maxw-1].dist;
-                       /* mind=ptn[0].dist; */ /* UNUSED */
-                       
-                       /* the weights here could be done better */
-                       for (w=0; w<maxw; w++) {
-                               parent[w]=ptn[w].index;
-                               pweight[w]=(float)pow(2.0,(double)(-6.0f*ptn[w].dist/maxd));
-                       }
-                       for (;w<10; w++) {
-                               parent[w]=-1;
-                               pweight[w]=0.0f;
-                       }
-
-                       for (w=0,i=0; w<maxw && i<4; w++) {
-                               if (parent[w]>=0) {
-                                       cpa->pa[i]=parent[w];
-                                       cpa->w[i]=pweight[w];
-                                       totw+=pweight[w];
-                                       i++;
-                               }
-                       }
-                       for (;i<4; i++) {
-                               cpa->pa[i]=-1;
-                               cpa->w[i]=0.0f;
-                       }
-
-                       if (totw>0.0f) for (w=0; w<4; w++)
-                               cpa->w[w]/=totw;
-
-                       cpa->parent=cpa->pa[0];
-               }
-       }
-
-       if (rng_skip_tot > 0) /* should never be below zero */
-               BLI_rng_skip(thread->rng, rng_skip_tot);
-}
-
-static void exec_distribute_parent(TaskPool *UNUSED(pool), void *taskdata, int UNUSED(threadid))
-{
-       ParticleTask *task = taskdata;
-       ParticleSystem *psys= task->ctx->sim.psys;
-       ParticleData *pa;
-       int p;
-       
-       pa= psys->particles + task->begin;
-       for (p = task->begin; p < task->end; ++p, ++pa)
-               distribute_threads_exec(task, pa, NULL, p);
-}
-
-static void exec_distribute_child(TaskPool *UNUSED(pool), void *taskdata, int UNUSED(threadid))
-{
-       ParticleTask *task = taskdata;
-       ParticleSystem *psys = task->ctx->sim.psys;
-       ChildParticle *cpa;
-       int p;
-       
-       /* RNG skipping at the beginning */
-       cpa = psys->child;
-       for (p = 0; p < task->begin; ++p, ++cpa) {
-               if (task->ctx->skip) /* simplification skip */
-                       BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP * task->ctx->skip[p]);
-               
-               BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP);
-       }
-               
-       for (; p < task->end; ++p, ++cpa) {
-               if (task->ctx->skip) /* simplification skip */
-                       BLI_rng_skip(task->rng, PSYS_RND_DIST_SKIP * task->ctx->skip[p]);
-               
-               distribute_threads_exec(task, NULL, cpa, p);
-       }
-}
-
-static int distribute_compare_orig_index(const void *p1, const void *p2, void *user_data)
-{
-       int *orig_index = (int *) user_data;
-       int index1 = orig_index[*(const int *)p1];
-       int index2 = orig_index[*(const int *)p2];
-
-       if (index1 < index2)
-               return -1;
-       else if (index1 == index2) {
-               /* this pointer comparison appears to make qsort stable for glibc,
-                * and apparently on solaris too, makes the renders reproducible */
-               if (p1 < p2)
-                       return -1;
-               else if (p1 == p2)
-                       return 0;
-               else
-                       return 1;
-       }
-       else
-               return 1;
-}
-
-static void distribute_invalid(Scene *scene, ParticleSystem *psys, int from)
-{
-       if (from == PART_FROM_CHILD) {
-               ChildParticle *cpa;
-               int p, totchild = get_psys_tot_child(scene, psys);
-
-               if (psys->child && totchild) {
-                       for (p=0,cpa=psys->child; p<totchild; p++,cpa++) {
-                               cpa->fuv[0]=cpa->fuv[1]=cpa->fuv[2]=cpa->fuv[3] = 0.0;
-                               cpa->foffset= 0.0f;
-                               cpa->parent=0;
-                               cpa->pa[0]=cpa->pa[1]=cpa->pa[2]=cpa->pa[3]=0;
-                               cpa->num= -1;
-                       }
-               }
-       }
-       else {
-               PARTICLE_P;
-               LOOP_PARTICLES {
-                       pa->fuv[0] = pa->fuv[1] = pa->fuv[2] = pa->fuv[3] = 0.0;
-                       pa->foffset= 0.0f;
-                       pa->num= -1;
-               }
-       }
-}
-
-/* Creates a distribution of coordinates on a DerivedMesh      */
-/* This is to denote functionality that does not yet work with mesh - only derived mesh */
-static int psys_thread_context_init_distribute(ParticleThreadContext *ctx, ParticleSimulationData *sim, int from)
-{
-       Scene *scene = sim->scene;
-       DerivedMesh *finaldm = sim->psmd->dm;
-       Object *ob = sim->ob;
-       ParticleSystem *psys= sim->psys;
-       ParticleData *pa=0, *tpars= 0;
-       ParticleSettings *part;
-       ParticleSeam *seams= 0;
-       KDTree *tree=0;
-       DerivedMesh *dm= NULL;
-       float *jit= NULL;
-       int i, p=0;
-       int cfrom=0;
-       int totelem=0, totpart, *particle_element=0, children=0, totseam=0;
-       int jitlevel= 1, distr;
-       float *element_weight=NULL,*element_sum=NULL,*jitter_offset=NULL, *vweight=NULL;
-       float cur, maxweight=0.0, tweight, totweight, inv_totweight, co[3], nor[3], orco[3];
-       
-       if (ELEM(NULL, ob, psys, psys->part))
-               return 0;
-       
-       part=psys->part;
-       totpart=psys->totpart;
-       if (totpart==0)
-               return 0;
-       
-       if (!finaldm->deformedOnly && !finaldm->getTessFaceDataArray(finaldm, CD_ORIGINDEX)) {
-               printf("Can't create particles with the current modifier stack, disable destructive modifiers\n");
-// XXX         error("Can't paint with the current modifier stack, disable destructive modifiers");
-               return 0;
-       }
-       
-       psys_thread_context_init(ctx, sim);
-       
-       /* First handle special cases */
-       if (from == PART_FROM_CHILD) {
-               /* Simple children */
-               if (part->childtype != PART_CHILD_FACES) {
-                       BLI_srandom(31415926 + psys->seed + psys->child_seed);
-                       distribute_simple_children(scene, ob, finaldm, psys);
-                       return 0;
-               }
-       }
-       else {
-               /* Grid distribution */
-               if (part->distr==PART_DISTR_GRID && from != PART_FROM_VERT) {
-                       BLI_srandom(31415926 + psys->seed);
-                       dm= CDDM_from_mesh((Mesh*)ob->data);
-                       DM_ensure_tessface(dm);
-                       distribute_grid(dm,psys);
-                       dm->release(dm);
-                       return 0;
-               }
-       }
-       
-       /* Create trees and original coordinates if needed */
-       if (from == PART_FROM_CHILD) {
-               distr=PART_DISTR_RAND;
-               BLI_srandom(31415926 + psys->seed + psys->child_seed);
-               dm= finaldm;
-
-               /* BMESH ONLY */
-               DM_ensure_tessface(dm);
-
-               children=1;
-
-               tree=BLI_kdtree_new(totpart);
-
-               for (p=0,pa=psys->particles; p<totpart; p++,pa++) {
-                       psys_particle_on_dm(dm,part->from,pa->num,pa->num_dmcache,pa->fuv,pa->foffset,co,nor,0,0,orco,NULL);
-                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &orco, 1, 1);
-                       BLI_kdtree_insert(tree, p, orco);
-               }
-
-               BLI_kdtree_balance(tree);
-
-               totpart = get_psys_tot_child(scene, psys);
-               cfrom = from = PART_FROM_FACE;
-       }
-       else {
-               distr = part->distr;
-               BLI_srandom(31415926 + psys->seed);
-               
-               if (psys->part->use_modifier_stack)
-                       dm = finaldm;
-               else
-                       dm= CDDM_from_mesh((Mesh*)ob->data);
-
-               /* BMESH ONLY, for verts we don't care about tessfaces */
-               if (from != PART_FROM_VERT) {
-                       DM_ensure_tessface(dm);
-               }
-
-               /* we need orco for consistent distributions */
-               if (!CustomData_has_layer(&dm->vertData, CD_ORCO))
-                       DM_add_vert_layer(dm, CD_ORCO, CD_ASSIGN, BKE_mesh_orco_verts_get(ob));
-
-               if (from == PART_FROM_VERT) {
-                       MVert *mv= dm->getVertDataArray(dm, CD_MVERT);
-                       float (*orcodata)[3] = dm->getVertDataArray(dm, CD_ORCO);
-                       int totvert = dm->getNumVerts(dm);
-
-                       tree=BLI_kdtree_new(totvert);
-
-                       for (p=0; p<totvert; p++) {
-                               if (orcodata) {
-                                       copy_v3_v3(co,orcodata[p]);
-                                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co, 1, 1);
-                               }
-                               else
-                                       copy_v3_v3(co,mv[p].co);
-                               BLI_kdtree_insert(tree, p, co);
-                       }
-
-                       BLI_kdtree_balance(tree);
-               }
-       }
-
-       /* Get total number of emission elements and allocate needed arrays */
-       totelem = (from == PART_FROM_VERT) ? dm->getNumVerts(dm) : dm->getNumTessFaces(dm);
-
-       if (totelem == 0) {
-               distribute_invalid(scene, psys, children ? PART_FROM_CHILD : 0);
-
-               if (G.debug & G_DEBUG)
-                       fprintf(stderr,"Particle distribution error: Nothing to emit from!\n");
-
-               if (dm != finaldm) dm->release(dm);
-
-               BLI_kdtree_free(tree);
-
-               return 0;
-       }
-
-       element_weight  = MEM_callocN(sizeof(float)*totelem, "particle_distribution_weights");
-       particle_element= MEM_callocN(sizeof(int)*totpart, "particle_distribution_indexes");
-       element_sum             = MEM_callocN(sizeof(float)*(totelem+1), "particle_distribution_sum");
-       jitter_offset   = MEM_callocN(sizeof(float)*totelem, "particle_distribution_jitoff");
-
-       /* Calculate weights from face areas */
-       if ((part->flag&PART_EDISTR || children) && from != PART_FROM_VERT) {
-               MVert *v1, *v2, *v3, *v4;
-               float totarea=0.f, co1[3], co2[3], co3[3], co4[3];
-               float (*orcodata)[3];
-               
-               orcodata= dm->getVertDataArray(dm, CD_ORCO);
-
-               for (i=0; i<totelem; i++) {
-                       MFace *mf=dm->getTessFaceData(dm,i,CD_MFACE);
-
-                       if (orcodata) {
-                               copy_v3_v3(co1, orcodata[mf->v1]);
-                               copy_v3_v3(co2, orcodata[mf->v2]);
-                               copy_v3_v3(co3, orcodata[mf->v3]);
-                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co1, 1, 1);
-                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co2, 1, 1);
-                               BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co3, 1, 1);
-                               if (mf->v4) {
-                                       copy_v3_v3(co4, orcodata[mf->v4]);
-                                       BKE_mesh_orco_verts_transform((Mesh*)ob->data, &co4, 1, 1);
-                               }
-                       }
-                       else {
-                               v1= (MVert*)dm->getVertData(dm,mf->v1,CD_MVERT);
-                               v2= (MVert*)dm->getVertData(dm,mf->v2,CD_MVERT);
-                               v3= (MVert*)dm->getVertData(dm,mf->v3,CD_MVERT);
-                               copy_v3_v3(co1, v1->co);
-                               copy_v3_v3(co2, v2->co);
-                               copy_v3_v3(co3, v3->co);
-                               if (mf->v4) {
-                                       v4= (MVert*)dm->getVertData(dm,mf->v4,CD_MVERT);
-                                       copy_v3_v3(co4, v4->co);
-                               }
-                       }
-
-                       cur = mf->v4 ? area_quad_v3(co1, co2, co3, co4) : area_tri_v3(co1, co2, co3);
-                       
-                       if (cur > maxweight)
-                               maxweight = cur;
-
-                       element_weight[i] = cur;
-                       totarea += cur;
-               }
-
-               for (i=0; i<totelem; i++)
-                       element_weight[i] /= totarea;
-
-               maxweight /= totarea;
-       }
-       else {
-               float min=1.0f/(float)(MIN2(totelem,totpart));
-               for (i=0; i<totelem; i++)
-                       element_weight[i]=min;
-               maxweight=min;
-       }
-
-       /* Calculate weights from vgroup */
-       vweight = psys_cache_vgroup(dm,psys,PSYS_VG_DENSITY);
-
-       if (vweight) {
-               if (from==PART_FROM_VERT) {
-                       for (i=0;i<totelem; i++)
-                               element_weight[i]*=vweight[i];
-               }
-               else { /* PART_FROM_FACE / PART_FROM_VOLUME */
-                       for (i=0;i<totelem; i++) {
-                               MFace *mf=dm->getTessFaceData(dm,i,CD_MFACE);
-                               tweight = vweight[mf->v1] + vweight[mf->v2] + vweight[mf->v3];
-                               
-                               if (mf->v4) {
-                                       tweight += vweight[mf->v4];
-                                       tweight /= 4.0f;
-                               }
-                               else {
-                                       tweight /= 3.0f;
-                               }
-
-                               element_weight[i]*=tweight;
-                       }
-               }
-               MEM_freeN(vweight);
-       }
-
-       /* Calculate total weight of all elements */
-       totweight= 0.0f;
-       for (i=0;i<totelem; i++)
-               totweight += element_weight[i];
-
-       inv_totweight = (totweight > 0.f ? 1.f/totweight : 0.f);
-
-       /* Calculate cumulative weights */
-       element_sum[0] = 0.0f;
-       for (i=0; i<totelem; i++)
-               element_sum[i+1] = element_sum[i] + element_weight[i] * inv_totweight;
-       
-       /* Finally assign elements to particles */
-       if ((part->flag&PART_TRAND) || (part->simplify_flag&PART_SIMPLIFY_ENABLE)) {
-               float pos;
-
-               for (p=0; p<totpart; p++) {
-                       /* In theory element_sum[totelem] should be 1.0, but due to float errors this is not necessarily always true, so scale pos accordingly. */
-                       pos= BLI_frand() * element_sum[totelem];
-                       particle_element[p] = distribute_binary_search(element_sum, totelem, pos);
-                       particle_element[p] = MIN2(totelem-1, particle_element[p]);
-                       jitter_offset[particle_element[p]] = pos;
-               }
-       }
-       else {
-               double step, pos;
-               
-               step= (totpart < 2) ? 0.5 : 1.0/(double)totpart;
-               pos= 1e-6; /* tiny offset to avoid zero weight face */
-               i= 0;
-
-               for (p=0; p<totpart; p++, pos+=step) {
-                       while ((i < totelem) && (pos > (double)element_sum[i + 1]))
-                               i++;
-
-                       particle_element[p] = MIN2(totelem-1, i);
-
-                       /* avoid zero weight face */
-                       if (p == totpart-1 && element_weight[particle_element[p]] == 0.0f)
-                               particle_element[p] = particle_element[p-1];
-
-                       jitter_offset[particle_element[p]] = pos;
-               }
-       }
-
-       MEM_freeN(element_sum);
-
-       /* For hair, sort by origindex (allows optimization's in rendering), */
-       /* however with virtual parents the children need to be in random order. */
-       if (part->type == PART_HAIR && !(part->childtype==PART_CHILD_FACES && part->parents!=0.0f)) {
-               int *orig_index = NULL;
-
-               if (from == PART_FROM_VERT) {
-                       if (dm->numVertData)
-                               orig_index = dm->getVertDataArray(dm, CD_ORIGINDEX);
-               }
-               else {
-                       if (dm->numTessFaceData)
-                               orig_index = dm->getTessFaceDataArray(dm, CD_ORIGINDEX);
-               }
-
-               if (orig_index) {
-                       BLI_qsort_r(particle_element, totpart, sizeof(int), distribute_compare_orig_index, orig_index);
-               }
-       }
-
-       /* Create jittering if needed */
-       if (distr==PART_DISTR_JIT && ELEM(from,PART_FROM_FACE,PART_FROM_VOLUME)) {
-               jitlevel= part->userjit;
-               
-               if (jitlevel == 0) {
-                       jitlevel= totpart/totelem;
-                       if (part->flag & PART_EDISTR) jitlevel*= 2;     /* looks better in general, not very scietific */
-                       if (jitlevel<3) jitlevel= 3;
-               }
-               
-               jit= MEM_callocN((2+ jitlevel*2)*sizeof(float), "jit");
-
-               /* for small amounts of particles we use regular jitter since it looks
-                * a bit better, for larger amounts we switch to hammersley sequence 
-                * because it is much faster */
-               if (jitlevel < 25)
-                       init_mv_jit(jit, jitlevel, psys->seed, part->jitfac);
-               else
-                       hammersley_create(jit, jitlevel+1, psys->seed, part->jitfac);
-               BLI_array_randomize(jit, 2*sizeof(float), jitlevel, psys->seed); /* for custom jit or even distribution */
-       }
-
-       /* Setup things for threaded distribution */
-       ctx->tree= tree;
-       ctx->seams= seams;
-       ctx->totseam= totseam;
-       ctx->sim.psys= psys;
-       ctx->index= particle_element;
-       ctx->jit= jit;
-       ctx->jitlevel= jitlevel;
-       ctx->jitoff= jitter_offset;
-       ctx->weight= element_weight;
-       ctx->maxweight= maxweight;
-       ctx->from= (children) ? PART_FROM_CHILD : from;
-       ctx->cfrom= cfrom;
-       ctx->distr= distr;
-       ctx->dm= dm;
-       ctx->tpars= tpars;
-
-       if (children) {
-               totpart= psys_render_simplify_distribution(ctx, totpart);
-               alloc_child_particles(psys, totpart);
-       }
-
-       return 1;
-}
-
-static void psys_task_init_distribute(ParticleTask *task, ParticleSimulationData *sim)
-{
-       /* init random number generator */
-       int seed = 31415926 + sim->psys->seed;
-       
-       task->rng = BLI_rng_new(seed);
-}
-
-static void distribute_particles_on_dm(ParticleSimulationData *sim, int from)
-{
-       TaskScheduler *task_scheduler;
-       TaskPool *task_pool;
-       ParticleThreadContext ctx;
-       ParticleTask *tasks;
-       DerivedMesh *finaldm = sim->psmd->dm;
-       int i, totpart, numtasks;
-       
-       /* create a task pool for distribution tasks */
-       if (!psys_thread_context_init_distribute(&ctx, sim, from))
-               return;
-       
-       task_scheduler = BLI_task_scheduler_get();
-       task_pool = BLI_task_pool_create(task_scheduler, &ctx);
-       
-       totpart = (from == PART_FROM_CHILD ? sim->psys->totchild : sim->psys->totpart);
-       psys_tasks_create(&ctx, totpart, &tasks, &numtasks);
-       for (i = 0; i < numtasks; ++i) {
-               ParticleTask *task = &tasks[i];
-               
-               psys_task_init_distribute(task, sim);
-               if (from == PART_FROM_CHILD)
-                       BLI_task_pool_push(task_pool, exec_distribute_child, task, false, TASK_PRIORITY_LOW);
-               else
-                       BLI_task_pool_push(task_pool, exec_distribute_parent, task, false, TASK_PRIORITY_LOW);
-       }
-       BLI_task_pool_work_and_wait(task_pool);
-       
-       BLI_task_pool_free(task_pool);
-       
-       psys_calc_dmcache(sim->ob, finaldm, sim->psys);
-       
-       if (ctx.dm != finaldm)
-               ctx.dm->release(ctx.dm);
-       
-       psys_tasks_free(tasks, numtasks);
-}
-
-/* ready for future use, to emit particles without geometry */
-static void distribute_particles_on_shape(ParticleSimulationData *sim, int UNUSED(from))
-{
-       distribute_invalid(sim->scene, sim->psys, 0);
-
-       fprintf(stderr,"Shape emission not yet possible!\n");
-}
-
-static void distribute_particles(ParticleSimulationData *sim, int from)
-{
-       PARTICLE_PSMD;
-       int distr_error=0;
-
-       if (psmd) {
-               if (psmd->dm)
-                       distribute_particles_on_dm(sim, from);
-               else
-                       distr_error=1;
-       }
-       else
-               distribute_particles_on_shape(sim, from);
-
-       if (distr_error) {
-               distribute_invalid(sim->scene, sim->psys, from);
-
-               fprintf(stderr,"Particle distribution error!\n");
-       }
-}
-
 /* threaded child particle distribution and path caching */
 void psys_thread_context_init(ParticleThreadContext *ctx, ParticleSimulationData *sim)
 {
@@ -3930,7 +2883,7 @@ static void psys_update_path_cache(ParticleSimulationData *sim, float cfra)
        Base *base;
        int distr=0, alloc=0, skip=0;
 
-       if ((psys->part->childtype && psys->totchild != get_psys_tot_child(sim->scene, psys)) || psys->recalc&PSYS_RECALC_RESET)
+       if ((psys->part->childtype && psys->totchild != psys_get_tot_child(sim->scene, psys)) || psys->recalc&PSYS_RECALC_RESET)
                alloc=1;
 
        if (alloc || psys->recalc&PSYS_RECALC_CHILD || (psys->vgroup[PSYS_VG_DENSITY] && (sim->ob && sim->ob->mode & OB_MODE_WEIGHT_PAINT)))
@@ -3940,7 +2893,7 @@ static void psys_update_path_cache(ParticleSimulationData *sim, float cfra)
                if (alloc)
                        realloc_particles(sim, sim->psys->totpart);
 
-               if (get_psys_tot_child(sim->scene, psys)) {
+               if (psys_get_tot_child(sim->scene, psys)) {
                        /* don't generate children while computing the hair keys */
                        if (!(psys->part->type == PART_HAIR) || (psys->flag & PSYS_HAIR_DONE)) {
                                distribute_particles(sim, PART_FROM_CHILD);
@@ -4600,7 +3553,7 @@ static void update_children(ParticleSimulationData *sim)
        /* don't generate children while growing hair - waste of time */
                psys_free_children(sim->psys);
        else if (sim->psys->part->childtype) {
-               if (sim->psys->totchild != get_psys_tot_child(sim->scene, sim->psys))
+               if (sim->psys->totchild != psys_get_tot_child(sim->scene, sim->psys))
                        distribute_particles(sim, PART_FROM_CHILD);
                else {
                        /* Children are up to date, nothing to do. */