KDTree: deprecate 'normal' argument
[blender-staging.git] / source / blender / blenkernel / intern / particle.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) 2007 by Janne Karhu.
19  * All rights reserved.
20  *
21  * The Original Code is: all of this file.
22  *
23  * Contributor(s): none yet.
24  *
25  * ***** END GPL LICENSE BLOCK *****
26  */
27
28 /** \file blender/blenkernel/intern/particle.c
29  *  \ingroup bke
30  */
31
32
33 #include <stdlib.h>
34 #include <math.h>
35 #include <string.h>
36
37 #include "MEM_guardedalloc.h"
38
39 #include "DNA_curve_types.h"
40 #include "DNA_group_types.h"
41 #include "DNA_key_types.h"
42 #include "DNA_material_types.h"
43 #include "DNA_mesh_types.h"
44 #include "DNA_meshdata_types.h"
45 #include "DNA_particle_types.h"
46 #include "DNA_smoke_types.h"
47 #include "DNA_scene_types.h"
48 #include "DNA_dynamicpaint_types.h"
49
50 #include "BLI_blenlib.h"
51 #include "BLI_noise.h"
52 #include "BLI_math.h"
53 #include "BLI_utildefines.h"
54 #include "BLI_kdtree.h"
55 #include "BLI_rand.h"
56 #include "BLI_threads.h"
57 #include "BLI_linklist.h"
58
59 #include "BLF_translation.h"
60
61 #include "BKE_anim.h"
62 #include "BKE_animsys.h"
63
64 #include "BKE_boids.h"
65 #include "BKE_cloth.h"
66 #include "BKE_effect.h"
67 #include "BKE_global.h"
68 #include "BKE_group.h"
69 #include "BKE_main.h"
70 #include "BKE_lattice.h"
71
72 #include "BKE_displist.h"
73 #include "BKE_particle.h"
74 #include "BKE_object.h"
75 #include "BKE_material.h"
76 #include "BKE_key.h"
77 #include "BKE_library.h"
78 #include "BKE_depsgraph.h"
79 #include "BKE_modifier.h"
80 #include "BKE_mesh.h"
81 #include "BKE_cdderivedmesh.h"
82 #include "BKE_pointcache.h"
83 #include "BKE_scene.h"
84 #include "BKE_deform.h"
85
86 #include "RE_render_ext.h"
87
88 unsigned int PSYS_FRAND_SEED_OFFSET[PSYS_FRAND_COUNT];
89 unsigned int PSYS_FRAND_SEED_MULTIPLIER[PSYS_FRAND_COUNT];
90 float PSYS_FRAND_BASE[PSYS_FRAND_COUNT];
91
92 void psys_init_rng(void)
93 {
94         int i;
95         BLI_srandom(5831); /* arbitrary */
96         for (i = 0; i < PSYS_FRAND_COUNT; ++i) {
97                 PSYS_FRAND_BASE[i] = BLI_frand();
98                 PSYS_FRAND_SEED_OFFSET[i] = (unsigned int)BLI_rand();
99                 PSYS_FRAND_SEED_MULTIPLIER[i] = (unsigned int)BLI_rand();
100         }
101 }
102
103 static void get_child_modifier_parameters(ParticleSettings *part, ParticleThreadContext *ctx,
104                                           ChildParticle *cpa, short cpa_from, int cpa_num, float *cpa_fuv, float *orco, ParticleTexture *ptex);
105 static void do_child_modifiers(ParticleSimulationData *sim,
106                                ParticleTexture *ptex, ParticleKey *par, float *par_rot, ChildParticle *cpa,
107                                float *orco, float mat[4][4], ParticleKey *state, float t);
108
109 /* few helpers for countall etc. */
110 int count_particles(ParticleSystem *psys)
111 {
112         ParticleSettings *part = psys->part;
113         PARTICLE_P;
114         int tot = 0;
115
116         LOOP_SHOWN_PARTICLES {
117                 if (pa->alive == PARS_UNBORN && (part->flag & PART_UNBORN) == 0) {}
118                 else if (pa->alive == PARS_DEAD && (part->flag & PART_DIED) == 0) {}
119                 else tot++;
120         }
121         return tot;
122 }
123 int count_particles_mod(ParticleSystem *psys, int totgr, int cur)
124 {
125         ParticleSettings *part = psys->part;
126         PARTICLE_P;
127         int tot = 0;
128
129         LOOP_SHOWN_PARTICLES {
130                 if (pa->alive == PARS_UNBORN && (part->flag & PART_UNBORN) == 0) {}
131                 else if (pa->alive == PARS_DEAD && (part->flag & PART_DIED) == 0) {}
132                 else if (p % totgr == cur) tot++;
133         }
134         return tot;
135 }
136 /* we allocate path cache memory in chunks instead of a big contiguous
137  * chunk, windows' memory allocater fails to find big blocks of memory often */
138
139 #define PATH_CACHE_BUF_SIZE 1024
140
141 static ParticleCacheKey **psys_alloc_path_cache_buffers(ListBase *bufs, int tot, int steps)
142 {
143         LinkData *buf;
144         ParticleCacheKey **cache;
145         int i, totkey, totbufkey;
146
147         tot = MAX2(tot, 1);
148         totkey = 0;
149         cache = MEM_callocN(tot * sizeof(void *), "PathCacheArray");
150
151         while (totkey < tot) {
152                 totbufkey = MIN2(tot - totkey, PATH_CACHE_BUF_SIZE);
153                 buf = MEM_callocN(sizeof(LinkData), "PathCacheLinkData");
154                 buf->data = MEM_callocN(sizeof(ParticleCacheKey) * totbufkey * steps, "ParticleCacheKey");
155
156                 for (i = 0; i < totbufkey; i++)
157                         cache[totkey + i] = ((ParticleCacheKey *)buf->data) + i * steps;
158
159                 totkey += totbufkey;
160                 BLI_addtail(bufs, buf);
161         }
162
163         return cache;
164 }
165
166 static void psys_free_path_cache_buffers(ParticleCacheKey **cache, ListBase *bufs)
167 {
168         LinkData *buf;
169
170         if (cache)
171                 MEM_freeN(cache);
172
173         for (buf = bufs->first; buf; buf = buf->next)
174                 MEM_freeN(buf->data);
175         BLI_freelistN(bufs);
176 }
177
178 /************************************************/
179 /*                      Getting stuff                                           */
180 /************************************************/
181 /* get object's active particle system safely */
182 ParticleSystem *psys_get_current(Object *ob)
183 {
184         ParticleSystem *psys;
185         if (ob == NULL) return NULL;
186
187         for (psys = ob->particlesystem.first; psys; psys = psys->next) {
188                 if (psys->flag & PSYS_CURRENT)
189                         return psys;
190         }
191         
192         return NULL;
193 }
194 short psys_get_current_num(Object *ob)
195 {
196         ParticleSystem *psys;
197         short i;
198
199         if (ob == NULL) return 0;
200
201         for (psys = ob->particlesystem.first, i = 0; psys; psys = psys->next, i++)
202                 if (psys->flag & PSYS_CURRENT)
203                         return i;
204         
205         return i;
206 }
207 void psys_set_current_num(Object *ob, int index)
208 {
209         ParticleSystem *psys;
210         short i;
211
212         if (ob == NULL) return;
213
214         for (psys = ob->particlesystem.first, i = 0; psys; psys = psys->next, i++) {
215                 if (i == index)
216                         psys->flag |= PSYS_CURRENT;
217                 else
218                         psys->flag &= ~PSYS_CURRENT;
219         }
220 }
221
222 #if 0 /* UNUSED */
223 Object *psys_find_object(Scene *scene, ParticleSystem *psys)
224 {
225         Base *base;
226         ParticleSystem *tpsys;
227
228         for (base = scene->base.first; base; base = base->next) {
229                 for (tpsys = base->object->particlesystem.first; psys; psys = psys->next) {
230                         if (tpsys == psys)
231                                 return base->object;
232                 }
233         }
234
235         return NULL;
236 }
237 #endif
238
239 struct LatticeDeformData *psys_create_lattice_deform_data(ParticleSimulationData *sim)
240 {
241         struct LatticeDeformData *lattice_deform_data = NULL;
242
243         if (psys_in_edit_mode(sim->scene, sim->psys) == 0) {
244                 Object *lattice = NULL;
245                 ModifierData *md = (ModifierData *)psys_get_modifier(sim->ob, sim->psys);
246
247                 for (; md; md = md->next) {
248                         if (md->type == eModifierType_Lattice) {
249                                 LatticeModifierData *lmd = (LatticeModifierData *)md;
250                                 lattice = lmd->object;
251                                 break;
252                         }
253                 }
254                 if (lattice)
255                         lattice_deform_data = init_latt_deform(lattice, NULL);
256         }
257
258         return lattice_deform_data;
259 }
260 void psys_disable_all(Object *ob)
261 {
262         ParticleSystem *psys = ob->particlesystem.first;
263
264         for (; psys; psys = psys->next)
265                 psys->flag |= PSYS_DISABLED;
266 }
267 void psys_enable_all(Object *ob)
268 {
269         ParticleSystem *psys = ob->particlesystem.first;
270
271         for (; psys; psys = psys->next)
272                 psys->flag &= ~PSYS_DISABLED;
273 }
274 int psys_in_edit_mode(Scene *scene, ParticleSystem *psys)
275 {
276         return (scene->basact && (scene->basact->object->mode & OB_MODE_PARTICLE_EDIT) && psys == psys_get_current((scene->basact)->object) && (psys->edit || psys->pointcache->edit) && !psys->renderdata);
277 }
278 int psys_check_enabled(Object *ob, ParticleSystem *psys)
279 {
280         ParticleSystemModifierData *psmd;
281
282         if (psys->flag & PSYS_DISABLED || psys->flag & PSYS_DELETE || !psys->part)
283                 return 0;
284
285         psmd = psys_get_modifier(ob, psys);
286         if (psys->renderdata || G.is_rendering) {
287                 if (!(psmd->modifier.mode & eModifierMode_Render))
288                         return 0;
289         }
290         else if (!(psmd->modifier.mode & eModifierMode_Realtime))
291                 return 0;
292
293         return 1;
294 }
295
296 int psys_check_edited(ParticleSystem *psys)
297 {
298         if (psys->part && psys->part->type == PART_HAIR)
299                 return (psys->flag & PSYS_EDITED || (psys->edit && psys->edit->edited));
300         else
301                 return (psys->pointcache->edit && psys->pointcache->edit->edited);
302 }
303
304 void psys_check_group_weights(ParticleSettings *part)
305 {
306         ParticleDupliWeight *dw, *tdw;
307         GroupObject *go;
308         int current = 0;
309
310         if (part->ren_as == PART_DRAW_GR && part->dup_group && part->dup_group->gobject.first) {
311                 /* first remove all weights that don't have an object in the group */
312                 dw = part->dupliweights.first;
313                 while (dw) {
314                         if (!BKE_group_object_exists(part->dup_group, dw->ob)) {
315                                 tdw = dw->next;
316                                 BLI_freelinkN(&part->dupliweights, dw);
317                                 dw = tdw;
318                         }
319                         else
320                                 dw = dw->next;
321                 }
322
323                 /* then add objects in the group to new list */
324                 go = part->dup_group->gobject.first;
325                 while (go) {
326                         dw = part->dupliweights.first;
327                         while (dw && dw->ob != go->ob)
328                                 dw = dw->next;
329                         
330                         if (!dw) {
331                                 dw = MEM_callocN(sizeof(ParticleDupliWeight), "ParticleDupliWeight");
332                                 dw->ob = go->ob;
333                                 dw->count = 1;
334                                 BLI_addtail(&part->dupliweights, dw);
335                         }
336
337                         go = go->next;
338                 }
339
340                 dw = part->dupliweights.first;
341                 for (; dw; dw = dw->next) {
342                         if (dw->flag & PART_DUPLIW_CURRENT) {
343                                 current = 1;
344                                 break;
345                         }
346                 }
347
348                 if (!current) {
349                         dw = part->dupliweights.first;
350                         if (dw)
351                                 dw->flag |= PART_DUPLIW_CURRENT;
352                 }
353         }
354         else {
355                 BLI_freelistN(&part->dupliweights);
356         }
357 }
358 int psys_uses_gravity(ParticleSimulationData *sim)
359 {
360         return sim->scene->physics_settings.flag & PHYS_GLOBAL_GRAVITY && sim->psys->part && sim->psys->part->effector_weights->global_gravity != 0.0f;
361 }
362 /************************************************/
363 /*                      Freeing stuff                                           */
364 /************************************************/
365 static void fluid_free_settings(SPHFluidSettings *fluid)
366 {
367         if (fluid)
368                 MEM_freeN(fluid); 
369 }
370
371 void BKE_particlesettings_free(ParticleSettings *part)
372 {
373         MTex *mtex;
374         int a;
375         BKE_free_animdata(&part->id);
376         free_partdeflect(part->pd);
377         free_partdeflect(part->pd2);
378
379         if (part->effector_weights)
380                 MEM_freeN(part->effector_weights);
381
382         BLI_freelistN(&part->dupliweights);
383
384         boid_free_settings(part->boids);
385         fluid_free_settings(part->fluid);
386
387         for (a = 0; a < MAX_MTEX; a++) {
388                 mtex = part->mtex[a];
389                 if (mtex && mtex->tex) mtex->tex->id.us--;
390                 if (mtex) MEM_freeN(mtex);
391         }
392 }
393
394 void free_hair(Object *UNUSED(ob), ParticleSystem *psys, int dynamics)
395 {
396         PARTICLE_P;
397
398         LOOP_PARTICLES {
399                 if (pa->hair)
400                         MEM_freeN(pa->hair);
401                 pa->hair = NULL;
402                 pa->totkey = 0;
403         }
404
405         psys->flag &= ~PSYS_HAIR_DONE;
406
407         if (psys->clmd) {
408                 if (dynamics) {
409                         BKE_ptcache_free_list(&psys->ptcaches);
410                         psys->clmd->point_cache = psys->pointcache = NULL;
411                         BLI_listbase_clear(&psys->clmd->ptcaches);
412
413                         modifier_free((ModifierData *)psys->clmd);
414                         
415                         psys->clmd = NULL;
416                         psys->pointcache = BKE_ptcache_add(&psys->ptcaches);
417                 }
418                 else {
419                         cloth_free_modifier(psys->clmd);
420                 }
421         }
422
423         if (psys->hair_in_dm)
424                 psys->hair_in_dm->release(psys->hair_in_dm);
425         psys->hair_in_dm = NULL;
426
427         if (psys->hair_out_dm)
428                 psys->hair_out_dm->release(psys->hair_out_dm);
429         psys->hair_out_dm = NULL;
430 }
431 void free_keyed_keys(ParticleSystem *psys)
432 {
433         PARTICLE_P;
434
435         if (psys->part->type == PART_HAIR)
436                 return;
437
438         if (psys->particles && psys->particles->keys) {
439                 MEM_freeN(psys->particles->keys);
440
441                 LOOP_PARTICLES {
442                         if (pa->keys) {
443                                 pa->keys = NULL;
444                                 pa->totkey = 0;
445                         }
446                 }
447         }
448 }
449 static void free_child_path_cache(ParticleSystem *psys)
450 {
451         psys_free_path_cache_buffers(psys->childcache, &psys->childcachebufs);
452         psys->childcache = NULL;
453         psys->totchildcache = 0;
454 }
455 void psys_free_path_cache(ParticleSystem *psys, PTCacheEdit *edit)
456 {
457         if (edit) {
458                 psys_free_path_cache_buffers(edit->pathcache, &edit->pathcachebufs);
459                 edit->pathcache = NULL;
460                 edit->totcached = 0;
461         }
462         if (psys) {
463                 psys_free_path_cache_buffers(psys->pathcache, &psys->pathcachebufs);
464                 psys->pathcache = NULL;
465                 psys->totcached = 0;
466
467                 free_child_path_cache(psys);
468         }
469 }
470 void psys_free_children(ParticleSystem *psys)
471 {
472         if (psys->child) {
473                 MEM_freeN(psys->child);
474                 psys->child = NULL;
475                 psys->totchild = 0;
476         }
477
478         free_child_path_cache(psys);
479 }
480 void psys_free_particles(ParticleSystem *psys)
481 {
482         PARTICLE_P;
483
484         if (psys->particles) {
485                 if (psys->part->type == PART_HAIR) {
486                         LOOP_PARTICLES {
487                                 if (pa->hair)
488                                         MEM_freeN(pa->hair);
489                         }
490                 }
491                 
492                 if (psys->particles->keys)
493                         MEM_freeN(psys->particles->keys);
494                 
495                 if (psys->particles->boid)
496                         MEM_freeN(psys->particles->boid);
497
498                 MEM_freeN(psys->particles);
499                 psys->particles = NULL;
500                 psys->totpart = 0;
501         }
502 }
503 void psys_free_pdd(ParticleSystem *psys)
504 {
505         if (psys->pdd) {
506                 if (psys->pdd->cdata)
507                         MEM_freeN(psys->pdd->cdata);
508                 psys->pdd->cdata = NULL;
509
510                 if (psys->pdd->vdata)
511                         MEM_freeN(psys->pdd->vdata);
512                 psys->pdd->vdata = NULL;
513
514                 if (psys->pdd->ndata)
515                         MEM_freeN(psys->pdd->ndata);
516                 psys->pdd->ndata = NULL;
517
518                 if (psys->pdd->vedata)
519                         MEM_freeN(psys->pdd->vedata);
520                 psys->pdd->vedata = NULL;
521
522                 psys->pdd->totpoint = 0;
523                 psys->pdd->tot_vec_size = 0;
524         }
525 }
526 /* free everything */
527 void psys_free(Object *ob, ParticleSystem *psys)
528 {       
529         if (psys) {
530                 int nr = 0;
531                 ParticleSystem *tpsys;
532                 
533                 psys_free_path_cache(psys, NULL);
534
535                 free_hair(ob, psys, 1);
536
537                 psys_free_particles(psys);
538
539                 if (psys->edit && psys->free_edit)
540                         psys->free_edit(psys->edit);
541
542                 if (psys->child) {
543                         MEM_freeN(psys->child);
544                         psys->child = NULL;
545                         psys->totchild = 0;
546                 }
547                 
548                 /* check if we are last non-visible particle system */
549                 for (tpsys = ob->particlesystem.first; tpsys; tpsys = tpsys->next) {
550                         if (tpsys->part) {
551                                 if (ELEM(tpsys->part->ren_as, PART_DRAW_OB, PART_DRAW_GR)) {
552                                         nr++;
553                                         break;
554                                 }
555                         }
556                 }
557                 /* clear do-not-draw-flag */
558                 if (!nr)
559                         ob->transflag &= ~OB_DUPLIPARTS;
560
561                 if (psys->part) {
562                         psys->part->id.us--;
563                         psys->part = NULL;
564                 }
565
566                 BKE_ptcache_free_list(&psys->ptcaches);
567                 psys->pointcache = NULL;
568                 
569                 BLI_freelistN(&psys->targets);
570
571                 BLI_bvhtree_free(psys->bvhtree);
572                 BLI_kdtree_free(psys->tree);
573  
574                 if (psys->fluid_springs)
575                         MEM_freeN(psys->fluid_springs);
576
577                 pdEndEffectors(&psys->effectors);
578
579                 if (psys->pdd) {
580                         psys_free_pdd(psys);
581                         MEM_freeN(psys->pdd);
582                 }
583
584                 MEM_freeN(psys);
585         }
586 }
587
588 /************************************************/
589 /*                      Rendering                                                       */
590 /************************************************/
591 /* these functions move away particle data and bring it back after
592  * rendering, to make different render settings possible without
593  * removing the previous data. this should be solved properly once */
594
595 typedef struct ParticleRenderElem {
596         int curchild, totchild, reduce;
597         float lambda, t, scalemin, scalemax;
598 } ParticleRenderElem;
599
600 typedef struct ParticleRenderData {
601         ChildParticle *child;
602         ParticleCacheKey **pathcache;
603         ParticleCacheKey **childcache;
604         ListBase pathcachebufs, childcachebufs;
605         int totchild, totcached, totchildcache;
606         DerivedMesh *dm;
607         int totdmvert, totdmedge, totdmface;
608
609         float mat[4][4];
610         float viewmat[4][4], winmat[4][4];
611         int winx, winy;
612
613         int do_simplify;
614         int timeoffset;
615         ParticleRenderElem *elems;
616
617         /* ORIGINDEX */
618         const int *index_mf_to_mpoly;
619         const int *index_mp_to_orig;
620 } ParticleRenderData;
621
622 static float psys_render_viewport_falloff(double rate, float dist, float width)
623 {
624         return pow(rate, dist / width);
625 }
626
627 static float psys_render_projected_area(ParticleSystem *psys, const float center[3], float area, double vprate, float *viewport)
628 {
629         ParticleRenderData *data = psys->renderdata;
630         float co[4], view[3], ortho1[3], ortho2[3], w, dx, dy, radius;
631         
632         /* transform to view space */
633         copy_v3_v3(co, center);
634         co[3] = 1.0f;
635         mul_m4_v4(data->viewmat, co);
636         
637         /* compute two vectors orthogonal to view vector */
638         normalize_v3_v3(view, co);
639         ortho_basis_v3v3_v3(ortho1, ortho2, view);
640
641         /* compute on screen minification */
642         w = co[2] * data->winmat[2][3] + data->winmat[3][3];
643         dx = data->winx * ortho2[0] * data->winmat[0][0];
644         dy = data->winy * ortho2[1] * data->winmat[1][1];
645         w = sqrtf(dx * dx + dy * dy) / w;
646
647         /* w squared because we are working with area */
648         area = area * w * w;
649
650         /* viewport of the screen test */
651
652         /* project point on screen */
653         mul_m4_v4(data->winmat, co);
654         if (co[3] != 0.0f) {
655                 co[0] = 0.5f * data->winx * (1.0f + co[0] / co[3]);
656                 co[1] = 0.5f * data->winy * (1.0f + co[1] / co[3]);
657         }
658
659         /* screen space radius */
660         radius = sqrt(area / (float)M_PI);
661
662         /* make smaller using fallof once over screen edge */
663         *viewport = 1.0f;
664
665         if (co[0] + radius < 0.0f)
666                 *viewport *= psys_render_viewport_falloff(vprate, -(co[0] + radius), data->winx);
667         else if (co[0] - radius > data->winx)
668                 *viewport *= psys_render_viewport_falloff(vprate, (co[0] - radius) - data->winx, data->winx);
669
670         if (co[1] + radius < 0.0f)
671                 *viewport *= psys_render_viewport_falloff(vprate, -(co[1] + radius), data->winy);
672         else if (co[1] - radius > data->winy)
673                 *viewport *= psys_render_viewport_falloff(vprate, (co[1] - radius) - data->winy, data->winy);
674         
675         return area;
676 }
677
678 void psys_render_set(Object *ob, ParticleSystem *psys, float viewmat[4][4], float winmat[4][4], int winx, int winy, int timeoffset)
679 {
680         ParticleRenderData *data;
681         ParticleSystemModifierData *psmd = psys_get_modifier(ob, psys);
682
683         if (psys->renderdata)
684                 return;
685
686         data = MEM_callocN(sizeof(ParticleRenderData), "ParticleRenderData");
687
688         data->child = psys->child;
689         data->totchild = psys->totchild;
690         data->pathcache = psys->pathcache;
691         data->pathcachebufs.first = psys->pathcachebufs.first;
692         data->pathcachebufs.last = psys->pathcachebufs.last;
693         data->totcached = psys->totcached;
694         data->childcache = psys->childcache;
695         data->childcachebufs.first = psys->childcachebufs.first;
696         data->childcachebufs.last = psys->childcachebufs.last;
697         data->totchildcache = psys->totchildcache;
698
699         if (psmd->dm)
700                 data->dm = CDDM_copy(psmd->dm);
701         data->totdmvert = psmd->totdmvert;
702         data->totdmedge = psmd->totdmedge;
703         data->totdmface = psmd->totdmface;
704
705         psys->child = NULL;
706         psys->pathcache = NULL;
707         psys->childcache = NULL;
708         psys->totchild = psys->totcached = psys->totchildcache = 0;
709         BLI_listbase_clear(&psys->pathcachebufs);
710         BLI_listbase_clear(&psys->childcachebufs);
711
712         copy_m4_m4(data->winmat, winmat);
713         mul_m4_m4m4(data->viewmat, viewmat, ob->obmat);
714         mul_m4_m4m4(data->mat, winmat, data->viewmat);
715         data->winx = winx;
716         data->winy = winy;
717
718         data->timeoffset = timeoffset;
719
720         psys->renderdata = data;
721
722         /* Hair can and has to be recalculated if everything isn't displayed. */
723         if (psys->part->disp != 100 && psys->part->type == PART_HAIR)
724                 psys->recalc |= PSYS_RECALC_RESET;
725 }
726
727 void psys_render_restore(Object *ob, ParticleSystem *psys)
728 {
729         ParticleRenderData *data;
730         ParticleSystemModifierData *psmd = psys_get_modifier(ob, psys);
731         float render_disp = psys_get_current_display_percentage(psys);
732         float disp;
733
734         data = psys->renderdata;
735         if (!data)
736                 return;
737         
738         if (data->elems)
739                 MEM_freeN(data->elems);
740
741         if (psmd->dm) {
742                 psmd->dm->needsFree = 1;
743                 psmd->dm->release(psmd->dm);
744         }
745
746         psys_free_path_cache(psys, NULL);
747
748         if (psys->child) {
749                 MEM_freeN(psys->child);
750                 psys->child = 0;
751                 psys->totchild = 0;
752         }
753
754         psys->child = data->child;
755         psys->totchild = data->totchild;
756         psys->pathcache = data->pathcache;
757         psys->pathcachebufs.first = data->pathcachebufs.first;
758         psys->pathcachebufs.last = data->pathcachebufs.last;
759         psys->totcached = data->totcached;
760         psys->childcache = data->childcache;
761         psys->childcachebufs.first = data->childcachebufs.first;
762         psys->childcachebufs.last = data->childcachebufs.last;
763         psys->totchildcache = data->totchildcache;
764
765         psmd->dm = data->dm;
766         psmd->totdmvert = data->totdmvert;
767         psmd->totdmedge = data->totdmedge;
768         psmd->totdmface = data->totdmface;
769         psmd->flag &= ~eParticleSystemFlag_psys_updated;
770
771         if (psmd->dm)
772                 psys_calc_dmcache(ob, psmd->dm, psys);
773
774         MEM_freeN(data);
775         psys->renderdata = NULL;
776
777         /* restore particle display percentage */
778         disp = psys_get_current_display_percentage(psys);
779
780         if (disp != render_disp) {
781                 PARTICLE_P;
782
783                 LOOP_PARTICLES {
784                         if (psys_frand(psys, p) > disp)
785                                 pa->flag |= PARS_NO_DISP;
786                         else
787                                 pa->flag &= ~PARS_NO_DISP;
788                 }
789         }
790 }
791
792 /* BMESH_TODO, for orig face data, we need to use MPoly */
793
794 int psys_render_simplify_distribution(ParticleThreadContext *ctx, int tot)
795 {
796         DerivedMesh *dm = ctx->dm;
797         Mesh *me = (Mesh *)(ctx->sim.ob->data);
798         MFace *mf, *mface;
799         MVert *mvert;
800         ParticleRenderData *data;
801         ParticleRenderElem *elems, *elem;
802         ParticleSettings *part = ctx->sim.psys->part;
803         float *facearea, (*facecenter)[3], size[3], fac, powrate, scaleclamp;
804         float co1[3], co2[3], co3[3], co4[3], lambda, arearatio, t, area, viewport;
805         double vprate;
806         int *facetotvert;
807         int a, b, totorigface, totface, newtot, skipped;
808
809         /* double lookup */
810         const int *index_mf_to_mpoly;
811         const int *index_mp_to_orig;
812
813         if (part->ren_as != PART_DRAW_PATH || !(part->draw & PART_DRAW_REN_STRAND))
814                 return tot;
815         if (!ctx->sim.psys->renderdata)
816                 return tot;
817
818         data = ctx->sim.psys->renderdata;
819         if (data->timeoffset)
820                 return 0;
821         if (!(part->simplify_flag & PART_SIMPLIFY_ENABLE))
822                 return tot;
823
824         mvert = dm->getVertArray(dm);
825         mface = dm->getTessFaceArray(dm);
826         totface = dm->getNumTessFaces(dm);
827         totorigface = me->totpoly;
828
829         if (totface == 0 || totorigface == 0)
830                 return tot;
831
832         index_mf_to_mpoly = dm->getTessFaceDataArray(dm, CD_ORIGINDEX);
833         index_mp_to_orig  = dm->getPolyDataArray(dm, CD_ORIGINDEX);
834         if (index_mf_to_mpoly == NULL) {
835                 index_mp_to_orig = NULL;
836         }
837
838         facearea = MEM_callocN(sizeof(float) * totorigface, "SimplifyFaceArea");
839         facecenter = MEM_callocN(sizeof(float[3]) * totorigface, "SimplifyFaceCenter");
840         facetotvert = MEM_callocN(sizeof(int) * totorigface, "SimplifyFaceArea");
841         elems = MEM_callocN(sizeof(ParticleRenderElem) * totorigface, "SimplifyFaceElem");
842
843         if (data->elems)
844                 MEM_freeN(data->elems);
845
846         data->do_simplify = TRUE;
847         data->elems = elems;
848         data->index_mf_to_mpoly = index_mf_to_mpoly;
849         data->index_mp_to_orig  = index_mp_to_orig;
850
851         /* compute number of children per original face */
852         for (a = 0; a < tot; a++) {
853                 b = (index_mf_to_mpoly) ? DM_origindex_mface_mpoly(index_mf_to_mpoly, index_mp_to_orig, ctx->index[a]) : ctx->index[a];
854                 if (b != ORIGINDEX_NONE) {
855                         elems[b].totchild++;
856                 }
857         }
858
859         /* compute areas and centers of original faces */
860         for (mf = mface, a = 0; a < totface; a++, mf++) {
861                 b = (index_mf_to_mpoly) ? DM_origindex_mface_mpoly(index_mf_to_mpoly, index_mp_to_orig, a) : a;
862
863                 if (b != ORIGINDEX_NONE) {
864                         copy_v3_v3(co1, mvert[mf->v1].co);
865                         copy_v3_v3(co2, mvert[mf->v2].co);
866                         copy_v3_v3(co3, mvert[mf->v3].co);
867
868                         add_v3_v3(facecenter[b], co1);
869                         add_v3_v3(facecenter[b], co2);
870                         add_v3_v3(facecenter[b], co3);
871
872                         if (mf->v4) {
873                                 copy_v3_v3(co4, mvert[mf->v4].co);
874                                 add_v3_v3(facecenter[b], co4);
875                                 facearea[b] += area_quad_v3(co1, co2, co3, co4);
876                                 facetotvert[b] += 4;
877                         }
878                         else {
879                                 facearea[b] += area_tri_v3(co1, co2, co3);
880                                 facetotvert[b] += 3;
881                         }
882                 }
883         }
884
885         for (a = 0; a < totorigface; a++)
886                 if (facetotvert[a] > 0)
887                         mul_v3_fl(facecenter[a], 1.0f / facetotvert[a]);
888
889         /* for conversion from BU area / pixel area to reference screen size */
890         BKE_mesh_texspace_get(me, 0, 0, size);
891         fac = ((size[0] + size[1] + size[2]) / 3.0f) / part->simplify_refsize;
892         fac = fac * fac;
893
894         powrate = log(0.5f) / log(part->simplify_rate * 0.5f);
895         if (part->simplify_flag & PART_SIMPLIFY_VIEWPORT)
896                 vprate = pow(1.0f - part->simplify_viewport, 5.0);
897         else
898                 vprate = 1.0;
899
900         /* set simplification parameters per original face */
901         for (a = 0, elem = elems; a < totorigface; a++, elem++) {
902                 area = psys_render_projected_area(ctx->sim.psys, facecenter[a], facearea[a], vprate, &viewport);
903                 arearatio = fac * area / facearea[a];
904
905                 if ((arearatio < 1.0f || viewport < 1.0f) && elem->totchild) {
906                         /* lambda is percentage of elements to keep */
907                         lambda = (arearatio < 1.0f) ? powf(arearatio, powrate) : 1.0f;
908                         lambda *= viewport;
909
910                         lambda = MAX2(lambda, 1.0f / elem->totchild);
911
912                         /* compute transition region */
913                         t = part->simplify_transition;
914                         elem->t = (lambda - t < 0.0f) ? lambda : (lambda + t > 1.0f) ? 1.0f - lambda : t;
915                         elem->reduce = 1;
916
917                         /* scale at end and beginning of the transition region */
918                         elem->scalemax = (lambda + t < 1.0f) ? 1.0f / lambda : 1.0f / (1.0f - elem->t * elem->t / t);
919                         elem->scalemin = (lambda + t < 1.0f) ? 0.0f : elem->scalemax * (1.0f - elem->t / t);
920
921                         elem->scalemin = sqrt(elem->scalemin);
922                         elem->scalemax = sqrt(elem->scalemax);
923
924                         /* clamp scaling */
925                         scaleclamp = (float)min_ii(elem->totchild, 10);
926                         elem->scalemin = MIN2(scaleclamp, elem->scalemin);
927                         elem->scalemax = MIN2(scaleclamp, elem->scalemax);
928
929                         /* extend lambda to include transition */
930                         lambda = lambda + elem->t;
931                         if (lambda > 1.0f)
932                                 lambda = 1.0f;
933                 }
934                 else {
935                         lambda = arearatio;
936
937                         elem->scalemax = 1.0f; //sqrt(lambda);
938                         elem->scalemin = 1.0f; //sqrt(lambda);
939                         elem->reduce = 0;
940                 }
941
942                 elem->lambda = lambda;
943                 elem->scalemin = sqrt(elem->scalemin);
944                 elem->scalemax = sqrt(elem->scalemax);
945                 elem->curchild = 0;
946         }
947
948         MEM_freeN(facearea);
949         MEM_freeN(facecenter);
950         MEM_freeN(facetotvert);
951
952         /* move indices and set random number skipping */
953         ctx->skip = MEM_callocN(sizeof(int) * tot, "SimplificationSkip");
954
955         skipped = 0;
956         for (a = 0, newtot = 0; a < tot; a++) {
957                 b = (index_mf_to_mpoly) ? DM_origindex_mface_mpoly(index_mf_to_mpoly, index_mp_to_orig, ctx->index[a]) : ctx->index[a];
958
959                 if (b != ORIGINDEX_NONE) {
960                         if (elems[b].curchild++ < ceil(elems[b].lambda * elems[b].totchild)) {
961                                 ctx->index[newtot] = ctx->index[a];
962                                 ctx->skip[newtot] = skipped;
963                                 skipped = 0;
964                                 newtot++;
965                         }
966                         else skipped++;
967                 }
968                 else skipped++;
969         }
970
971         for (a = 0, elem = elems; a < totorigface; a++, elem++)
972                 elem->curchild = 0;
973
974         return newtot;
975 }
976
977 int psys_render_simplify_params(ParticleSystem *psys, ChildParticle *cpa, float *params)
978 {
979         ParticleRenderData *data;
980         ParticleRenderElem *elem;
981         float x, w, scale, alpha, lambda, t, scalemin, scalemax;
982         int b;
983
984         if (!(psys->renderdata && (psys->part->simplify_flag & PART_SIMPLIFY_ENABLE)))
985                 return 0;
986         
987         data = psys->renderdata;
988         if (!data->do_simplify)
989                 return 0;
990         b = (data->index_mf_to_mpoly) ? DM_origindex_mface_mpoly(data->index_mf_to_mpoly, data->index_mp_to_orig, cpa->num) : cpa->num;
991         if (b == ORIGINDEX_NONE) {
992                 return 0;
993         }
994
995         elem = &data->elems[b];
996
997         lambda = elem->lambda;
998         t = elem->t;
999         scalemin = elem->scalemin;
1000         scalemax = elem->scalemax;
1001
1002         if (!elem->reduce) {
1003                 scale = scalemin;
1004                 alpha = 1.0f;
1005         }
1006         else {
1007                 x = (elem->curchild + 0.5f) / elem->totchild;
1008                 if (x < lambda - t) {
1009                         scale = scalemax;
1010                         alpha = 1.0f;
1011                 }
1012                 else if (x >= lambda + t) {
1013                         scale = scalemin;
1014                         alpha = 0.0f;
1015                 }
1016                 else {
1017                         w = (lambda + t - x) / (2.0f * t);
1018                         scale = scalemin + (scalemax - scalemin) * w;
1019                         alpha = w;
1020                 }
1021         }
1022
1023         params[0] = scale;
1024         params[1] = alpha;
1025
1026         elem->curchild++;
1027
1028         return 1;
1029 }
1030
1031 /************************************************/
1032 /*                      Interpolation                                           */
1033 /************************************************/
1034 static float interpolate_particle_value(float v1, float v2, float v3, float v4, const float w[4], int four)
1035 {
1036         float value;
1037
1038         value = w[0] * v1 + w[1] * v2 + w[2] * v3;
1039         if (four)
1040                 value += w[3] * v4;
1041
1042         CLAMP(value, 0.f, 1.f);
1043         
1044         return value;
1045 }
1046
1047 void psys_interpolate_particle(short type, ParticleKey keys[4], float dt, ParticleKey *result, int velocity)
1048 {
1049         float t[4];
1050
1051         if (type < 0) {
1052                 interp_cubic_v3(result->co, result->vel, keys[1].co, keys[1].vel, keys[2].co, keys[2].vel, dt);
1053         }
1054         else {
1055                 key_curve_position_weights(dt, t, type);
1056
1057                 interp_v3_v3v3v3v3(result->co, keys[0].co, keys[1].co, keys[2].co, keys[3].co, t);
1058
1059                 if (velocity) {
1060                         float temp[3];
1061
1062                         if (dt > 0.999f) {
1063                                 key_curve_position_weights(dt - 0.001f, t, type);
1064                                 interp_v3_v3v3v3v3(temp, keys[0].co, keys[1].co, keys[2].co, keys[3].co, t);
1065                                 sub_v3_v3v3(result->vel, result->co, temp);
1066                         }
1067                         else {
1068                                 key_curve_position_weights(dt + 0.001f, t, type);
1069                                 interp_v3_v3v3v3v3(temp, keys[0].co, keys[1].co, keys[2].co, keys[3].co, t);
1070                                 sub_v3_v3v3(result->vel, temp, result->co);
1071                         }
1072                 }
1073         }
1074 }
1075
1076
1077
1078 typedef struct ParticleInterpolationData {
1079         HairKey *hkey[2];
1080
1081         DerivedMesh *dm;
1082         MVert *mvert[2];
1083
1084         int keyed;
1085         ParticleKey *kkey[2];
1086
1087         PointCache *cache;
1088         PTCacheMem *pm;
1089
1090         PTCacheEditPoint *epoint;
1091         PTCacheEditKey *ekey[2];
1092
1093         float birthtime, dietime;
1094         int bspline;
1095 } ParticleInterpolationData;
1096 /* Assumes pointcache->mem_cache exists, so for disk cached particles call psys_make_temp_pointcache() before use */
1097 /* It uses ParticleInterpolationData->pm to store the current memory cache frame so it's thread safe. */
1098 static void get_pointcache_keys_for_time(Object *UNUSED(ob), PointCache *cache, PTCacheMem **cur, int index, float t, ParticleKey *key1, ParticleKey *key2)
1099 {
1100         static PTCacheMem *pm = NULL;
1101         int index1, index2;
1102
1103         if (index < 0) { /* initialize */
1104                 *cur = cache->mem_cache.first;
1105
1106                 if (*cur)
1107                         *cur = (*cur)->next;
1108         }
1109         else {
1110                 if (*cur) {
1111                         while (*cur && (*cur)->next && (float)(*cur)->frame < t)
1112                                 *cur = (*cur)->next;
1113
1114                         pm = *cur;
1115
1116                         index2 = BKE_ptcache_mem_index_find(pm, index);
1117                         index1 = BKE_ptcache_mem_index_find(pm->prev, index);
1118
1119                         BKE_ptcache_make_particle_key(key2, index2, pm->data, (float)pm->frame);
1120                         if (index1 < 0)
1121                                 copy_particle_key(key1, key2, 1);
1122                         else
1123                                 BKE_ptcache_make_particle_key(key1, index1, pm->prev->data, (float)pm->prev->frame);
1124                 }
1125                 else if (cache->mem_cache.first) {
1126                         pm = cache->mem_cache.first;
1127                         index2 = BKE_ptcache_mem_index_find(pm, index);
1128                         BKE_ptcache_make_particle_key(key2, index2, pm->data, (float)pm->frame);
1129                         copy_particle_key(key1, key2, 1);
1130                 }
1131         }
1132 }
1133 static int get_pointcache_times_for_particle(PointCache *cache, int index, float *start, float *end)
1134 {
1135         PTCacheMem *pm;
1136         int ret = 0;
1137
1138         for (pm = cache->mem_cache.first; pm; pm = pm->next) {
1139                 if (BKE_ptcache_mem_index_find(pm, index) >= 0) {
1140                         *start = pm->frame;
1141                         ret++;
1142                         break;
1143                 }
1144         }
1145
1146         for (pm = cache->mem_cache.last; pm; pm = pm->prev) {
1147                 if (BKE_ptcache_mem_index_find(pm, index) >= 0) {
1148                         *end = pm->frame;
1149                         ret++;
1150                         break;
1151                 }
1152         }
1153
1154         return ret == 2;
1155 }
1156
1157 float psys_get_dietime_from_cache(PointCache *cache, int index)
1158 {
1159         PTCacheMem *pm;
1160         int dietime = 10000000; /* some max value so that we can default to pa->time+lifetime */
1161
1162         for (pm = cache->mem_cache.last; pm; pm = pm->prev) {
1163                 if (BKE_ptcache_mem_index_find(pm, index) >= 0)
1164                         return (float)pm->frame;
1165         }
1166
1167         return (float)dietime;
1168 }
1169
1170 static void init_particle_interpolation(Object *ob, ParticleSystem *psys, ParticleData *pa, ParticleInterpolationData *pind)
1171 {
1172
1173         if (pind->epoint) {
1174                 PTCacheEditPoint *point = pind->epoint;
1175
1176                 pind->ekey[0] = point->keys;
1177                 pind->ekey[1] = point->totkey > 1 ? point->keys + 1 : NULL;
1178
1179                 pind->birthtime = *(point->keys->time);
1180                 pind->dietime = *((point->keys + point->totkey - 1)->time);
1181         }
1182         else if (pind->keyed) {
1183                 ParticleKey *key = pa->keys;
1184                 pind->kkey[0] = key;
1185                 pind->kkey[1] = pa->totkey > 1 ? key + 1 : NULL;
1186
1187                 pind->birthtime = key->time;
1188                 pind->dietime = (key + pa->totkey - 1)->time;
1189         }
1190         else if (pind->cache) {
1191                 float start = 0.0f, end = 0.0f;
1192                 get_pointcache_keys_for_time(ob, pind->cache, &pind->pm, -1, 0.0f, NULL, NULL);
1193                 pind->birthtime = pa ? pa->time : pind->cache->startframe;
1194                 pind->dietime = pa ? pa->dietime : pind->cache->endframe;
1195
1196                 if (get_pointcache_times_for_particle(pind->cache, pa - psys->particles, &start, &end)) {
1197                         pind->birthtime = MAX2(pind->birthtime, start);
1198                         pind->dietime = MIN2(pind->dietime, end);
1199                 }
1200         }
1201         else {
1202                 HairKey *key = pa->hair;
1203                 pind->hkey[0] = key;
1204                 pind->hkey[1] = key + 1;
1205
1206                 pind->birthtime = key->time;
1207                 pind->dietime = (key + pa->totkey - 1)->time;
1208
1209                 if (pind->dm) {
1210                         pind->mvert[0] = CDDM_get_vert(pind->dm, pa->hair_index);
1211                         pind->mvert[1] = pind->mvert[0] + 1;
1212                 }
1213         }
1214 }
1215 static void edit_to_particle(ParticleKey *key, PTCacheEditKey *ekey)
1216 {
1217         copy_v3_v3(key->co, ekey->co);
1218         if (ekey->vel) {
1219                 copy_v3_v3(key->vel, ekey->vel);
1220         }
1221         key->time = *(ekey->time);
1222 }
1223 static void hair_to_particle(ParticleKey *key, HairKey *hkey)
1224 {
1225         copy_v3_v3(key->co, hkey->co);
1226         key->time = hkey->time;
1227 }
1228
1229 static void mvert_to_particle(ParticleKey *key, MVert *mvert, HairKey *hkey)
1230 {
1231         copy_v3_v3(key->co, mvert->co);
1232         key->time = hkey->time;
1233 }
1234
1235 static void do_particle_interpolation(ParticleSystem *psys, int p, ParticleData *pa, float t, ParticleInterpolationData *pind, ParticleKey *result)
1236 {
1237         PTCacheEditPoint *point = pind->epoint;
1238         ParticleKey keys[4];
1239         int point_vel = (point && point->keys->vel);
1240         float real_t, dfra, keytime, invdt = 1.f;
1241
1242         /* billboards wont fill in all of these, so start cleared */
1243         memset(keys, 0, sizeof(keys));
1244
1245         /* interpret timing and find keys */
1246         if (point) {
1247                 if (result->time < 0.0f)
1248                         real_t = -result->time;
1249                 else
1250                         real_t = *(pind->ekey[0]->time) + t * (*(pind->ekey[0][point->totkey - 1].time) - *(pind->ekey[0]->time));
1251
1252                 while (*(pind->ekey[1]->time) < real_t)
1253                         pind->ekey[1]++;
1254
1255                 pind->ekey[0] = pind->ekey[1] - 1;
1256         }
1257         else if (pind->keyed) {
1258                 /* we have only one key, so let's use that */
1259                 if (pind->kkey[1] == NULL) {
1260                         copy_particle_key(result, pind->kkey[0], 1);
1261                         return;
1262                 }
1263
1264                 if (result->time < 0.0f)
1265                         real_t = -result->time;
1266                 else
1267                         real_t = pind->kkey[0]->time + t * (pind->kkey[0][pa->totkey - 1].time - pind->kkey[0]->time);
1268
1269                 if (psys->part->phystype == PART_PHYS_KEYED && psys->flag & PSYS_KEYED_TIMING) {
1270                         ParticleTarget *pt = psys->targets.first;
1271
1272                         pt = pt->next;
1273
1274                         while (pt && pa->time + pt->time < real_t)
1275                                 pt = pt->next;
1276
1277                         if (pt) {
1278                                 pt = pt->prev;
1279
1280                                 if (pa->time + pt->time + pt->duration > real_t)
1281                                         real_t = pa->time + pt->time;
1282                         }
1283                         else
1284                                 real_t = pa->time + ((ParticleTarget *)psys->targets.last)->time;
1285                 }
1286
1287                 CLAMP(real_t, pa->time, pa->dietime);
1288
1289                 while (pind->kkey[1]->time < real_t)
1290                         pind->kkey[1]++;
1291                 
1292                 pind->kkey[0] = pind->kkey[1] - 1;
1293         }
1294         else if (pind->cache) {
1295                 if (result->time < 0.0f) /* flag for time in frames */
1296                         real_t = -result->time;
1297                 else
1298                         real_t = pa->time + t * (pa->dietime - pa->time);
1299         }
1300         else {
1301                 if (result->time < 0.0f)
1302                         real_t = -result->time;
1303                 else
1304                         real_t = pind->hkey[0]->time + t * (pind->hkey[0][pa->totkey - 1].time - pind->hkey[0]->time);
1305
1306                 while (pind->hkey[1]->time < real_t) {
1307                         pind->hkey[1]++;
1308                         pind->mvert[1]++;
1309                 }
1310
1311                 pind->hkey[0] = pind->hkey[1] - 1;
1312         }
1313
1314         /* set actual interpolation keys */
1315         if (point) {
1316                 edit_to_particle(keys + 1, pind->ekey[0]);
1317                 edit_to_particle(keys + 2, pind->ekey[1]);
1318         }
1319         else if (pind->dm) {
1320                 pind->mvert[0] = pind->mvert[1] - 1;
1321                 mvert_to_particle(keys + 1, pind->mvert[0], pind->hkey[0]);
1322                 mvert_to_particle(keys + 2, pind->mvert[1], pind->hkey[1]);
1323         }
1324         else if (pind->keyed) {
1325                 memcpy(keys + 1, pind->kkey[0], sizeof(ParticleKey));
1326                 memcpy(keys + 2, pind->kkey[1], sizeof(ParticleKey));
1327         }
1328         else if (pind->cache) {
1329                 get_pointcache_keys_for_time(NULL, pind->cache, &pind->pm, p, real_t, keys + 1, keys + 2);
1330         }
1331         else {
1332                 hair_to_particle(keys + 1, pind->hkey[0]);
1333                 hair_to_particle(keys + 2, pind->hkey[1]);
1334         }
1335
1336         /* set secondary interpolation keys for hair */
1337         if (!pind->keyed && !pind->cache && !point_vel) {
1338                 if (point) {
1339                         if (pind->ekey[0] != point->keys)
1340                                 edit_to_particle(keys, pind->ekey[0] - 1);
1341                         else
1342                                 edit_to_particle(keys, pind->ekey[0]);
1343                 }
1344                 else if (pind->dm) {
1345                         if (pind->hkey[0] != pa->hair)
1346                                 mvert_to_particle(keys, pind->mvert[0] - 1, pind->hkey[0] - 1);
1347                         else
1348                                 mvert_to_particle(keys, pind->mvert[0], pind->hkey[0]);
1349                 }
1350                 else {
1351                         if (pind->hkey[0] != pa->hair)
1352                                 hair_to_particle(keys, pind->hkey[0] - 1);
1353                         else
1354                                 hair_to_particle(keys, pind->hkey[0]);
1355                 }
1356
1357                 if (point) {
1358                         if (pind->ekey[1] != point->keys + point->totkey - 1)
1359                                 edit_to_particle(keys + 3, pind->ekey[1] + 1);
1360                         else
1361                                 edit_to_particle(keys + 3, pind->ekey[1]);
1362                 }
1363                 else if (pind->dm) {
1364                         if (pind->hkey[1] != pa->hair + pa->totkey - 1)
1365                                 mvert_to_particle(keys + 3, pind->mvert[1] + 1, pind->hkey[1] + 1);
1366                         else
1367                                 mvert_to_particle(keys + 3, pind->mvert[1], pind->hkey[1]);
1368                 }
1369                 else {
1370                         if (pind->hkey[1] != pa->hair + pa->totkey - 1)
1371                                 hair_to_particle(keys + 3, pind->hkey[1] + 1);
1372                         else
1373                                 hair_to_particle(keys + 3, pind->hkey[1]);
1374                 }
1375         }
1376
1377         dfra = keys[2].time - keys[1].time;
1378         keytime = (real_t - keys[1].time) / dfra;
1379
1380         /* convert velocity to timestep size */
1381         if (pind->keyed || pind->cache || point_vel) {
1382                 invdt = dfra * 0.04f * (psys ? psys->part->timetweak : 1.f);
1383                 mul_v3_fl(keys[1].vel, invdt);
1384                 mul_v3_fl(keys[2].vel, invdt);
1385                 interp_qt_qtqt(result->rot, keys[1].rot, keys[2].rot, keytime);
1386         }
1387
1388         /* now we should have in chronologiacl order k1<=k2<=t<=k3<=k4 with keytime between [0, 1]->[k2, k3] (k1 & k4 used for cardinal & bspline interpolation)*/
1389         psys_interpolate_particle((pind->keyed || pind->cache || point_vel) ? -1 /* signal for cubic interpolation */
1390                                   : (pind->bspline ? KEY_BSPLINE : KEY_CARDINAL),
1391                                   keys, keytime, result, 1);
1392
1393         /* the velocity needs to be converted back from cubic interpolation */
1394         if (pind->keyed || pind->cache || point_vel)
1395                 mul_v3_fl(result->vel, 1.f / invdt);
1396 }
1397
1398 static void interpolate_pathcache(ParticleCacheKey *first, float t, ParticleCacheKey *result)
1399 {
1400         int i = 0;
1401         ParticleCacheKey *cur = first;
1402
1403         /* scale the requested time to fit the entire path even if the path is cut early */
1404         t *= (first + first->steps)->time;
1405
1406         while (i < first->steps && cur->time < t)
1407                 cur++;
1408
1409         if (cur->time == t)
1410                 *result = *cur;
1411         else {
1412                 float dt = (t - (cur - 1)->time) / (cur->time - (cur - 1)->time);
1413                 interp_v3_v3v3(result->co, (cur - 1)->co, cur->co, dt);
1414                 interp_v3_v3v3(result->vel, (cur - 1)->vel, cur->vel, dt);
1415                 interp_qt_qtqt(result->rot, (cur - 1)->rot, cur->rot, dt);
1416                 result->time = t;
1417         }
1418
1419         /* first is actual base rotation, others are incremental from first */
1420         if (cur == first || cur - 1 == first)
1421                 copy_qt_qt(result->rot, first->rot);
1422         else
1423                 mul_qt_qtqt(result->rot, first->rot, result->rot);
1424 }
1425
1426 /************************************************/
1427 /*                      Particles on a dm                                       */
1428 /************************************************/
1429 /* interpolate a location on a face based on face coordinates */
1430 void psys_interpolate_face(MVert *mvert, MFace *mface, MTFace *tface, float (*orcodata)[3],
1431                            float w[4], float vec[3], float nor[3], float utan[3], float vtan[3],
1432                            float orco[3], float ornor[3])
1433 {
1434         float *v1 = 0, *v2 = 0, *v3 = 0, *v4 = 0;
1435         float e1[3], e2[3], s1, s2, t1, t2;
1436         float *uv1, *uv2, *uv3, *uv4;
1437         float n1[3], n2[3], n3[3], n4[3];
1438         float tuv[4][2];
1439         float *o1, *o2, *o3, *o4;
1440
1441         v1 = mvert[mface->v1].co;
1442         v2 = mvert[mface->v2].co;
1443         v3 = mvert[mface->v3].co;
1444
1445         normal_short_to_float_v3(n1, mvert[mface->v1].no);
1446         normal_short_to_float_v3(n2, mvert[mface->v2].no);
1447         normal_short_to_float_v3(n3, mvert[mface->v3].no);
1448
1449         if (mface->v4) {
1450                 v4 = mvert[mface->v4].co;
1451                 normal_short_to_float_v3(n4, mvert[mface->v4].no);
1452                 
1453                 interp_v3_v3v3v3v3(vec, v1, v2, v3, v4, w);
1454
1455                 if (nor) {
1456                         if (mface->flag & ME_SMOOTH)
1457                                 interp_v3_v3v3v3v3(nor, n1, n2, n3, n4, w);
1458                         else
1459                                 normal_quad_v3(nor, v1, v2, v3, v4);
1460                 }
1461         }
1462         else {
1463                 interp_v3_v3v3v3(vec, v1, v2, v3, w);
1464
1465                 if (nor) {
1466                         if (mface->flag & ME_SMOOTH)
1467                                 interp_v3_v3v3v3(nor, n1, n2, n3, w);
1468                         else
1469                                 normal_tri_v3(nor, v1, v2, v3);
1470                 }
1471         }
1472         
1473         /* calculate tangent vectors */
1474         if (utan && vtan) {
1475                 if (tface) {
1476                         uv1 = tface->uv[0];
1477                         uv2 = tface->uv[1];
1478                         uv3 = tface->uv[2];
1479                         uv4 = tface->uv[3];
1480                 }
1481                 else {
1482                         uv1 = tuv[0]; uv2 = tuv[1]; uv3 = tuv[2]; uv4 = tuv[3];
1483                         map_to_sphere(uv1, uv1 + 1, v1[0], v1[1], v1[2]);
1484                         map_to_sphere(uv2, uv2 + 1, v2[0], v2[1], v2[2]);
1485                         map_to_sphere(uv3, uv3 + 1, v3[0], v3[1], v3[2]);
1486                         if (v4)
1487                                 map_to_sphere(uv4, uv4 + 1, v4[0], v4[1], v4[2]);
1488                 }
1489
1490                 if (v4) {
1491                         s1 = uv3[0] - uv1[0];
1492                         s2 = uv4[0] - uv1[0];
1493
1494                         t1 = uv3[1] - uv1[1];
1495                         t2 = uv4[1] - uv1[1];
1496
1497                         sub_v3_v3v3(e1, v3, v1);
1498                         sub_v3_v3v3(e2, v4, v1);
1499                 }
1500                 else {
1501                         s1 = uv2[0] - uv1[0];
1502                         s2 = uv3[0] - uv1[0];
1503
1504                         t1 = uv2[1] - uv1[1];
1505                         t2 = uv3[1] - uv1[1];
1506
1507                         sub_v3_v3v3(e1, v2, v1);
1508                         sub_v3_v3v3(e2, v3, v1);
1509                 }
1510
1511                 vtan[0] = (s1 * e2[0] - s2 * e1[0]);
1512                 vtan[1] = (s1 * e2[1] - s2 * e1[1]);
1513                 vtan[2] = (s1 * e2[2] - s2 * e1[2]);
1514
1515                 utan[0] = (t1 * e2[0] - t2 * e1[0]);
1516                 utan[1] = (t1 * e2[1] - t2 * e1[1]);
1517                 utan[2] = (t1 * e2[2] - t2 * e1[2]);
1518         }
1519
1520         if (orco) {
1521                 if (orcodata) {
1522                         o1 = orcodata[mface->v1];
1523                         o2 = orcodata[mface->v2];
1524                         o3 = orcodata[mface->v3];
1525
1526                         if (mface->v4) {
1527                                 o4 = orcodata[mface->v4];
1528
1529                                 interp_v3_v3v3v3v3(orco, o1, o2, o3, o4, w);
1530
1531                                 if (ornor)
1532                                         normal_quad_v3(ornor, o1, o2, o3, o4);
1533                         }
1534                         else {
1535                                 interp_v3_v3v3v3(orco, o1, o2, o3, w);
1536
1537                                 if (ornor)
1538                                         normal_tri_v3(ornor, o1, o2, o3);
1539                         }
1540                 }
1541                 else {
1542                         copy_v3_v3(orco, vec);
1543                         if (ornor && nor)
1544                                 copy_v3_v3(ornor, nor);
1545                 }
1546         }
1547 }
1548 void psys_interpolate_uvs(const MTFace *tface, int quad, const float w[4], float uvco[2])
1549 {
1550         float v10 = tface->uv[0][0];
1551         float v11 = tface->uv[0][1];
1552         float v20 = tface->uv[1][0];
1553         float v21 = tface->uv[1][1];
1554         float v30 = tface->uv[2][0];
1555         float v31 = tface->uv[2][1];
1556         float v40, v41;
1557
1558         if (quad) {
1559                 v40 = tface->uv[3][0];
1560                 v41 = tface->uv[3][1];
1561
1562                 uvco[0] = w[0] * v10 + w[1] * v20 + w[2] * v30 + w[3] * v40;
1563                 uvco[1] = w[0] * v11 + w[1] * v21 + w[2] * v31 + w[3] * v41;
1564         }
1565         else {
1566                 uvco[0] = w[0] * v10 + w[1] * v20 + w[2] * v30;
1567                 uvco[1] = w[0] * v11 + w[1] * v21 + w[2] * v31;
1568         }
1569 }
1570
1571 void psys_interpolate_mcol(const MCol *mcol, int quad, const float w[4], MCol *mc)
1572 {
1573         char *cp, *cp1, *cp2, *cp3, *cp4;
1574
1575         cp = (char *)mc;
1576         cp1 = (char *)&mcol[0];
1577         cp2 = (char *)&mcol[1];
1578         cp3 = (char *)&mcol[2];
1579         
1580         if (quad) {
1581                 cp4 = (char *)&mcol[3];
1582
1583                 cp[0] = (int)(w[0] * cp1[0] + w[1] * cp2[0] + w[2] * cp3[0] + w[3] * cp4[0]);
1584                 cp[1] = (int)(w[0] * cp1[1] + w[1] * cp2[1] + w[2] * cp3[1] + w[3] * cp4[1]);
1585                 cp[2] = (int)(w[0] * cp1[2] + w[1] * cp2[2] + w[2] * cp3[2] + w[3] * cp4[2]);
1586                 cp[3] = (int)(w[0] * cp1[3] + w[1] * cp2[3] + w[2] * cp3[3] + w[3] * cp4[3]);
1587         }
1588         else {
1589                 cp[0] = (int)(w[0] * cp1[0] + w[1] * cp2[0] + w[2] * cp3[0]);
1590                 cp[1] = (int)(w[0] * cp1[1] + w[1] * cp2[1] + w[2] * cp3[1]);
1591                 cp[2] = (int)(w[0] * cp1[2] + w[1] * cp2[2] + w[2] * cp3[2]);
1592                 cp[3] = (int)(w[0] * cp1[3] + w[1] * cp2[3] + w[2] * cp3[3]);
1593         }
1594 }
1595
1596 static float psys_interpolate_value_from_verts(DerivedMesh *dm, short from, int index, const float fw[4], const float *values)
1597 {
1598         if (values == 0 || index == -1)
1599                 return 0.0;
1600
1601         switch (from) {
1602                 case PART_FROM_VERT:
1603                         return values[index];
1604                 case PART_FROM_FACE:
1605                 case PART_FROM_VOLUME:
1606                 {
1607                         MFace *mf = dm->getTessFaceData(dm, index, CD_MFACE);
1608                         return interpolate_particle_value(values[mf->v1], values[mf->v2], values[mf->v3], values[mf->v4], fw, mf->v4);
1609                 }
1610                         
1611         }
1612         return 0.0f;
1613 }
1614
1615 /* conversion of pa->fw to origspace layer coordinates */
1616 static void psys_w_to_origspace(const float w[4], float uv[2])
1617 {
1618         uv[0] = w[1] + w[2];
1619         uv[1] = w[2] + w[3];
1620 }
1621
1622 /* conversion of pa->fw to weights in face from origspace */
1623 static void psys_origspace_to_w(OrigSpaceFace *osface, int quad, const float w[4], float neww[4])
1624 {
1625         float v[4][3], co[3];
1626
1627         v[0][0] = osface->uv[0][0]; v[0][1] = osface->uv[0][1]; v[0][2] = 0.0f;
1628         v[1][0] = osface->uv[1][0]; v[1][1] = osface->uv[1][1]; v[1][2] = 0.0f;
1629         v[2][0] = osface->uv[2][0]; v[2][1] = osface->uv[2][1]; v[2][2] = 0.0f;
1630
1631         psys_w_to_origspace(w, co);
1632         co[2] = 0.0f;
1633         
1634         if (quad) {
1635                 v[3][0] = osface->uv[3][0]; v[3][1] = osface->uv[3][1]; v[3][2] = 0.0f;
1636                 interp_weights_poly_v3(neww, v, 4, co);
1637         }
1638         else {
1639                 interp_weights_poly_v3(neww, v, 3, co);
1640                 neww[3] = 0.0f;
1641         }
1642 }
1643
1644 /* find the derived mesh face for a particle, set the mf passed. this is slow
1645  * and can be optimized but only for many lookups. returns the face index. */
1646 int psys_particle_dm_face_lookup(Object *ob, DerivedMesh *dm, int index, const float fw[4], struct LinkNode *node)
1647 {
1648         Mesh *me = (Mesh *)ob->data;
1649         MPoly *mpoly;
1650         OrigSpaceFace *osface;
1651         int quad, findex, totface;
1652         float uv[2], (*faceuv)[2];
1653
1654         /* double lookup */
1655         const int *index_mf_to_mpoly = dm->getTessFaceDataArray(dm, CD_ORIGINDEX);
1656         const int *index_mp_to_orig  = dm->getPolyDataArray(dm, CD_ORIGINDEX);
1657         if (index_mf_to_mpoly == NULL) {
1658                 index_mp_to_orig = NULL;
1659         }
1660
1661         mpoly = dm->getPolyArray(dm);
1662         osface = dm->getTessFaceDataArray(dm, CD_ORIGSPACE);
1663
1664         totface = dm->getNumTessFaces(dm);
1665         
1666         if (osface == NULL || index_mf_to_mpoly == NULL) {
1667                 /* Assume we don't need osface data */
1668                 if (index < totface) {
1669                         //printf("\tNO CD_ORIGSPACE, assuming not needed\n");
1670                         return index;
1671                 }
1672                 else {
1673                         printf("\tNO CD_ORIGSPACE, error out of range\n");
1674                         return DMCACHE_NOTFOUND;
1675                 }
1676         }
1677         else if (index >= me->totpoly)
1678                 return DMCACHE_NOTFOUND;  /* index not in the original mesh */
1679
1680         psys_w_to_origspace(fw, uv);
1681         
1682         if (node) { /* we have a linked list of faces that we use, faster! */
1683                 for (; node; node = node->next) {
1684                         findex = GET_INT_FROM_POINTER(node->link);
1685                         faceuv = osface[findex].uv;
1686                         quad = (mpoly[findex].totloop == 4);
1687
1688                         /* check that this intersects - Its possible this misses :/ -
1689                          * could also check its not between */
1690                         if (quad) {
1691                                 if (isect_point_quad_v2(uv, faceuv[0], faceuv[1], faceuv[2], faceuv[3]))
1692                                         return findex;
1693                         }
1694                         else if (isect_point_tri_v2(uv, faceuv[0], faceuv[1], faceuv[2]))
1695                                 return findex;
1696                 }
1697         }
1698         else { /* if we have no node, try every face */
1699                 for (findex = 0; findex < totface; findex++) {
1700                         const int findex_orig = DM_origindex_mface_mpoly(index_mf_to_mpoly, index_mp_to_orig, findex);
1701                         if (findex_orig == index) {
1702                                 faceuv = osface[findex].uv;
1703                                 quad = (mpoly[findex].totloop == 4);
1704
1705                                 /* check that this intersects - Its possible this misses :/ -
1706                                  * could also check its not between */
1707                                 if (quad) {
1708                                         if (isect_point_quad_v2(uv, faceuv[0], faceuv[1], faceuv[2], faceuv[3]))
1709                                                 return findex;
1710                                 }
1711                                 else if (isect_point_tri_v2(uv, faceuv[0], faceuv[1], faceuv[2]))
1712                                         return findex;
1713                         }
1714                 }
1715         }
1716
1717         return DMCACHE_NOTFOUND;
1718 }
1719
1720 static int psys_map_index_on_dm(DerivedMesh *dm, int from, int index, int index_dmcache, const float fw[4], float UNUSED(foffset), int *mapindex, float mapfw[4])
1721 {
1722         if (index < 0)
1723                 return 0;
1724
1725         if (dm->deformedOnly || index_dmcache == DMCACHE_ISCHILD) {
1726                 /* for meshes that are either only deformed or for child particles, the
1727                  * index and fw do not require any mapping, so we can directly use it */
1728                 if (from == PART_FROM_VERT) {
1729                         if (index >= dm->getNumVerts(dm))
1730                                 return 0;
1731
1732                         *mapindex = index;
1733                 }
1734                 else { /* FROM_FACE/FROM_VOLUME */
1735                         if (index >= dm->getNumTessFaces(dm))
1736                                 return 0;
1737
1738                         *mapindex = index;
1739                         copy_v4_v4(mapfw, fw);
1740                 }
1741         }
1742         else {
1743                 /* for other meshes that have been modified, we try to map the particle
1744                  * to their new location, which means a different index, and for faces
1745                  * also a new face interpolation weights */
1746                 if (from == PART_FROM_VERT) {
1747                         if (index_dmcache == DMCACHE_NOTFOUND || index_dmcache > dm->getNumVerts(dm))
1748                                 return 0;
1749
1750                         *mapindex = index_dmcache;
1751                 }
1752                 else { /* FROM_FACE/FROM_VOLUME */
1753                            /* find a face on the derived mesh that uses this face */
1754                         MFace *mface;
1755                         OrigSpaceFace *osface;
1756                         int i;
1757
1758                         i = index_dmcache;
1759
1760                         if (i == DMCACHE_NOTFOUND || i >= dm->getNumTessFaces(dm))
1761                                 return 0;
1762
1763                         *mapindex = i;
1764
1765                         /* modify the original weights to become
1766                          * weights for the derived mesh face */
1767                         osface = dm->getTessFaceDataArray(dm, CD_ORIGSPACE);
1768                         mface = dm->getTessFaceData(dm, i, CD_MFACE);
1769
1770                         if (osface == NULL)
1771                                 mapfw[0] = mapfw[1] = mapfw[2] = mapfw[3] = 0.0f;
1772                         else
1773                                 psys_origspace_to_w(&osface[i], mface->v4, fw, mapfw);
1774                 }
1775         }
1776
1777         return 1;
1778 }
1779
1780 /* interprets particle data to get a point on a mesh in object space */
1781 void psys_particle_on_dm(DerivedMesh *dm, int from, int index, int index_dmcache,
1782                          const float fw[4], float foffset, float vec[3], float nor[3], float utan[3], float vtan[3],
1783                          float orco[3], float ornor[3])
1784 {
1785         float tmpnor[3], mapfw[4];
1786         float (*orcodata)[3];
1787         int mapindex;
1788
1789         if (!psys_map_index_on_dm(dm, from, index, index_dmcache, fw, foffset, &mapindex, mapfw)) {
1790                 if (vec) { vec[0] = vec[1] = vec[2] = 0.0; }
1791                 if (nor) { nor[0] = nor[1] = 0.0; nor[2] = 1.0; }
1792                 if (orco) { orco[0] = orco[1] = orco[2] = 0.0; }
1793                 if (ornor) { ornor[0] = ornor[1] = 0.0; ornor[2] = 1.0; }
1794                 if (utan) { utan[0] = utan[1] = utan[2] = 0.0; }
1795                 if (vtan) { vtan[0] = vtan[1] = vtan[2] = 0.0; }
1796
1797                 return;
1798         }
1799
1800         orcodata = dm->getVertDataArray(dm, CD_ORCO);
1801
1802         if (from == PART_FROM_VERT) {
1803                 dm->getVertCo(dm, mapindex, vec);
1804
1805                 if (nor) {
1806                         dm->getVertNo(dm, mapindex, nor);
1807                         normalize_v3(nor);
1808                 }
1809
1810                 if (orco)
1811                         copy_v3_v3(orco, orcodata[mapindex]);
1812
1813                 if (ornor) {
1814                         dm->getVertNo(dm, mapindex, ornor);
1815                         normalize_v3(ornor);
1816                 }
1817
1818                 if (utan && vtan) {
1819                         utan[0] = utan[1] = utan[2] = 0.0f;
1820                         vtan[0] = vtan[1] = vtan[2] = 0.0f;
1821                 }
1822         }
1823         else { /* PART_FROM_FACE / PART_FROM_VOLUME */
1824                 MFace *mface;
1825                 MTFace *mtface;
1826                 MVert *mvert;
1827
1828                 mface = dm->getTessFaceData(dm, mapindex, CD_MFACE);
1829                 mvert = dm->getVertDataArray(dm, CD_MVERT);
1830                 mtface = CustomData_get_layer(&dm->faceData, CD_MTFACE);
1831
1832                 if (mtface)
1833                         mtface += mapindex;
1834
1835                 if (from == PART_FROM_VOLUME) {
1836                         psys_interpolate_face(mvert, mface, mtface, orcodata, mapfw, vec, tmpnor, utan, vtan, orco, ornor);
1837                         if (nor)
1838                                 copy_v3_v3(nor, tmpnor);
1839
1840                         normalize_v3(tmpnor);  /* XXX Why not normalize tmpnor before copying it into nor??? -- mont29 */
1841                         mul_v3_fl(tmpnor, -foffset);
1842                         add_v3_v3(vec, tmpnor);
1843                 }
1844                 else
1845                         psys_interpolate_face(mvert, mface, mtface, orcodata, mapfw, vec, nor, utan, vtan, orco, ornor);
1846         }
1847 }
1848
1849 float psys_particle_value_from_verts(DerivedMesh *dm, short from, ParticleData *pa, float *values)
1850 {
1851         float mapfw[4];
1852         int mapindex;
1853
1854         if (!psys_map_index_on_dm(dm, from, pa->num, pa->num_dmcache, pa->fuv, pa->foffset, &mapindex, mapfw))
1855                 return 0.0f;
1856         
1857         return psys_interpolate_value_from_verts(dm, from, mapindex, mapfw, values);
1858 }
1859
1860 ParticleSystemModifierData *psys_get_modifier(Object *ob, ParticleSystem *psys)
1861 {
1862         ModifierData *md;
1863         ParticleSystemModifierData *psmd;
1864
1865         for (md = ob->modifiers.first; md; md = md->next) {
1866                 if (md->type == eModifierType_ParticleSystem) {
1867                         psmd = (ParticleSystemModifierData *) md;
1868                         if (psmd->psys == psys) {
1869                                 return psmd;
1870                         }
1871                 }
1872         }
1873         return NULL;
1874 }
1875 /************************************************/
1876 /*                      Particles on a shape                            */
1877 /************************************************/
1878 /* ready for future use */
1879 static void psys_particle_on_shape(int UNUSED(distr), int UNUSED(index),
1880                                    float *UNUSED(fuv), float vec[3], float nor[3], float utan[3], float vtan[3],
1881                                    float orco[3], float ornor[3])
1882 {
1883         /* TODO */
1884         float zerovec[3] = {0.0f, 0.0f, 0.0f};
1885         if (vec) {
1886                 copy_v3_v3(vec, zerovec);
1887         }
1888         if (nor) {
1889                 copy_v3_v3(nor, zerovec);
1890         }
1891         if (utan) {
1892                 copy_v3_v3(utan, zerovec);
1893         }
1894         if (vtan) {
1895                 copy_v3_v3(vtan, zerovec);
1896         }
1897         if (orco) {
1898                 copy_v3_v3(orco, zerovec);
1899         }
1900         if (ornor) {
1901                 copy_v3_v3(ornor, zerovec);
1902         }
1903 }
1904 /************************************************/
1905 /*                      Particles on emitter                            */
1906 /************************************************/
1907 void psys_particle_on_emitter(ParticleSystemModifierData *psmd, int from, int index, int index_dmcache,
1908                               float fuv[4], float foffset, float vec[3], float nor[3], float utan[3], float vtan[3],
1909                               float orco[3], float ornor[3])
1910 {
1911         if (psmd) {
1912                 if (psmd->psys->part->distr == PART_DISTR_GRID && psmd->psys->part->from != PART_FROM_VERT) {
1913                         if (vec)
1914                                 copy_v3_v3(vec, fuv);
1915
1916                         if (orco)
1917                                 copy_v3_v3(orco, fuv);
1918                         return;
1919                 }
1920                 /* we cant use the num_dmcache */
1921                 psys_particle_on_dm(psmd->dm, from, index, index_dmcache, fuv, foffset, vec, nor, utan, vtan, orco, ornor);
1922         }
1923         else
1924                 psys_particle_on_shape(from, index, fuv, vec, nor, utan, vtan, orco, ornor);
1925
1926 }
1927 /************************************************/
1928 /*                      Path Cache                                                      */
1929 /************************************************/
1930
1931 static void do_kink(ParticleKey *state, ParticleKey *par, float *par_rot, float time, float freq, float shape, float amplitude, float flat, short type, short axis, float obmat[4][4], int smooth_start)
1932 {
1933         float kink[3] = {1.f, 0.f, 0.f}, par_vec[3], q1[4] = {1.f, 0.f, 0.f, 0.f};
1934         float t, dt = 1.f, result[3];
1935
1936         if (par == NULL || type == PART_KINK_NO)
1937                 return;
1938
1939         CLAMP(time, 0.f, 1.f);
1940
1941         if (shape != 0.0f && type != PART_KINK_BRAID) {
1942                 if (shape < 0.0f)
1943                         time = (float)pow(time, 1.f + shape);
1944                 else
1945                         time = (float)pow(time, 1.f / (1.f - shape));
1946         }
1947
1948         t = time * freq * (float)M_PI;
1949         
1950         if (smooth_start) {
1951                 dt = fabs(t);
1952                 /* smooth the beginning of kink */
1953                 CLAMP(dt, 0.f, (float)M_PI);
1954                 dt = sin(dt / 2.f);
1955         }
1956
1957         if (type != PART_KINK_RADIAL) {
1958                 float temp[3];
1959
1960                 kink[axis] = 1.f;
1961
1962                 if (obmat)
1963                         mul_mat3_m4_v3(obmat, kink);
1964                 
1965                 if (par_rot)
1966                         mul_qt_v3(par_rot, kink);
1967
1968                 /* make sure kink is normal to strand */
1969                 project_v3_v3v3(temp, kink, par->vel);
1970                 sub_v3_v3(kink, temp);
1971                 normalize_v3(kink);
1972         }
1973
1974         copy_v3_v3(result, state->co);
1975         sub_v3_v3v3(par_vec, par->co, state->co);
1976
1977         switch (type) {
1978                 case PART_KINK_CURL:
1979                 {
1980                         negate_v3(par_vec);
1981
1982                         if (flat > 0.f) {
1983                                 float proj[3];
1984                                 project_v3_v3v3(proj, par_vec, par->vel);
1985                                 madd_v3_v3fl(par_vec, proj, -flat);
1986
1987                                 project_v3_v3v3(proj, par_vec, kink);
1988                                 madd_v3_v3fl(par_vec, proj, -flat);
1989                         }
1990
1991                         axis_angle_to_quat(q1, kink, (float)M_PI / 2.f);
1992
1993                         mul_qt_v3(q1, par_vec);
1994
1995                         madd_v3_v3fl(par_vec, kink, amplitude);
1996
1997                         /* rotate kink vector around strand tangent */
1998                         if (t != 0.f) {
1999                                 axis_angle_to_quat(q1, par->vel, t);
2000                                 mul_qt_v3(q1, par_vec);
2001                         }
2002
2003                         add_v3_v3v3(result, par->co, par_vec);
2004                         break;
2005                 }
2006                 case PART_KINK_RADIAL:
2007                 {
2008                         if (flat > 0.f) {
2009                                 float proj[3];
2010                                 /* flatten along strand */
2011                                 project_v3_v3v3(proj, par_vec, par->vel);
2012                                 madd_v3_v3fl(result, proj, flat);
2013                         }
2014
2015                         madd_v3_v3fl(result, par_vec, -amplitude * (float)sin(t));
2016                         break;
2017                 }
2018                 case PART_KINK_WAVE:
2019                 {
2020                         madd_v3_v3fl(result, kink, amplitude * (float)sin(t));
2021
2022                         if (flat > 0.f) {
2023                                 float proj[3];
2024                                 /* flatten along wave */
2025                                 project_v3_v3v3(proj, par_vec, kink);
2026                                 madd_v3_v3fl(result, proj, flat);
2027
2028                                 /* flatten along strand */
2029                                 project_v3_v3v3(proj, par_vec, par->vel);
2030                                 madd_v3_v3fl(result, proj, flat);
2031                         }
2032                         break;
2033                 }
2034                 case PART_KINK_BRAID:
2035                 {
2036                         float y_vec[3] = {0.f, 1.f, 0.f};
2037                         float z_vec[3] = {0.f, 0.f, 1.f};
2038                         float vec_one[3], state_co[3];
2039                         float inp_y, inp_z, length;
2040                 
2041                         if (par_rot) {
2042                                 mul_qt_v3(par_rot, y_vec);
2043                                 mul_qt_v3(par_rot, z_vec);
2044                         }
2045
2046                         negate_v3(par_vec);
2047                         normalize_v3_v3(vec_one, par_vec);
2048
2049                         inp_y = dot_v3v3(y_vec, vec_one);
2050                         inp_z = dot_v3v3(z_vec, vec_one);
2051
2052                         if (inp_y > 0.5f) {
2053                                 copy_v3_v3(state_co, y_vec);
2054
2055                                 mul_v3_fl(y_vec, amplitude * (float)cos(t));
2056                                 mul_v3_fl(z_vec, amplitude / 2.f * (float)sin(2.f * t));
2057                         }
2058                         else if (inp_z > 0.0f) {
2059                                 mul_v3_v3fl(state_co, z_vec, (float)sin((float)M_PI / 3.f));
2060                                 madd_v3_v3fl(state_co, y_vec, -0.5f);
2061
2062                                 mul_v3_fl(y_vec, -amplitude * (float)cos(t + (float)M_PI / 3.f));
2063                                 mul_v3_fl(z_vec, amplitude / 2.f * (float)cos(2.f * t + (float)M_PI / 6.f));
2064                         }
2065                         else {
2066                                 mul_v3_v3fl(state_co, z_vec, -(float)sin((float)M_PI / 3.f));
2067                                 madd_v3_v3fl(state_co, y_vec, -0.5f);
2068
2069                                 mul_v3_fl(y_vec, amplitude * (float)-sin(t + (float)M_PI / 6.f));
2070                                 mul_v3_fl(z_vec, amplitude / 2.f * (float)-sin(2.f * t + (float)M_PI / 3.f));
2071                         }
2072
2073                         mul_v3_fl(state_co, amplitude);
2074                         add_v3_v3(state_co, par->co);
2075                         sub_v3_v3v3(par_vec, state->co, state_co);
2076
2077                         length = normalize_v3(par_vec);
2078                         mul_v3_fl(par_vec, MIN2(length, amplitude / 2.f));
2079
2080                         add_v3_v3v3(state_co, par->co, y_vec);
2081                         add_v3_v3(state_co, z_vec);
2082                         add_v3_v3(state_co, par_vec);
2083
2084                         shape = 2.f * (float)M_PI * (1.f + shape);
2085
2086                         if (t < shape) {
2087                                 shape = t / shape;
2088                                 shape = (float)sqrt((double)shape);
2089                                 interp_v3_v3v3(result, result, state_co, shape);
2090                         }
2091                         else {
2092                                 copy_v3_v3(result, state_co);
2093                         }
2094                         break;
2095                 }
2096         }
2097
2098         /* blend the start of the kink */
2099         if (dt < 1.f)
2100                 interp_v3_v3v3(state->co, state->co, result, dt);
2101         else
2102                 copy_v3_v3(state->co, result);
2103 }
2104
2105 static float do_clump(ParticleKey *state, ParticleKey *par, float time, float clumpfac, float clumppow, float pa_clump)
2106 {
2107         float clump = 0.f;
2108
2109         if (par && clumpfac != 0.0f) {
2110                 float cpow;
2111
2112                 if (clumppow < 0.0f)
2113                         cpow = 1.0f + clumppow;
2114                 else
2115                         cpow = 1.0f + 9.0f * clumppow;
2116
2117                 if (clumpfac < 0.0f) /* clump roots instead of tips */
2118                         clump = -clumpfac * pa_clump * (float)pow(1.0 - (double)time, (double)cpow);
2119                 else
2120                         clump = clumpfac * pa_clump * (float)pow((double)time, (double)cpow);
2121
2122                 interp_v3_v3v3(state->co, state->co, par->co, clump);
2123         }
2124
2125         return clump;
2126 }
2127 void precalc_guides(ParticleSimulationData *sim, ListBase *effectors)
2128 {
2129         EffectedPoint point;
2130         ParticleKey state;
2131         EffectorData efd;
2132         EffectorCache *eff;
2133         ParticleSystem *psys = sim->psys;
2134         EffectorWeights *weights = sim->psys->part->effector_weights;
2135         GuideEffectorData *data;
2136         PARTICLE_P;
2137
2138         if (!effectors)
2139                 return;
2140
2141         LOOP_PARTICLES {
2142                 psys_particle_on_emitter(sim->psmd, sim->psys->part->from, pa->num, pa->num_dmcache, pa->fuv, pa->foffset, state.co, 0, 0, 0, 0, 0);
2143                 
2144                 mul_m4_v3(sim->ob->obmat, state.co);
2145                 mul_mat3_m4_v3(sim->ob->obmat, state.vel);
2146                 
2147                 pd_point_from_particle(sim, pa, &state, &point);
2148
2149                 for (eff = effectors->first; eff; eff = eff->next) {
2150                         if (eff->pd->forcefield != PFIELD_GUIDE)
2151                                 continue;
2152
2153                         if (!eff->guide_data)
2154                                 eff->guide_data = MEM_callocN(sizeof(GuideEffectorData) * psys->totpart, "GuideEffectorData");
2155
2156                         data = eff->guide_data + p;
2157
2158                         sub_v3_v3v3(efd.vec_to_point, state.co, eff->guide_loc);
2159                         copy_v3_v3(efd.nor, eff->guide_dir);
2160                         efd.distance = len_v3(efd.vec_to_point);
2161
2162                         copy_v3_v3(data->vec_to_point, efd.vec_to_point);
2163                         data->strength = effector_falloff(eff, &efd, &point, weights);
2164                 }
2165         }
2166 }
2167 int do_guides(ListBase *effectors, ParticleKey *state, int index, float time)
2168 {
2169         EffectorCache *eff;
2170         PartDeflect *pd;
2171         Curve *cu;
2172         ParticleKey key, par;
2173         GuideEffectorData *data;
2174
2175         float effect[3] = {0.0f, 0.0f, 0.0f}, veffect[3] = {0.0f, 0.0f, 0.0f};
2176         float guidevec[4], guidedir[3], rot2[4], temp[3];
2177         float guidetime, radius, weight, angle, totstrength = 0.0f;
2178         float vec_to_point[3];
2179
2180         if (effectors) for (eff = effectors->first; eff; eff = eff->next) {
2181                         pd = eff->pd;
2182
2183                         if (pd->forcefield != PFIELD_GUIDE)
2184                                 continue;
2185
2186                         data = eff->guide_data + index;
2187
2188                         if (data->strength <= 0.0f)
2189                                 continue;
2190
2191                         guidetime = time / (1.0f - pd->free_end);
2192
2193                         if (guidetime > 1.0f)
2194                                 continue;
2195
2196                         cu = (Curve *)eff->ob->data;
2197
2198                         if (pd->flag & PFIELD_GUIDE_PATH_ADD) {
2199                                 if (where_on_path(eff->ob, data->strength * guidetime, guidevec, guidedir, NULL, &radius, &weight) == 0)
2200                                         return 0;
2201                         }
2202                         else {
2203                                 if (where_on_path(eff->ob, guidetime, guidevec, guidedir, NULL, &radius, &weight) == 0)
2204                                         return 0;
2205                         }
2206
2207                         mul_m4_v3(eff->ob->obmat, guidevec);
2208                         mul_mat3_m4_v3(eff->ob->obmat, guidedir);
2209
2210                         normalize_v3(guidedir);
2211
2212                         copy_v3_v3(vec_to_point, data->vec_to_point);
2213
2214                         if (guidetime != 0.0f) {
2215                                 /* curve direction */
2216                                 cross_v3_v3v3(temp, eff->guide_dir, guidedir);
2217                                 angle = dot_v3v3(eff->guide_dir, guidedir) / (len_v3(eff->guide_dir));
2218                                 angle = saacos(angle);
2219                                 axis_angle_to_quat(rot2, temp, angle);
2220                                 mul_qt_v3(rot2, vec_to_point);
2221
2222                                 /* curve tilt */
2223                                 axis_angle_to_quat(rot2, guidedir, guidevec[3] - eff->guide_loc[3]);
2224                                 mul_qt_v3(rot2, vec_to_point);
2225                         }
2226
2227                         /* curve taper */
2228                         if (cu->taperobj)
2229                                 mul_v3_fl(vec_to_point, BKE_displist_calc_taper(eff->scene, cu->taperobj, (int)(data->strength * guidetime * 100.0f), 100));
2230
2231                         else { /* curve size*/
2232                                 if (cu->flag & CU_PATH_RADIUS) {
2233                                         mul_v3_fl(vec_to_point, radius);
2234                                 }
2235                         }
2236                         par.co[0] = par.co[1] = par.co[2] = 0.0f;
2237                         copy_v3_v3(key.co, vec_to_point);
2238                         do_kink(&key, &par, 0, guidetime, pd->kink_freq, pd->kink_shape, pd->kink_amp, 0.f, pd->kink, pd->kink_axis, 0, 0);
2239                         do_clump(&key, &par, guidetime, pd->clump_fac, pd->clump_pow, 1.0f);
2240                         copy_v3_v3(vec_to_point, key.co);
2241
2242                         add_v3_v3(vec_to_point, guidevec);
2243
2244                         //sub_v3_v3v3(pa_loc, pa_loc, pa_zero);
2245                         madd_v3_v3fl(effect, vec_to_point, data->strength);
2246                         madd_v3_v3fl(veffect, guidedir, data->strength);
2247                         totstrength += data->strength;
2248
2249                         if (pd->flag & PFIELD_GUIDE_PATH_WEIGHT)
2250                                 totstrength *= weight;
2251                 }
2252
2253         if (totstrength != 0.0f) {
2254                 if (totstrength > 1.0f)
2255                         mul_v3_fl(effect, 1.0f / totstrength);
2256                 CLAMP(totstrength, 0.0f, 1.0f);
2257                 //add_v3_v3(effect, pa_zero);
2258                 interp_v3_v3v3(state->co, state->co, effect, totstrength);
2259
2260                 normalize_v3(veffect);
2261                 mul_v3_fl(veffect, len_v3(state->vel));
2262                 copy_v3_v3(state->vel, veffect);
2263                 return 1;
2264         }
2265         return 0;
2266 }
2267 static void do_rough(float *loc, float mat[4][4], float t, float fac, float size, float thres, ParticleKey *state)
2268 {
2269         float rough[3];
2270         float rco[3];
2271
2272         if (thres != 0.0f)
2273                 if ((float)fabs((float)(-1.5f + loc[0] + loc[1] + loc[2])) < 1.5f * thres) return;
2274
2275         copy_v3_v3(rco, loc);
2276         mul_v3_fl(rco, t);
2277         rough[0] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[0], rco[1], rco[2], 2, 0, 2);
2278         rough[1] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[1], rco[2], rco[0], 2, 0, 2);
2279         rough[2] = -1.0f + 2.0f * BLI_gTurbulence(size, rco[2], rco[0], rco[1], 2, 0, 2);
2280
2281         madd_v3_v3fl(state->co, mat[0], fac * rough[0]);
2282         madd_v3_v3fl(state->co, mat[1], fac * rough[1]);
2283         madd_v3_v3fl(state->co, mat[2], fac * rough[2]);
2284 }
2285 static void do_rough_end(float *loc, float mat[4][4], float t, float fac, float shape, ParticleKey *state)
2286 {
2287         float rough[2];
2288         float roughfac;
2289
2290         roughfac = fac * (float)pow((double)t, shape);
2291         copy_v2_v2(rough, loc);
2292         rough[0] = -1.0f + 2.0f * rough[0];
2293         rough[1] = -1.0f + 2.0f * rough[1];
2294         mul_v2_fl(rough, roughfac);
2295
2296         madd_v3_v3fl(state->co, mat[0], rough[0]);
2297         madd_v3_v3fl(state->co, mat[1], rough[1]);
2298 }
2299 static void do_path_effectors(ParticleSimulationData *sim, int i, ParticleCacheKey *ca, int k, int steps, float *UNUSED(rootco), float effector, float UNUSED(dfra), float UNUSED(cfra), float *length, float *vec)
2300 {
2301         float force[3] = {0.0f, 0.0f, 0.0f};
2302         ParticleKey eff_key;
2303         EffectedPoint epoint;
2304
2305         /* Don't apply effectors for dynamic hair, otherwise the effectors don't get applied twice. */
2306         if (sim->psys->flag & PSYS_HAIR_DYNAMICS)
2307                 return;
2308
2309         copy_v3_v3(eff_key.co, (ca - 1)->co);
2310         copy_v3_v3(eff_key.vel, (ca - 1)->vel);
2311         copy_qt_qt(eff_key.rot, (ca - 1)->rot);
2312
2313         pd_point_from_particle(sim, sim->psys->particles + i, &eff_key, &epoint);
2314         pdDoEffectors(sim->psys->effectors, sim->colliders, sim->psys->part->effector_weights, &epoint, force, NULL);
2315
2316         mul_v3_fl(force, effector * powf((float)k / (float)steps, 100.0f * sim->psys->part->eff_hair) / (float)steps);
2317
2318         add_v3_v3(force, vec);
2319
2320         normalize_v3(force);
2321
2322         if (k < steps)
2323                 sub_v3_v3v3(vec, (ca + 1)->co, ca->co);
2324
2325         madd_v3_v3v3fl(ca->co, (ca - 1)->co, force, *length);
2326
2327         if (k < steps)
2328                 *length = len_v3(vec);
2329 }
2330 static int check_path_length(int k, ParticleCacheKey *keys, ParticleCacheKey *state, float max_length, float *cur_length, float length, float *dvec)
2331 {
2332         if (*cur_length + length > max_length) {
2333                 mul_v3_fl(dvec, (max_length - *cur_length) / length);
2334                 add_v3_v3v3(state->co, (state - 1)->co, dvec);
2335                 keys->steps = k;
2336                 /* something over the maximum step value */
2337                 return k = 100000;
2338         }
2339         else {
2340                 *cur_length += length;
2341                 return k;
2342         }
2343 }
2344 static void offset_child(ChildParticle *cpa, ParticleKey *par, float *par_rot, ParticleKey *child, float flat, float radius)
2345 {
2346         copy_v3_v3(child->co, cpa->fuv);
2347         mul_v3_fl(child->co, radius);
2348
2349         child->co[0] *= flat;
2350
2351         copy_v3_v3(child->vel, par->vel);
2352
2353         if (par_rot) {
2354                 mul_qt_v3(par_rot, child->co);
2355                 copy_qt_qt(child->rot, par_rot);
2356         }
2357         else
2358                 unit_qt(child->rot);
2359
2360         add_v3_v3(child->co, par->co);
2361 }
2362 float *psys_cache_vgroup(DerivedMesh *dm, ParticleSystem *psys, int vgroup)
2363 {
2364         float *vg = 0;
2365
2366         if (vgroup < 0) {
2367                 /* hair dynamics pinning vgroup */
2368
2369         }
2370         else if (psys->vgroup[vgroup]) {
2371                 MDeformVert *dvert = dm->getVertDataArray(dm, CD_MDEFORMVERT);
2372                 if (dvert) {
2373                         int totvert = dm->getNumVerts(dm), i;
2374                         vg = MEM_callocN(sizeof(float) * totvert, "vg_cache");
2375                         if (psys->vg_neg & (1 << vgroup)) {
2376                                 for (i = 0; i < totvert; i++)
2377                                         vg[i] = 1.0f - defvert_find_weight(&dvert[i], psys->vgroup[vgroup] - 1);
2378                         }
2379                         else {
2380                                 for (i = 0; i < totvert; i++)
2381                                         vg[i] =  defvert_find_weight(&dvert[i], psys->vgroup[vgroup] - 1);
2382                         }
2383                 }
2384         }
2385         return vg;
2386 }
2387 void psys_find_parents(ParticleSimulationData *sim)
2388 {
2389         ParticleSettings *part = sim->psys->part;
2390         KDTree *tree;
2391         ChildParticle *cpa;
2392         int p, totparent, totchild = sim->psys->totchild;
2393         float co[3], orco[3];
2394         int from = PART_FROM_FACE;
2395         totparent = (int)(totchild * part->parents * 0.3f);
2396
2397         if ((sim->psys->renderdata || G.is_rendering) && part->child_nbr && part->ren_child_nbr)
2398                 totparent *= (float)part->child_nbr / (float)part->ren_child_nbr;
2399
2400         /* hard limit, workaround for it being ignored above */
2401         if (sim->psys->totpart < totparent) {
2402                 totparent = sim->psys->totpart;
2403         }
2404
2405         tree = BLI_kdtree_new(totparent);
2406
2407         for (p = 0, cpa = sim->psys->child; p < totparent; p++, cpa++) {
2408                 psys_particle_on_emitter(sim->psmd, from, cpa->num, DMCACHE_ISCHILD, cpa->fuv, cpa->foffset, co, 0, 0, 0, orco, 0);
2409                 BLI_kdtree_insert(tree, p, orco);
2410         }
2411
2412         BLI_kdtree_balance(tree);
2413
2414         for (; p < totchild; p++, cpa++) {
2415                 psys_particle_on_emitter(sim->psmd, from, cpa->num, DMCACHE_ISCHILD, cpa->fuv, cpa->foffset, co, 0, 0, 0, orco, 0);
2416                 cpa->parent = BLI_kdtree_find_nearest(tree, orco, NULL);
2417         }
2418
2419         BLI_kdtree_free(tree);
2420 }
2421
2422 static void get_strand_normal(Material *ma, const float surfnor[3], float surfdist, float nor[3])
2423 {
2424         float cross[3], nstrand[3], vnor[3], blend;
2425
2426         if (!((ma->mode & MA_STR_SURFDIFF) || (ma->strand_surfnor > 0.0f)))
2427                 return;
2428
2429         if (ma->mode & MA_STR_SURFDIFF) {
2430                 cross_v3_v3v3(cross, surfnor, nor);
2431                 cross_v3_v3v3(nstrand, nor, cross);
2432
2433                 blend = dot_v3v3(nstrand, surfnor);
2434                 CLAMP(blend, 0.0f, 1.0f);
2435
2436                 interp_v3_v3v3(vnor, nstrand, surfnor, blend);
2437                 normalize_v3(vnor);
2438         }
2439         else {
2440                 copy_v3_v3(vnor, nor);
2441         }
2442         
2443         if (ma->strand_surfnor > 0.0f) {
2444                 if (ma->strand_surfnor > surfdist) {
2445                         blend = (ma->strand_surfnor - surfdist) / ma->strand_surfnor;
2446                         interp_v3_v3v3(vnor, vnor, surfnor, blend);
2447                         normalize_v3(vnor);
2448                 }
2449         }
2450
2451         copy_v3_v3(nor, vnor);
2452 }
2453
2454 static int psys_threads_init_path(ParticleThread *threads, Scene *scene, float cfra, int editupdate)
2455 {
2456         ParticleThreadContext *ctx = threads[0].ctx;
2457 /*      Object *ob = ctx->sim.ob; */
2458         ParticleSystem *psys = ctx->sim.psys;
2459         ParticleSettings *part = psys->part;
2460 /*      ParticleEditSettings *pset = &scene->toolsettings->particle; */
2461         int totparent = 0, between = 0;
2462         int steps = (int)pow(2.0, (double)part->draw_step);
2463         int totchild = psys->totchild;
2464         int i, seed, totthread = threads[0].tot;
2465
2466         /*---start figuring out what is actually wanted---*/
2467         if (psys_in_edit_mode(scene, psys)) {
2468                 ParticleEditSettings *pset = &scene->toolsettings->particle;
2469
2470                 if (psys->renderdata == 0 && (psys->edit == NULL || pset->flag & PE_DRAW_PART) == 0)
2471                         totchild = 0;
2472
2473                 steps = (int)pow(2.0, (double)pset->draw_step);
2474         }
2475
2476         if (totchild && part->childtype == PART_CHILD_FACES) {
2477                 totparent = (int)(totchild * part->parents * 0.3f);
2478                 
2479                 if ((psys->renderdata || G.is_rendering) && part->child_nbr && part->ren_child_nbr)
2480                         totparent *= (float)part->child_nbr / (float)part->ren_child_nbr;
2481
2482                 /* part->parents could still be 0 so we can't test with totparent */
2483                 between = 1;
2484         }
2485
2486         if (psys->renderdata)
2487                 steps = (int)pow(2.0, (double)part->ren_step);
2488         else {
2489                 totchild = (int)((float)totchild * (float)part->disp / 100.0f);
2490                 totparent = MIN2(totparent, totchild);
2491         }
2492
2493         if (totchild == 0) return 0;
2494
2495         /* init random number generator */
2496         seed = 31415926 + ctx->sim.psys->seed;
2497         
2498         if (ctx->editupdate || totchild < 10000)
2499                 totthread = 1;
2500         
2501         for (i = 0; i < totthread; i++) {
2502                 threads[i].rng_path = BLI_rng_new(seed);
2503                 threads[i].tot = totthread;
2504         }
2505
2506         /* fill context values */
2507         ctx->between = between;
2508         ctx->steps = steps;
2509         ctx->totchild = totchild;
2510         ctx->totparent = totparent;
2511         ctx->parent_pass = 0;
2512         ctx->cfra = cfra;
2513         ctx->editupdate = editupdate;
2514
2515         psys->lattice_deform_data = psys_create_lattice_deform_data(&ctx->sim);
2516
2517         /* cache all relevant vertex groups if they exist */
2518         ctx->vg_length = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_LENGTH);
2519         ctx->vg_clump = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_CLUMP);
2520         ctx->vg_kink = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_KINK);
2521         ctx->vg_rough1 = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_ROUGH1);
2522         ctx->vg_rough2 = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_ROUGH2);
2523         ctx->vg_roughe = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_ROUGHE);
2524         if (psys->part->flag & PART_CHILD_EFFECT)
2525                 ctx->vg_effector = psys_cache_vgroup(ctx->dm, psys, PSYS_VG_EFFECTOR);
2526
2527         /* set correct ipo timing */
2528 #if 0 // XXX old animation system
2529         if (part->flag & PART_ABS_TIME && part->ipo) {
2530                 calc_ipo(part->ipo, cfra);
2531                 execute_ipo((ID *)part, part->ipo);
2532         }
2533 #endif // XXX old animation system
2534
2535         return 1;
2536 }
2537
2538 /* note: this function must be thread safe, except for branching! */
2539 static void psys_thread_create_path(ParticleThread *thread, struct ChildParticle *cpa, ParticleCacheKey *child_keys, int i)
2540 {
2541         ParticleThreadContext *ctx = thread->ctx;
2542         Object *ob = ctx->sim.ob;
2543         ParticleSystem *psys = ctx->sim.psys;
2544         ParticleSettings *part = psys->part;
2545         ParticleCacheKey **cache = psys->childcache;
2546         ParticleCacheKey **pcache = psys_in_edit_mode(ctx->sim.scene, psys) ? psys->edit->pathcache : psys->pathcache;
2547         ParticleCacheKey *child, *par = NULL, *key[4];
2548         ParticleTexture ptex;
2549         float *cpa_fuv = 0, *par_rot = 0, rot[4];
2550         float orco[3], ornor[3], hairmat[4][4], t, dvec[3], off1[4][3], off2[4][3];
2551         float length, max_length = 1.0f, cur_length = 0.0f;
2552         float eff_length, eff_vec[3], weight[4];
2553         int k, cpa_num;
2554         short cpa_from;
2555
2556         if (!pcache)
2557                 return;
2558
2559         if (ctx->between) {
2560                 ParticleData *pa = psys->particles + cpa->pa[0];
2561                 int w, needupdate;
2562                 float foffset, wsum = 0.f;
2563                 float co[3];
2564                 float p_min = part->parting_min;
2565                 float p_max = part->parting_max;
2566                 /* Virtual parents don't work nicely with parting. */
2567                 float p_fac = part->parents > 0.f ? 0.f : part->parting_fac;
2568
2569                 if (ctx->editupdate) {
2570                         needupdate = 0;
2571                         w = 0;
2572                         while (w < 4 && cpa->pa[w] >= 0) {
2573                                 if (psys->edit->points[cpa->pa[w]].flag & PEP_EDIT_RECALC) {
2574                                         needupdate = 1;
2575                                         break;
2576                                 }
2577                                 w++;
2578                         }
2579
2580                         if (!needupdate)
2581                                 return;
2582                         else
2583                                 memset(child_keys, 0, sizeof(*child_keys) * (ctx->steps + 1));
2584                 }
2585
2586                 /* get parent paths */
2587                 for (w = 0; w < 4; w++) {
2588                         if (cpa->pa[w] >= 0) {
2589                                 key[w] = pcache[cpa->pa[w]];
2590                                 weight[w] = cpa->w[w];
2591                         }
2592                         else {
2593                                 key[w] = pcache[0];
2594                                 weight[w] = 0.f;
2595                         }
2596                 }
2597
2598                 /* modify weights to create parting */
2599                 if (p_fac > 0.f) {
2600                         for (w = 0; w < 4; w++) {
2601                                 if (w && weight[w] > 0.f) {
2602                                         float d;
2603                                         if (part->flag & PART_CHILD_LONG_HAIR) {
2604                                                 /* For long hair use tip distance/root distance as parting factor instead of root to tip angle. */
2605                                                 float d1 = len_v3v3(key[0]->co, key[w]->co);
2606                                                 float d2 = len_v3v3((key[0] + key[0]->steps - 1)->co, (key[w] + key[w]->steps - 1)->co);
2607
2608                                                 d = d1 > 0.f ? d2 / d1 - 1.f : 10000.f;
2609                                         }
2610                                         else {
2611                                                 float v1[3], v2[3];
2612                                                 sub_v3_v3v3(v1, (key[0] + key[0]->steps - 1)->co, key[0]->co);
2613                                                 sub_v3_v3v3(v2, (key[w] + key[w]->steps - 1)->co, key[w]->co);
2614                                                 normalize_v3(v1);
2615                                                 normalize_v3(v2);
2616
2617                                                 d = RAD2DEGF(saacos(dot_v3v3(v1, v2)));
2618                                         }
2619
2620                                         if (p_max > p_min)
2621                                                 d = (d - p_min) / (p_max - p_min);
2622                                         else
2623                                                 d = (d - p_min) <= 0.f ? 0.f : 1.f;
2624
2625                                         CLAMP(d, 0.f, 1.f);
2626
2627                                         if (d > 0.f)
2628                                                 weight[w] *= (1.f - d);
2629                                 }
2630                                 wsum += weight[w];
2631                         }
2632                         for (w = 0; w < 4; w++)
2633                                 weight[w] /= wsum;
2634
2635                         interp_v4_v4v4(weight, cpa->w, weight, p_fac);
2636                 }
2637
2638                 /* get the original coordinates (orco) for texture usage */
2639                 cpa_num = cpa->num;
2640                 
2641                 foffset = cpa->foffset;
2642                 cpa_fuv = cpa->fuv;
2643                 cpa_from = PART_FROM_FACE;
2644
2645                 psys_particle_on_emitter(ctx->sim.psmd, cpa_from, cpa_num, DMCACHE_ISCHILD, cpa->fuv, foffset, co, ornor, 0, 0, orco, 0);
2646
2647                 mul_m4_v3(ob->obmat, co);
2648
2649                 for (w = 0; w < 4; w++)
2650                         sub_v3_v3v3(off1[w], co, key[w]->co);
2651
2652                 psys_mat_hair_to_global(ob, ctx->sim.psmd->dm, psys->part->from, pa, hairmat);
2653         }
2654         else {
2655                 ParticleData *pa = psys->particles + cpa->parent;
2656                 float co[3];
2657                 if (ctx->editupdate) {
2658                         if (!(psys->edit->points[cpa->parent].flag & PEP_EDIT_RECALC))
2659                                 return;
2660
2661                         memset(child_keys, 0, sizeof(*child_keys) * (ctx->steps + 1));
2662                 }
2663
2664                 /* get the parent path */
2665                 key[0] = pcache[cpa->parent];
2666
2667                 /* get the original coordinates (orco) for texture usage */
2668                 cpa_from = part->from;
2669                 cpa_num = pa->num;
2670                 cpa_fuv = pa->fuv;
2671
2672                 psys_particle_on_emitter(ctx->sim.psmd, cpa_from, cpa_num, DMCACHE_ISCHILD, cpa_fuv, pa->foffset, co, ornor, 0, 0, orco, 0);
2673
2674                 psys_mat_hair_to_global(ob, ctx->sim.psmd->dm, psys->part->from, pa, hairmat);
2675         }
2676
2677         child_keys->steps = ctx->steps;
2678
2679         /* get different child parameters from textures & vgroups */
2680         get_child_modifier_parameters(part, ctx, cpa, cpa_from, cpa_num, cpa_fuv, orco, &ptex);
2681
2682         if (ptex.exist < psys_frand(psys, i + 24)) {
2683                 child_keys->steps = -1;
2684                 return;
2685         }
2686
2687         /* create the child path */
2688         for (k = 0, child = child_keys; k <= ctx->steps; k++, child++) {
2689                 if (ctx->between) {
2690                         int w = 0;
2691
2692                         zero_v3(child->co);
2693                         zero_v3(child->vel);
2694                         unit_qt(child->rot);
2695
2696                         for (w = 0; w < 4; w++) {
2697                                 copy_v3_v3(off2[w], off1[w]);
2698
2699                                 if (part->flag & PART_CHILD_LONG_HAIR) {
2700                                         /* Use parent rotation (in addition to emission location) to determine child offset. */
2701                                         if (k)
2702                                                 mul_qt_v3((key[w] + k)->rot, off2[w]);
2703
2704                                         /* Fade the effect of rotation for even lengths in the end */
2705                                         project_v3_v3v3(dvec, off2[w], (key[w] + k)->vel);
2706                                         madd_v3_v3fl(off2[w], dvec, -(float)k / (float)ctx->steps);
2707                                 }
2708
2709                                 add_v3_v3(off2[w], (key[w] + k)->co);
2710                         }
2711
2712                         /* child position is the weighted sum of parent positions */
2713                         interp_v3_v3v3v3v3(child->co, off2[0], off2[1], off2[2], off2[3], weight);
2714                         interp_v3_v3v3v3v3(child->vel, (key[0] + k)->vel, (key[1] + k)->vel, (key[2] + k)->vel, (key[3] + k)->vel, weight);
2715
2716                         copy_qt_qt(child->rot, (key[0] + k)->rot);
2717                 }
2718                 else {
2719                         if (k) {
2720                                 mul_qt_qtqt(rot, (key[0] + k)->rot, key[0]->rot);
2721                                 par_rot = rot;
2722                         }
2723                         else {
2724                                 par_rot = key[0]->rot;
2725                         }
2726                         /* offset the child from the parent position */
2727                         offset_child(cpa, (ParticleKey *)(key[0] + k), par_rot, (ParticleKey *)child, part->childflat, part->childrad);
2728                 }
2729
2730                 child->time = (float)k / (float)ctx->steps;
2731         }
2732
2733         /* apply effectors */
2734         if (part->flag & PART_CHILD_EFFECT) {
2735                 for (k = 0, child = child_keys; k <= ctx->steps; k++, child++) {
2736                         if (k) {
2737                                 do_path_effectors(&ctx->sim, cpa->pa[0], child, k, ctx->steps, child_keys->co, ptex.effector, 0.0f, ctx->cfra, &eff_length, eff_vec);
2738                         }
2739                         else {
2740                                 sub_v3_v3v3(eff_vec, (child + 1)->co, child->co);
2741                                 eff_length = len_v3(eff_vec);
2742                         }
2743                 }
2744         }
2745
2746         for (k = 0, child = child_keys; k <= ctx->steps; k++, child++) {
2747                 t = (float)k / (float)ctx->steps;
2748
2749                 if (ctx->totparent)
2750                         /* this is now threadsafe, virtual parents are calculated before rest of children */
2751                         par = (i >= ctx->totparent) ? cache[cpa->parent] : NULL;
2752                 else if (cpa->parent >= 0)
2753                         par = pcache[cpa->parent];
2754
2755                 if (par) {
2756                         if (k) {
2757                                 mul_qt_qtqt(rot, (par + k)->rot, par->rot);
2758                                 par_rot = rot;
2759                         }
2760                         else {
2761                                 par_rot = par->rot;
2762                         }
2763                         par += k;
2764                 }
2765
2766                 /* apply different deformations to the child path */
2767                 do_child_modifiers(&ctx->sim, &ptex, (ParticleKey *)par, par_rot, cpa, orco, hairmat, (ParticleKey *)child, t);
2768
2769                 /* we have to correct velocity because of kink & clump */
2770                 if (k > 1) {
2771                         sub_v3_v3v3((child - 1)->vel, child->co, (child - 2)->co);
2772                         mul_v3_fl((child - 1)->vel, 0.5);
2773
2774                         if (ctx->ma && (part->draw_col == PART_DRAW_COL_MAT))
2775                                 get_strand_normal(ctx->ma, ornor, cur_length, (child - 1)->vel);
2776                 }
2777
2778                 if (k == ctx->steps)
2779                         sub_v3_v3v3(child->vel, child->co, (child - 1)->co);
2780
2781                 /* check if path needs to be cut before actual end of data points */
2782                 if (k) {
2783                         sub_v3_v3v3(dvec, child->co, (child - 1)->co);
2784                         length = 1.0f / (float)ctx->steps;
2785                         k = check_path_length(k, child_keys, child, max_length, &cur_length, length, dvec);
2786                 }
2787                 else {
2788                         /* initialize length calculation */
2789                         max_length = ptex.length;
2790                         cur_length = 0.0f;
2791                 }
2792
2793                 if (ctx->ma && (part->draw_col == PART_DRAW_COL_MAT)) {
2794                         copy_v3_v3(child->col, &ctx->ma->r);
2795                         get_strand_normal(ctx->ma, ornor, cur_length, child->vel);
2796                 }
2797         }
2798
2799         /* Hide virtual parents */
2800         if (i < ctx->totparent)
2801                 child_keys->steps = -1;
2802 }
2803
2804 static void *exec_child_path_cache(void *data)
2805 {
2806         ParticleThread *thread = (ParticleThread *)data;
2807         ParticleThreadContext *ctx = thread->ctx;
2808         ParticleSystem *psys = ctx->sim.psys;
2809         ParticleCacheKey **cache = psys->childcache;
2810         ChildParticle *cpa;
2811         int i, totchild = ctx->totchild, first = 0;
2812
2813         if (thread->tot > 1) {
2814                 first = ctx->parent_pass ? 0 : ctx->totparent;
2815                 totchild = ctx->parent_pass ? ctx->totparent : ctx->totchild;
2816         }
2817         
2818         cpa = psys->child + first + thread->num;
2819         for (i = first + thread->num; i < totchild; i += thread->tot, cpa += thread->tot)
2820                 psys_thread_create_path(thread, cpa, cache[i], i);
2821
2822         return 0;
2823 }
2824
2825 void psys_cache_child_paths(ParticleSimulationData *sim, float cfra, int editupdate)
2826 {
2827         ParticleThread *pthreads;
2828         ParticleThreadContext *ctx;
2829         ListBase threads;
2830         int i, totchild, totparent, totthread;
2831
2832         if (sim->psys->flag & PSYS_GLOBAL_HAIR)
2833                 return;
2834
2835         pthreads = psys_threads_create(sim);
2836
2837         if (!psys_threads_init_path(pthreads, sim->scene, cfra, editupdate)) {
2838                 psys_threads_free(pthreads);
2839                 return;
2840         }
2841
2842         ctx = pthreads[0].ctx;
2843         totchild = ctx->totchild;
2844         totparent = ctx->totparent;
2845
2846         if (editupdate && sim->psys->childcache && totchild == sim->psys->totchildcache) {
2847                 ; /* just overwrite the existing cache */
2848         }
2849         else {
2850                 /* clear out old and create new empty path cache */
2851                 free_child_path_cache(sim->psys);
2852                 sim->psys->childcache = psys_alloc_path_cache_buffers(&sim->psys->childcachebufs, totchild, ctx->steps + 1);
2853                 sim->psys->totchildcache = totchild;
2854         }
2855
2856         totthread = pthreads[0].tot;
2857
2858         if (totthread > 1) {
2859
2860                 /* make virtual child parents thread safe by calculating them first */
2861                 if (totparent) {
2862                         BLI_init_threads(&threads, exec_child_path_cache, totthread);
2863                         
2864                         for (i = 0; i < totthread; i++) {
2865                                 pthreads[i].ctx->parent_pass = 1;
2866                                 BLI_insert_thread(&threads, &pthreads[i]);
2867                         }
2868
2869                         BLI_end_threads(&threads);
2870
2871                         for (i = 0; i < totthread; i++)
2872                                 pthreads[i].ctx->parent_pass = 0;
2873                 }
2874
2875                 BLI_init_threads(&threads, exec_child_path_cache, totthread);
2876
2877                 for (i = 0; i < totthread; i++)
2878                         BLI_insert_thread(&threads, &pthreads[i]);
2879
2880                 BLI_end_threads(&threads);
2881         }
2882         else
2883                 exec_child_path_cache(&pthreads[0]);
2884
2885         psys_threads_free(pthreads);
2886 }
2887 /* figure out incremental rotations along path starting from unit quat */
2888 static void cache_key_incremental_rotation(ParticleCacheKey *key0, ParticleCacheKey *key1, ParticleCacheKey *key2, float *prev_tangent, int i)
2889 {
2890         float cosangle, angle, tangent[3], normal[3], q[4];
2891
2892         switch (i) {
2893                 case 0:
2894                         /* start from second key */
2895                         break;
2896                 case 1:
2897                         /* calculate initial tangent for incremental rotations */
2898                         sub_v3_v3v3(prev_tangent, key0->co, key1->co);
2899                         normalize_v3(prev_tangent);
2900                         unit_qt(key1->rot);
2901                         break;
2902                 default:
2903                         sub_v3_v3v3(tangent, key0->co, key1->co);
2904                         normalize_v3(tangent);
2905
2906                         cosangle = dot_v3v3(tangent, prev_tangent);
2907
2908                         /* note we do the comparison on cosangle instead of
2909                          * angle, since floating point accuracy makes it give
2910                          * different results across platforms */
2911                         if (cosangle > 0.999999f) {
2912                                 copy_v4_v4(key1->rot, key2->rot);
2913                         }
2914                         else {
2915                                 angle = saacos(cosangle);
2916                                 cross_v3_v3v3(normal, prev_tangent, tangent);
2917                                 axis_angle_to_quat(q, normal, angle);
2918                                 mul_qt_qtqt(key1->rot, q, key2->rot);
2919                         }
2920
2921                         copy_v3_v3(prev_tangent, tangent);
2922         }
2923 }
2924
2925 /**
2926  * Calculates paths ready for drawing/rendering
2927  * - Useful for making use of opengl vertex arrays for super fast strand drawing.
2928  * - Makes child strands possible and creates them too into the cache.
2929  * - Cached path data is also used to determine cut position for the editmode tool. */
2930 void psys_cache_paths(ParticleSimulationData *sim, float cfra)
2931 {
2932         PARTICLE_PSMD;
2933         ParticleEditSettings *pset = &sim->scene->toolsettings->particle;
2934         ParticleSystem *psys = sim->psys;
2935         ParticleSettings *part = psys->part;
2936         ParticleCacheKey *ca, **cache;
2937
2938         DerivedMesh *hair_dm = (psys->part->type == PART_HAIR && psys->flag & PSYS_HAIR_DYNAMICS) ? psys->hair_out_dm : NULL;
2939         
2940         ParticleKey result;
2941         
2942         Material *ma;
2943         ParticleInterpolationData pind;
2944         ParticleTexture ptex;
2945
2946         PARTICLE_P;
2947         
2948         float birthtime = 0.0, dietime = 0.0;
2949         float t, time = 0.0, dfra = 1.0 /* , frs_sec = sim->scene->r.frs_sec*/ /*UNUSED*/;
2950         float col[4] = {0.5f, 0.5f, 0.5f, 1.0f};
2951         float prev_tangent[3] = {0.0f, 0.0f, 0.0f}, hairmat[4][4];
2952         float rotmat[3][3];
2953         int k;
2954         int steps = (int)pow(2.0, (double)(psys->renderdata ? part->ren_step : part->draw_step));
2955         int totpart = psys->totpart;
2956         float length, vec[3];
2957         float *vg_effector = NULL;
2958         float *vg_length = NULL, pa_length = 1.0f;
2959         int keyed, baked;
2960
2961         /* we don't have anything valid to create paths from so let's quit here */
2962         if ((psys->flag & PSYS_HAIR_DONE || psys->flag & PSYS_KEYED || psys->pointcache) == 0)
2963                 return;
2964
2965         if (psys_in_edit_mode(sim->scene, psys))
2966                 if (psys->renderdata == 0 && (psys->edit == NULL || pset->flag & PE_DRAW_PART) == 0)
2967                         return;
2968
2969         keyed = psys->flag & PSYS_KEYED;
2970         baked = psys->pointcache->mem_cache.first && psys->part->type != PART_HAIR;
2971
2972         /* clear out old and create new empty path cache */
2973         psys_free_path_cache(psys, psys->edit);
2974         cache = psys->pathcache = psys_alloc_path_cache_buffers(&psys->pathcachebufs, totpart, steps + 1);
2975
2976         psys->lattice_deform_data = psys_create_lattice_deform_data(sim);
2977         ma = give_current_material(sim->ob, psys->part->omat);
2978         if (ma && (psys->part->draw_col == PART_DRAW_COL_MAT))
2979                 copy_v3_v3(col, &ma->r);
2980
2981         if ((psys->flag & PSYS_GLOBAL_HAIR) == 0) {
2982                 if ((psys->part->flag & PART_CHILD_EFFECT) == 0)
2983                         vg_effector = psys_cache_vgroup(psmd->dm, psys, PSYS_VG_EFFECTOR);
2984                 
2985                 if (!psys->totchild)
2986                         vg_length = psys_cache_vgroup(psmd->dm, psys, PSYS_VG_LENGTH);
2987         }
2988
2989         /* ensure we have tessfaces to be used for mapping */
2990         if (part->from != PART_FROM_VERT) {
2991                 DM_ensure_tessface(psmd->dm);
2992         }
2993
2994         /*---first main loop: create all actual particles' paths---*/
2995         LOOP_SHOWN_PARTICLES {
2996                 if (!psys->totchild) {
2997                         psys_get_texture(sim, pa, &ptex, PAMAP_LENGTH, 0.f);
2998                         pa_length = ptex.length * (1.0f - part->randlength * psys_frand(psys, psys->seed + p));
2999                         if (vg_length)
3000                                 pa_length *= psys_particle_value_from_verts(psmd->dm, part->from, pa, vg_length);
3001                 }
3002
3003                 pind.keyed = keyed;
3004                 pind.cache = baked ? psys->pointcache : NULL;
3005                 pind.epoint = NULL;
3006                 pind.bspline = (psys->part->flag & PART_HAIR_BSPLINE);
3007                 pind.dm = hair_dm;
3008
3009                 memset(cache[p], 0, sizeof(*cache[p]) * (steps + 1));
3010
3011                 cache[p]->steps = steps;
3012
3013                 /*--get the first data points--*/
3014                 init_particle_interpolation(sim->ob, sim->psys, pa, &pind);
3015
3016                 /* hairmat is needed for for non-hair particle too so we get proper rotations */
3017                 psys_mat_hair_to_global(sim->ob, psmd->dm, psys->part->from, pa, hairmat);
3018                 copy_v3_v3(rotmat[0], hairmat[2]);
3019                 copy_v3_v3(rotmat[1], hairmat[1]);
3020                 copy_v3_v3(rotmat[2], hairmat[0]);
3021
3022                 if (part->draw & PART_ABS_PATH_TIME) {
3023                         birthtime = MAX2(pind.birthtime, part->path_start);
3024                         dietime = MIN2(pind.dietime, part->path_end);
3025                 }
3026                 else {
3027                         float tb = pind.birthtime;
3028                         birthtime = tb + part->path_start * (pind.dietime - tb);
3029                         dietime = tb + part->path_end * (pind.dietime - tb);
3030                 }
3031
3032                 if (birthtime >= dietime) {
3033                         cache[p]->steps = -1;
3034                         continue;
3035                 }
3036
3037                 dietime = birthtime + pa_length * (dietime - birthtime);
3038
3039                 /*--interpolate actual path from data points--*/
3040                 for (k = 0, ca = cache[p]; k <= steps; k++, ca++) {
3041                         time = (float)k / (float)steps;
3042                         t = birthtime + time * (dietime - birthtime);
3043                         result.time = -t;
3044                         do_particle_interpolation(psys, p, pa, t, &pind, &result);
3045                         copy_v3_v3(ca->co, result.co);
3046
3047                         /* dynamic hair is in object space */
3048                         /* keyed and baked are already in global space */
3049                         if (hair_dm)
3050                                 mul_m4_v3(sim->ob->obmat, ca->co);
3051                         else if (!keyed && !baked && !(psys->flag & PSYS_GLOBAL_HAIR))
3052                                 mul_m4_v3(hairmat, ca->co);
3053
3054                         copy_v3_v3(ca->col, col);
3055                 }
3056                 
3057                 /*--modify paths and calculate rotation & velocity--*/
3058
3059                 if (!(psys->flag & PSYS_GLOBAL_HAIR)) {
3060                         /* apply effectors */
3061                         if ((psys->part->flag & PART_CHILD_EFFECT) == 0) {
3062                                 float effector = 1.0f;
3063                                 if (vg_effector)
3064                                         effector *= psys_particle_value_from_verts(psmd->dm, psys->part->from, pa, vg_effector);
3065
3066                                 sub_v3_v3v3(vec, (cache[p] + 1)->co, cache[p]->co);
3067                                 length = len_v3(vec);
3068
3069                                 for (k = 1, ca = cache[p] + 1; k <= steps; k++, ca++)
3070                                         do_path_effectors(sim, p, ca, k, steps, cache[p]->co, effector, dfra, cfra, &length, vec);
3071                         }
3072
3073                         /* apply guide curves to path data */
3074                         if (sim->psys->effectors && (psys->part->flag & PART_CHILD_EFFECT) == 0) {