Orange: tuesday sync with bf-blender
[blender.git] / source / blender / blenkernel / intern / effect.c
1 /*  effect.c
2  * 
3  * 
4  * $Id$
5  *
6  * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
7  *
8  * This program is free software; you can redistribute it and/or
9  * modify it under the terms of the GNU General Public License
10  * as published by the Free Software Foundation; either version 2
11  * of the License, or (at your option) any later version. The Blender
12  * Foundation also sells licenses for use in proprietary software under
13  * the Blender License.  See http://www.blender.org/BL/ for information
14  * about this.
15  *
16  * This program is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
19  * GNU General Public License for more details.
20  *
21  * You should have received a copy of the GNU General Public License
22  * along with this program; if not, write to the Free Software Foundation,
23  * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
24  *
25  * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
26  * All rights reserved.
27  *
28  * The Original Code is: all of this file.
29  *
30  * Contributor(s): none yet.
31  *
32  * ***** END GPL/BL DUAL LICENSE BLOCK *****
33  */
34
35 #include <math.h>
36 #include <stdlib.h>
37
38 #include "MEM_guardedalloc.h"
39
40 #include "DNA_curve_types.h"
41 #include "DNA_effect_types.h"
42 #include "DNA_group_types.h"
43 #include "DNA_ipo_types.h"
44 #include "DNA_key_types.h"
45 #include "DNA_lattice_types.h"
46 #include "DNA_listBase.h"
47 #include "DNA_mesh_types.h"
48 #include "DNA_meshdata_types.h"
49 #include "DNA_material_types.h"
50 #include "DNA_object_types.h"
51 #include "DNA_object_force.h"
52 #include "DNA_texture_types.h"
53 #include "DNA_scene_types.h"
54
55 #include "BLI_blenlib.h"
56 #include "BLI_arithb.h"
57 #include "BLI_rand.h"
58
59 #include "BKE_action.h"
60 #include "BKE_anim.h"           /* needed for where_on_path */
61 #include "BKE_armature.h"
62 #include "BKE_bad_level_calls.h"
63 #include "BKE_blender.h"
64 #include "BKE_constraint.h"
65 #include "BKE_deform.h"
66 #include "BKE_displist.h"
67 #include "BKE_DerivedMesh.h"
68 #include "BKE_effect.h"
69 #include "BKE_global.h"
70 #include "BKE_ipo.h"
71 #include "BKE_key.h"
72 #include "BKE_lattice.h"
73 #include "BKE_mesh.h"
74 #include "BKE_material.h"
75 #include "BKE_main.h"
76 #include "BKE_object.h"
77 #include "BKE_screen.h"
78 #include "BKE_utildefines.h"
79
80 #include "render.h"             // externtex, bad level call (ton)
81 #include "PIL_time.h"
82
83 /* temporal struct, used for reading return of mesh_get_mapped_verts_nors() */
84 typedef struct VeNoCo {
85         float co[3], no[3];
86 } VeNoCo;
87
88 Effect *add_effect(int type)
89 {
90         Effect *eff=0;
91         PartEff *paf;
92         int a;
93         
94         switch(type) {
95         case EFF_PARTICLE:
96                 paf= MEM_callocN(sizeof(PartEff), "neweff");
97                 eff= (Effect *)paf;
98                 
99                 paf->sta= 1.0;
100                 paf->end= 100.0;
101                 paf->lifetime= 50.0;
102                 for(a=0; a<PAF_MAXMULT; a++) {
103                         paf->life[a]= 50.0;
104                         paf->child[a]= 4;
105                         paf->mat[a]= 1;
106                 }
107                 
108                 paf->totpart= 1000;
109                 paf->totkey= 8;
110                 paf->staticstep= 5;
111                 paf->defvec[2]= 1.0f;
112                 paf->nabla= 0.05f;
113                 paf->disp = 100;
114                 paf->speedtex = 8;
115                 paf->omat = 1;
116                 paf->flag= PAF_FACE;
117
118                 break;
119         }
120         
121         eff->type= eff->buttype= type;
122         eff->flag |= SELECT;
123         
124         return eff;
125 }
126
127 void free_effect(Effect *eff)
128 {
129         PartEff *paf;
130         
131         if(eff->type==EFF_PARTICLE) {
132                 paf= (PartEff *)eff;
133                 if(paf->keys) MEM_freeN(paf->keys);
134         }
135         MEM_freeN(eff); 
136 }
137
138
139 void free_effects(ListBase *lb)
140 {
141         Effect *eff;
142         
143         eff= lb->first;
144         while(eff) {
145                 BLI_remlink(lb, eff);
146                 free_effect(eff);
147                 eff= lb->first;
148         }
149 }
150
151 Effect *copy_effect(Effect *eff) 
152 {
153         Effect *effn;
154
155         effn= MEM_dupallocN(eff);
156         if(effn->type==EFF_PARTICLE) ((PartEff *)effn)->keys= 0;
157
158         return effn;    
159 }
160
161 void copy_act_effect(Object *ob)
162 {
163         /* return a copy of the active effect */
164         Effect *effn, *eff;
165         
166         eff= ob->effect.first;
167         while(eff) {
168                 if(eff->flag & SELECT) {
169                         
170                         effn= copy_effect(eff);
171                         BLI_addtail(&ob->effect, effn);
172                         
173                         eff->flag &= ~SELECT;
174                         return;
175                         
176                 }
177                 eff= eff->next;
178         }
179         
180         /* when it comes here: add new effect */
181         eff= add_effect(EFF_PARTICLE);
182         BLI_addtail(&ob->effect, eff);
183                         
184 }
185
186 void copy_effects(ListBase *lbn, ListBase *lb)
187 {
188         Effect *eff, *effn;
189
190         lbn->first= lbn->last= 0;
191
192         eff= lb->first;
193         while(eff) {
194                 effn= copy_effect(eff);
195                 BLI_addtail(lbn, effn);
196                 
197                 eff= eff->next;
198         }
199         
200 }
201
202 void deselectall_eff(Object *ob)
203 {
204         Effect *eff= ob->effect.first;
205         
206         while(eff) {
207                 eff->flag &= ~SELECT;
208                 eff= eff->next;
209         }
210 }
211
212 /* ***************** PARTICLES ***************** */
213
214 static Particle *new_particle(PartEff *paf)
215 {
216         static Particle *pa;
217         static int cur;
218
219         /* we agree: when paf->keys==0: alloc */        
220         if(paf->keys==NULL) {
221                 pa= paf->keys= MEM_callocN( paf->totkey*paf->totpart*sizeof(Particle), "particlekeys" );
222                 cur= 0;
223         }
224         else {
225                 if(cur && cur<paf->totpart) pa+=paf->totkey;
226                 cur++;
227         }
228         return pa;
229 }
230
231 PartEff *give_parteff(Object *ob)
232 {
233         PartEff *paf;
234         
235         paf= ob->effect.first;
236         while(paf) {
237                 if(paf->type==EFF_PARTICLE) return paf;
238                 paf= paf->next;
239         }
240         return 0;
241 }
242
243 void where_is_particle(PartEff *paf, Particle *pa, float ctime, float *vec)
244 {
245         Particle *p[4];
246         float dt, t[4];
247         int a;
248         
249         if(paf->totkey==1 || ctime < pa->time) {
250                 VECCOPY(vec, pa->co);
251                 return;
252         }
253         
254         /* first find the first particlekey */
255         a= (int)((paf->totkey-1)*(ctime-pa->time)/pa->lifetime);
256         if(a>=paf->totkey) a= paf->totkey-1;
257         else if(a<0) a= 0;
258         
259         pa+= a;
260         
261         if(a>0) p[0]= pa-1; else p[0]= pa;
262         p[1]= pa;
263         
264         if(a+1<paf->totkey) p[2]= pa+1; else p[2]= pa;
265         if(a+2<paf->totkey) p[3]= pa+2; else p[3]= p[2];
266         
267         if(p[1]==p[2] || p[2]->time == p[1]->time) dt= 0.0;
268         else dt= (ctime-p[1]->time)/(p[2]->time - p[1]->time);
269
270         if(paf->flag & PAF_BSPLINE) set_four_ipo(dt, t, KEY_BSPLINE);
271         else set_four_ipo(dt, t, KEY_CARDINAL);
272
273         vec[0]= t[0]*p[0]->co[0] + t[1]*p[1]->co[0] + t[2]*p[2]->co[0] + t[3]*p[3]->co[0];
274         vec[1]= t[0]*p[0]->co[1] + t[1]*p[1]->co[1] + t[2]*p[2]->co[1] + t[3]*p[3]->co[1];
275         vec[2]= t[0]*p[0]->co[2] + t[1]*p[1]->co[2] + t[2]*p[2]->co[2] + t[3]*p[3]->co[2];
276
277 }
278
279 static void particle_tex(MTex *mtex, PartEff *paf, float *co, float *no)
280 {                               
281         float tin, tr, tg, tb, ta;
282         float old;
283         
284         externtex(mtex, co, &tin, &tr, &tg, &tb, &ta);
285
286         if(paf->texmap==PAF_TEXINT) {
287                 tin*= paf->texfac;
288                 no[0]+= tin*paf->defvec[0];
289                 no[1]+= tin*paf->defvec[1];
290                 no[2]+= tin*paf->defvec[2];
291         }
292         else if(paf->texmap==PAF_TEXRGB) {
293                 no[0]+= (tr-0.5f)*paf->texfac;
294                 no[1]+= (tg-0.5f)*paf->texfac;
295                 no[2]+= (tb-0.5f)*paf->texfac;
296         }
297         else {  /* PAF_TEXGRAD */
298                 
299                 old= tin;
300                 co[0]+= paf->nabla;
301                 externtex(mtex, co, &tin, &tr, &tg, &tb, &ta);
302                 no[0]+= (old-tin)*paf->texfac;
303                 
304                 co[0]-= paf->nabla;
305                 co[1]+= paf->nabla;
306                 externtex(mtex, co, &tin, &tr, &tg, &tb, &ta);
307                 no[1]+= (old-tin)*paf->texfac;
308                 
309                 co[1]-= paf->nabla;
310                 co[2]+= paf->nabla;
311                 externtex(mtex, co, &tin, &tr, &tg, &tb, &ta);
312                 no[2]+= (old-tin)*paf->texfac;
313                 
314         }
315 }
316
317 /* -------------------------- Effectors ------------------ */
318
319 typedef struct pEffectorCache {
320         struct pEffectorCache *next, *prev;
321         struct Object *ob;
322         
323         /* precalculated variables */
324         float oldloc[3], oldspeed[3];
325         float scale, time_scale;
326         float guide_dist;
327         
328         Object obcopy;  /* for restoring transformation data */
329 } pEffectorCache;
330
331 static void add_to_effectorcache(ListBase *lb, Object *ob, Object *obsrc)
332 {
333         pEffectorCache *ec;
334         PartDeflect *pd= ob->pd;
335                         
336         if(pd->forcefield == PFIELD_GUIDE) {
337                 if(ob->type==OB_CURVE && obsrc->type==OB_MESH) {        /* guides only do mesh particles */
338                         Curve *cu= ob->data;
339                         if(cu->flag & CU_PATH) {
340                                 if(cu->path==NULL || cu->path->data==NULL)
341                                         makeDispListCurveTypes(ob, 0);
342                                 if(cu->path && cu->path->data) {
343                                         ec= MEM_callocN(sizeof(pEffectorCache), "effector cache");
344                                         ec->ob= ob;
345                                         BLI_addtail(lb, ec);
346                                 }
347                         }
348                 }
349         }
350         else if(pd->forcefield) {
351                 ec= MEM_callocN(sizeof(pEffectorCache), "effector cache");
352                 ec->ob= ob;
353                 BLI_addtail(lb, ec);
354         }
355 }
356
357 /* returns ListBase handle with objects taking part in the effecting */
358 ListBase *pdInitEffectors(Object *obsrc, Group *group)
359 {
360         static ListBase listb={NULL, NULL};
361         pEffectorCache *ec;
362         Base *base;
363         unsigned int layer= obsrc->lay;
364         
365         if(group) {
366                 GroupObject *go;
367                 
368                 for(go= group->gobject.first; go; go= go->next) {
369                         if( (go->ob->lay & layer) && go->ob->pd && go->ob!=obsrc) {
370                                 add_to_effectorcache(&listb, go->ob, obsrc);
371                         }
372                 }
373         }
374         else {
375                 for(base = G.scene->base.first; base; base= base->next) {
376                         if( (base->lay & layer) && base->object->pd && base->object!=obsrc) {
377                                 add_to_effectorcache(&listb, base->object, obsrc);
378                         }
379                 }
380         }
381         
382         /* make a full copy */
383         for(ec= listb.first; ec; ec= ec->next) {
384                 ec->obcopy= *(ec->ob);
385         }
386
387         if(listb.first)
388                 return &listb;
389         
390         return NULL;
391 }
392
393 void pdEndEffectors(ListBase *lb)
394 {
395         if(lb) {
396                 pEffectorCache *ec;
397                 /* restore full copy */
398                 for(ec= lb->first; ec; ec= ec->next)
399                         *(ec->ob)= ec->obcopy;
400
401                 BLI_freelistN(lb);
402         }
403 }
404
405 /* local for this c file, only for guides now */
406 static void precalc_effectors(Object *ob, PartEff *paf, Particle *pa, ListBase *lb)
407 {
408         pEffectorCache *ec;
409         
410         for(ec= lb->first; ec; ec= ec->next) {
411                 PartDeflect *pd= ec->ob->pd;
412                 
413                 ec->oldspeed[0]= ec->oldspeed[1]= ec->oldspeed[2]= 0.0f;
414                 
415                 if(pd->forcefield==PFIELD_GUIDE && ec->ob->type==OB_CURVE) {
416                         float vec[4], dir[3];
417                                 
418                         /* scale corrects speed vector to curve size */
419                         if(paf->totkey>1) ec->scale= (paf->totkey-1)/pa->lifetime;
420                         else ec->scale= 1.0f;
421                         
422                         /* time_scale is for random life */
423                         if(pa->lifetime>paf->lifetime)
424                                 ec->time_scale= paf->lifetime/pa->lifetime;
425                         else
426                                 ec->time_scale= pa->lifetime/paf->lifetime;
427
428                         /* distance of first path point to particle origin */
429                         where_on_path(ec->ob, 0.0f, vec, dir);
430                         VECCOPY(ec->oldloc, vec);       /* store local coord for differences */
431                         Mat4MulVecfl(ec->ob->obmat, vec);
432                         
433                         /* for static we need to move to global space */
434                         if(paf->flag & PAF_STATIC) {
435                                 VECCOPY(dir, pa->co);
436                                 Mat4MulVecfl(ob->obmat, dir);
437                                 ec->guide_dist= VecLenf(vec, dir);
438                         }
439                         else 
440                                 ec->guide_dist= VecLenf(vec, pa->co);
441                 }
442         }
443 }
444
445
446 /*  -------- pdDoEffectors() --------
447     generic force/speed system, now used for particles and softbodies
448         lb                      = listbase with objects that take part in effecting
449         opco            = global coord, as input
450     force               = force accumulator
451     speed               = actual current speed which can be altered
452         cur_time        = "external" time in frames, is constant for static particles
453         loc_time        = "local" time in frames, range <0-1> for the lifetime of particle
454     par_layer   = layer the caller is in
455         flags           = only used for softbody wind now
456         guide           = old speed of particle
457
458 */
459 void pdDoEffectors(ListBase *lb, float *opco, float *force, float *speed, float cur_time, float loc_time, unsigned int flags)
460 {
461 /*
462         Modifies the force on a particle according to its
463         relation with the effector object
464         Different kind of effectors include:
465                 Forcefields: Gravity-like attractor
466                 (force power is related to the inverse of distance to the power of a falloff value)
467                 Vortex fields: swirling effectors
468                 (particles rotate around Z-axis of the object. otherwise, same relation as)
469                 (Forcefields, but this is not done through a force/acceleration)
470                 Guide: particles on a path
471                 (particles are guided along a curve bezier or old nurbs)
472                 (is independent of other effectors)
473 */
474         Object *ob;
475         pEffectorCache *ec;
476         PartDeflect *pd;
477         float vect_to_vert[3];
478         float f_force, force_vec[3];
479         float *obloc;
480         float distance, force_val, ffall_val;
481         float guidecollect[3], guidedist= 0.0f;
482         short cur_frame;
483         
484         guidecollect[0]= guidecollect[1]= guidecollect[2]=0.0f;
485
486         /* Cycle through collected objects, get total of (1/(gravity_strength * dist^gravity_power)) */
487         /* Check for min distance here? (yes would be cool to add that, ton) */
488         
489         for(ec = lb->first; ec; ec= ec->next) {
490                 /* object effectors were fully checked to be OK to evaluate! */
491                 ob= ec->ob;
492                 pd= ob->pd;
493                         
494                 /* Get IPO force strength and fall off values here */
495                 if (has_ipo_code(ob->ipo, OB_PD_FSTR))
496                         force_val = IPO_GetFloatValue(ob->ipo, OB_PD_FSTR, cur_time);
497                 else 
498                         force_val = pd->f_strength;
499                 
500                 if (has_ipo_code(ob->ipo, OB_PD_FFALL)) 
501                         ffall_val = IPO_GetFloatValue(ob->ipo, OB_PD_FFALL, cur_time);
502                 else 
503                         ffall_val = pd->f_power;
504                         
505                 /* Need to set r.cfra for paths (investigate, ton) (uses ob->ctime now, ton) */
506                 if(ob->ctime!=cur_time) {
507                         cur_frame = G.scene->r.cfra;
508                         G.scene->r.cfra = (short)cur_time;
509                         where_is_object_time(ob, cur_time);
510                         G.scene->r.cfra = cur_frame;
511                 }
512                         
513                 /* use center of object for distance calculus */
514                 obloc= ob->obmat[3];
515                 VECSUB(vect_to_vert, obloc, opco);
516                 distance = VecLength(vect_to_vert);
517                         
518                 if((pd->flag & PFIELD_USEMAX) && distance>pd->maxdist && pd->forcefield != PFIELD_GUIDE)
519                         ;       /* don't do anything */
520                 else if(pd->forcefield == PFIELD_WIND) {
521                         VECCOPY(force_vec, ob->obmat[2]);
522                         
523                         /* wind works harder perpendicular to normal, would be nice for softbody later (ton) */
524                         
525                         /* Limit minimum distance to vertex so that */
526                         /* the force is not too big */
527                         if (distance < 0.001) distance = 0.001f;
528                         f_force = (force_val)*(1/(1000 * (float)pow((double)distance, (double)ffall_val)));
529                         /* this option for softbody only */
530                         if(flags && PE_WIND_AS_SPEED){
531                                 speed[0] -= (force_vec[0] * f_force );
532                                 speed[1] -= (force_vec[1] * f_force );
533                                 speed[2] -= (force_vec[2] * f_force );
534                         }
535                         else{
536                                 force[0] += force_vec[0]*f_force;
537                                 force[1] += force_vec[1]*f_force;
538                                 force[2] += force_vec[2]*f_force;
539                         }
540                 }
541                 else if(pd->forcefield == PFIELD_FORCE) {
542                         
543                         /* only use center of object */
544                         obloc= ob->obmat[3];
545
546                         /* Now calculate the gravitational force */
547                         VECSUB(vect_to_vert, obloc, opco);
548                         distance = VecLength(vect_to_vert);
549
550                         /* Limit minimum distance to vertex so that */
551                         /* the force is not too big */
552                         if (distance < 0.001) distance = 0.001f;
553                         f_force = (force_val)*(1.0/(1000.0 * (float)pow((double)distance, (double)ffall_val)));
554                         force[0] += (vect_to_vert[0] * f_force );
555                         force[1] += (vect_to_vert[1] * f_force );
556                         force[2] += (vect_to_vert[2] * f_force );
557                 }
558                 else if(pd->forcefield == PFIELD_VORTEX) {
559                         float vortexvec[3];
560                         
561                         /* only use center of object */
562                         obloc= ob->obmat[3];
563
564                         /* Now calculate the vortex force */
565                         VECSUB(vect_to_vert, obloc, opco);
566                         distance = VecLength(vect_to_vert);
567
568                         Crossf(force_vec, ob->obmat[2], vect_to_vert);
569                         Normalise(force_vec);
570
571                         /* Limit minimum distance to vertex so that */
572                         /* the force is not too big */
573                         if (distance < 0.001) distance = 0.001f;
574                         f_force = (force_val)*(1.0/(100.0 * (float)pow((double)distance, (double)ffall_val)));
575                         vortexvec[0]= -(force_vec[0] * f_force );
576                         vortexvec[1]= -(force_vec[1] * f_force );
577                         vortexvec[2]= -(force_vec[2] * f_force );
578                         
579                         /* this option for softbody only */
580                         if(flags &&PE_WIND_AS_SPEED) {
581                                 speed[0]+= vortexvec[0];
582                                 speed[1]+= vortexvec[1];
583                                 speed[2]+= vortexvec[2];
584                         }
585                         else {
586                                 /* since vortex alters the speed, we have to correct for the previous vortex result */
587                                 speed[0]+= vortexvec[0] - ec->oldspeed[0];
588                                 speed[1]+= vortexvec[1] - ec->oldspeed[1];
589                                 speed[2]+= vortexvec[2] - ec->oldspeed[2];
590                                 
591                                 VECCOPY(ec->oldspeed, vortexvec);
592                         }
593                 }
594                 else if(pd->forcefield == PFIELD_GUIDE) {
595                         float guidevec[4], guidedir[3];
596                         float mindist= force_val; /* force_val is actually mindist in the UI */
597                         
598                         distance= ec->guide_dist;
599                         
600                         /* WARNING: bails out with continue here */
601                         if((pd->flag & PFIELD_USEMAX) && distance>pd->maxdist) continue;
602                         
603                         /* calculate contribution factor for this guide */
604                         if(distance<=mindist) f_force= 1.0f;
605                         else if(pd->flag & PFIELD_USEMAX) {
606                                 if(distance>pd->maxdist || mindist>=pd->maxdist) f_force= 0.0f;
607                                 else {
608                                         f_force= 1.0f - (distance-mindist)/(pd->maxdist - mindist);
609                                         if(ffall_val!=0.0f)
610                                                 f_force = (float)pow(f_force, ffall_val+1.0);
611                                 }
612                         }
613                         else {
614                                 f_force= 1.0f/(1.0f + distance-mindist);
615                                 if(ffall_val!=0.0f)
616                                         f_force = (float)pow(f_force, ffall_val+1.0);
617                         }
618                         
619                         /* now derive path point from loc_time */
620                         if(pd->flag & PFIELD_GUIDE_PATH_ADD)
621                                 where_on_path(ob, f_force*loc_time*ec->time_scale, guidevec, guidedir);
622                         else
623                                 where_on_path(ob, loc_time*ec->time_scale, guidevec, guidedir);
624                         
625                         VECSUB(guidedir, guidevec, ec->oldloc);
626                         VECCOPY(ec->oldloc, guidevec);
627                         
628                         Mat4Mul3Vecfl(ob->obmat, guidedir);
629                         VecMulf(guidedir, ec->scale);   /* correction for lifetime and speed */
630                         
631                         /* we subtract the speed we gave it previous step */
632                         VECCOPY(guidevec, guidedir);
633                         VECSUB(guidedir, guidedir, ec->oldspeed);
634                         VECCOPY(ec->oldspeed, guidevec);
635                         
636                         /* if it fully contributes, we stop */
637                         if(f_force==1.0) {
638                                 VECCOPY(guidecollect, guidedir);
639                                 guidedist= 1.0f;
640                                 break;
641                         }
642                         else if(guidedist<1.0f) {
643                                 VecMulf(guidedir, f_force);
644                                 VECADD(guidecollect, guidecollect, guidedir);
645                                 guidedist += f_force;
646                         }                                       
647                 }
648         }
649
650         /* all guides are accumulated here */
651         if(guidedist!=0.0f) {
652                 if(guidedist!=1.0f) VecMulf(guidecollect, 1.0f/guidedist);
653                 VECADD(speed, speed, guidecollect);
654         }
655 }
656
657 static void cache_object_vertices(Object *ob)
658 {
659         Mesh *me;
660         MVert *mvert;
661         float *fp;
662         int a;
663         
664         me= ob->data;
665         if(me->totvert==0) return;
666
667         fp= ob->sumohandle= MEM_mallocN(3*sizeof(float)*me->totvert, "cache particles");
668         mvert= me->mvert;
669         a= me->totvert;
670         while(a--) {
671                 VECCOPY(fp, mvert->co);
672                 Mat4MulVecfl(ob->obmat, fp);
673                 mvert++;
674                 fp+= 3;
675         }
676 }
677
678 static int pdDoDeflection(RNG *rng, float opco[3], float npco[3], float opno[3],
679         float npno[3], float life, float force[3], int def_depth,
680         float cur_time, unsigned int par_layer, int *last_object,
681                 int *last_face, int *same_face)
682 {
683         /* Particle deflection code */
684         /* The code is in two sections: the first part checks whether a particle has            */
685         /* intersected a face of a deflector mesh, given its old and new co-ords, opco and npco */
686         /* and which face it hit first                                                          */
687         /* The second part calculates the new co-ordinates given that collision and updates     */
688         /* the new co-ordinates accordingly */
689         Base *base;
690         Object *ob, *deflection_object = NULL;
691         Mesh *def_mesh;
692         MFace *mface, *deflection_face = NULL;
693         float *v1, *v2, *v3, *v4, *vcache=NULL;
694         float nv1[3], nv2[3], nv3[3], nv4[3], edge1[3], edge2[3];
695         float dv1[3], dv2[3], dv3[3];
696         float vect_to_int[3], refl_vel[3];
697         float d_intersect_co[3], d_intersect_vect[3], d_nvect[3], d_i_co_above[3];
698         float forcec[3];
699         float k_point3, dist_to_plane;
700         float first_dist, ref_plane_mag;
701         float dk_plane=0, dk_point1=0;
702         float icalctop, icalcbot, n_mag;
703         float mag_iv, x_m,y_m,z_m;
704         float damping, perm_thresh;
705         float perm_val, rdamp_val;
706         int a, deflected=0, deflected_now=0;
707         float t,t2, min_t;
708         float mat[3][3], obloc[3];
709         short cur_frame;
710         float time_before, time_after;
711         float force_mag_norm;
712         int d_object=0, d_face=0, ds_object=0, ds_face=0;
713
714         first_dist = 200000;
715         min_t = 200000;
716
717         /* The first part of the code, finding the first intersected face*/
718         base= G.scene->base.first;
719         while (base) {
720                 /*Only proceed for mesh object in same layer */
721                 if(base->object->type==OB_MESH && (base->lay & par_layer)) {
722                         ob= base->object;
723                         /* only with deflecting set */
724                         if(ob->pd && ob->pd->deflect) {
725                                 def_mesh= ob->data;
726                         
727                                 d_object = d_object + 1;
728
729                                 d_face = d_face + 1;
730                                 mface= def_mesh->mface;
731                                 a = def_mesh->totface;
732                                 
733                                 
734                                 if(ob->parent==NULL && ob->ipo==NULL) { // static
735                                         if(ob->sumohandle==NULL) cache_object_vertices(ob);
736                                         vcache= ob->sumohandle;
737                                 }
738                                 else {
739                                         /*Find out where the object is at this time*/
740                                         cur_frame = G.scene->r.cfra;
741                                         G.scene->r.cfra = (short)cur_time;
742                                         where_is_object_time(ob, cur_time);
743                                         G.scene->r.cfra = cur_frame;
744                                         
745                                         /*Pass the values from ob->obmat to mat*/
746                                         /*and the location values to obloc           */
747                                         Mat3CpyMat4(mat,ob->obmat);
748                                         obloc[0] = ob->obmat[3][0];
749                                         obloc[1] = ob->obmat[3][1];
750                                         obloc[2] = ob->obmat[3][2];
751                                         vcache= NULL;
752
753                                 }
754                                 
755                                 while (a--) {
756
757                                         if(vcache) {
758                                                 v1= vcache+ 3*(mface->v1);
759                                                 VECCOPY(nv1, v1);
760                                                 v1= vcache+ 3*(mface->v2);
761                                                 VECCOPY(nv2, v1);
762                                                 v1= vcache+ 3*(mface->v3);
763                                                 VECCOPY(nv3, v1);
764                                                 v1= vcache+ 3*(mface->v4);
765                                                 VECCOPY(nv4, v1);
766                                         }
767                                         else {
768                                                 /* Calculate the global co-ordinates of the vertices*/
769                                                 v1= (def_mesh->mvert+(mface->v1))->co;
770                                                 v2= (def_mesh->mvert+(mface->v2))->co;
771                                                 v3= (def_mesh->mvert+(mface->v3))->co;
772                                                 v4= (def_mesh->mvert+(mface->v4))->co;
773         
774                                                 VECCOPY(nv1, v1);
775                                                 VECCOPY(nv2, v2);
776                                                 VECCOPY(nv3, v3);
777                                                 VECCOPY(nv4, v4);
778         
779                                                 /*Apply the objects deformation matrix*/
780                                                 Mat3MulVecfl(mat, nv1);
781                                                 Mat3MulVecfl(mat, nv2);
782                                                 Mat3MulVecfl(mat, nv3);
783                                                 Mat3MulVecfl(mat, nv4);
784         
785                                                 VECADD(nv1, nv1, obloc);
786                                                 VECADD(nv2, nv2, obloc);
787                                                 VECADD(nv3, nv3, obloc);
788                                                 VECADD(nv4, nv4, obloc);
789                                         }
790                                         
791                                         deflected_now = 0;
792
793                                                 
794                                                 
795 //                                      t= 0.5; // this is labda of line, can use it optimize quad intersection
796 // sorry but no .. see below (BM)                                       
797                                         if( LineIntersectsTriangle(opco, npco, nv1, nv2, nv3, &t) ) {
798                                                 if (t < min_t) {
799                                                         deflected = 1;
800                                                         deflected_now = 1;
801                                                 }
802                                         }
803 //                                      else if (mface->v4 && (t>=0.0 && t<=1.0)) {
804 // no, you can't skip testing the other triangle
805 // it might give a smaller t on (close to) the edge .. this is numerics not esoteric maths :)
806 // note: the 2 triangles don't need to share a plane ! (BM)
807                                         if (mface->v4) {
808                                                 if( LineIntersectsTriangle(opco, npco, nv1, nv3, nv4, &t2) ) {
809                                                         if (t2 < min_t) {
810                                                                 deflected = 1;
811                                                                 deflected_now = 2;
812                                                         }
813                                                 }
814                                         }
815                                         
816                                         if ((deflected_now > 0) && ((t < min_t) ||(t2 < min_t))) {
817                         min_t = t;
818                         ds_object = d_object;
819                                                 ds_face = d_face;
820                                                 deflection_object = ob;
821                                                 deflection_face = mface;
822                                                 if (deflected_now==1) {
823                         min_t = t;
824                                                         VECCOPY(dv1, nv1);
825                                                         VECCOPY(dv2, nv2);
826                                                         VECCOPY(dv3, nv3);
827                                                 }
828                                                 else {
829                         min_t = t2;
830                                                         VECCOPY(dv1, nv1);
831                                                         VECCOPY(dv2, nv3);
832                                                         VECCOPY(dv3, nv4);
833                                                 }
834                                         }
835                                         mface++;
836                                 }
837                         }
838                 }
839                 base = base->next;
840         }
841
842
843         /* Here's the point to do the permeability calculation */
844         /* Set deflected to 0 if a random number is below the value */
845         /* Get the permeability IPO here*/
846         if (deflected) {
847                 
848                 if (has_ipo_code(deflection_object->ipo, OB_PD_PERM)) 
849                         perm_val = IPO_GetFloatValue(deflection_object->ipo, OB_PD_PERM, cur_time);
850                 else 
851                         perm_val = deflection_object->pd->pdef_perm;
852
853                 perm_thresh = rng_getFloat(rng) - perm_val;
854                 if (perm_thresh < 0 ) {
855                         deflected = 0;
856                 }
857         }
858
859         /* Now for the second part of the deflection code - work out the new speed */
860         /* and position of the particle if a collision occurred */
861         if (deflected) {
862         VECSUB(edge1, dv1, dv2);
863                 VECSUB(edge2, dv3, dv2);
864                 Crossf(d_nvect, edge2, edge1);
865                 n_mag = Normalise(d_nvect);
866                 dk_plane = INPR(d_nvect, nv1);
867                 dk_point1 = INPR(d_nvect,opco);
868
869                 VECSUB(d_intersect_vect, npco, opco);
870
871                 d_intersect_co[0] = opco[0] + (min_t * (npco[0] - opco[0]));
872                 d_intersect_co[1] = opco[1] + (min_t * (npco[1] - opco[1]));
873                 d_intersect_co[2] = opco[2] + (min_t * (npco[2] - opco[2]));
874                 
875                 d_i_co_above[0] = (d_intersect_co[0] + (0.001f * d_nvect[0]));
876                 d_i_co_above[1] = (d_intersect_co[1] + (0.001f * d_nvect[1]));
877                 d_i_co_above[2] = (d_intersect_co[2] + (0.001f * d_nvect[2]));
878                 mag_iv = Normalise(d_intersect_vect);
879                 VECCOPY(npco, d_intersect_co);
880                 
881                 VECSUB(vect_to_int, opco, d_intersect_co);
882                 first_dist = Normalise(vect_to_int);
883
884                 /* Work out the lengths of time before and after collision*/
885                 time_before = (life*(first_dist / (mag_iv)));
886                 time_after =  (life*((mag_iv - first_dist) / (mag_iv)));
887
888                 /* We have to recalculate what the speed would have been at the */
889                 /* point of collision, not the key frame time */
890                 npno[0]= opno[0] + time_before*force[0];
891                 npno[1]= opno[1] + time_before*force[1];
892                 npno[2]= opno[2] + time_before*force[2];
893
894
895                 /* Reflect the speed vector in the face */
896                 x_m = (2 * npno[0] * d_nvect[0]);
897                 y_m = (2 * npno[1] * d_nvect[1]);
898                 z_m = (2 * npno[2] * d_nvect[2]);
899                 refl_vel[0] = npno[0] - (d_nvect[0] * (x_m + y_m + z_m));
900                 refl_vel[1] = npno[1] - (d_nvect[1] * (x_m + y_m + z_m));
901                 refl_vel[2] = npno[2] - (d_nvect[2] * (x_m + y_m + z_m));
902
903                 /*A random variation in the damping factor........ */
904                 /*Get the IPO values for damping here*/
905                 
906                 if (has_ipo_code(deflection_object->ipo, OB_PD_SDAMP)) 
907                         damping = IPO_GetFloatValue(deflection_object->ipo, OB_PD_SDAMP, cur_time);
908                 else 
909                         damping = deflection_object->pd->pdef_damp;
910                 
911                 if (has_ipo_code(deflection_object->ipo, OB_PD_RDAMP)) 
912                         rdamp_val = IPO_GetFloatValue(deflection_object->ipo, OB_PD_RDAMP, cur_time);
913                 else 
914                         rdamp_val = deflection_object->pd->pdef_rdamp;
915
916                 damping = damping + ((1.0f - damping) * rng_getFloat(rng) *rdamp_val);
917                 damping = damping * damping;
918         ref_plane_mag = INPR(refl_vel,d_nvect);
919
920                 if (damping > 0.999) damping = 0.999f;
921
922                 /* Now add in the damping force - only damp in the direction of */
923                 /* the faces normal vector */
924                 npno[0] = (refl_vel[0] - (d_nvect[0] * ref_plane_mag * damping));
925                 npno[1] = (refl_vel[1] - (d_nvect[1] * ref_plane_mag * damping));
926                 npno[2] = (refl_vel[2] - (d_nvect[2] * ref_plane_mag * damping));
927
928                 /* Now reset opno */
929                 VECCOPY(opno,npno);
930                 VECCOPY(forcec, force);
931
932                 /* If the particle has bounced more than four times on the same */
933                 /* face within this cycle (depth > 4, same face > 4 )           */
934                 /* Then set the force to be only that component of the force    */
935                 /* in the same direction as the face normal                     */
936                 /* i.e. subtract the component of the force in the direction    */
937                 /* of the face normal from the actual force                     */
938                 if ((ds_object == *last_object) && (ds_face == *last_face)) {
939                         /* Increment same_face */
940                         *same_face = *same_face + 1;
941                         if ((*same_face > 3) && (def_depth > 3)) {
942                 force_mag_norm = INPR(forcec, d_nvect);
943                 forcec[0] = forcec[0] - (d_nvect[0] * force_mag_norm);
944                 forcec[1] = forcec[1] - (d_nvect[1] * force_mag_norm);
945                 forcec[2] = forcec[2] - (d_nvect[2] * force_mag_norm);
946                         }
947                 }
948                 else *same_face = 1;
949
950                 *last_object = ds_object;
951                 *last_face = ds_face;
952
953                 /* We have the particles speed at the point of collision    */
954                 /* Now we want the particles speed at the current key frame */
955
956                 npno[0]= npno[0] + time_after*forcec[0];
957                 npno[1]= npno[1] + time_after*forcec[1];
958                 npno[2]= npno[2] + time_after*forcec[2];
959
960                 /* Now we have to recalculate pa->co for the remainder*/
961                 /* of the time since the intersect*/
962                 npco[0]= npco[0] + time_after*npno[0];
963                 npco[1]= npco[1] + time_after*npno[1];
964                 npco[2]= npco[2] + time_after*npno[2];
965
966                 /* And set the old co-ordinates back to the point just above the intersection */
967                 VECCOPY(opco, d_i_co_above);
968
969                 /* Finally update the time */
970                 life = time_after;
971                 cur_time += time_before;
972
973                 /* The particle may have fallen through the face again by now!!*/
974                 /* So check if the particle has changed sides of the plane compared*/
975                 /* the co-ordinates at the last keyframe*/
976                 /* But only do this as a last resort, if we've got to the end of the */
977                 /* number of collisions allowed */
978                 if (def_depth==9) {
979                         k_point3 = INPR(d_nvect,npco);
980                         if (((dk_plane > k_point3) && (dk_plane < dk_point1))||((dk_plane < k_point3) && (dk_plane > dk_point1))) {
981
982                                 /* Yup, the pesky particle may have fallen through a hole!!! */
983                 /* So we'll cheat a bit and move the particle along the normal vector */
984                 /* until it's just the other side of the plane */
985                 icalctop = (dk_plane - d_nvect[0]*npco[0] - d_nvect[1]*npco[1] - d_nvect[2]*npco[2]);
986                 icalcbot = (d_nvect[0]*d_nvect[0] + d_nvect[1]*d_nvect[1] + d_nvect[2]*d_nvect[2]);
987                 dist_to_plane = icalctop / icalcbot;
988
989                 /*  Now just increase the distance a little to place */
990                 /* the point the other side of the plane */
991                 dist_to_plane *= 1.1f;
992                 npco[0]= npco[0] + (dist_to_plane * d_nvect[0]);
993                 npco[1]= npco[1] + (dist_to_plane * d_nvect[1]);
994                 npco[2]= npco[2] + (dist_to_plane * d_nvect[2]);
995
996                         }
997                 }
998         }
999         return deflected;
1000 }
1001
1002 /*
1003         rng= random number generator
1004         ob = object that spawns the particles
1005         depth = for fireworks
1006         nr = index nr of current particle
1007         paf = the particle system
1008         part = current particle
1009         force = force vector
1010         deform = flag to indicate lattice deform
1011  */
1012 static void make_particle_keys(RNG *rng, Object *ob, int depth, int nr, PartEff *paf, Particle *part, float *force, int deform, MTex *mtex, ListBase *effectorbase)
1013 {
1014         Particle *pa, *opa = NULL;
1015         float damp, deltalife, life;
1016         float cur_time;
1017         float opco[3], opno[3], npco[3], npno[3], new_force[3], new_speed[3];
1018         int b, rt1, rt2, deflected, deflection, finish_defs, def_count;
1019         int last_ob, last_fc, same_fc;
1020
1021         damp= 1.0f-paf->damp;
1022         pa= part;
1023
1024         /* start speed: random */
1025         if(paf->randfac!=0.0) {
1026                 pa->no[0]+= paf->randfac*(rng_getFloat(rng) - 0.5f);
1027                 pa->no[1]+= paf->randfac*(rng_getFloat(rng) - 0.5f);
1028                 pa->no[2]+= paf->randfac*(rng_getFloat(rng) - 0.5f);
1029         }
1030
1031         /* start speed: texture */
1032         if(mtex && paf->texfac!=0.0) {
1033                 particle_tex(mtex, paf, pa->co, pa->no);
1034         }
1035
1036         /* effectors here? */
1037         if(effectorbase)
1038                 precalc_effectors(ob, paf, pa, effectorbase);
1039         
1040         if(paf->totkey>1) deltalife= pa->lifetime/(paf->totkey-1);
1041         else deltalife= pa->lifetime;
1042
1043         /* longer lifetime results in longer distance covered */
1044         VecMulf(pa->no, deltalife);
1045         
1046         opa= pa;
1047         pa++;
1048
1049         for(b=1; b<paf->totkey; b++) {
1050
1051                 /* new time */
1052                 pa->time= opa->time+deltalife;
1053                 cur_time = pa->time;
1054
1055                 /* set initial variables                                */
1056                 VECCOPY(opco, opa->co);
1057                 VECCOPY(new_force, force);
1058                 VECCOPY(new_speed, opa->no);
1059                 VecMulf(new_speed, 1.0f/deltalife);
1060                 //new_speed[0] = new_speed[1] = new_speed[2] = 0.0f;
1061
1062                 /* handle differences between static (local coords, fixed frame) and dynamic */
1063                 if(effectorbase) {
1064                         float loc_time= ((float)b)/(float)(paf->totkey-1);
1065                         
1066                         if(paf->flag & PAF_STATIC) {
1067                                 float opco1[3], new_force1[3];
1068                                 
1069                                 /* move co and force to global coords */
1070                                 VECCOPY(opco1, opco);
1071                                 Mat4MulVecfl(ob->obmat, opco1);
1072                                 VECCOPY(new_force1, new_force);
1073                                 Mat4Mul3Vecfl(ob->obmat, new_force1);
1074                                 Mat4Mul3Vecfl(ob->obmat, new_speed);
1075                                 
1076                                 cur_time = G.scene->r.cfra;
1077                                 
1078                                 /* force fields */
1079                                 pdDoEffectors(effectorbase, opco1, new_force1, new_speed, cur_time, loc_time, 0);
1080                                 
1081                                 /* move co, force and newspeed back to local */
1082                                 VECCOPY(opco, opco1);
1083                                 Mat4MulVecfl(ob->imat, opco);
1084                                 VECCOPY(new_force, new_force1);
1085                                 Mat4Mul3Vecfl(ob->imat, new_force);
1086                                 Mat4Mul3Vecfl(ob->imat, new_speed);
1087                         }
1088                         else {
1089                                  /* force fields */
1090                                 pdDoEffectors(effectorbase, opco, new_force, new_speed, cur_time, loc_time, 0);
1091                         }
1092                 }
1093                 
1094                 /* new speed */
1095                 pa->no[0]= deltalife * (new_speed[0] + new_force[0]);
1096                 pa->no[1]= deltalife * (new_speed[1] + new_force[1]);
1097                 pa->no[2]= deltalife * (new_speed[2] + new_force[2]);
1098                 
1099                 /* new location */
1100                 pa->co[0]= opa->co[0] + pa->no[0];
1101                 pa->co[1]= opa->co[1] + pa->no[1];
1102                 pa->co[2]= opa->co[2] + pa->no[2];
1103
1104                 /* Particle deflection code */
1105                 if((paf->flag & PAF_STATIC)==0) {
1106                         deflection = 0;
1107                         finish_defs = 1;
1108                         def_count = 0;
1109
1110                         VECCOPY(opno, opa->no);
1111                         VECCOPY(npco, pa->co);
1112                         VECCOPY(npno, pa->no);
1113
1114                         life = deltalife;
1115                         cur_time -= deltalife;
1116
1117                         last_ob = -1;
1118                         last_fc = -1;
1119                         same_fc = 0;
1120
1121                         /* First call the particle deflection check for the particle moving   */
1122                         /* between the old co-ordinates and the new co-ordinates              */
1123                         /* If a deflection occurs, call the code again, this time between the */
1124                         /* intersection point and the updated new co-ordinates                */
1125                         /* Bail out if we've done the calculation 10 times - this seems ok     */
1126                         /* for most scenes I've tested */
1127                         while (finish_defs) {
1128                                 deflected =  pdDoDeflection(rng, opco, npco, opno, npno, life, new_force,
1129                                                                 def_count, cur_time, ob->lay,
1130                                                                 &last_ob, &last_fc, &same_fc);
1131                                 if (deflected) {
1132                                         def_count = def_count + 1;
1133                                         deflection = 1;
1134                                         if (def_count==10) finish_defs = 0;
1135                                 }
1136                                 else {
1137                                         finish_defs = 0;
1138                                 }
1139                         }
1140
1141                         /* Only update the particle positions and speed if we had a deflection */
1142                         if (deflection) {
1143                                 pa->co[0] = npco[0];
1144                                 pa->co[1] = npco[1];
1145                                 pa->co[2] = npco[2];
1146                                 pa->no[0] = npno[0];
1147                                 pa->no[1] = npno[1];
1148                                 pa->no[2] = npno[2];
1149                         }
1150                 }
1151                 
1152                 /* speed: texture */
1153                 if(mtex && paf->texfac!=0.0) {
1154                         particle_tex(mtex, paf, pa->co, pa->no);
1155                 }
1156                 if(damp!=1.0) {
1157                         pa->no[0]*= damp;
1158                         pa->no[1]*= damp;
1159                         pa->no[2]*= damp;
1160                 }
1161                 
1162                 opa= pa;
1163                 pa++;
1164                 /* opa is used later on too! */
1165         }
1166
1167         if(deform) {
1168                 /* deform all keys */
1169                 pa= part;
1170                 b= paf->totkey;
1171                 while(b--) {
1172                         calc_latt_deform(pa->co, 1.0f);
1173                         pa++;
1174                 }
1175         }
1176         
1177         /* the big multiplication */
1178         if(depth<PAF_MAXMULT && paf->mult[depth]!=0.0) {
1179                 
1180                 /* new 'child' emerges from an average 'mult' part from 
1181                         the particles */
1182                 damp = (float)nr;
1183                 rt1= (int)(damp*paf->mult[depth]);
1184                 rt2= (int)((damp+1.0)*paf->mult[depth]);
1185                 if(rt1!=rt2) {
1186                         
1187                         for(b=0; b<paf->child[depth]; b++) {
1188                                 pa= new_particle(paf);
1189                                 *pa= *opa;
1190                                 pa->lifetime= paf->life[depth];
1191                                 if(paf->randlife!=0.0) {
1192                                         pa->lifetime*= 1.0f + paf->randlife*(rng_getFloat(rng) - 0.5f);
1193                                 }
1194                                 pa->mat_nr= paf->mat[depth];
1195
1196                                 make_particle_keys(rng, ob, depth+1, b, paf, pa, force, deform, mtex, effectorbase);
1197                         }
1198                 }
1199         }
1200 }
1201
1202 static void init_mv_jit(float *jit, int num, int seed2)
1203 {
1204         RNG *rng;
1205         float *jit2, x, rad1, rad2, rad3;
1206         int i, num2;
1207
1208         if(num==0) return;
1209
1210         rad1= (float)(1.0/sqrt((float)num));
1211         rad2= (float)(1.0/((float)num));
1212         rad3= (float)sqrt((float)num)/((float)num);
1213
1214         rng = rng_new(31415926 + num + seed2);
1215         x= 0;
1216         num2 = 2 * num;
1217         for(i=0; i<num2; i+=2) {
1218         
1219                 jit[i]= x + rad1*(0.5f - rng_getFloat(rng));
1220                 jit[i+1]= i/(2.0f*num) + rad1*(0.5f - rng_getFloat(rng));
1221                 
1222                 jit[i]-= (float)floor(jit[i]);
1223                 jit[i+1]-= (float)floor(jit[i+1]);
1224                 
1225                 x+= rad3;
1226                 x -= (float)floor(x);
1227         }
1228
1229         jit2= MEM_mallocN(12 + 2*sizeof(float)*num, "initjit");
1230
1231         for (i=0 ; i<4 ; i++) {
1232                 RE_jitterate1(jit, jit2, num, rad1);
1233                 RE_jitterate1(jit, jit2, num, rad1);
1234                 RE_jitterate2(jit, jit2, num, rad2);
1235         }
1236         MEM_freeN(jit2);
1237         rng_free(rng);
1238 }
1239
1240 #define JIT_RAND        32
1241
1242 /* for a position within a face, tot is total amount of faces */
1243 static void give_mesh_particle_coord(PartEff *paf, VeNoCo *noco, MFace *mface, int partnr, int subnr, float *co, float *no)
1244 {
1245         static float *jit= NULL;
1246         static float *trands= NULL;
1247         static int jitlevel= 1;
1248         float *v1, *v2, *v3, *v4;
1249         float u, v;
1250         float *n1, *n2, *n3, *n4;
1251         
1252         /* free signal */
1253         if(paf==NULL) {
1254                 if(jit) MEM_freeN(jit);
1255                 jit= NULL;
1256                 if(trands) MEM_freeN(trands);
1257                 trands= NULL;
1258                 return;
1259         }
1260         
1261         /* first time initialize jitter or trand, partnr then is total amount of particles, subnr total amount of faces */
1262         if(trands==NULL && jit==NULL) {
1263                 RNG *rng = rng_new(31415926 + paf->seed);
1264                 int i, tot;
1265
1266                 if(paf->flag & PAF_TRAND)
1267                         tot= partnr;
1268                 else 
1269                         tot= JIT_RAND;  /* arbitrary... allows JIT_RAND times more particles in a face for jittered distro */
1270                         
1271                 trands= MEM_callocN(2+2*tot*sizeof(float), "trands");
1272                 for(i=0; i<tot; i++) {
1273                         trands[2*i]= rng_getFloat(rng);
1274                         trands[2*i+1]= rng_getFloat(rng);
1275                 }
1276                 rng_free(rng);
1277
1278                 if((paf->flag & PAF_TRAND)==0) {
1279                         jitlevel= paf->userjit;
1280                         
1281                         if(jitlevel == 0) {
1282                                 jitlevel= partnr/subnr;
1283                                 if(paf->flag & PAF_EDISTR) jitlevel*= 2;        /* looks better in general, not very scietific */
1284                                 if(jitlevel<3) jitlevel= 3;
1285                                 if(jitlevel>100) jitlevel= 100;
1286                         }
1287                         
1288                         jit= MEM_callocN(2+ jitlevel*2*sizeof(float), "jit");
1289                         init_mv_jit(jit, jitlevel, paf->seed);
1290                         BLI_array_randomize(jit, 2*sizeof(float), jitlevel, paf->seed); /* for custom jit or even distribution */
1291                 }
1292                 return;
1293         }
1294         
1295         if(paf->flag & PAF_TRAND) {
1296                 u= trands[2*partnr];
1297                 v= trands[2*partnr+1];
1298         }
1299         else {
1300                 /* jittered distribution gets fixed random offset */
1301                 if(subnr>=jitlevel) {
1302                         int jitrand= (subnr/jitlevel) % JIT_RAND;
1303                 
1304                         subnr %= jitlevel;
1305                         u= jit[2*subnr] + trands[2*jitrand];
1306                         v= jit[2*subnr+1] + trands[2*jitrand+1];
1307                         if(u > 1.0f) u-= 1.0f;
1308                         if(v > 1.0f) v-= 1.0f;
1309                 }
1310                 else {
1311                         u= jit[2*subnr];
1312                         v= jit[2*subnr+1];
1313                 }
1314         }
1315         
1316         v1= (noco+(mface->v1))->co;
1317         v2= (noco+(mface->v2))->co;
1318         v3= (noco+(mface->v3))->co;
1319         n1= (noco+(mface->v1))->no;
1320         n2= (noco+(mface->v2))->no;
1321         n3= (noco+(mface->v3))->no;
1322         
1323         if(mface->v4) {
1324                 float uv= u*v;
1325                 float muv= (1.0f-u)*(v);
1326                 float umv= (u)*(1.0f-v);
1327                 float mumv= (1.0f-u)*(1.0f-v);
1328                 
1329                 v4= (noco+(mface->v4))->co;
1330                 n4= (noco+(mface->v4))->no;
1331                 
1332                 co[0]= mumv*v1[0] + muv*v2[0] + uv*v3[0] + umv*v4[0];
1333                 co[1]= mumv*v1[1] + muv*v2[1] + uv*v3[1] + umv*v4[1];
1334                 co[2]= mumv*v1[2] + muv*v2[2] + uv*v3[2] + umv*v4[2];
1335
1336                 no[0]= mumv*n1[0] + muv*n2[0] + uv*n3[0] + umv*n4[0];
1337                 no[1]= mumv*n1[1] + muv*n2[1] + uv*n3[1] + umv*n4[1];
1338                 no[2]= mumv*n1[2] + muv*n2[2] + uv*n3[2] + umv*n4[2];
1339         }
1340         else {
1341                 /* mirror triangle uv coordinates when on other side */
1342                 if(u + v > 1.0f) {
1343                         u= 1.0f-u;
1344                         v= 1.0f-v;
1345                 }
1346                 co[0]= v1[0] + u*(v3[0]-v1[0]) + v*(v2[0]-v1[0]);
1347                 co[1]= v1[1] + u*(v3[1]-v1[1]) + v*(v2[1]-v1[1]);
1348                 co[2]= v1[2] + u*(v3[2]-v1[2]) + v*(v2[2]-v1[2]);
1349                 
1350                 no[0]= n1[0] + u*(n3[0]-n1[0]) + v*(n2[0]-n1[0]);
1351                 no[1]= n1[1] + u*(n3[1]-n1[1]) + v*(n2[1]-n1[1]);
1352                 no[2]= n1[2] + u*(n3[2]-n1[2]) + v*(n2[2]-n1[2]);
1353         }
1354 }
1355
1356
1357 /* Gets a MDeformVert's weight in group (0 if not in group) */
1358 /* note; this call could be in mesh.c or deform.c, but OK... it's in armature.c too! (ton) */
1359 static float vert_weight(MDeformVert *dvert, int group)
1360 {
1361         MDeformWeight *dw;
1362         int i;
1363         
1364         if(dvert) {
1365                 dw= dvert->dw;
1366                 for(i= dvert->totweight; i>0; i--, dw++) {
1367                         if(dw->def_nr == group) return dw->weight;
1368                         if(i==1) break; /*otherwise dw will point to somewhere it shouldn't*/
1369                 }
1370         }
1371         return 0.0;
1372 }
1373
1374 /* Gets a faces average weight in a group, helper for below, face and weights are always set */
1375 static float face_weight(MFace *face, float *weights)
1376 {
1377         float tweight;
1378         
1379         tweight = weights[face->v1] + weights[face->v2] + weights[face->v3];
1380         
1381         if(face->v4) {
1382                 tweight += weights[face->v4];
1383                 tweight /= 4.0;
1384         }
1385         else {
1386                 tweight /= 3.0;
1387         }
1388
1389         return tweight;
1390 }
1391
1392 /* helper function for build_particle_system() */
1393 static void make_weight_tables(PartEff *paf, Mesh *me, int totpart, VeNoCo *vertlist, int totvert, MFace *facelist, int totface, float **vweights, float **fweights)
1394 {
1395         MFace *mface;
1396         float *foweights=NULL, *voweights=NULL;
1397         float totvweight=0.0f, totfweight=0.0f;
1398         int a;
1399         
1400         if((paf->flag & PAF_FACE)==0) totface= 0;
1401         
1402         /* collect emitting vertices & faces if vert groups used */
1403         if(paf->vertgroup && me->dvert) {
1404                 
1405                 /* allocate weights array for all vertices, also for lookup of faces later on. note it's a malloc */
1406                 *vweights= voweights= MEM_mallocN( totvert*sizeof(float), "pafvoweights" );
1407                 totvweight= 0.0f;
1408                 for(a=0; a<totvert; a++) {
1409                         voweights[a]= vert_weight(me->dvert+a, paf->vertgroup-1);
1410                         totvweight+= voweights[a];
1411                 }
1412                 
1413                 if(totface) {
1414                         /* allocate weights array for faces, note it's a malloc */
1415                         *fweights= foweights= MEM_mallocN(totface*sizeof(float), "paffoweights" );
1416                         for(a=0, mface=facelist; a<totface; a++, mface++) {
1417                                 foweights[a] = face_weight(mface, voweights);
1418                         }
1419                 }
1420         }
1421         
1422         /* make weights for faces or for even area distribution */
1423         if(totface && (paf->flag & PAF_EDISTR)) {
1424                 float maxfarea= 0.0f, curfarea;
1425                 
1426                 /* two cases for area distro, second case we already have group weights */
1427                 if(foweights==NULL) {
1428                         /* allocate weights array for faces, note it's a malloc */
1429                         *fweights= foweights= MEM_mallocN(totface*sizeof(float), "paffoweights" );
1430                         
1431                         for(a=0, mface=facelist; a<totface; a++, mface++) {
1432                                 if (mface->v4)
1433                                         curfarea= AreaQ3Dfl(vertlist[mface->v1].co, vertlist[mface->v2].co, vertlist[mface->v3].co, vertlist[mface->v4].co);
1434                                 else
1435                                         curfarea= AreaT3Dfl(vertlist[mface->v1].co,  vertlist[mface->v2].co, vertlist[mface->v3].co);
1436                                 if(curfarea>maxfarea)
1437                                         maxfarea = curfarea;
1438                                 foweights[a]= curfarea;
1439                         }
1440                 }
1441                 else {
1442                         for(a=0, mface=facelist; a<totface; a++, mface++) {
1443                                 if(foweights[a]!=0.0f) {
1444                                         if (mface->v4)
1445                                                 curfarea= AreaQ3Dfl(vertlist[mface->v1].co, vertlist[mface->v2].co, vertlist[mface->v3].co, vertlist[mface->v4].co);
1446                                         else
1447                                                 curfarea= AreaT3Dfl(vertlist[mface->v1].co,  vertlist[mface->v2].co, vertlist[mface->v3].co);
1448                                         if(curfarea>maxfarea)
1449                                                 maxfarea = curfarea;
1450                                         foweights[a]*= curfarea;
1451                                 }
1452                         }
1453                 }
1454                 
1455                 /* normalize weights for max face area, calculate tot */
1456                 if(maxfarea!=0.0f) {
1457                         maxfarea= 1.0f/maxfarea;
1458                         for(a=0; a< totface; a++) {
1459                                 if(foweights[a]!=0.0) {
1460                                         foweights[a] *= maxfarea;
1461                                         totfweight+= foweights[a];
1462                                 }
1463                         }
1464                 }
1465         }
1466         else if(foweights) {
1467                 /* only add totfweight value */
1468                 for(a=0; a< totface; a++) {
1469                         if(foweights[a]!=0.0) {
1470                                 totfweight+= foweights[a];
1471                         }
1472                 }
1473         }
1474         
1475         /* if weight arrays, we turn these arrays into the amount of particles */
1476         if(totvert && voweights) {
1477                 float mult= (float)totpart/totvweight;
1478                 
1479                 for(a=0; a< totvert; a++) {
1480                         if(voweights[a]!=0.0)
1481                                 voweights[a] *= mult;
1482                 }
1483         }
1484         
1485         if(totface && foweights) {
1486                 float mult= (float)totpart/totfweight;
1487                 
1488                 for(a=0; a< totface; a++) {
1489                         if(foweights[a]!=0.0)
1490                                 foweights[a] *= mult;
1491                 }
1492         }
1493 }
1494
1495 /* helper function for build_particle_system() */
1496 static void make_length_tables(PartEff *paf, Mesh *me, int totvert, MFace *facelist, int totface, float **vlengths, float **flengths)
1497 {
1498         MFace *mface;
1499         float *folengths=NULL, *volengths=NULL;
1500         int a;
1501         
1502         if((paf->flag & PAF_FACE)==0) totface= 0;
1503         
1504         /* collect emitting vertices & faces if vert groups used */
1505         if(paf->vertgroup_v && me->dvert) {
1506                 
1507                 /* allocate lengths array for all vertices, also for lookup of faces later on. note it's a malloc */
1508                 *vlengths= volengths= MEM_mallocN( totvert*sizeof(float), "pafvolengths" );
1509                 for(a=0; a<totvert; a++) {
1510                         volengths[a]= vert_weight(me->dvert+a, paf->vertgroup_v-1);
1511                 }
1512                 
1513                 if(totface) {
1514                         /* allocate lengths array for faces, note it's a malloc */
1515                         *flengths= folengths= MEM_mallocN(totface*sizeof(float), "paffolengths" );
1516                         for(a=0, mface=facelist; a<totface; a++, mface++) {
1517                                 folengths[a] = face_weight(mface, volengths);
1518                         }
1519                 }
1520         }
1521 }
1522
1523
1524 /* for paf start to end, store all matrices for objects */
1525 typedef struct pMatrixCache {
1526         float obmat[4][4];
1527         float imat[3][3];
1528 } pMatrixCache;
1529
1530 static pMatrixCache *cache_object_matrices(Object *ob, int start, int end)
1531 {
1532         pMatrixCache *mcache, *mc;
1533         Object *par, ob_store;
1534         float framelenold, cfrao;
1535         
1536         mcache= mc= MEM_mallocN( (end-start+1)*sizeof(pMatrixCache), "ob matrix cache");
1537         
1538         framelenold= G.scene->r.framelen;
1539         G.scene->r.framelen= 1.0f;
1540         cfrao= G.scene->r.cfra;
1541         ob_store= *ob;  /* quick copy of all settings */
1542         ob->sf= 0.0f;
1543         
1544         for(G.scene->r.cfra= start; G.scene->r.cfra<=end; G.scene->r.cfra++, mc++) {
1545                 
1546                 par= ob;
1547                 while(par) {
1548                         par->ctime= -1234567.0;         /* hrms? */
1549                         do_ob_key(par);
1550                         if(par->type==OB_ARMATURE) {
1551                                 do_all_pose_actions(par);       // only does this object actions
1552                                 where_is_pose(par);
1553                         }
1554                         par= par->parent;
1555                 }
1556
1557                 where_is_object(ob);
1558                 Mat4CpyMat4(mc->obmat, ob->obmat);
1559                 Mat4Invert(ob->imat, ob->obmat);
1560                 Mat3CpyMat4(mc->imat, ob->imat);
1561                 Mat3Transp(mc->imat);
1562         }
1563         
1564         
1565         /* restore */
1566         G.scene->r.cfra= cfrao;
1567         G.scene->r.framelen= framelenold;
1568         *ob= ob_store;
1569         
1570         /* restore hierarchy, weak code destroying potential depgraph stuff... */
1571         par= ob;
1572         while(par) {
1573                 /* do not do ob->ipo: keep insertkey */
1574                 do_ob_key(par);
1575                 
1576                 if(par->type==OB_ARMATURE) {
1577                         do_all_pose_actions(par);       // only does this object actions
1578                         where_is_pose(par);
1579                 }
1580                 par= par->parent;
1581         }
1582         
1583         return mcache;
1584 }
1585
1586 /* main particle building function 
1587    one day particles should become dynamic (realtime) with the current method as a 'bake' (ton) */
1588 void build_particle_system(Object *ob)
1589 {
1590         RNG *rng;
1591         PartEff *paf;
1592         Particle *pa;
1593         Mesh *me;
1594         Base *base;
1595         MTex *mtexmove=0, *mtextime=0;
1596         Material *ma;
1597         MFace *facelist= NULL;
1598         pMatrixCache *mcache=NULL, *mcnow, *mcprev;
1599         ListBase *effectorbase;
1600         VeNoCo *vertexcosnos;
1601         double startseconds= PIL_check_seconds_timer();
1602         float ftime, dtime, force[3], vec[3], fac, co[3], no[3];
1603         float *voweights= NULL, *foweights= NULL, maxw=1.0f;
1604         float *volengths= NULL, *folengths= NULL;
1605         int deform=0, a, totpart, paf_sta, paf_end;
1606         int waitcursor_set= 0, totvert, totface, curface, curvert;
1607         
1608         /* return conditions */
1609         if(ob->type!=OB_MESH) return;
1610         me= ob->data;
1611
1612         paf= give_parteff(ob);
1613         if(paf==NULL) return;
1614         if(paf->keys) MEM_freeN(paf->keys);     /* free as early as possible, for returns */
1615         paf->keys= NULL;
1616         
1617         if(paf->end < paf->sta) return;
1618         
1619         if( (paf->flag & PAF_OFACE) && (paf->flag & PAF_FACE)==0) return;
1620         
1621         if(me->totvert==0) return;
1622         
1623         if(ob==G.obedit) return;
1624         
1625         totpart= (R.flag & R_RENDERING)?paf->totpart:(paf->disp*paf->totpart)/100;
1626         if(totpart==0) return;
1627         
1628         /* No returns after this line! */
1629         
1630         /* material */
1631         ma= give_current_material(ob, paf->omat);
1632         if(ma) {
1633                 if(paf->speedtex)
1634                         mtexmove= ma->mtex[paf->speedtex-1];
1635                 mtextime= ma->mtex[paf->timetex-1];
1636         }
1637
1638         disable_speed_curve(1); /* check this... */
1639
1640         /* initialize particles */
1641         new_particle(paf);
1642
1643         /* reset deflector cache, sumohandle is free, but its still sorta abuse... (ton) */
1644         for(base= G.scene->base.first; base; base= base->next)
1645                 base->object->sumohandle= NULL;
1646
1647         /* all object positions from start to end */
1648         paf_sta= (int)floor(paf->sta);
1649         paf_end= (int)ceil(paf->end);
1650         if((paf->flag & PAF_STATIC)==0)
1651                 mcache= cache_object_matrices(ob, paf_sta, paf_end);
1652         
1653         /* mult generations? */
1654         for(a=0; a<PAF_MAXMULT; a++) {
1655                 if(paf->mult[a]!=0.0) {
1656                         /* interesting formula! this way after 'x' generations the total is paf->totpart */
1657                         totpart= (int)(totpart / (1.0+paf->mult[a]*paf->child[a]));
1658                 }
1659                 else break;
1660         }
1661
1662         /* for static particles, calculate system on current frame (? ton) */
1663         if(ma) do_mat_ipo(ma);
1664         
1665         /* matrix invert for static too */
1666         Mat4Invert(ob->imat, ob->obmat);
1667         
1668         /* new random generator */
1669         rng = rng_new(paf->seed);
1670         
1671         /* otherwise it goes way too fast */
1672         force[0]= paf->force[0]*0.05f;
1673         force[1]= paf->force[1]*0.05f;
1674         force[2]= paf->force[2]*0.05f;
1675         
1676         if( paf->flag & PAF_STATIC ) deform= 0;
1677         else {
1678                 deform= (ob->parent && ob->parent->type==OB_LATTICE && ob->partype==PARSKEL);
1679                 if(deform) init_latt_deform(ob->parent, 0);
1680         }
1681         
1682         /* get the effectors */
1683         effectorbase= pdInitEffectors(ob, paf->group);
1684         
1685         /* init geometry, return is 6 x float * me->totvert in size */
1686         vertexcosnos= (VeNoCo *)mesh_get_mapped_verts_nors(ob);
1687         facelist= me->mface;
1688         totvert= me->totvert;
1689         totface= me->totface;
1690         
1691         /* if vertexweights or even distribution, it makes weight tables, also checks where it emits from */
1692         make_weight_tables(paf, me, totpart, vertexcosnos, totvert, facelist, totface, &voweights, &foweights);
1693         
1694         /* vertexweights can define lengths too */
1695         make_length_tables(paf, me, totvert, facelist, totface, &volengths, &folengths);
1696         
1697         /* now define where to emit from, if there are face weights we skip vertices */
1698         if(paf->flag & PAF_OFACE) totvert= 0;
1699         if((paf->flag & PAF_FACE)==0) totface= 0;
1700         if(foweights) totvert= 0;
1701         
1702         /* initialize give_mesh_particle_coord */
1703         if(totface)
1704                 give_mesh_particle_coord(paf, vertexcosnos, facelist, totpart, totface, NULL, NULL);
1705         
1706         /* correction for face timing when using weighted average */
1707         if(totface && foweights) {
1708                 maxw= (paf->end-paf->sta)/foweights[0];
1709         }
1710         else if(totvert && voweights) {
1711                 maxw= (paf->end-paf->sta)/voweights[0];
1712         }
1713         
1714         /* for loop below */
1715         if (paf->flag & PAF_STATIC) {
1716                 ftime = G.scene->r.cfra;
1717                 dtime= 0.0f;
1718         } else {
1719                 ftime= paf->sta;
1720                 dtime= (paf->end - paf->sta)/(float)totpart;
1721         }
1722         
1723         curface= curvert= 0;
1724         for(a=0; a<totpart; a++, ftime+=dtime) {
1725                 
1726                 /* we set waitcursor only when a half second expired, particles now are realtime updated */
1727                 if(waitcursor_set==0 && (a % 256)==255) {
1728                         double seconds= PIL_check_seconds_timer();
1729                         if(seconds - startseconds > 0.5) {
1730                                 waitcursor(1);
1731                                 waitcursor_set= 1;
1732                         }
1733                 }
1734                 
1735                 pa= new_particle(paf);
1736                 pa->time= ftime;
1737
1738                 /* get coordinates from faces, only when vertices set to zero */
1739                 if(totvert==0 && totface) {
1740                         int curjit;
1741                         
1742                         /* use weight table, we have to do faces in order to be able to use jitter table... */
1743                         if(foweights) {
1744                                 
1745                                 if(foweights[curface] < 1.0f) {
1746                                         float remainder= 0.0f;
1747                                         
1748                                         while(remainder + foweights[curface] < 1.0f && curface<totface-1) {
1749                                                 remainder += foweights[curface];
1750                                                 curface++;
1751                                         }
1752                                         foweights[curface] += remainder;
1753                                         maxw= (paf->end-paf->sta)/foweights[curface];
1754                                 }
1755                                 
1756                                 if(foweights[curface]==0.0f)
1757                                         break;  /* WARN skips here out of particle generating */
1758                                 else {
1759                                         if(foweights[curface] > 1.0f)
1760                                                 foweights[curface] -= 1.0f;
1761                                         
1762                                         curjit= (int) foweights[curface];
1763                                         give_mesh_particle_coord(paf, vertexcosnos, facelist+curface, a, curjit, co, no);
1764                                         
1765                                         /* time correction to make particles appear evenly, maxw does interframe (0-1) */
1766                                         pa->time= paf->sta + maxw*foweights[curface];
1767                                 }
1768                         }
1769                         else {
1770                                 curface= a % totface;
1771                                 curjit= a/totface;
1772                                 give_mesh_particle_coord(paf, vertexcosnos, facelist+curface, a, curjit, co, no);
1773                         }
1774                 }
1775                 /* get coordinates from vertices */
1776                 if(totvert) {
1777                         /* use weight table */
1778                         if(voweights) {
1779                                 
1780                                 if(voweights[curvert] < 1.0f) {
1781                                         float remainder= 0.0f;
1782                                         
1783                                         while(remainder + voweights[curvert] < 1.0f && curvert<totvert-1) {
1784                                                 remainder += voweights[curvert];
1785                                                 curvert++;
1786                                         }
1787                                         voweights[curvert] += remainder;
1788                                         maxw= (paf->end-paf->sta)/voweights[curvert];
1789                                 }
1790                                 
1791                                 if(voweights[curvert]==0.0f)
1792                                         break;  /* WARN skips here out of particle generating */
1793                                 else {
1794                                         if(voweights[curvert] > 1.0f)
1795                                                 voweights[curvert] -= 1.0f;
1796                                         
1797                                         /* time correction to make particles appear evenly */
1798                                         pa->time= paf->sta + maxw*voweights[curvert];
1799                                 }
1800                         }
1801                         else {
1802                                 curvert= a % totvert;
1803                                 if(a >= totvert && totface)
1804                                         totvert= 0;
1805                         }
1806                         
1807                         VECCOPY(co, vertexcosnos[curvert].co);
1808                         VECCOPY(no, vertexcosnos[curvert].no);
1809                 }
1810                 
1811                 VECCOPY(pa->co, co);
1812                 
1813                 /* dynamic options */
1814                 if((paf->flag & PAF_STATIC)==0) {
1815                         int cur;
1816                         
1817                         /* particle retiming with texture */
1818                         if(mtextime && (paf->flag2 & PAF_TEXTIME)) {
1819                                 float tin, tr, tg, tb, ta, orco[3];
1820                                 
1821                                 /* calculate normalized orco */
1822                                 orco[0] = (co[0]-me->loc[0])/me->size[0];
1823                                 orco[1] = (co[1]-me->loc[1])/me->size[1];
1824                                 orco[2] = (co[2]-me->loc[2])/me->size[2];
1825                                 externtex(mtextime, orco, &tin, &tr, &tg, &tb, &ta);
1826                                 
1827                                 if(paf->flag2neg & PAF_TEXTIME)
1828                                         pa->time = paf->sta + (paf->end - paf->sta)*tin;
1829                                 else
1830                                         pa->time = paf->sta + (paf->end - paf->sta)*(1.0f-tin);
1831                         }
1832
1833                         /* set ob at correct time, we use cached matrices */
1834                         cur= (int)floor(pa->time) + 1 ;         /* + 1 has a reason: (obmat/prevobmat) otherwise comet-tails start too late */
1835                         
1836                         if(cur <= paf_end) mcnow= mcache + cur - paf_sta;
1837                         else mcnow= mcache + paf_end - paf_sta + 1;
1838                         
1839                         if(cur > paf_sta) mcprev= mcnow-1;
1840                         else mcprev= mcache;
1841                         
1842                         /* move to global space */
1843                         Mat4MulVecfl(mcnow->obmat, pa->co);
1844                 
1845                         VECCOPY(vec, co);
1846                         Mat4MulVecfl(mcprev->obmat, vec);
1847                         
1848                         /* first start speed: object */
1849                         VECSUB(pa->no, pa->co, vec);
1850
1851                         VecMulf(pa->no, paf->obfac);
1852                         
1853                         /* calculate the correct inter-frame */ 
1854                         fac= (pa->time- (float)floor(pa->time));
1855                         pa->co[0]= fac*pa->co[0] + (1.0f-fac)*vec[0];
1856                         pa->co[1]= fac*pa->co[1] + (1.0f-fac)*vec[1];
1857                         pa->co[2]= fac*pa->co[2] + (1.0f-fac)*vec[2];
1858
1859                         /* start speed: normal */
1860                         if(paf->normfac!=0.0) {
1861                                 /* imat is transpose ! */
1862                                 VECCOPY(vec, no);
1863                                 Mat3MulVecfl(mcnow->imat, vec);
1864                         
1865                                 Normalise(vec);
1866                                 VecMulf(vec, paf->normfac);
1867                                 VECADD(pa->no, pa->no, vec);
1868                         }
1869                 }
1870                 else {
1871                         if(paf->normfac!=0.0) {
1872                                 VECCOPY(pa->no, no);
1873                                 Normalise(pa->no);
1874                                 VecMulf(pa->no, paf->normfac);
1875                         }
1876                 }
1877                 
1878                 pa->lifetime= paf->lifetime;
1879                 if(paf->randlife!=0.0) {
1880                         pa->lifetime*= 1.0f + paf->randlife*(rng_getFloat(rng) - 0.5f);
1881                 }
1882                 pa->mat_nr= paf->omat;
1883                 
1884                 if(folengths)
1885                         pa->lifetime*= folengths[curface];
1886
1887                 make_particle_keys(rng, ob, 0, a, paf, pa, force, deform, mtexmove, effectorbase);
1888         }
1889         
1890         /* free stuff */
1891         give_mesh_particle_coord(NULL, NULL, NULL, 0, 0, NULL, NULL);
1892         MEM_freeN(vertexcosnos);
1893         if(voweights) MEM_freeN(voweights);
1894         if(foweights) MEM_freeN(foweights);
1895         if(volengths) MEM_freeN(volengths);
1896         if(folengths) MEM_freeN(folengths);
1897         if(mcache) MEM_freeN(mcache);
1898         rng_free(rng);
1899
1900         if(deform) end_latt_deform();
1901         
1902         if(effectorbase)
1903                 pdEndEffectors(effectorbase);   
1904         
1905         /* reset deflector cache */
1906         for(base= G.scene->base.first; base; base= base->next) {
1907                 if(base->object->sumohandle) {
1908                         MEM_freeN(base->object->sumohandle);
1909                         base->object->sumohandle= NULL;
1910                 }
1911         }
1912
1913         disable_speed_curve(0);
1914         
1915         if(waitcursor_set) waitcursor(0);
1916 }
1917