New hair child options:

* Renamed children to "simple" and "interpolated" as this is
  easier to explain and more descriptive than "from particles"
  and "from faces".
* Also shuffled the child ui around a bit to make it clearer.
* Child seed parameter allows to change the seed for children
  independent of the main seed value.
* Long hair mode for interpolated children:
  - Making even haircuts was impossible before as the child
    strand lengths were even, but their root coordinates were
    not similar in relation to the parent strands.
  - The "long hair" option uses the tips of the parent strands
    to calculate the child strand tips.
* Hair parting options:
  - Hair parting can now be calculated dynamically on the fly
    when in 2.49 there was a cumbersome way of using emitter mesh
    seams to define parting lines.
  - For long hair parting can be created by a tip distance/root
    distance threshold. For example setting the minimum threshold
    to 2.0 creates partings between children belonging to parents
    with tip distance of three times the root distance
    ((1+2)*root distance).
  - For short hair the parting thresholds are used as angles
    between the root directions.
* New kink parameters:
  - Kink flatness calculates kink into a shape that would have
    been achieved with an actual curling iron.
  - Kink amplitude clump determines how much the main clump value
    effects the kink amplitude.
  - The beginning of kink is now smoothed to make the hair look
    more natural close to the roots.
* Some bugs fixed along the way too:
  - Child parent's were not determined correctly in some cases.
  - Children didn't always look correct in particle mode.
  - Changing child parameters caused actual particles to be
    recalculated.
* Also cleaned up some deprecated code.

All in all there should be no real changes to how old files look
(except perhaps a bit better!), but the new options should make
hair/fur creation a bit more enjoyable. I'll try to make a video
demonstrating the new stuff shortly.

1 files changed, 401 insertions, 391 deletions

diff --git a/source/blender/blenkernel/intern/particle.c b/source/blender/blenkernel/intern/particle.c
index 95e14542e28..d3a0cd5df7c 100644
--- a/source/blender/blenkernel/intern/particle.c
+++ b/source/blender/blenkernel/intern/particle.c
@@ -1862,144 +1862,184 @@ static float vert_weight(MDeformVert *dvert, int group)
 	return 0.0;
 }
 
-static void do_prekink(ParticleKey *state, ParticleKey *par, float *par_rot, float time, float freq, float shape, float amplitude, short type, short axis, float obmat[][4])
+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], int smooth_start)
 {
-	float vec[3]={0.0,0.0,0.0}, q1[4]={1,0,0,0},q2[4];
-	float t;
+	float kink[3]={1.f,0.f,0.f}, par_vec[3], q1[4]={1.f,0.f,0.f,0.f};
+	float t, dt=1.f, result[3];
 
-	CLAMP(time,0.0,1.0);
+	if(par == NULL || type == PART_KINK_NO)
+		return;
+
+	CLAMP(time, 0.f, 1.f);
 
 	if(shape!=0.0f && type!=PART_KINK_BRAID) {
 		if(shape<0.0f)
-			time= (float)pow(time, 1.0+shape);
+			time= (float)pow(time, 1.f+shape);
 		else
-			time= (float)pow(time, 1.0/(1.0-shape));
+			time= (float)pow(time, 1.f/(1.f-shape));
 	}
 
-	t=time;
+	t = time * freq *(float)M_PI;
+	
+	if(smooth_start) {
+		dt = fabs(t);
+		/* smooth the beginning of kink */
+		CLAMP(dt, 0.f, (float)M_PI);
+		dt = sin(dt/2.f);
+	}
 
-	t*=(float)M_PI*freq;
+	if(type != PART_KINK_RADIAL) {
+		float temp[3];
 
-	if(par==0) return;
+		kink[axis]=1.f;
 
-	switch(type){
-		case PART_KINK_CURL:
-			vec[axis]=1.0;
-			if(par_rot)
-				QUATCOPY(q2,par_rot)
-			else
-				vec_to_quat( q2,par->vel,axis,(axis+1)%3);
-			mul_qt_v3(q2,vec);
-			mul_v3_fl(vec,amplitude);
-			VECADD(state->co,state->co,vec);
+		if(obmat)
+			mul_mat3_m4_v3(obmat, kink);
+		
+		if(par_rot)
+			mul_qt_v3(par_rot, kink);
 
-			VECSUB(vec,state->co,par->co);
+		/* make sure kink is normal to strand */
+		project_v3_v3v3(temp, kink, par->vel);
+		sub_v3_v3(kink, temp);
+		normalize_v3(kink);
+	}
 
-			if(t!=0.0)
-				axis_angle_to_quat(q1,par->vel,t);
-			
-			mul_qt_v3(q1,vec);
-			
-			VECADD(state->co,par->co,vec);
-			break;
-		case PART_KINK_RADIAL:
-			VECSUB(vec,state->co,par->co);
+	copy_v3_v3(result, state->co);
+	sub_v3_v3v3(par_vec, par->co, state->co);
 
-			normalize_v3(vec);
-			mul_v3_fl(vec,amplitude*(float)sin(t));
+	switch(type) {
+	case PART_KINK_CURL:
+	{
+		mul_v3_fl(par_vec, -1.f);
 
-			VECADD(state->co,state->co,vec);
-			break;
-		case PART_KINK_WAVE:
-			vec[axis]=1.0;
-			if(obmat)
-				mul_mat3_m4_v3(obmat,vec);
+		if(flat > 0.f) {
+			float proj[3];
+			project_v3_v3v3(proj, par_vec, par->vel);
+			madd_v3_v3fl(par_vec, proj, -flat);
 
-			if(par_rot)
-				mul_qt_v3(par_rot,vec);
+			project_v3_v3v3(proj, par_vec, kink);
+			madd_v3_v3fl(par_vec, proj, -flat);
+		}
 
-			project_v3_v3v3(q1,vec,par->vel);
-			
-			VECSUB(vec,vec,q1);
-			normalize_v3(vec);
+		axis_angle_to_quat(q1, kink, (float)M_PI/2.f);
 
-			mul_v3_fl(vec,amplitude*(float)sin(t));
+		mul_qt_v3(q1, par_vec);
 
-			VECADD(state->co,state->co,vec);
-			break;
-		case PART_KINK_BRAID:
-			if(par){
-				float y_vec[3]={0.0,1.0,0.0};
-				float z_vec[3]={0.0,0.0,1.0};
-				float vec_from_par[3], vec_one[3], radius, state_co[3];
-				float inp_y,inp_z,length;
-				
-				if(par_rot)
-					QUATCOPY(q2,par_rot)
-				else
-					vec_to_quat(q2,par->vel,axis,(axis+1)%3);
-				mul_qt_v3(q2,y_vec);
-				mul_qt_v3(q2,z_vec);
-				
-				VECSUB(vec_from_par,state->co,par->co);
-				radius= normalize_v3_v3(vec_one, vec_from_par);
+		madd_v3_v3fl(par_vec, kink, amplitude);
 
-				inp_y=dot_v3v3(y_vec,vec_one);
-				inp_z=dot_v3v3(z_vec,vec_one);
+		/* rotate kink vector around strand tangent */
+		if(t!=0.f) {
+			axis_angle_to_quat(q1, par->vel, t);
+			mul_qt_v3(q1, par_vec);
+		}
 
-				if(inp_y>0.5){
-					VECCOPY(state_co,y_vec);
+		add_v3_v3v3(result, par->co, par_vec);
+		break;
+	}
+	case PART_KINK_RADIAL:
+	{
+		if(flat > 0.f) {
+			float proj[3];
+			/* flatten along strand */
+			project_v3_v3v3(proj, par_vec, par->vel);
+			madd_v3_v3fl(result, proj, flat);
+		}
 
-					mul_v3_fl(y_vec,amplitude*(float)cos(t));
-					mul_v3_fl(z_vec,amplitude/2.0f*(float)sin(2.0f*t));
-				}
-				else if(inp_z>0.0){
-					VECCOPY(state_co,z_vec);
-					mul_v3_fl(state_co,(float)sin(M_PI/3.0f));
-					VECADDFAC(state_co,state_co,y_vec,-0.5f);
+		madd_v3_v3fl(result, par_vec, -amplitude*(float)sin(t));
+		break;
+	}
+	case PART_KINK_WAVE:
+	{
+		madd_v3_v3fl(result, kink, amplitude*(float)sin(t));
 
-					mul_v3_fl(y_vec,-amplitude*(float)cos(t + M_PI/3.0f));
-					mul_v3_fl(z_vec,amplitude/2.0f*(float)cos(2.0f*t + M_PI/6.0f));
-				}
-				else{
-					VECCOPY(state_co,z_vec);
-					mul_v3_fl(state_co,-(float)sin(M_PI/3.0f));
-					VECADDFAC(state_co,state_co,y_vec,-0.5f);
+		if(flat > 0.f) {
+			float proj[3];
+			/* flatten along wave */
+			project_v3_v3v3(proj, par_vec, kink);
+			madd_v3_v3fl(result, proj, flat);
 
-					mul_v3_fl(y_vec,amplitude*(float)-sin(t+M_PI/6.0f));
-					mul_v3_fl(z_vec,amplitude/2.0f*(float)-sin(2.0f*t+M_PI/3.0f));
-				}
+			/* flatten along strand */
+			project_v3_v3v3(proj, par_vec, par->vel);
+			madd_v3_v3fl(result, proj, flat);
+		}
+		break;
+	}
+	case PART_KINK_BRAID:
+	{
+		float y_vec[3]={0.f,1.f,0.f};
+		float z_vec[3]={0.f,0.f,1.f};
+		float vec_one[3], radius, state_co[3];
+		float inp_y, inp_z, length;
 
-				mul_v3_fl(state_co,amplitude);
-				VECADD(state_co,state_co,par->co);
-				VECSUB(vec_from_par,state->co,state_co);
+		mul_qt_v3(par_rot, y_vec);
+		mul_qt_v3(par_rot, z_vec);
+		
+		mul_v3_fl(par_vec, -1.f);
+		radius= normalize_v3_v3(vec_one, par_vec);
 
-				length=normalize_v3(vec_from_par);
-				mul_v3_fl(vec_from_par,MIN2(length,amplitude/2.0f));
+		inp_y=dot_v3v3(y_vec, vec_one);
+		inp_z=dot_v3v3(z_vec, vec_one);
 
-				VECADD(state_co,par->co,y_vec);
-				VECADD(state_co,state_co,z_vec);
-				VECADD(state_co,state_co,vec_from_par);
+		if(inp_y>0.5){
+			copy_v3_v3(state_co, y_vec);
 
-				shape=(2.0f*(float)M_PI)*(1.0f+shape);
+			mul_v3_fl(y_vec, amplitude*(float)cos(t));
+			mul_v3_fl(z_vec, amplitude/2.f*(float)sin(2.f*t));
+		}
+		else if(inp_z>0.0){
+			mul_v3_v3fl(state_co, z_vec, (float)sin(M_PI/3.f));
+			VECADDFAC(state_co,state_co,y_vec,-0.5f);
 
-				if(t<shape){
-					shape=t/shape;
-					shape=(float)sqrt((double)shape);
-					interp_v3_v3v3(state->co,state->co,state_co,shape);
-				}
-				else{
-					VECCOPY(state->co,state_co);
-				}
-			}
-			break;
+			mul_v3_fl(y_vec, -amplitude * (float)cos(t + M_PI/3.f));
+			mul_v3_fl(z_vec, amplitude/2.f * (float)cos(2.f*t + M_PI/6.f));
+		}
+		else{
+			mul_v3_v3fl(state_co, z_vec, -(float)sin(M_PI/3.f));
+			madd_v3_v3fl(state_co, y_vec, -0.5f);
+
+			mul_v3_fl(y_vec, amplitude * (float)-sin(t + M_PI/6.f));
+			mul_v3_fl(z_vec, amplitude/2.f * (float)-sin(2.f*t + M_PI/3.f));
+		}
+
+		mul_v3_fl(state_co, amplitude);
+		add_v3_v3(state_co, par->co);
+		sub_v3_v3v3(par_vec, state->co, state_co);
+
+		length = normalize_v3(par_vec);
+		mul_v3_fl(par_vec, MIN2(length, amplitude/2.f));
+
+		add_v3_v3v3(state_co, par->co, y_vec);
+		add_v3_v3(state_co, z_vec);
+		add_v3_v3(state_co, par_vec);
+
+		shape = 2.f*(float)M_PI * (1.f+shape);
+
+		if(t<shape){
+			shape = t/shape;
+			shape = (float)sqrt((double)shape);
+			interp_v3_v3v3(result, result, state_co, shape);
+		}
+		else{
+			copy_v3_v3(result, state_co);
+		}
+		break;
 	}
+	}
+
+	/* blend the start of the kink */
+	if(dt < 1.f)
+		interp_v3_v3v3(state->co, state->co, result, dt);
+	else
+		copy_v3_v3(state->co, result);
 }
 
-static void do_clump(ParticleKey *state, ParticleKey *par, float time, float clumpfac, float clumppow, float pa_clump)
+static float do_clump(ParticleKey *state, ParticleKey *par, float time, float clumpfac, float clumppow, float pa_clump)
 {
+	float clump = 0.f;
+
 	if(par && clumpfac!=0.0){
-		float clump, cpow;
+		float cpow;
 
 		if(clumppow<0.0)
 			cpow=1.0f+clumppow;
@@ -2010,8 +2050,11 @@ static void do_clump(ParticleKey *state, ParticleKey *par, float time, float clu
 			clump = -clumpfac*pa_clump*(float)pow(1.0-(double)time,(double)cpow);
 		else
 			clump = clumpfac*pa_clump*(float)pow((double)time,(double)cpow);
+
 		interp_v3_v3v3(state->co,state->co,par->co,clump);
 	}
+
+	return clump;
 }
 void precalc_guides(ParticleSimulationData *sim, ListBase *effectors)
 {
@@ -2124,7 +2167,7 @@ int do_guides(ListBase *effectors, ParticleKey *state, int index, float time)
 		}
 		par.co[0] = par.co[1] = par.co[2] = 0.0f;
 		VECCOPY(key.co, vec_to_point);
-		do_prekink(&key, &par, 0, guidetime, pd->kink_freq, pd->kink_shape, pd->kink_amp, pd->kink, pd->kink_axis, 0);
+		do_kink(&key, &par, 0, guidetime, pd->kink_freq, pd->kink_shape, pd->kink_amp, 0.f, pd->kink, pd->kink_axis, 0, 0);
 		do_clump(&key, &par, guidetime, pd->clump_fac, pd->clump_pow, 1.0f);
 		VECCOPY(vec_to_point, key.co);
 
@@ -2229,20 +2272,23 @@ static int check_path_length(int k, ParticleCacheKey *keys, ParticleCacheKey *st
 		return k;
 	}
 }
-static void offset_child(ChildParticle *cpa, ParticleKey *par, ParticleKey *child, float flat, float radius)
+static void offset_child(ChildParticle *cpa, ParticleKey *par, float *par_rot, ParticleKey *child, float flat, float radius)
 {
-	VECCOPY(child->co,cpa->fuv);
-	mul_v3_fl(child->co,radius);
+	copy_v3_v3(child->co, cpa->fuv);
+	mul_v3_fl(child->co, radius);
 
 	child->co[0]*=flat;
 
-	VECCOPY(child->vel,par->vel);
-
-	mul_qt_v3(par->rot,child->co);
+	copy_v3_v3(child->vel, par->vel);
 
-	QUATCOPY(child->rot,par->rot);
+	if(par_rot) {
+		mul_qt_v3(par_rot, child->co);
+		copy_qt_qt(child->rot, par_rot);
+	}
+	else
+		unit_qt(child->rot);
 
-	VECADD(child->co,child->co,par->co);
+	add_v3_v3(child->co, par->co);
 }
 float *psys_cache_vgroup(DerivedMesh *dm, ParticleSystem *psys, int vgroup)
 {
@@ -2372,12 +2418,9 @@ static int psys_threads_init_path(ParticleThread *threads, Scene *scene, float c
 	if(totchild==0) return 0;
 
 	/* init random number generator */
-	if(ctx->sim.psys->part->flag & PART_ANIM_BRANCHING)
-		seed= 31415926 + ctx->sim.psys->seed + (int)cfra;
-	else
-		seed= 31415926 + ctx->sim.psys->seed;
+	seed= 31415926 + ctx->sim.psys->seed;
 	
-	if(part->flag & PART_BRANCHING || ctx->editupdate || totchild < 10000)
+	if(ctx->editupdate || totchild < 10000)
 		totthread= 1;
 	
 	for(i=0; i<totthread; i++) {
@@ -2420,7 +2463,7 @@ static int psys_threads_init_path(ParticleThread *threads, Scene *scene, float c
 }
 
 /* note: this function must be thread safe, except for branching! */
-static void psys_thread_create_path(ParticleThread *thread, struct ChildParticle *cpa, ParticleCacheKey *keys, int i)
+static void psys_thread_create_path(ParticleThread *thread, struct ChildParticle *cpa, ParticleCacheKey *child_keys, int i)
 {
 	ParticleThreadContext *ctx= thread->ctx;
 	Object *ob= ctx->sim.ob;
@@ -2428,44 +2471,29 @@ static void psys_thread_create_path(ParticleThread *thread, struct ChildParticle
 	ParticleSettings *part = psys->part;
 	ParticleCacheKey **cache= psys->childcache;
 	ParticleCacheKey **pcache= psys_in_edit_mode(ctx->sim.scene, psys) ? psys->edit->pathcache : psys->pathcache;
-	ParticleCacheKey *state, *par = NULL, *key[4];
-	ParticleData *pa=NULL;
+	ParticleCacheKey *child, *par = NULL, *key[4];
 	ParticleTexture ptex;
-	float *cpa_fuv=0, *par_rot=0;
-	float co[3], orco[3], ornor[3], hairmat[4][4], t, cpa_1st[3], dvec[3];
-	float branch_begin, branch_end, branch_prob, rough_rand;
+	float *cpa_fuv=0, *par_rot=0, rot[4];
+	float orco[3], ornor[3], hairmat[4][4], t, dvec[3], off1[4][3], off2[4][3];
 	float length, max_length = 1.0f, cur_length = 0.0f;
-	float eff_length, eff_vec[3];
-	int k, cpa_num;
+	float eff_length, eff_vec[3], weight[4];
+	int k, cpa_num, maxw=0;
 	short cpa_from;
 
 	if(!pcache)
 		return;
 
-	if(part->flag & PART_BRANCHING) {
-		branch_begin=rng_getFloat(thread->rng_path);
-		branch_end=branch_begin+(1.0f-branch_begin)*rng_getFloat(thread->rng_path);
-		branch_prob=rng_getFloat(thread->rng_path);
-		rough_rand=rng_getFloat(thread->rng_path);
-	}
-	else {
-		branch_begin= 0.0f;
-		branch_end= 0.0f;
-		branch_prob= 0.0f;
-		rough_rand= 0.0f;
-	}
-
-	if(i<psys->totpart){
-		branch_begin=0.0f;
-		branch_end=1.0f;
-		branch_prob=0.0f;
-	}
-
 	if(ctx->between){
+		ParticleData *pa = psys->particles + cpa->pa[0];
 		int w, needupdate;
-		float foffset;
-
-		if(ctx->editupdate && !(part->flag & PART_BRANCHING)) {
+		float foffset, wsum=0.f;
+		float co[3];
+		float p_min = part->parting_min;
+		float p_max = part->parting_max;
+		/* Virtual parents don't work nicely with parting. */
+		float p_fac = part->parents > 0.f ? 0.f : part->parting_fac;
+
+		if(ctx->editupdate) {
 			needupdate= 0;
 			w= 0;
 			while(w<4 && cpa->pa[w]>=0) {
@@ -2479,223 +2507,223 @@ static void psys_thread_create_path(ParticleThread *thread, struct ChildParticle
 			if(!needupdate)
 				return;
 			else
-				memset(keys, 0, sizeof(*keys)*(ctx->steps+1));
+				memset(child_keys, 0, sizeof(*child_keys)*(ctx->steps+1));
 		}
 
 		/* get parent paths */
-		w= 0;
-		while(w<4 && cpa->pa[w]>=0){
-			key[w] = pcache[cpa->pa[w]];
-			w++;
+		for(w=0; w<4; w++) {
+			if(cpa->pa[w] >= 0) {
+				key[w] = pcache[cpa->pa[w]];
+				weight[w] = cpa->w[w];
+			}
+			else {
+				key[w] = pcache[0];
+				weight[w] = 0.f;
+			}
+		}
+
+		/* modify weights to create parting */
+		if(p_fac > 0.f) {
+			for(w=0; w<4; w++) {
+				if(w && weight[w] > 0.f) {
+					float d;
+					if(part->flag & PART_CHILD_LONG_HAIR) {
+						/* For long hair use tip distance/root distance as parting factor instead of root to tip angle. */
+						float d1 = len_v3v3(key[0]->co, key[w]->co);
+						float d2 = len_v3v3((key[0]+key[0]->steps-1)->co, (key[w]+key[w]->steps-1)->co);
+
+						d = d1 > 0.f ? d2/d1 - 1.f : 10000.f;
+					}
+					else {
+						float v1[3], v2[3];
+						sub_v3_v3v3(v1, (key[0]+key[0]->steps-1)->co, key[0]->co);
+						sub_v3_v3v3(v2, (key[w]+key[w]->steps-1)->co, key[w]->co);
+						normalize_v3(v1);
+						normalize_v3(v2);
+
+						d = saacos(dot_v3v3(v1, v2)) * 180.f / M_PI;
+					}
+
+					if(p_max > p_min)
+						d = (d - p_min)/(p_max - p_min);
+					else
+						d = (d - p_min) <= 0.f ? 0.f : 1.f;
+
+					CLAMP(d, 0.f, 1.f);
+
+					if(d > 0.f)
+						weight[w] *= (1.f - d);
+				}
+				wsum += weight[w];
+			}
+			for(w=0; w<4; w++)
+				weight[w] /= wsum;
+
+			interp_v4_v4v4(weight, cpa->w, weight, p_fac);
 		}
 
 		/* get the original coordinates (orco) for texture usage */
 		cpa_num = cpa->num;
 		
-		foffset= cpa->foffset;
+		foffset = cpa->foffset;
 		cpa_fuv = cpa->fuv;
 		cpa_from = PART_FROM_FACE;
 
 		psys_particle_on_emitter(ctx->sim.psmd,cpa_from,cpa_num,DMCACHE_ISCHILD,cpa->fuv,foffset,co,ornor,0,0,orco,0);
 
-		if(part->path_start==0.0f) {
-			/* we need to save the actual root position of the child for positioning it accurately to the surface of the emitter */
-			VECCOPY(cpa_1st,co);
-			mul_m4_v3(ob->obmat,cpa_1st);
-		}
+		mul_m4_v3(ob->obmat, co);
 
-		pa = psys->particles + cpa->pa[0];
+		for(w=0; w<4; w++)
+			sub_v3_v3v3(off1[w], co, key[w]->co);
 
 		psys_mat_hair_to_global(ob, ctx->sim.psmd->dm, psys->part->from, pa, hairmat);
-
-		pa=0;
 	}
 	else{
-		if(ctx->editupdate && !(part->flag & PART_BRANCHING)) {
+		ParticleData *pa = psys->particles + cpa->parent;
+		float co[3];
+		if(ctx->editupdate) {
 			if(!(psys->edit->points[cpa->parent].flag & PEP_EDIT_RECALC))
 				return;
 
-			memset(keys, 0, sizeof(*keys)*(ctx->steps+1));
+			memset(child_keys, 0, sizeof(*child_keys)*(ctx->steps+1));
 		}
 
 		/* get the parent path */
-		key[0]=pcache[cpa->parent];
+		key[0] = pcache[cpa->parent];
 
 		/* get the original coordinates (orco) for texture usage */
-		pa=psys->particles+cpa->parent;
-
-		cpa_from=part->from;
-		cpa_num=pa->num;
-		cpa_fuv=pa->fuv;
+		cpa_from = part->from;
+		cpa_num = pa->num;
+		cpa_fuv = pa->fuv;
 
 		psys_particle_on_emitter(ctx->sim.psmd,cpa_from,cpa_num,DMCACHE_ISCHILD,cpa_fuv,pa->foffset,co,ornor,0,0,orco,0);
 
 		psys_mat_hair_to_global(ob, ctx->sim.psmd->dm, psys->part->from, pa, hairmat);
 	}
 
-	keys->steps = ctx->steps;
-
-	/* correct child ipo timing */
-#if 0 // XXX old animation system
-	if((part->flag&PART_ABS_TIME)==0 && part->ipo){
-		float dsta=part->end-part->sta;
-		calc_ipo(part->ipo, 100.0f*(ctx->cfra-(part->sta+dsta*cpa->rand[1]))/(part->lifetime*(1.0f - part->randlife*cpa->rand[0])));
-		execute_ipo((ID *)part, part->ipo);
-	}
-#endif // XXX old animation system
+	child_keys->steps = ctx->steps;
 
 	/* get different child parameters from textures & vgroups */
 	get_child_modifier_parameters(part, ctx, cpa, cpa_from, cpa_num, cpa_fuv, orco, &ptex);
 
 	if(ptex.exist < PSYS_FRAND(i + 24)) {
-		keys->steps = -1;
+		child_keys->steps = -1;
 		return;
 	}
 
 	/* create the child path */
-	for(k=0,state=keys; k<=ctx->steps; k++,state++){
+	for(k=0,child=child_keys; k<=ctx->steps; k++,child++){
 		if(ctx->between){
 			int w=0;
 
-			state->co[0] = state->co[1] = state->co[2] = 0.0f;
-			state->vel[0] = state->vel[1] = state->vel[2] = 0.0f;
-			state->rot[0] = state->rot[1] = state->rot[2] = state->rot[3] = 0.0f;
+			zero_v3(child->co);
+			zero_v3(child->vel);
+			unit_qt(child->rot);
 
-			//QUATCOPY(state->rot,key[0]->rot);
+			for(w=0; w<4; w++) {
+				copy_v3_v3(off2[w], off1[w]);
 
-			/* child position is the weighted sum of parent positions */
-			while(w<4 && cpa->pa[w]>=0){
-				state->co[0] += cpa->w[w] * key[w]->co[0];
-				state->co[1] += cpa->w[w] * key[w]->co[1];
-				state->co[2] += cpa->w[w] * key[w]->co[2];
-
-				state->vel[0] += cpa->w[w] * key[w]->vel[0];
-				state->vel[1] += cpa->w[w] * key[w]->vel[1];
-				state->vel[2] += cpa->w[w] * key[w]->vel[2];
-				key[w]++;
-				w++;
-			}
-			if(part->path_start==0.0f) {
-				if(k==0){
-					/* calculate the offset between actual child root position and first position interpolated from parents */
-					VECSUB(cpa_1st,cpa_1st,state->co);
+				if(part->flag & PART_CHILD_LONG_HAIR) {
+					/* Use parent rotation (in addition to emission location) to determine child offset. */
+					if(k)
+						mul_qt_v3((key[w]+k)->rot, off2[w]);
+
+					/* Fade the effect of rotation for even lengths in the end */
+					project_v3_v3v3(dvec, off2[w], (key[w]+k)->vel);
+					madd_v3_v3fl(off2[w], dvec, -(float)k/(float)ctx->steps);
 				}
-				/* apply offset for correct positioning */
-				VECADD(state->co,state->co,cpa_1st);
+
+				add_v3_v3(off2[w], (key[w]+k)->co);
 			}
+
+			/* child position is the weighted sum of parent positions */
+			interp_v3_v3v3v3v3(child->co, off2[0], off2[1], off2[2], off2[3], weight);
+			interp_v3_v3v3v3v3(child->vel, (key[0]+k)->vel, (key[1]+k)->vel, (key[2]+k)->vel, (key[3]+k)->vel, weight);
+
+			copy_qt_qt(child->rot, (key[0]+k)->rot);
 		}
 		else{
+			if(k) {
+				mul_qt_qtqt(rot, (key[0]+k)->rot, key[0]->rot);
+				par_rot = rot;
+			}
+			else {
+				par_rot = key[0]->rot;
+			}
 			/* offset the child from the parent position */
-			offset_child(cpa, (ParticleKey*)key[0], (ParticleKey*)state, part->childflat, part->childrad);
-
-			key[0]++;
+			offset_child(cpa, (ParticleKey*)(key[0]+k), par_rot, (ParticleKey*)child, part->childflat, part->childrad);
 		}
 	}
 
 	/* apply effectors */
 	if(part->flag & PART_CHILD_EFFECT) {
-		for(k=0,state=keys; k<=ctx->steps; k++,state++) {
+		for(k=0,child=child_keys; k<=ctx->steps; k++,child++) {
 			if(k) {
-				do_path_effectors(&ctx->sim, cpa->pa[0], state, k, ctx->steps, keys->co, ptex.effector, 0.0f, ctx->cfra, &eff_length, eff_vec);
+				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);
 			}
 			else {
-				sub_v3_v3v3(eff_vec,(state+1)->co,state->co);
-				eff_length= len_v3(eff_vec);
+				sub_v3_v3v3(eff_vec, (child+1)->co, child->co);
+				eff_length = len_v3(eff_vec);
 			}
 		}
 	}
 
-	for(k=0,state=keys; k<=ctx->steps; k++,state++){
-		t=(float)k/(float)ctx->steps;
+	for(k=0,child=child_keys; k<=ctx->steps; k++,child++){
+		t = (float)k/(float)ctx->steps;
+
+		if(ctx->totparent)
+			/* this is now threadsafe, virtual parents are calculated before rest of children */
+			par = (i >= ctx->totparent) ? cache[cpa->parent] : NULL;
+		else if(cpa->parent >= 0)
+			par = pcache[cpa->parent];
 
-		if(ctx->totparent){
-			if(i>=ctx->totparent) {
-				/* this is now threadsafe, virtual parents are calculated before rest of children */
-				par = cache[cpa->parent] + k;
+		if(par) {
+			if(k) {
+				mul_qt_qtqt(rot, (par+k)->rot, par->rot);
+				par_rot = rot;
 			}
-			else
-				par=0;
-		}
-		else if(cpa->parent>=0){
-			par=pcache[cpa->parent]+k;
-			par_rot = par->rot;
+			else {
+				par_rot = par->rot;
+			}
+			par += k;
 		}
 
 		/* apply different deformations to the child path */
-		do_child_modifiers(&ctx->sim, &ptex, (ParticleKey *)par, par_rot, cpa, orco, hairmat, (ParticleKey *)state, t);
-
-		/* TODO: better branching */
-		//if(part->flag & PART_BRANCHING && ctx->between == 0 && part->flag & PART_ANIM_BRANCHING)
-		//	rough_t = t * rough_rand;
-		//else
-		//	rough_t = t;
-
-		/* TODO: better branching */
-		//if(part->flag & PART_BRANCHING && ctx->between==0){
-		//	if(branch_prob > part->branch_thres){
-		//		branchfac=0.0f;
-		//	}
-		//	else{
-		//		if(part->flag & PART_SYMM_BRANCHING){
-		//			if(t < branch_begin || t > branch_end)
-		//				branchfac=0.0f;
-		//			else{
-		//				if((t-branch_begin)/(branch_end-branch_begin)<0.5)
-		//					branchfac=2.0f*(t-branch_begin)/(branch_end-branch_begin);
-		//				else
-		//					branchfac=2.0f*(branch_end-t)/(branch_end-branch_begin);
-
-		//				CLAMP(branchfac,0.0f,1.0f);
-		//			}
-		//		}
-		//		else{
-		//			if(t < branch_begin){
-		//				branchfac=0.0f;
-		//			}
-		//			else{
-		//				branchfac=(t-branch_begin)/((1.0f-branch_begin)*0.5f);
-		//				CLAMP(branchfac,0.0f,1.0f);
-		//			}
-		//		}
-		//	}
-
-		//	if(i<psys->totpart)
-		//		interp_v3_v3v3(state->co, (pcache[i] + k)->co, state->co, branchfac);
-		//	else
-		//		/* this is not threadsafe, but should only happen for
-		//		 * branching particles particles, which are not threaded */
-		//		interp_v3_v3v3(state->co, (cache[i - psys->totpart] + k)->co, state->co, branchfac);
-		//}
+		do_child_modifiers(&ctx->sim, &ptex, (ParticleKey *)par, par_rot, cpa, orco, hairmat, (ParticleKey *)child, t);
 
 		/* we have to correct velocity because of kink & clump */
 		if(k>1){
-			VECSUB((state-1)->vel,state->co,(state-2)->co);
-			mul_v3_fl((state-1)->vel,0.5);
+			sub_v3_v3v3((child-1)->vel, child->co, (child-2)->co);
+			mul_v3_fl((child-1)->vel, 0.5);
 
 			if(ctx->ma && (part->draw & PART_DRAW_MAT_COL))
-				get_strand_normal(ctx->ma, ornor, cur_length, (state-1)->vel);
+				get_strand_normal(ctx->ma, ornor, cur_length, (child-1)->vel);
 		}
 
 		if(k == ctx->steps)
-			VECSUB(state->vel,state->co,(state-1)->co);
+			sub_v3_v3v3(child->vel, child->co, (child-1)->co);
 
 		/* check if path needs to be cut before actual end of data points */
 		if(k){
-			VECSUB(dvec,state->co,(state-1)->co);
-			length=1.0f/(float)ctx->steps;
-			k=check_path_length(k,keys,state,max_length,&cur_length,length,dvec);
+			sub_v3_v3v3(dvec, child->co, (child-1)->co);
+			length = 1.0f/(float)ctx->steps;
+			k = check_path_length(k, child_keys, child, max_length, &cur_length, length, dvec);
 		}
 		else{
 			/* initialize length calculation */
-			max_length= ptex.length;
-			cur_length= 0.0f;
+			max_length = ptex.length;
+			cur_length = 0.0f;
 		}
 
 		if(ctx->ma && (part->draw & PART_DRAW_MAT_COL)) {
-			VECCOPY(state->col, &ctx->ma->r)
-			get_strand_normal(ctx->ma, ornor, cur_length, state->vel);
+			VECCOPY(child->col, &ctx->ma->r)
+			get_strand_normal(ctx->ma, ornor, cur_length, child->vel);
 		}
 	}
+
+	/* Hide virtual parents */
+	if(i < ctx->totparent)
+		child_keys->steps = -1;
 }
 
 static void *exec_child_path_cache(void *data)
@@ -2742,7 +2770,7 @@ void psys_cache_child_paths(ParticleSimulationData *sim, float cfra, int editupd
 	totchild= ctx->totchild;
 	totparent= ctx->totparent;
 
-	if(editupdate && sim->psys->childcache && !(part->flag & PART_BRANCHING) && totchild == sim->psys->totchildcache) {
+	if(editupdate && sim->psys->childcache && totchild == sim->psys->totchildcache) {
 		cache = sim->psys->childcache;
 	}
 	else {
@@ -2783,6 +2811,43 @@ void psys_cache_child_paths(ParticleSimulationData *sim, float cfra, int editupd
 
 	psys_threads_free(pthreads);
 }
+/* figure out incremental rotations along path starting from unit quat */
+static void cache_key_incremental_rotation(ParticleCacheKey *key0, ParticleCacheKey *key1, ParticleCacheKey *key2, float *prev_tangent, int i)
+{
+	float cosangle, angle, tangent[3], normal[3], q[4];
+
+	switch(i) {
+	case 0:
+		/* start from second key */
+		break;
+	case 1:
+		/* calculate initial tangent for incremental rotations */
+		sub_v3_v3v3(prev_tangent, key0->co, key1->co);
+		normalize_v3(prev_tangent);
+		unit_qt(key1->rot);
+		break;
+	default:
+		sub_v3_v3v3(tangent, key0->co, key1->co);
+		normalize_v3(tangent);
+
+		cosangle= dot_v3v3(tangent, prev_tangent);
+
+		/* note we do the comparison on cosangle instead of
+		* angle, since floating point accuracy makes it give
+		* different results across platforms */
+		if(cosangle > 0.999999f) {
+			QUATCOPY(key1->rot, key2->rot);
+		}
+		else {
+			angle= saacos(cosangle);
+			cross_v3_v3v3(normal, prev_tangent, tangent);
+			axis_angle_to_quat( q,normal, angle);
+			mul_qt_qtqt(key1->rot, q, key2->rot);
+		}
+
+		copy_v3_v3(prev_tangent, tangent);
+	}
+}
 /* Calculates paths ready for drawing/rendering.									*/
 /* -Usefull for making use of opengl vertex arrays for super fast strand drawing.	*/
 /* -Makes child strands possible and creates them too into the cache.				*/
@@ -2895,22 +2960,19 @@ void psys_cache_paths(ParticleSimulationData *sim, float cfra)
 		/*--interpolate actual path from data points--*/
 		for(k=0, ca=cache[p]; k<=steps; k++, ca++){
 			time = (float)k / (float)steps;
-
 			t = birthtime + time * (dietime - birthtime);
-
 			result.time = -t;
-
 			do_particle_interpolation(psys, p, pa, t, frs_sec, &pind, &result);
+			copy_v3_v3(ca->co, result.co);
 
 			/* dynamic hair is in object space */
 			/* keyed and baked are already in global space */
 			if(hair_dm)
-				mul_m4_v3(sim->ob->obmat, result.co);
+				mul_m4_v3(sim->ob->obmat, ca->co);
 			else if(!keyed && !baked && !(psys->flag & PSYS_GLOBAL_HAIR))
-				mul_m4_v3(hairmat, result.co);
+				mul_m4_v3(hairmat, ca->co);
 
-			VECCOPY(ca->co, result.co);
-			VECCOPY(ca->col, col);
+			copy_v3_v3(ca->col, col);
 		}
 		
 		/*--modify paths and calculate rotation & velocity--*/
@@ -2945,54 +3007,25 @@ void psys_cache_paths(ParticleSimulationData *sim, float cfra)
 
 		/* finally do rotation & velocity */
 		for(k=1, ca=cache[p]+1; k<=steps; k++, ca++) {
-			/* figure out rotation */
-			float cosangle, angle, tangent[3], normal[3], q[4];
-
-			if(k == 1) {
-				/* calculate initial tangent for incremental rotations */
-				VECSUB(tangent, ca->co, (ca - 1)->co);
-				normalize_v3_v3(prev_tangent, tangent);
-
-				/* First rotation is based on emitting face orientation.		*/
-				/* This is way better than having flipping rotations resulting	*/
-				/* from using a global axis as a rotation pole (vec_to_quat()). */
-				/* It's not an ideal solution though since it disregards the	*/
-				/* initial tangent, but taking that in to account will allow	*/
-				/* the possibility of flipping again. -jahka					*/
-				mat3_to_quat_is_ok( (ca-1)->rot,rotmat);
-			}
-			else {
-				VECSUB(tangent, ca->co, (ca - 1)->co);
-				normalize_v3(tangent);
-
-				cosangle= dot_v3v3(tangent, prev_tangent);
-
-				/* note we do the comparison on cosangle instead of
-				* angle, since floating point accuracy makes it give
-				* different results across platforms */
-				if(cosangle > 0.999999f) {
-					QUATCOPY((ca - 1)->rot, (ca - 2)->rot);
-				}
-				else {
-					angle= saacos(cosangle);
-					cross_v3_v3v3(normal, prev_tangent, tangent);
-					axis_angle_to_quat( q,normal, angle);
-					mul_qt_qtqt((ca - 1)->rot, q, (ca - 2)->rot);
-				}
-
-				VECCOPY(prev_tangent, tangent);
-			}
+			cache_key_incremental_rotation(ca, ca - 1, ca - 2, prev_tangent, k);
 
 			if(k == steps)
-				QUATCOPY(ca->rot, (ca - 1)->rot);
-			
+				copy_qt_qt(ca->rot, (ca - 1)->rot);
 
 			/* set velocity */
-			VECSUB(ca->vel, ca->co, (ca-1)->co);
+			sub_v3_v3v3(ca->vel, ca->co, (ca-1)->co);
 
 			if(k==1)
-				VECCOPY((ca-1)->vel, ca->vel);
+				copy_v3_v3((ca-1)->vel, ca->vel);
 		}
+		/* First rotation is based on emitting face orientation.
+		 * This is way better than having flipping rotations resulting
+		 * from using a global axis as a rotation pole (vec_to_quat()).
+		 * It's not an ideal solution though since it disregards the
+		 * initial tangent, but taking that in to account will allow
+		 * the possibility of flipping again. -jahka
+		 */
+		mat3_to_quat_is_ok(cache[p]->rot, rotmat);
 	}
 
 	psys->totcached = totpart;
@@ -3047,12 +3080,8 @@ void psys_cache_edit_paths(Scene *scene, Object *ob, PTCacheEdit *edit, float cf
 
 	frs_sec = (psys || edit->pid.flag & PTCACHE_VEL_PER_SEC) ? 25.0f : 1.0f;
 
-	if(pset->brushtype == PE_BRUSH_WEIGHT){
-		/* use weight painting colors now... */
-#if 0
-		sel_col[0] = sel_col[1] = sel_col[2] = 1.0f;
-		nosel_col[0] = nosel_col[1] = nosel_col[2] = 0.0f;
-#endif
+	if(pset->brushtype == PE_BRUSH_WEIGHT) {
+		;/* use weight painting colors now... */
 	}
 	else{
 		sel_col[0] = (float)edit->sel_col[0] / 255.0f;
@@ -3093,9 +3122,9 @@ void psys_cache_edit_paths(Scene *scene, Object *ob, PTCacheEdit *edit, float cf
 
 		if(psys) {
 			psys_mat_hair_to_global(ob, psmd->dm, psys->part->from, pa, hairmat);
-			VECCOPY(rotmat[0], hairmat[2]);
-			VECCOPY(rotmat[1], hairmat[1]);
-			VECCOPY(rotmat[2], hairmat[0]);
+			copy_v3_v3(rotmat[0], hairmat[2]);
+			copy_v3_v3(rotmat[1], hairmat[1]);
+			copy_v3_v3(rotmat[2], hairmat[0]);
 		}
 
 		birthtime = pind.birthtime;
@@ -3109,66 +3138,32 @@ void psys_cache_edit_paths(Scene *scene, Object *ob, PTCacheEdit *edit, float cf
 		/*--interpolate actual path from data points--*/
 		for(k=0, ca=cache[i]; k<=steps; k++, ca++){
 			time = (float)k / (float)steps;
-
 			t = birthtime + time * (dietime - birthtime);
-
 			result.time = -t;
-
 			do_particle_interpolation(psys, i, pa, t, frs_sec, &pind, &result);
+			copy_v3_v3(ca->co, result.co);
 
 			 /* non-hair points are already in global space */
 			if(psys && !(psys->flag & PSYS_GLOBAL_HAIR)) {
-				mul_m4_v3(hairmat, result.co);
+				mul_m4_v3(hairmat, ca->co);
 
-				/* create rotations for proper creation of children */
 				if(k) {
-					float cosangle, angle, tangent[3], normal[3], q[4];
-
-					if(k == 1) {
-						/* calculate initial tangent for incremental rotations */
-						VECSUB(tangent, ca->co, (ca - 1)->co);
-						normalize_v3_v3(prev_tangent, tangent);
-
-						/* First rotation is based on emitting face orientation.		*/
-						/* This is way better than having flipping rotations resulting	*/
-						/* from using a global axis as a rotation pole (vec_to_quat()). */
-						/* It's not an ideal solution though since it disregards the	*/
-						/* initial tangent, but taking that in to account will allow	*/
-						/* the possibility of flipping again. -jahka					*/
-						mat3_to_quat_is_ok( (ca-1)->rot,rotmat);
-					}
-					else {
-						VECSUB(tangent, ca->co, (ca - 1)->co);
-						normalize_v3(tangent);
-
-						cosangle= dot_v3v3(tangent, prev_tangent);
-
-						/* note we do the comparison on cosangle instead of
-						* angle, since floating point accuracy makes it give
-						* different results across platforms */
-						if(cosangle > 0.999999f) {
-							QUATCOPY((ca - 1)->rot, (ca - 2)->rot);
-						}
-						else {
-							angle= saacos(cosangle);
-							cross_v3_v3v3(normal, prev_tangent, tangent);
-							axis_angle_to_quat( q,normal, angle);
-							mul_qt_qtqt((ca - 1)->rot, q, (ca - 2)->rot);
-						}
-
-						VECCOPY(prev_tangent, tangent);
-					}
+					cache_key_incremental_rotation(ca, ca - 1, ca - 2, prev_tangent, k);
 
 					if(k == steps)
-						QUATCOPY(ca->rot, (ca - 1)->rot);
-				}
+						copy_qt_qt(ca->rot, (ca - 1)->rot);
 
-			}
+					/* set velocity */
+					sub_v3_v3v3(ca->vel, ca->co, (ca - 1)->co);
 
-			VECCOPY(ca->co, result.co);
-
-			ca->vel[0] = ca->vel[1] = 0.0f;
-			ca->vel[1] = 1.0f;
+					if(k==1)
+						copy_v3_v3((ca - 1)->vel, ca->vel);
+				}
+			}
+			else {
+				ca->vel[0] = ca->vel[1] = 0.0f;
+				ca->vel[1] = 1.0f;
+			}
 
 			/* selection coloring in edit mode */
 			if(pset->brushtype==PE_BRUSH_WEIGHT){
@@ -3216,6 +3211,16 @@ void psys_cache_edit_paths(Scene *scene, Object *ob, PTCacheEdit *edit, float cf
 
 			ca->time = t;
 		}
+		if(psys && !(psys->flag & PSYS_GLOBAL_HAIR)) {
+			/* First rotation is based on emitting face orientation.
+			 * This is way better than having flipping rotations resulting
+			 * from using a global axis as a rotation pole (vec_to_quat()).
+			 * It's not an ideal solution though since it disregards the
+			 * initial tangent, but taking that in to account will allow
+			 * the possibility of flipping again. -jahka
+			 */
+			mat3_to_quat_is_ok(cache[i]->rot, rotmat);
+		}
 	}
 
 	edit->totcached = totpart;
@@ -3492,6 +3497,7 @@ static void default_particle_settings(ParticleSettings *part)
 	part->adapt_angle= 5;
 	part->adapt_pix= 3;
 	part->kink_axis= 2;
+	part->kink_amp_clump= 1.f;
 	part->reactevent= PART_EVENT_DEATH;
 	part->disp=100;
 	part->from= PART_FROM_FACE;
@@ -3921,11 +3927,15 @@ static void do_child_modifiers(ParticleSimulationData *sim, ParticleTexture *pte
 		guided = do_guides(sim->psys->effectors, (ParticleKey*)state, cpa->parent, t);
 
 	if(guided==0){
-		if(kink_freq > 0.f)
-			do_prekink(state, par, par_rot, t, kink_freq, part->kink_shape,
-			part->kink_amp, part->kink, part->kink_axis, sim->ob->obmat);
-				
-		do_clump(state, par, t, part->clumpfac, part->clumppow, ptex ? ptex->clump : 1.f);
+		float clump = do_clump(state, par, t, part->clumpfac, part->clumppow, ptex ? ptex->clump : 1.f);
+
+		if(kink_freq != 0.f) {
+			float kink_amp = part->kink_amp * (1.f - part->kink_amp_clump * clump);
+
+			do_kink(state, par, par_rot, t, kink_freq, part->kink_shape,
+					kink_amp, part->kink_flat, part->kink, part->kink_axis,
+					sim->ob->obmat, sim->psys->part->childtype == PART_CHILD_FACES);
+		}
 	}
 
 	if(rough1 > 0.f)
@@ -4114,7 +4124,7 @@ void psys_get_particle_on_path(ParticleSimulationData *sim, int p, ParticleKey *
 		}
 		else{
 			/* offset the child from the parent position */
-			offset_child(cpa, keys, state, part->childflat, part->childrad);
+			offset_child(cpa, keys, keys->rot, state, part->childflat, part->childrad);
 		}
 
 		par = keys;
@@ -4212,7 +4222,7 @@ int psys_get_particle_state(ParticleSimulationData *sim, int p, ParticleKey *sta
 			float t = (cfra - pa->time) / pa->lifetime;
 
 			key1=&pa->state;
-			offset_child(cpa, key1, state, part->childflat, part->childrad);
+			offset_child(cpa, key1, key1->rot, state, part->childflat, part->childrad);
 
 			CLAMP(t,0.0,1.0);