extracting a quadrilateral image to a rectangle

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.


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.]]
"""