| Andrew Cooke | Contents | Latest | RSS | Twitter | Previous | Next


Welcome to my blog, which was once a mailing list of the same name and is still generated by mail. Please reply via the "comment" links.

Always interested in offers/projects/new ideas. Eclectic experience in fields like: numerical computing; Python web; Java enterprise; functional languages; GPGPU; SQL databases; etc. Based in Santiago, Chile; telecommute worldwide. CV; email.

Personal Projects

Lepl parser for Python.

Colorless Green.

Photography around Santiago.

SVG experiment.

Professional Portfolio

Calibration of seismometers.

Data access via web services.

Cache rewrite.

Extending OpenSSH.

Last 100 entries

Chinese Govt Intercepts External Web To DDOS github; Numbering Permutations; Teenage Engineering - Low Price Synths; GCHQ Can Do Whatever It Wants; Dublinesque; A Cryptographic SAT Solver; Security Challenges; Word Lists for Crosswords; 3D Printing and Speaker Design; Searchable Snowden Archive; XCode Backdoored; Derived Apps Have Malware (CIA); Rowhammer - Hacking Software Via Hardware (DRAM) Bugs; Immutable SQL Database (Kinda); Tor GPS Tracker; That PyCon Dongle Mess...; ASCII Fluid Dynamics; Brandalism; Table of Shifter, Cassette and Derailleur Compatability; Lenovo Demonstrates How Bad HTTPS Is; Telegraph Owned by HSBC; Smaptop - Sunrise (Music); Equation Group (NSA); UK Torture in NI; And - A Natural Extension To Regexps; This Is The Future Of Religion; The Shazam (Music Matching) Algorithm; Tributes To Lesbian Community From AIDS Survivors; Nice Rust Summary; List of Good Fiction Books; Constructing JSON From Postgres (Part 2); Constructing JSON From Postgres (Part 1); Postgres in Docker; Why Poor Places Are More Diverse; Smart Writing on Graceland; Satire in France; Free Speech in France; MTB Cornering - Where Should We Point Our Thrusters?; Secure Secure Shell; Java Generics over Primitives; 2014 (Charlie Brooker); How I am 7; Neural Nets Applied to Go; Programming, Business, Social Contracts; Distributed Systems for Fun and Profit; XML and Scheme; Internet Radio Stations (Curated List); Solid Data About Placebos; Half of Americans Think Climate Change Is a Sign of the Apocalypse; Saturday Surf Sessions With Juvenile Delinquents; Ssh, tty, stdout and stderr; Feathers falling in a vacuum; Santiago 30m Bike Route; Mapa de Ciclovias en Santiago; How Unreliable is UDP?; SE Santiago 20m Bike Route; Cameron's Rap; Configuring libxml with Eclipse; Reducing Combinatorial Complexity With Occam - AI; Sentidos Comunes (Chilean Online Magazine); Hilary Mantel: The Assassination of Margaret Thatcher - August 6th 1983; NSA Interceptng Gmail During Delivery; General IIR Filters; What's happening with Scala?; Interesting (But Largely Illegible) Typeface; Retiring Essentialism; Poorest in UK, Poorest in N Europe; I Want To Be A Redneck!; Reverse Racism; The Lost Art Of Nomography; IBM Data Center (Photo); Interesting Account Of Gamma Hack; The Most Interesting Audiophile In The World; How did the first world war actually end?; Ky - Restaurant Santiago; The Black Dork Lives!; The UN Requires Unaninmous Decisions; LPIR - Steganography in Practice; How I Am 6; Clear Explanation of Verizon / Level 3 / Netflix; Teenage Girls; Formalising NSA Attacks; Switching Brakes (Tektro Hydraulic); Naim NAP 100 (Power Amp); AKG 550 First Impressions; Facebook manipulates emotions (no really); Map Reduce "No Longer Used" At Google; Removing RAID metadata; New Bike (Good Bike Shop, Santiago Chile); Removing APE Tags in Linux; Compiling Python 3.0 With GCC 4.8; Maven is Amazing; Generating Docs from a GitHub Wiki; Modular Shelves; Bash Best Practices; Good Emergency Gasfiter (Santiago, Chile); Readings in Recent Architecture; Roger Casement; Integrated Information Theory (Or Not); Possibly undefined macro AC_ENABLE_SHARED; Update on Charges; Sunburst Visualisation

© 2006-2013 Andrew Cooke (site) / post authors (content).

Cube - Series of Images for Laser Printer

From: "andrew cooke" <andrew@...>

Date: Fri, 26 Dec 2008 17:29:04 -0300 (CLST)

I've been messing around with Python (using NumPy for the calculations and
PyScript for the final Postscript rendering), generating a series of
images based on a cube that is "wrapped" by one or more lines.

Initially points on the surface of a cube are generated in 3D - this is
implemented using a Python generator.  Then they are filtered by a
"pipeline" of transforms that filter, project to 2D, and scale/rotate.

The initial generator takes as parameters further generators, one for each
line).  So, for example, the length of each line might be given by
"constant(10)", or "repeated_uniform(1, 5)" (repeated sampling from a
uniform distribution from 1 to 5), or [1, 2, 3] (lengths for 3 lines).

The transforms are not necessarily generators - they are functions that
are composed.  This is because some transforms require the entire set of
points to function correctly.  But in practice most, those that operate on
individual points, are implemented as generators.


I'm wondering about printing some out on "nice" paper and seeing if I can
get a gallery to stock them.  Is that crazy?  I have no idea...



From: "andrew cooke" <andrew@...>

Date: Sat, 27 Dec 2008 08:43:40 -0300 (CLST)

The conversion from 3D (spatial model) to 2D points (paper plot) is worth
recording.  This is a summary of what I learnt from googling around - it
may not be correct (the info on the web about projective geometry isn't
completely consistent).

There are two common kinds of projections - orthographic and perspective.

Orthographic projections are calculated using an affine transformation - 
a linear calculation using a matrix multiplication (with a shift that can
be implemented using an additional vector addition, or by using an extra
dimension with "1" extending the vectors).

Orthographic projections can be further divided into three groups.  First,
those perpendicular to "principal" planes (like plans), which are not very
interesting here.  Second, axonometric projections which show a "skew"
view.  These are equivalent to a perspective projection "from infinity"
(so parallel lines remain parallel in projection) and the most common is
probably the isometric view.  Third, oblique projections where one plane
is shown "face on" while others are skewed.  The most common here is
"cavalier perspective" with an angle of 45 degrees and no scaling of

My code for the orthographic projections is:

  def _orthographic(transform, offset):
    t = np.asarray(transform)
    o = np.asarray(offset)
    return lambda p: np.dot(t, p) + o

  def cavalier(offset=(0,0)):
    rsq2 = 1 / sqrt(2)
    return per_point(
	transform=((-1,rsq2,0), (0,-rsq2,1)),

  def isometric(offset=(0,0)):
    s32 = sqrt(3) / 2
    return per_point(

where "per_point" is the function to convert from a pointwise transform to
one for the whole image (image is a list of paths; paths are lists of

I couldn't find a simple description  of how to calculate a perspective
projection, but it's simple to derive by hand if you construct the
projection by taking the figure generated by rays from eye to object that
intersect a plane (the "image plane").

I imagine the image plane as a sheet hanging in space, through which I can
vaguely see the object that is to be projected; the perspective projection
is just tracing the outline seen.

For an "undistorted" view the image plane is perpendicular to the line of
sight to the "centre" of the object (ie the point on the object that will
be at the centre of the image).  Also, vertical lines in the object will
be vertical in the image.  This is sufficient to define the geometry if we
place the eye at E (capital letters for vectors) and the centre of the
image at C (so the image plane contains C and is perpendicular to EC).

I will use, for example, EC as shorthand for C-E - the vector from E to C
(the reversal of order between "EC" and "C-E" is a real pain).

There are two coordinate systems in use here, confusingly: Z is the
spatial vertical direction, while X and Y are the horizontal and vertical
directions in the image plane.  From the discussion above
("..undistorted...") we want Y to appear (from E) parallel to Z.  We can
use this to define X and Y in terms of C, E and Z.

First, let us define unit(V) = V / |V| (where |V| is the length of V, an
arbitrary vector) and N = unit(CE) (the normal to the image plane,
pointing towards E).

Then X = unit(Z^N) since it is perpendicular to N (in the image plane) and
perpendicular to the vertical.  I am using "^" for the "cross" or "vector"

Similarly, Y = N^X

Now we can do the projection.  The eye looks towards a point P on the
object.  The line EP intersects the image plane at Q (there is no need for
Q to be between E and P - it may lie on the extended line, in which case
we get weird "fisheye" effects).

So the projection consists in finding Q for each P on the object and then
resolving Q along X and Y.

How to find Q?

We know that Q lies on EP, so for some parameter t, Q = E + t EP (a simple
parametrised line).  We also know that Q is in the image line and so is
perpendicular to N: 0 = N.Q (where "." is the "dot" or "scalar" product).

We can therefore solve for t:
Q = E + t EP
0 = N.Q
0 = N.(E + t EP)
-(N.E) = t (N.EP)
t = -(N.E) / (N.EP)

So Q = E - EP (N.E) / (N.EP)

The final projection to 2D coords is then just the dot product "along" X
and Y: (x, y) = (Q.X, Q.Y).

My code is a bit ugly, but seems to work:

  def perspective(eye, screen):
    def unit(v):
      return v / sqrt(np.dot(v, v))
    eye = np.asarray(eye)
    screen = np.asarray(screen)
    normal = eye - screen
    normal2 = np.dot(normal, normal)
    n = normal / sqrt(normal2)
    # screen axes
    x = unit(np.cross(Z, n))
    y = np.cross(n, x)
    def process(p1):
      los = p1 - eye
      t = -normal2 / np.dot(los, normal)
      p2 = eye + t * los # point in the plane
      return (np.dot(p2, x), np.dot(p2, y))
    return per_point(process)


Comment on this post