W06: Linear algebra becomes graphics and spatial computation

Prerequisites and observable goals#

The central question: which space contains the quantity?#

A coordinate tuple is a representation relative to a basis, origin, and units. The same physical point may have model, world, camera, clip, normalized-device, and screen coordinates. Applying a matrix without naming its source and destination spaces can produce a plausible image that is nevertheless mathematically inconsistent. Treat a transform as a typed function from one representation to another.

Use column vectors in this workshop. In two dimensions, represent a point as $(x,y,1)^T$ and a direction as $(d_x,d_y,0)^T$. A homogeneous translation matrix

moves points but not directions. The final component encodes this distinction; it is not decorative padding. A general projective output $(X,Y,W)$ with $W\ne0$ becomes Cartesian $(X/W,Y/W)$. If $W=0$, the same division is unavailable and the representation has a different projective meaning.

Order matters because the output of one map is the next input#

Let scaling be $S=\diag(2,3,1)$ and translation $(5,-1)$. For point $(1,2,1)^T$, $TS$ first scales to $(2,6,1)$ and then translates to $(7,5,1)$. In the opposite order, $ST$ first translates to $(6,1,1)$ and then scales to $(12,3,1)$. Matrices generally do not commute. Row-vector conventions reverse how composed matrices are written; copying a formula from another convention without changing the interpretation is a representation error, not a mysterious graphics bug.

The same idea explains normals. If tangent vectors transform by an invertible linear map $A$, a transformed normal $n'$ must satisfy $(n')^T(Av)=0$ whenever $n^Tv=0$. Choosing $n'=A^{-T}n$ gives $(A^{-T}n)^TAv=n^Tv$. Normalize afterward when a unit normal is required. Transforming normals with $A$ works in some special cases, such as rotations, but not arbitrary nonuniform scaling. Singular transforms require separate treatment because the inverse does not exist.

Derive barycentric coordinates rather than guess weights#

For a nondegenerate triangle with vertices $a,b,c\in\R^2$, write a point as

Subtract $a$ and eliminate $\alpha$:

The determinant is twice the signed triangle area. A zero determinant means a degenerate triangle and no unique two-dimensional barycentric coordinates. A very small determinant relative to coordinate scales indicates sensitivity; a floating-point nonzero check is not a robust geometric predicate for every adversarial mesh.

For $a=(0,0),b=(4,0),c=(0,4)$ and $p=(1,1)$, the equations are $4\beta=1,4\gamma=1$, so weights are $(1/2,1/4,1/4)$. Their sum is one and the reconstructed position is exactly $(1,1)$. A point inside the closed triangle has nonnegative weights. An outside point may have a negative weight even though the sum and reconstruction still hold. Interpolation and containment are different uses of the same representation.

Why screen-linear attributes can be wrong under perspective#

Perspective division is nonlinear. Suppose the vertices’ pre-division homogeneous components are $w_i>0$, and screen-space barycentric weights are $\lambda_i$. An attribute that is affine in the corresponding pre-projection geometry must be interpolated as

To see the normalization, let pre-projection weights be $\mu_i$. Projection gives screen weights proportional to $\mu_iw_i$, so $\lambda_i=\mu_iw_i/\sum_j\mu_jw_j$. Solving for $\mu_i$ and imposing $\sum_i\mu_i=1$ yields $\mu_i=(\lambda_i/w_i)/\sum_j\lambda_j/w_j$. Substituting $a=\sum_i\mu_i a_i$ gives the formula. This derivation assumes a consistent projective representation and a valid visible triangle; clipping and sign-changing $w$ require additional rules.

Use the weights $(1/2,1/4,1/4)$, $w=(1,2,4)$, and scalar attributes $(0,8,16)$. The numerator is $0+1+1=2$ and denominator $1/2+1/8+1/16=11/16$, giving $32/11\approx2.9091$. Naive screen-linear interpolation gives $0+2+4=6$, a different result. The issue is not that one is more aesthetically pleasing; they interpolate under different coordinate relationships.

Geometry, spatial indexing, and verification#

A spatial index usually produces candidates, not a geometric theorem for free. A grid bucket overlapping a query can contain points outside the exact region; candidate filtering must be followed by the final predicate. J05 makes this distinction explicit for axis-aligned boxes. In triangle work, a bounding box supplies a broad filter and barycentric coordinates supply a more precise test, with boundary and numerical policies documented.

Parallel execution adds another layer. Independent pixel or particle outputs can often be partitioned among workers, but shared reductions change accumulation order. Floating-point addition is not associative, so bitwise agreement with a sequential order may be unavailable even when both approximate the same real sum. Choose a tolerance or reproducibility contract appropriate to the use, test against an analytically justified or higher-accuracy reference, and include transfer/setup costs when claiming acceleration. A faster image with different geometry or different precision is not automatically a like-for-like speedup.

Independent practice#

F foundation C standard course work T transfer S optional stretch L runnable lab or project. Attempt the complete set before opening its solutions.

Complete solutions: W06#

Compare the first step where your work differs. A corrected answer is not yet an independent success.

Open every solution

"""W06: triangle coordinates and perspective interpolation, CPU-only."""
import math


def barycentric(point,a,b,c):
    if any(len(p)!=2 or any(not math.isfinite(x) for x in p) for p in (point,a,b,c)):
        raise ValueError("finite two-dimensional points required")
    ux,uy=b[0]-a[0],b[1]-a[1]
    vx,vy=c[0]-a[0],c[1]-a[1]
    px,py=point[0]-a[0],point[1]-a[1]
    determinant=ux*vy-uy*vx
    if not math.isfinite(determinant) or determinant==0:
        raise ValueError("degenerate or unrepresentable triangle")
    beta=(px*vy-py*vx)/determinant
    gamma=(ux*py-uy*px)/determinant
    result=(1-beta-gamma,beta,gamma)
    if not all(math.isfinite(x) for x in result):raise ValueError("nonfinite weights")
    return result


def perspective(weights,w,attributes):
    if len(weights)!=3 or len(w)!=3 or len(attributes)!=3 or not attributes[0]:
        raise ValueError("three vertices and nonempty attributes required")
    width=len(attributes[0])
    if width>64 or any(len(row)!=width for row in attributes):raise ValueError("attribute shape")
    if any(not math.isfinite(x) for x in [*weights,*w,*[v for row in attributes for v in row]]):
        raise ValueError("finite inputs required")
    if any(x<0 for x in weights) or not math.isclose(math.fsum(weights),1.,rel_tol=0.,abs_tol=1e-12):
        raise ValueError("inside-triangle weights summing to one required")
    if any(x<=0 for x in w):raise ValueError("positive clip w required in this bounded lab")
    coefficients=[a/b for a,b in zip(weights,w)]
    denominator=math.fsum(coefficients)
    if not math.isfinite(denominator) or denominator<=0:raise ValueError("invalid denominator")
    result=[math.fsum(c*row[k] for c,row in zip(coefficients,attributes))/denominator
            for k in range(width)]
    if not all(math.isfinite(x) for x in result):raise ValueError("nonfinite attribute")
    return result


if __name__=='__main__':
    print(barycentric((1,1),(0,0),(4,0),(0,4)))
    print(perspective([.5,.25,.25],[1,2,4],[[0],[8],[16]]))

Retention and readiness#

Label the coordinate space and units on every vector in a graphics derivation. Reconstruct a point from barycentric weights, and explain why perspective division changes interpolation. In a performance report, state which geometric and precision contract stayed fixed.