From: horst.kraemer@t-online.de (Horst Kraemer) Subject: Re: Tetrahedron Interpolation Date: Sat, 06 Nov 1999 03:48:19 GMT Newsgroups: sci.math On Thu, 04 Nov 1999 13:48:48 -0500, Jonathan Hoyle wrote: > >> Could you elaborate on what you mean by "tetrahedron interpolation" ? > > Sure. This is the name I've heard it referred to, but there if there is > a better term let me know... > > Given four points in three dimensional space (xi,yi,zi) for i=0,1,2, and > a mapping to f(xi,yi,zi) for each xi,yi,zi, this algorithm interpolates > values for f on any (x,y,z) within the tetrahedron defined by the four > corners (xi,yi,zi). > > Although in my case my function is f:R^3->R^3, but I believe I can > generalize from an algorithm handling f:R^3->R. This appears to be a simple linear interpolation in barycentric coordinates (c1,c2,...,c_n,c_{n+1}) within the simplex defined by N+1 points and N+1 values in R^n. In one dimension it would correspond to interpolate x1 -> f1 x2 -> f2 f(c1*x1+c2*x2) = c1*f1 + c2*f2 (c1+c2=1), all c_i>=0 In two dimensions P1=(x1,y1) -> f1 P2=(x2,y2) -> f2 P3=(x3,y3) -> f3 f(c1*P1+c2*P2+c3*P3) = c1*f1+c2*f2+c3*f3, c1+c2+c3=1, all c_i>=0 In three dimensions P1=(x1,y1,z1) -> f1 P2=(x2,y2,z2) -> f2 P3=(x3,y3,z3) -> f3 P4=(x4,y4,z4) -> f4 f(c1*P1+c2*P2+c3*P3+c4*P4) := c1*f1+c2*f2+c3*f3+c3*f4 , c1+c2+c3+c4=1, all c_i>=0. 2 points which are not identical For any set of 3 points which are not on the same straigth line 4 points which are not in the same plane R etc. there is one and only one affine mapping from R^2 etc. R^3 into a linear space E such that f(P_i) = f_i, and that's probably the mapping you want because this is the "natural" mapping. Any set of 4 nonnegative real numbers (c1,c2,c3,c4) with sum 1 - in the so-called barycentric coordinate system relative to the tetrahedron P1,P2,P3,P4 - will yield some (x,y,z) = c1*P1+c2*P2+c3*P3+c4*P4 inside the tetrahedron and some f by f=c1*f1+c2*f2+c3*f3+c4*f4 where f fulfills the four initial conditions f(P_i) = f_i if c_i=1 and the rest of the c_j = 0 and where f is the unique affine mapping defined by these conditions. If one or more of the c_i are negative then the point is outside the tetrahedron but the mapping is still a continuation of the same mapping. You may rewrite the above formula in a asymetric way by choosing any P_i as the "anchor point", for instance P1 f((1-c2-c3-c4)P1 + c2P2 + c3P3 + c4 ) = f(P1 + c2(P2-P1) + c3(P3-P1) + c3(P4-P1) ) = f1 + c2*(f2-f1) + c3*(f3-f1) + c4*(f4-f1) where P is in the tetrahedron if and only if c2,c3,c4 are >=0 and c2+c3+c4<=1, i.e. c1=1-c2-c3-c4 >=0. No problem so far. What is your problem ? If you start from the c_i your problem is solved. Do you want to get the corresponding c_i for an arbitrary (x,y,z) ? This may be written in an elegant way. The c2,c3,c4 are the solutions of the linear equation system x2-x1 x3-x1 x4-x1 c2 x-x1 y2-y1 y3-y1 y4-y1 * c3 = y-y1 z2-z1 z3-z1 z4-z1 c4 z-z1 for c2,c3,c4 where x,y,z are the input parameters. If you substitute Q2:=P2-P1,Q3:=P3-P1,Q4=P4-P1 then x-x1 c2 = < y-y1 , [Q3.Q4]/det(Q2.Q3.Q4) > z-y1 x-x1 c3 = < y-y1 , [Q4.Q2]/det(Q2,Q3,Q4) > z-y1 x-x1 c4 = < y-y1 , [Q2.Q3]/det(Q2,Q3,Q4) > z-y1 Where <*.*> denotes the dot product and [*,*] the cross product and det(*,*,*) the determinant which may be evaluated as det(Q2,Q3,Q4) = = = f(x,y,z) = f1 + c2*(f2-f1) + c3*(f3-f1) + c4*(f4-f1) = (1-c2-c3-c4)*f1 + c2*f2 + c3*f3 *c4*f4 Hope it helps Horst