Fix for bug: cross platform strand render differences with kink/branch.

This time is was due to different accuracy of floating point computation,
now it uses does a comparison a bit different to avoid this.

Also changed the vectoquat function to be threadsafe.

11 files changed, 63 insertions, 71 deletions

diff --git a/source/blender/blenkernel/intern/anim.c b/source/blender/blenkernel/intern/anim.c
index 21b50cd01b2..764b8aee032 100644
--- a/source/blender/blenkernel/intern/anim.c
+++ b/source/blender/blenkernel/intern/anim.c
@@ -387,7 +387,7 @@ static void vertex_dupli__mapFunc(void *userData, int index, float *co, float *n
 {
 	DupliObject *dob;
 	struct vertexDupliData *vdd= userData;
-	float vec[3], *q2, mat[3][3], tmat[4][4], obmat[4][4];
+	float vec[3], q2[4], mat[3][3], tmat[4][4], obmat[4][4];
 	
 	VECCOPY(vec, co);
 	Mat4MulVecfl(vdd->pmat, vec);
@@ -405,7 +405,7 @@ static void vertex_dupli__mapFunc(void *userData, int index, float *co, float *n
 			vec[0]= -no_s[0]; vec[1]= -no_s[1]; vec[2]= -no_s[2];
 		}
 		
-		q2= vectoquat(vec, vdd->ob->trackflag, vdd->ob->upflag);
+		vectoquat(vec, vdd->ob->trackflag, vdd->ob->upflag, q2);
 		
 		QuatToMat3(q2, mat);
 		Mat4CpyMat4(tmat, obmat);
@@ -737,7 +737,7 @@ static void new_particle_duplilist(ListBase *lb, ID *id, Object *par, float par_
 	float ctime, pa_time, scale = 1.0f;
 	float tmat[4][4], mat[4][4], obrotmat[4][4], pamat[4][4], size=0.0;
 	float obmat[4][4], (*obmatlist)[4][4]=0;
-	float xvec[3] = {-1.0, 0.0, 0.0}, *q;
+	float xvec[3] = {-1.0, 0.0, 0.0}, q[4];
 	int lay, a, b, k, step_nbr = 0, counter, hair = 0;
 	int totpart, totchild, totgroup=0, pa_num;
 
@@ -891,7 +891,7 @@ static void new_particle_duplilist(ListBase *lb, ID *id, Object *par, float par_
 					where_is_object_time(ob, ctime-pa_time);
 					
 					if(!hair) {
-						q = vectoquat(xvec, ob->trackflag, ob->upflag);
+						vectoquat(xvec, ob->trackflag, ob->upflag, q);
 						QuatToMat4(q, obrotmat);
 						obrotmat[3][3]= 1.0f;
 						Mat4MulMat4(mat, obrotmat, pamat);
diff --git a/source/blender/blenkernel/intern/constraint.c b/source/blender/blenkernel/intern/constraint.c
index 572c91de7f3..7c19563d932 100644
--- a/source/blender/blenkernel/intern/constraint.c
+++ b/source/blender/blenkernel/intern/constraint.c
@@ -1259,7 +1259,7 @@ static void followpath_get_tarmat (bConstraint *con, bConstraintOb *cob, bConstr
 	
 	if (VALID_CONS_TARGET(ct)) {
 		Curve *cu= ct->tar->data;
-		float q[4], vec[4], dir[3], *quat, x1;
+		float q[4], vec[4], dir[3], quat[4], x1;
 		float totmat[4][4];
 		float curvetime;
 		
@@ -1284,7 +1284,7 @@ static void followpath_get_tarmat (bConstraint *con, bConstraintOb *cob, bConstr
 			
 			if ( where_on_path(ct->tar, curvetime, vec, dir) ) {
 				if (data->followflag) {
-					quat= vectoquat(dir, (short) data->trackflag, (short) data->upflag);
+					vectoquat(dir, (short) data->trackflag, (short) data->upflag, quat);
 					
 					Normalize(dir);
 					q[0]= (float)cos(0.5*vec[3]);
diff --git a/source/blender/blenkernel/intern/curve.c b/source/blender/blenkernel/intern/curve.c
index c23cc90aeb8..dcd8bc6d9eb 100644
--- a/source/blender/blenkernel/intern/curve.c
+++ b/source/blender/blenkernel/intern/curve.c
@@ -1788,13 +1788,13 @@ void makeBevelList(Object *ob)
 			bevp2->cosa= bevp1->cosa;
 
 			if(cu->flag & CU_3D) {	/* 3D */
-				float *quat, q[4];
+				float quat[4], q[4];
 			
 				vec[0]= bevp1->x - bevp2->x;
 				vec[1]= bevp1->y - bevp2->y;
 				vec[2]= bevp1->z - bevp2->z;
 				
-				quat= vectoquat(vec, 5, 1);
+				vectoquat(vec, 5, 1, quat);
 				
 				Normalize(vec);
 				q[0]= (float)cos(0.5*bevp1->alfa);
@@ -1820,7 +1820,7 @@ void makeBevelList(Object *ob)
 			while(nr--) {
 	
 				if(cu->flag & CU_3D) {	/* 3D */
-					float *quat, q[4];
+					float quat[4], q[4];
 				
 					vec[0]= bevp2->x - bevp0->x;
 					vec[1]= bevp2->y - bevp0->y;
@@ -1828,7 +1828,7 @@ void makeBevelList(Object *ob)
 					
 					Normalize(vec);
 
-					quat= vectoquat(vec, 5, 1);
+					vectoquat(vec, 5, 1, quat);
 					
 					q[0]= (float)cos(0.5*bevp1->alfa);
 					x1= (float)sin(0.5*bevp1->alfa);
diff --git a/source/blender/blenkernel/intern/exotic.c b/source/blender/blenkernel/intern/exotic.c
index 0f1f8c6078a..b78efcd4ecf 100644
--- a/source/blender/blenkernel/intern/exotic.c
+++ b/source/blender/blenkernel/intern/exotic.c
@@ -878,7 +878,7 @@ static void read_videoscape_lamp(char *str)
 	Object *ob;
 	Lamp *la;
 	FILE *fp;
-	float vec[3], *q1;
+	float vec[3], q1[4];
 	int tot, val;
 	char s[50];
 	
@@ -906,7 +906,7 @@ static void read_videoscape_lamp(char *str)
 		
 		fscanf(fp, "%f %f %f\n", ob->loc, ob->loc+1, ob->loc+2);
 		val= fscanf(fp, "%f %f %f\n", vec, vec+1, vec+2);
-		q1= vectoquat(vec, 5, 2);
+		vectoquat(vec, 5, 2, q1);
 		QuatToEul(q1, ob->rot);
 		
 		if(val<=0) break;
diff --git a/source/blender/blenkernel/intern/lattice.c b/source/blender/blenkernel/intern/lattice.c
index 2b815c28cc5..8e89b0d7966 100644
--- a/source/blender/blenkernel/intern/lattice.c
+++ b/source/blender/blenkernel/intern/lattice.c
@@ -575,14 +575,13 @@ static float *calc_curve_deform(Object *par, float *co, short axis, CurveDeform
 	}
 	
 	if( where_on_path_deform(par, fac, loc, dir)) {	/* returns OK */
-		float q[4], mat[3][3];
-		float *quat;
+		float q[4], mat[3][3], quat[4];
 		
 		if(cd->no_rot_axis)	/* set by caller */
 			dir[cd->no_rot_axis-1]= 0.0f;
 		
 		/* -1 for compatibility with old track defines */
-		quat= vectoquat(dir, axis-1, upflag);	/* gives static quat */
+		vectoquat(dir, axis-1, upflag, quat);
 		
 		/* the tilt */
 		if(loc[3]!=0.0) {
diff --git a/source/blender/blenkernel/intern/object.c b/source/blender/blenkernel/intern/object.c
index e25c89c9f83..16934783a07 100644
--- a/source/blender/blenkernel/intern/object.c
+++ b/source/blender/blenkernel/intern/object.c
@@ -1438,7 +1438,7 @@ int enable_cu_speed= 1;
 static void ob_parcurve(Object *ob, Object *par, float mat[][4])
 {
 	Curve *cu;
-	float q[4], vec[4], dir[3], *quat, x1, ctime;
+	float q[4], vec[4], dir[3], quat[4], x1, ctime;
 	float timeoffs = 0.0, sf_orig = 0.0;
 	
 	Mat4One(mat);
@@ -1487,7 +1487,7 @@ static void ob_parcurve(Object *ob, Object *par, float mat[][4])
  	if( where_on_path(par, ctime, vec, dir) ) {
 
 		if(cu->flag & CU_FOLLOW) {
-			quat= vectoquat(dir, ob->trackflag, ob->upflag);
+			vectoquat(dir, ob->trackflag, ob->upflag, quat);
 			
 			/* the tilt */
 			Normalize(dir);
@@ -1887,13 +1887,13 @@ static void solve_parenting (Object *ob, Object *par, float obmat[][4], float sl
 }
 void solve_tracking (Object *ob, float targetmat[][4])
 {
-	float *quat;
+	float quat[4];
 	float vec[3];
 	float totmat[3][3];
 	float tmat[4][4];
 	
 	VecSubf(vec, ob->obmat[3], targetmat[3]);
-	quat= vectoquat(vec, ob->trackflag, ob->upflag);
+	vectoquat(vec, ob->trackflag, ob->upflag, quat);
 	QuatToMat3(quat, totmat);
 	
 	if(ob->parent && (ob->transflag & OB_POWERTRACK)) {
diff --git a/source/blender/blenkernel/intern/particle.c b/source/blender/blenkernel/intern/particle.c
index 93d9d972b80..0af959597a6 100644
--- a/source/blender/blenkernel/intern/particle.c
+++ b/source/blender/blenkernel/intern/particle.c
@@ -1369,7 +1369,7 @@ static float vert_weight(MDeformVert *dvert, int group)
 }
 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])
 {
-	float vec[3]={0.0,0.0,0.0}, q1[4]={1,0,0,0},*q2;
+	float vec[3]={0.0,0.0,0.0}, q1[4]={1,0,0,0},q2[4];
 	float t;
 
 	CLAMP(time,0.0,1.0);
@@ -1391,10 +1391,9 @@ static void do_prekink(ParticleKey *state, ParticleKey *par, float *par_rot, flo
 		case PART_KINK_CURL:
 			vec[axis]=1.0;
 			if(par_rot)
-				q2=par_rot;
-			else{
-				q2=vectoquat(par->vel,axis,(axis+1)%3);
-			}
+				QUATCOPY(q2,par_rot)
+			else
+				vectoquat(par->vel,axis,(axis+1)%3, q2);
 			QuatMulVecf(q2,vec);
 			VecMulf(vec,amplitude);
 			VECADD(state->co,state->co,vec);
@@ -1441,9 +1440,9 @@ static void do_prekink(ParticleKey *state, ParticleKey *par, float *par_rot, flo
 				float inp_y,inp_z,length;
 				
 				if(par_rot)
-					q2=par_rot;
+					QUATCOPY(q2,par_rot)
 				else
-					q2=vectoquat(par->vel,axis,(axis+1)%3);
+					vectoquat(par->vel,axis,(axis+1)%3,q2);
 				QuatMulVecf(q2,y_vec);
 				QuatMulVecf(q2,z_vec);
 				
@@ -1592,7 +1591,7 @@ int do_guide(ParticleKey *state, int pa_num, float time, ListBase *lb)
 						VecRotToQuat(guidedir,guidevec[3]-ec->firstloc[3],rot2);
 						QuatMulVecf(rot2,pa_loc);
 
-						//q=vectoquat(guidedir, pd->kink_axis, (pd->kink_axis+1)%3);
+						//vectoquat(guidedir, pd->kink_axis, (pd->kink_axis+1)%3, q);
 						//QuatMul(par.rot,rot2,q);
 					}
 					//else{
@@ -2536,16 +2535,12 @@ void psys_cache_paths(Object *ob, ParticleSystem *psys, float cfra, int editupda
 			/* figure out rotation */
 			
 			if(k) {
-				float angle, tangent[3], normal[3], q[4];
+				float cosangle, angle, tangent[3], normal[3], q[4];
 
 				if(k == 1) {
-					float *q2;
-
 					VECSUB(tangent, ca->co, (ca - 1)->co);
 
-					q2 = vectoquat(tangent, OB_POSX, OB_POSZ);
-
-					QUATCOPY((ca - 1)->rot, q2);
+					vectoquat(tangent, OB_POSX, OB_POSZ, (ca-1)->rot);
 
 					VECCOPY(prev_tangent, tangent);
 					Normalize(prev_tangent);
@@ -2553,12 +2548,17 @@ void psys_cache_paths(Object *ob, ParticleSystem *psys, float cfra, int editupda
 				else {
 					VECSUB(tangent, ca->co, (ca - 1)->co);
 					Normalize(tangent);
-					angle = saacos(Inpf(tangent, prev_tangent));
 
-					if((angle > -0.000001) && (angle < 0.000001)){
+					cosangle= Inpf(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{
+					else {
+						angle= saacos(cosangle);
 						Crossf(normal, prev_tangent, tangent);
 						VecRotToQuat(normal, angle, q);
 						QuatMul((ca - 1)->rot, q, (ca - 2)->rot);
@@ -2567,9 +2567,8 @@ void psys_cache_paths(Object *ob, ParticleSystem *psys, float cfra, int editupda
 					VECCOPY(prev_tangent, tangent);
 				}
 
-				if(k == steps) {
+				if(k == steps)
 					QUATCOPY(ca->rot, (ca - 1)->rot);
-				}
 			}
 
 			
@@ -3412,11 +3411,9 @@ void psys_get_particle_on_path(Object *ob, ParticleSystem *psys, int p, Particle
 		//else{
 		//	/* TODO: different rotations */
 		//	float nvel[3];
-		//	float *q2;
 		//	VECCOPY(nvel,state->vel);
 		//	VecMulf(nvel,-1.0f);
-		//	q2=vectoquat(nvel, OB_POSX, OB_POSZ);
-		//	QUATCOPY(state->rot,q2);
+		//	vectoquat(nvel, OB_POSX, OB_POSZ, state->rot);
 		//}
 
 		dfra = keys[2].time - keys[1].time;
@@ -3740,7 +3737,7 @@ void psys_get_dupli_texture(Object *ob, ParticleSettings *part, ParticleSystemMo
 void psys_get_dupli_path_transform(Object *ob, ParticleSystem *psys, ParticleSystemModifierData *psmd, ParticleData *pa, ChildParticle *cpa, ParticleCacheKey *cache, float mat[][4], float *scale)
 {
 	float loc[3], nor[3], vec[3], side[3], len, obrotmat[4][4], qmat[4][4];
-	float xvec[3] = {-1.0, 0.0, 0.0}, *q;
+	float xvec[3] = {-1.0, 0.0, 0.0}, q[4];
 
 	VecSubf(vec, (cache+cache->steps-1)->co, cache->co);
 	len= Normalize(vec);
@@ -3756,7 +3753,7 @@ void psys_get_dupli_path_transform(Object *ob, ParticleSystem *psys, ParticleSys
 		if(!pa)
 			pa= psys->particles+cpa->pa[0];
 
-		q = vectoquat(xvec, ob->trackflag, ob->upflag);
+		vectoquat(xvec, ob->trackflag, ob->upflag, q);
 		QuatToMat4(q, obrotmat);
 		obrotmat[3][3]= 1.0f;
 
diff --git a/source/blender/blenkernel/intern/particle_system.c b/source/blender/blenkernel/intern/particle_system.c
index 37a14543ee6..02f19597389 100644
--- a/source/blender/blenkernel/intern/particle_system.c
+++ b/source/blender/blenkernel/intern/particle_system.c
@@ -1573,7 +1573,7 @@ void reset_particle(ParticleData *pa, ParticleSystem *psys, ParticleSystemModifi
 	ParticleTexture ptex;
 	ParticleKey state;
 	IpoCurve *icu=0;
-	float fac, phasefac, nor[3]={0,0,0},loc[3],tloc[3],vel[3]={0.0,0.0,0.0},rot[4],*q2=0;
+	float fac, phasefac, nor[3]={0,0,0},loc[3],tloc[3],vel[3]={0.0,0.0,0.0},rot[4],q2[4];
 	float r_vel[3],r_ave[3],r_rot[4],p_vel[3]={0.0,0.0,0.0};
 	float x_vec[3]={1.0,0.0,0.0}, utan[3]={0.0,1.0,0.0}, vtan[3]={0.0,0.0,1.0}, rot_vec[3]={0.0,0.0,0.0};
 	float q_phase[4];
@@ -1749,7 +1749,7 @@ void reset_particle(ParticleData *pa, ParticleSystem *psys, ParticleSystemModifi
 		
 		/* create rotation quat */
 		VecMulf(rot_vec,-1.0);
-		q2= vectoquat(rot_vec, OB_POSX, OB_POSZ);
+		vectoquat(rot_vec, OB_POSX, OB_POSZ, q2);
 
 		/* randomize rotation quat */
 		if(part->randrotfac!=0.0f)
@@ -3085,7 +3085,7 @@ static void deflect_particle(Object *pob, ParticleSystemModifierData *psmd, Part
 	ParticleKey cstate;
 	float imat[4][4];
 	float co1[3],co2[3],def_loc[3],def_nor[3],unit_nor[3],def_tan[3],dvec[3],def_vel[3],dave[3],dvel[3];
-	float t_co1[3],t_co2[3];
+	float t_co1[3]={0.0,0.0,0.0},t_co2[3]={0.0,0.0,0.0};
 	float pa_minmax[6];
 	float min_w[4], zerovec[3]={0.0,0.0,0.0}, ipoint[3];
 	float min_d,dotprod,damp,frict,o_len,d_len,radius=-1.0f;
@@ -3811,7 +3811,7 @@ static void boid_brain(BoidVecFunc *bvf, ParticleData *pa, Object *ob, ParticleS
 static void boid_body(BoidVecFunc *bvf, ParticleData *pa, ParticleSystem *psys, ParticleSettings *part, float timestep, float *acc, ParticleKey *state)
 {
 	float dvec[3], bvec[3], length, max_vel=part->max_vel;
-	float *q2, q[4];
+	float q2[4], q[4];
 	float g=9.81f, pa_mass=part->mass;
 	float yvec[3]={0.0,1.0,0.0}, zvec[3]={0.0,0.0,-1.0}, bank;
 
@@ -3893,7 +3893,7 @@ static void boid_body(BoidVecFunc *bvf, ParticleData *pa, ParticleSystem *psys,
 
 	VECCOPY(dvec,state->vel);
 	VecMulf(dvec,-1.0f);
-	q2= vectoquat(dvec, OB_POSX, OB_POSZ);
+	vectoquat(dvec, OB_POSX, OB_POSZ, q2);
 
 	QuatMul(state->rot,q,q2);
 
diff --git a/source/blender/blenlib/BLI_arithb.h b/source/blender/blenlib/BLI_arithb.h
index b3c94ca79ba..84348d59448 100644
--- a/source/blender/blenlib/BLI_arithb.h
+++ b/source/blender/blenlib/BLI_arithb.h
@@ -261,7 +261,7 @@ void Vec2Addf(float *v, float *v1, float *v2);
 void Vec2Subf(float *v, float *v1, float *v2);
 void Vec2Copyf(float *v1, float *v2);
 
-float *vectoquat(float *vec, short axis, short upflag);
+void vectoquat(float *vec, short axis, short upflag, float *q);
 
 float VecAngle2(float *v1, float *v2);
 float VecAngle3(float *v1, float *v2, float *v3);
diff --git a/source/blender/blenlib/intern/arithb.c b/source/blender/blenlib/intern/arithb.c
index cf02fc67ec4..bbb8b418fc5 100644
--- a/source/blender/blenlib/intern/arithb.c
+++ b/source/blender/blenlib/intern/arithb.c
@@ -1337,9 +1337,8 @@ void NormalQuat(float *q)
 	}
 }
 
-float *vectoquat( float *vec, short axis, short upflag)
+void vectoquat(float *vec, short axis, short upflag, float *q)
 {
-	static float q1[4];
 	float q2[4], nor[3], *fp, mat[3][3], angle, si, co, x2, y2, z2, len1;
 	
 	/* first rotate to axis */
@@ -1351,11 +1350,11 @@ float *vectoquat( float *vec, short axis, short upflag)
 		x2= -vec[0] ; y2= -vec[1] ; z2= -vec[2];
 	}
 	
-	q1[0]=1.0; 
-	q1[1]=q1[2]=q1[3]= 0.0;
+	q[0]=1.0; 
+	q[1]=q[2]=q[3]= 0.0;
 
 	len1= (float)sqrt(x2*x2+y2*y2+z2*z2);
-	if(len1 == 0.0) return(q1);
+	if(len1 == 0.0) return;
 
 	/* nasty! I need a good routine for this...
 	 * problem is a rotation of an Y axis to the negative Y-axis for example.
@@ -1366,9 +1365,8 @@ float *vectoquat( float *vec, short axis, short upflag)
 		nor[1]= -z2;
 		nor[2]= y2;
 
-		if( fabs(y2)+fabs(z2)<0.0001 ) {
+		if(fabs(y2)+fabs(z2)<0.0001)
 			nor[1]= 1.0;
-		}
 
 		co= x2;
 	}
@@ -1377,9 +1375,8 @@ float *vectoquat( float *vec, short axis, short upflag)
 		nor[1]= 0.0;
 		nor[2]= -x2;
 		
-		if( fabs(x2)+fabs(z2)<0.0001 ) {
+		if(fabs(x2)+fabs(z2)<0.0001)
 			nor[2]= 1.0;
-		}
 		
 		co= y2;
 	}
@@ -1388,9 +1385,8 @@ float *vectoquat( float *vec, short axis, short upflag)
 		nor[1]= x2;
 		nor[2]= 0.0;
 
-		if( fabs(x2)+fabs(y2)<0.0001 ) {
+		if(fabs(x2)+fabs(y2)<0.0001)
 			nor[0]= 1.0;
-		}
 
 		co= z2;
 	}
@@ -1400,13 +1396,13 @@ float *vectoquat( float *vec, short axis, short upflag)
 	
 	angle= 0.5f*saacos(co);
 	si= (float)sin(angle);
-	q1[0]= (float)cos(angle);
-	q1[1]= nor[0]*si;
-	q1[2]= nor[1]*si;
-	q1[3]= nor[2]*si;
+	q[0]= (float)cos(angle);
+	q[1]= nor[0]*si;
+	q[2]= nor[1]*si;
+	q[3]= nor[2]*si;
 	
 	if(axis!=upflag) {
-		QuatToMat3(q1, mat);
+		QuatToMat3(q, mat);
 
 		fp= mat[2];
 		if(axis==0) {
@@ -1429,10 +1425,8 @@ float *vectoquat( float *vec, short axis, short upflag)
 		q2[2]= y2*si;
 		q2[3]= z2*si;
 			
-		QuatMul(q1,q2,q1);
+		QuatMul(q,q2,q);
 	}
-
-	return(q1);
 }
 
 void VecUpMat3old( float *vec, float mat[][3], short axis)
diff --git a/source/blender/python/api2_2x/vector.c b/source/blender/python/api2_2x/vector.c
index 45ad36d0093..b3bcdf17762 100644
--- a/source/blender/python/api2_2x/vector.c
+++ b/source/blender/python/api2_2x/vector.c
@@ -1,5 +1,5 @@
 /*
- * $Id: vector.c 12314 2007-10-20 20:24:09Z campbellbarton $
+ * $Id$
  * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
  *
  * This program is free software; you can redistribute it and/or
@@ -168,7 +168,7 @@ PyObject *Vector_Resize4D(VectorObject * self)
   extract a quaternion from the vector and the track and up axis */
 PyObject *Vector_ToTrackQuat( VectorObject * self, PyObject * args )
 {
-	float vec[3];
+	float vec[3], quat[4];
 	char *strack, *sup;
 	short track = 2, up = 1;
 
@@ -271,7 +271,9 @@ PyObject *Vector_ToTrackQuat( VectorObject * self, PyObject * args )
 	vec[1] = -self->vec[1];
 	vec[2] = -self->vec[2];
 
-	return newQuaternionObject(vectoquat(vec, track, up), Py_NEW);
+	vectoquat(vec, track, up, quat);
+
+	return newQuaternionObject(quat, Py_NEW);
 }