A: 

I think what you need is affine transformation which can be accomplished using matrix math.

Andy West
Actually, this is a projective transformation.
tzaman
You're right, I forgot about the perspective. An affine transformation would result in a distorted image. Sorry about that.
Andy West
+4  A: 

What you want is called planar rectification, and it's not all that simple, I'm afraid. What you need to do is recover the homography that maps this oblique view of the van side onto the front-facing view. Photoshop / etc. have tools to do this for you given some control points; if you want to implement it for yourself you'll have to start delving into computer vision.

Edit - OK, here you go: a Python script to do the warping, using the OpenCV library which has convenient functions to calculate the homography and warp the image for you:

import cv

def warpImage(image, corners, target):
    mat = cv.CreateMat(3, 3, cv.CV_32F)
    cv.GetPerspectiveTransform(corners, target, mat)
    out = cv.CreateMat(height, width, cv.CV_8UC3)
    cv.WarpPerspective(image, out, mat, cv.CV_INTER_CUBIC)
    return out

if __name__ == '__main__':
    width, height = 400, 250
    corners = [(171,72),(331,93),(333,188),(177,210)]
    target = [(0,0),(width,0),(width,height),(0,height)]
    image = cv.LoadImageM('fries.jpg')
    out = warpImage(image, corners, target)
    cv.SaveImage('fries_warped.jpg', out)

And the output:
Warped image

OpenCV also has C and C++ bindings, or you can use EmguCV for a .NET wrapper; the API is fairly consistent across all languages so you can replicate this in whichever language suits your fancy.

tzaman
I fear I phrased my question too generally - I know the dimensions of the rectangular face I'm trying to extract, and the coordinates of it's corners as in the image, leaving me wondering only how to get best-quality pixels out - how to calculate the coverage of each pixel, and how to best mix them
Will
@Will - selecting four quadrilateral corner points in the image, and four corresponding rectangle points is enough information to recover the homography (i.e. transformation matrix), but you need to do some math; see this pdf: http://www.cs.brown.edu/courses/csci1950-g/asgn/proj6/resources/ProjectiveMappings.pdf
tzaman
+2  A: 

Look up "quad to quad" transform, e.g. threeblindmiceandamonkey.
A 3x3 transform on 2d homogeneous coordinates can transform any 4 points (a quad) to any other quad; conversely, any fromquad and toquad, such as the corners of your truck and a target rectangle, give a 3 x 3 transform.

Qt has quadToQuad and can transform pixmaps with it, but I guess you don't have Qt ?
Added 10Jun: from labs.trolltech.com/page/Graphics/Examples there's a nice demo which quad-to-quads a pixmap as you move the corners:

alt text

Added 11Jun: @Will, here's translate.h in Python (which you know a bit ? """ ...""" are multiline comments.)
perstrans() is the key; hope that makes sense, if not ask.

Bytheway, you could map the pixels one by one, mapQuadToQuad( target rect, orig quad ), but without pixel interpolation it'll look terrible; OpenCV does it all.

#!/usr/bin/env python
""" square <-> quad maps
    from http://threeblindmiceandamonkey.com/?p=16 matrix.h
"""

from __future__ import division
import numpy as np

__date__ = "2010-06-11 jun denis"

def det2(a, b, c, d):
    return a*d - b*c

def mapSquareToQuad( quad ):  # [4][2]
    SQ = np.zeros((3,3))
    px = quad[0,0] - quad[1,0] + quad[2,0] - quad[3,0]
    py = quad[0,1] - quad[1,1] + quad[2,1] - quad[3,1]
    if abs(px) < 1e-10 and abs(py) < 1e-10:
        SQ[0,0] = quad[1,0] - quad[0,0]
        SQ[1,0] = quad[2,0] - quad[1,0]
        SQ[2,0] = quad[0,0]
        SQ[0,1] = quad[1,1] - quad[0,1]
        SQ[1,1] = quad[2,1] - quad[1,1]
        SQ[2,1] = quad[0,1]
        SQ[0,2] = 0.
        SQ[1,2] = 0.
        SQ[2,2] = 1.
        return SQ
    else:
        dx1 = quad[1,0] - quad[2,0]
        dx2 = quad[3,0] - quad[2,0]
        dy1 = quad[1,1] - quad[2,1]
        dy2 = quad[3,1] - quad[2,1]
        det = det2(dx1,dx2, dy1,dy2)
        if det == 0.:
            return None
        SQ[0,2] = det2(px,dx2, py,dy2) / det
        SQ[1,2] = det2(dx1,px, dy1,py) / det
        SQ[2,2] = 1.
        SQ[0,0] = quad[1,0] - quad[0,0] + SQ[0,2]*quad[1,0]
        SQ[1,0] = quad[3,0] - quad[0,0] + SQ[1,2]*quad[3,0]
        SQ[2,0] = quad[0,0]
        SQ[0,1] = quad[1,1] - quad[0,1] + SQ[0,2]*quad[1,1]
        SQ[1,1] = quad[3,1] - quad[0,1] + SQ[1,2]*quad[3,1]
        SQ[2,1] = quad[0,1]
        return SQ

#...............................................................................
def mapQuadToSquare( quad ):
    return np.linalg.inv( mapSquareToQuad( quad ))

def mapQuadToQuad( a, b ):
    return np.dot( mapQuadToSquare(a), mapSquareToQuad(b) )

def perstrans( X, t ):
    """ perspective transform X Nx2, t 3x3:
        [x0 y0 1] t = [a0 b0 w0] -> [a0/w0 b0/w0]
        [x1 y1 1] t = [a1 b1 w1] -> [a1/w1 b1/w1]
        ...
    """
    x1 = np.vstack(( X.T, np.ones(len(X)) ))
    y = np.dot( t.T, x1 )
    return (y[:-1] / y[-1]) .T

#...............................................................................
if __name__ == "__main__":
    np.set_printoptions( 2, threshold=100, suppress=True )  # .2f

    sq = np.array([[0,0], [1,0], [1,1], [0,1]])
    quad = np.array([[171, 72], [331, 93], [333, 188], [177, 210]])
    print "quad:", quad
    print "square to quad:", perstrans( sq, mapSquareToQuad(quad) )
    print "quad to square:", perstrans( quad, mapQuadToSquare(quad) )

    dw, dh = 300, 250
    rect = np.array([[0, 0], [dw, 0], [dw, dh], [0, dh]])
    quadquad = mapQuadToQuad( quad, rect )
    print "quad to quad transform:", quadquad
    print "quad to rect:", perstrans( quad, quadquad )
"""
quad: [[171  72]
 [331  93]
 [333 188]
 [177 210]]
square to quad: [[ 171.   72.]
 [ 331.   93.]
 [ 333.  188.]
 [ 177.  210.]]
quad to square: [[-0.  0.]
 [ 1.  0.]
 [ 1.  1.]
 [ 0.  1.]]
quad to quad transform: [[   1.29   -0.23   -0.  ]
 [  -0.06    1.79   -0.  ]
 [-217.24  -88.54    1.34]]
quad to rect: [[   0.    0.]
 [ 300.    0.]
 [ 300.  250.]
 [   0.  250.]]
"""
Denis
@Denis I'm struggling and I've updated the question, but thx for the excellent lead
Will
Thank you for showing the workings
Will
@Will, you're welcome. What did you end up using ?
Denis
@Denis the threeblindmiceandamonkey code with my own supersampling; the opencv bilinear-etc parameter gives me an uneasy feeling
Will