Update raymath.h

A spherical linear interpolation from one Vector3 point to another, without quaternions, and with selectable radius and orientation of arc.

These comments were in the function definition but to keep the header uncluttered I've put them here.  It explains the math.

    //Provides a constant-rate interpolation along an arc from v1 to v2 as t goes from 0 to 1.0, i.e. a spherical linear interpolation (slerp)
    //The basic formula is as follows:   (https://en.wikipedia.org/wiki/Slerp)
    //Define omega as the angle between the two points using omega=acos( a-dot-b / ||a||*||b|| ) (Law of cosines)
    //slerp(v1,v2,t) = ( sin((1-t)omega)/sin(omega) )*v1 + ( sin(t*omega)/sin(omega) )*v2
    //However, the basic formula has two drawbacks.  First, the center of the arc always lies on the midpoint, so the arc is a semi-circle.  One might want a flatter arc with an increased radius and thus a center located off the midpoint
    //Second, the arc is drawn with reference to the origin (i.e., it arcs away from wherever the origin is).  One might wish to specify a different direction of the arc.
    //This formula permits both parameters: iradius increases the radius from the minimum which is half the distance between the two points.  So the resulting radius of the arc is iradius+Length(midpoint-v1)
    //Theta (in radians) rotates the direction of the arc around the line connecting the two points (think of a jumprope)

    //Overview of formula below: We define L as the line between v1 and v2.  L=v2-v1.  Midpoint is the midpoint of this line: M=v1+0.5*L.  From the midpoint we will construct an arc centerpoint, on a radial line perpendicular to L and rotated around it by theta.
    //For the plane of rotation of the radial line, we construct a plane normal to L at M. The plane by definition has the form (L.x)*X+(L.y)*Y+(L.z)*Z+D=0.  To locate the plane at the midpoint, we substitute the midpoint into the equation and solve for D.  So D=-(L-dot-M).
    //Next we need to pick two orthogonal vectors in this plane as a basis for rotating the the radial line of the arc around it.  We call these vectors F and G.  The choice of lines on the plane is infinite and arbitrary.
    //So we choose an easy point to solve: the intersection of the plane and the x-axis (where y==0 && z==0).  This works in all cases except those where (1) the plane is parallel to the x-axis, which happens when L.x==0.  So in this case we choose y-axis or z-axis instead which are just as easy, and (2) the midpoint is very close to the origin, in which case we just move aside a small amount to avoid denominators close to zero
    //The x-intercept (called Incept below) is some vector {x,0,0} so we substitute this into the equation for the plane and obtain x=-D/L.x. Incept={-D/L.x,0,0} F is therefore the line between the midpoint and this intercept.  F=Incept-midpoint.  We normalize F because we will use it to rotate the radial line.
    //The second vector must be orthogonal to F and L.  It is easily obtained by LxF (cross product).  We normalize G as well.
    //Now F and G permit us to rotate the radial line by theta.
    //The centerpoint of the arc, called perp below, becomes midpoint+(F*cos(theta)+G*sin(theta))*radius).
    //Finally we translate our start and end point into a coordinate system with its center at perp: a=v1-perp and b=v2-perp.  We apply our basic slerp formula to points a and b.  We translate our result back into the original coordinate frame returning r+perp

src/raymath.h

@@ -621,6 +621,40 @@ RMDEF float3 Vector3ToFloatV(Vector3 v)
     return buffer;
 }
 
+RMDEF Vector3 Vector3Slerp(Vector3 v1, Vector3 v2, float t, float iradius, float theta) {
+
+    if (Vector3LengthSqr(Vector3Subtract(v2,v1))==0.0) return v1; //If our two points are the same, simply return v1
+
+    Vector3 r,a,b,incept,F,G;
+    Vector3 L=Vector3Subtract(v2,v1);  // L=v2-v1
+    Vector3 midpoint=Vector3Add(v1,Vector3Scale(L,0.5)); // midpoint=v1+L*0.5
+    float D=-1*Vector3DotProduct(L,midpoint);  //Solve plane equation for D
+    if (L.x!=0.0) {  //Find intercept of midpoint and x-axis as long as the plane isn't parallel to the x-axis
+        incept=(Vector3){-D/L.x,0.0,0.0};
+    }
+    else if (L.y!=0.0) { //Otherwise use y-axis
+        incept=(Vector3){0.0,-D/L.y,0.0};
+    }
+    else { //or z-axis.  All three cannot be zero because of the first if-statement in routine
+        incept=(Vector3){0.0,0.0,-D/L.z};
+    }
+
+    if (Vector3LengthSqr(Vector3Subtract(incept,midpoint))<0.1) {  //If our midpoint and axis intercept are very close together, displace the midpoint slightly so that normalization of F does not result in divide-by-zero or G=0 because F={1,0,0}
+        midpoint=Vector3Add(v1,Vector3Scale(L,0.49));
+        midpoint=Vector3Add(midpoint,Vector3CrossProduct(midpoint,{0.2,0.0,0.0}));
+    }
+    F=Vector3Normalize(Vector3Subtract(incept,midpoint)); // First orthogonal vector
+    G=Vector3Normalize(Vector3CrossProduct(L,F));   // Second orthogonal vector
+
+    Vector3 perp=Vector3Add(midpoint,Vector3Scale(Vector3Add( Vector3Scale(F,cos(theta)),Vector3Scale(G,sin(theta))),iradius)); //perp=midpoint+(F*cos(theta)+G*sin(theta))*iradius
+    a=Vector3Subtract(v1,perp); //translate coordinates to center at perp
+    b=Vector3Subtract(v2,perp);
+    float om=acos( Vector3DotProduct(a,b) / (Vector3Length(a)*Vector3Length(b)) ); //Apply slerp formula for points at origin
+    r= Vector3Add(Vector3Scale(a,sin((1.0F-t)*om)/sin(om)),Vector3Scale(b,sin(t*om)/sin(om)));
+return Vector3Add(r,perp); //translate coordinates back to original frame
+}
+
+
 //----------------------------------------------------------------------------------
 // Module Functions Definition - Matrix math
 //----------------------------------------------------------------------------------