You are currently browsing the monthly archive for December 2012.


I’m in Sweden for Christmas and new year, firefly escaping the miserably wet weather at home. Until normal service is resumed, here’s a snowy scene from Stockholm.

Happy Christmas!

With the end of the Mayan Long Count Calendar today, it seems only fair to join the melee and talk a little about the end of the world, oh – excuse me: the End of the World. It at least deserves to be capitalised. First of all let me say that if you genuinely believe the world is going to end I am not out to mock you (although why you would spend your last hours on the internet reading my blog is not entirely clear). As ever, I’m just going to talk about some thoughts I’ve had on the subject and explain what I believe is a rational response to prophecies of impending doom.

But more of that later. Prophecies of impending doom seem to come around every so often, my first encounter with the idea was in about 1984 when I was 9. It was round about the end of September, a documentary on BBC2 announced that one particular prophecy had it that the end of the world was nigh. “in fact,” claimed the presenter, “very nigh: next Thursday”. That, it is safe to say, scared the 9 year old me very much indeed – so much so that I had to have the day off school the following Thursday. As I recall, we had a class outing planned to Saddlers-Wells theatre that day which I missed as a result. I was kind of looking forward to it at the time.

Needless to say, the world did not end but I learned a lesson that day about uncertainty and not believing everything you see on TV and the fact that adults don’t always necessarily know what they’re talking about. Since then predictions of the end of the world have come and gone, and I’ve developed what I believe is a sane a sensible response to them. It is as follows: You think the world is going to end on Saturday? I’ll bet you £50 it won’t. Or an equivalent quantity in your local currency.

This is based on the fact that, in the absence of evidence that a doomsday scenario is on the cards, I don’t believe that it will. If the world is still there, I’m £50 up. If the world ends, well I’m not likely to have to make good on the debt and this fact may actually be of some small comfort as the meteor slams into the Earth’s crust / the tail of the comet seeds a deadly virus / the antichrist, harbinger of the unrighteous, punishes the unrighteous.

This is the essence of Game theory. I play a certain strategy, I get a certain payoff in a certain circumstance. In this case I choose to bet. If I win the bet I am richer, if I lose the bet everyone is dead. On the other had I could choose not to bet, in which case if the world survives I am no richer but if the world ends everyone still dies. A rational agent should therefore choose to make the bet. Of course, this implies that someone is willing to bet with me, but seeing as this hypothetical person must presumably believe that the world is going to end they have nothing to lose and are probably thinking about rewarding ways to spend their last hours, in which case a bit of a wager might be fun, and the prospect of me having to admit that all my fancy book learnin’ hasn’t helped me one bit might also soften the blow a little.

I could look at this another way: way do I believe that the world isn’t going to end? You could argue that none of the previous predictions have proved true so why should this one. This, however, is very faulty reasoning. The end of the world is something that that, by definition, will only happen once. So the fact that it is still here after all those false predictions is just confirmation bias – if any of them had proved true we would not be here to ask the question.

This is where we get into some interesting questions about statistics and how we interpret probability. Classically, probability was interpreted as a frequency of events: given a large number of trials you would expect one result x percent of the time, and other y percent of the time and so on. This is called a frequentist approach, and it works well when you are dealing with lots of repeatable events, like coin tosses. If I flip a coin, I expect heads 50% of the time and tails 50% of the time (in a vanishingly small number of cases it might land on its side, but we’ll ignore that), so the probability of heads is 0.5.

This idea works well for coins tosses, but less well for the end of the world. The end of the world won’t repeat, but I would still like to attach a probability to it happening: surely on a day when there is no prophecy or threat of nuclear destruction the probability is lower than a day when there is? Here we get into the realm of Bayesean statistics. Bayesean statistics is built on the idea of condition probabilities: “given that this has happened, the probability that that will happen is…” and allows you to include prior information, such as the amount of evidence for a particular outcome. This means that I can say things like

“given the recent observation of a huge fiery comet of doom headed directly for central London, the probability of the end of the world is…”


“given that major changes in various numbering systems have come and gone without incident and the evidence of any effect of overflow in numbering time based on arbitrary start dates and dynamical range is thin at best, the probability of the end of the world on 21st December 2012 is not statistically different than in was on the 20th or 22nd Decamber 2012”.

So a combination of Baysean statistics, with its ability to incorporate prior knowledge into estimates of probability, and Game Theory, with its notion of strategy and expectation, I will therefore bet you £50 that the world does not end today. In fact, I’m not the first one to come to this conclusion. See you tomorrow in the 13th Long Count period!

It’s been a while since I posted anything techie, so this week I thought I’d revisit the visualiser I mentioned a while back. Before we continue, however, a short disclaimer: this post contains maths and (a little) code. The maths involves nothing more strenuous than vector cross-products, and the code is a few lines of JoGL but if those things don’t float your boat then by all means feel free to ignore the rest of this post. I won’t hold it against you.

So, that said: I’ve been doing some work adding interactivity to Seer, the visualiser for the Camino Monte-Carlo simulation. The visualiser is primarily a public engagement and demonstration tool, but it’s also pretty handy for debugging. I’ve talked about it before here and here. What the Seer does is visualise a live diffusion MRI simulation. It shows the diffusion environment and the positions of diffusing particles. It currently looks like this:

visualising walkers on a mesh

What we’ve got here is the tissue mesh I’ve talked about before with spins executing random walks rendered in red. Their positions are update by a Camino simulation running on the same mesh in a separate thread. The smaller plot in the bottom-left is a 3D scatter plot of the net displacement of each spin, also updated live. This shows the sort of data that we would measure from the simulation using a diffusion MRI pulse sequence.

Point picking & JoGL code

What we decided, though, was that it needed a bit more interactivity. Specifically we wnt to be able to reset all the spins to a location that you select with the mouse. Since I’m already using a single click to move, the right mousebutton to reset the visualisation, and the mouse wheel to zoom, I decided to go with a double-click to set spin positions.

This presents an interesting challenge, though. How do you translate the mouse coordinates (which are 2D) into a position in 3D space? Setting out I had this plan about projecting a vector into the OpenGL viewport and checking for intersections with the mesh, sorting them along the arclength and then projecting into visualisation (modelview) coords and then on to simulation (substrate) coordinates. What was quite nice, though, was that it turns out thatOpenGL, or rather the GLU API, does quite a bit of this for you.

Point picking works by taking the coordinates of the pixel you click on, converting to a position in the plane at the frot of the view frustrum (this can be done in 2D), then projecting into the scene along the current z-axis until you hit something. You then use the z-coordinate of the object as your third coordinate. This gives you a 3D point that you then project into the model coordinates via the current projection and modelview matrices. GLU provides methods to do this, specifically they’re called glReadPixels() and gluUnProject(). There’s an excellent tutorial on NeHe’s website here.

Because this is a tutorial for OpenGL in C/C++, I’ll also add my code snippet in JoGL:

public final void resetWalkerPositions(GLAutoDrawable drawable, Mesh mesh){

GL gl= drawable.getGL();

IntBuffer viewport= BufferUtil.newIntBuffer(4);
DoubleBuffer modelview= BufferUtil.newDoubleBuffer(16);
DoubleBuffer projection= BufferUtil.newDoubleBuffer(16);

gl.glGetIntegerv(GL.GL_VIEWPORT, viewport);

int winx= walkerX;
int winy= viewport.get(3)-walkerY;

FloatBuffer posZ= BufferUtil.newFloatBuffer(1);
DoubleBuffer pos= BufferUtil.newDoubleBuffer(3);

gl.glReadPixels(winx, winy, 1, 1, GL.GL_DEPTH_COMPONENT, GL.GL_FLOAT, posZ);

gl.glGetDoublev(GL.GL_MODELVIEW_MATRIX, modelview);
gl.glGetDoublev(GL.GL_PROJECTION_MATRIX, projection);

glu.gluUnProject((double)winx, (double)winy, (double)posZ.get(0), modelview, projection, viewport, pos);

// transform into substrate coords
boolean onMesh= mesh.GLtoSubstrate(pos, subsCoords);

// if the coordinates are on the mesh, reset walker positions
// tell the simulation thread to reset walker positions
System.err.println("coords are off-mesh. walkers not reset.");


In addition to the use of gluUnProject(), there’s one additional JoGL-specific issue here: the GL object itself. The way I’d designed the code meant that the method that catches the double-click event was nowhere near the rendering code that does the unprojection and talks to the simulation. I spent a bit of time trying to get hold of a GL object and hand it over to the event handler, but nothing I tried worked so instead I realised that all the event handler actually needed to do was to provide the mouse coordinates and instruct the render method to do the rest. So all it does is set a flag and hand over the coords via a class-level variable. That’s a theme that’s emerged a little recently: making instructions to other parts of the system via global flags rather than method calls. It works pretty well when you’ve got functionality that’s spread across different parts f the code. (I suppose I could also have used static variables but the principle is the same and this way things are more self-contained).

Planes, projection and a bit of maths

So: sorted. Well, actually no. Unfortunately The meshes that I’m clocking on have a lot of holes in them, and sometimes I was to click on a hole instead of a triangle. In this case, glUnProject() gives a point at infinity, which isn’t what I want. I want a point half way across my substrate. This means there’s a special case to catch. Fortunately, points at infinity are easy enough to catch as the coordinate will be equal to 1, but what to do once you’ve caught it?

Firstly, we need to recognise that this is essentially a projection into a plane. The plane in question bisects the substrate half way along its z-axis ans so is easily defined but in viewport coords will depend on the current view of the mesh (the modelview matrix). Given a plane ax + by + cz + d =0 , and a point \left( X, Y, Z \right) we just choose a new z-coord such that

Z' = -\frac{aX + bY +d}{c}

The tricky part is knowing what a b , c and d are. My initial thought was to back rotate into substrate coords and project into the appropriate plane, but this requires you to invert the modelview matrix, which frankly a cannot be bothered to write code to do (and in any case is an expensive operation) so I need to be working in viewport coordinates, not modelview coordinates. So then I thought I’d use the modelview matrix to rotate the plane normal but it turns out that plane normals transform with the inverse of the rotation matrix so once gain we’re back to square one.

The answer is to define the plan using three points and use the cross product to get the plane normal. Any three non-collinear points define a plane. These points transform using the modelview matrix, not the inverse, and the the components of the normal to the plane are the coefficients we want. The algebr works out like this [Cracks knuckles], [flexes fingers in the manner of concert pianist]:

\hat{\mathbf{n}}= \frac{\mathbf{n}}{|\mathbf{n}|}

\mathbf{n} = \left(\mathbf{v}_3 \times \mathbf{v}_1\right)-\left(\mathbf{v}_2 \times \mathbf{v}_1\right)-\left(\mathbf{v}_3 \times \mathbf{v}_2\right)


(a - b) \times (c - b) = a \times (c-b) - b \times (c-b) = (c-b) \times b - (c-b) \times a

(a - b) \times (c - b) = c \times b - b \times b - c \times a + b \times a

(a - b) \times (c - b) = c \times b - c \times a + b \times a

and we’re away. Three cross products and no matrix inversion.

I’ll call it a day there. I’ll post some code snippets once they’re done.

I was in a tile shop the other day and it got me thinking about the fourth dimension. What do you mean “what are you talking about?”, I would have thought the connection was obvious. Oh, all right then…

Anyway, I was in a tile shop looking at different tiling patterns. Many of them were very lovely and its impressive to see how a good tile shop can embrace patterns form so many different countries and cultures. This shop had a fine line in British, Spanish, Italian, Moroccan, Arabic, Greek, you name it. Bathroom designs to satisfy even the most demanding and culturally promiscuous time-traveler.

Except after a while you start to notice that all these patterns have something in common: they are all periodic. That is to say, all the patterns repeat. In maths this is known as a translational symmetry – if you slide a repeating pattern around it will eventually line up with the tiles a few feet over. This could be as simple as a checkerboard pattern or something much more complex but typically a tiling that perfectly fills the floor space available repeats (i.e. it doesn’t have any gaps in it) will be repetitive.

There’s actually a whole branch of maths that investigates tilings and symmetries. It’s provable that certain patterns will fill the plane perfectly (like squares or hexagons) and other (like pentagons or octagons) won’t, and that certain combinations of shapes (like octagons and squares) can be used together to perfectly tile the plane. It’s related to the wider study of symmetry and a field called Group Theory, which turns out to be hugely important in physics, cropping up in everything from molecular chemistry and spectroscopy to particle physics. Heady stuff.

Anyway, what’s all this got to do with the 4th dimension. Well, aside from culturally promiscuous time-travel (which I’ll cover in a separate post). The fourth dimension is link to tiling in a rather surprising way. I mentioned that tilings tend to be periodic, but does that mean that EVERY tiling of the plane has to be periodic? The idea is related to a thing called the Domino problem in maths. Specifically, given a set of shapes, is it possible to design an algorithm that will decide if they can tile the plane or not. During the 1960s, a mathematician called Hao Wang suggested that this problem was solvable if all tilings of the plane were in periodic. You’d just need to decide if your set of tiles could be arranged periodically or not and you’d know the answer. Neat, huh?

Except there’s a snag. What if there are tilings of the plane that aren’t periodic? If those existed the test would fail. In 1966 a non-periodic tiling was found: it used over 20000 different tiles (pretty hard to construct, and all the more impressive for being found before powerful computers were widespread). Then anther with 104 tiles was found. Then another with 40 and another with just 13. Finally in 1974 Roger Penrose found an aperiodic tiling that required only 2 tiles. These are pretty interesting patterns: patterns that cover the plane by don’t (quite) repeat. You can’t slide the pattern around and match it to itself, and it only contains two types of tile! It looks like this:

penrose tilingAt first it looks regular, but stare at it for a while and you can see that it doesn’t quite repeat. Also notice that this piece doesn’t form a repeating unit on its own either. It’s a pretty cool thing. What turns out to be extra interesting, though, is that tiling patterns are related to packing problems in 3D, and in particular to the arrangement of atoms in crystals. You could ask the same question here: if I want to pack a load of (spherical) atoms into a 3D space, does the pattern have to be repetitive?

Well, no. crystals are periodic, glasses are disordered arrangements of atoms. More organic forms like wood are the result of nonlinear growth processes and many rocks are made up of mixtures of crystals, glasses and dusts all mixed up together but none of these have the same properties as the Penrose Tiling. There are things that do, though: they’re called quasicrystals, and they have a whole field of research attached to them.

A quasi-crystal is a regular but non-periodic packing of atoms. Using electronic force microscopy you can actually make images of the individual atoms and the surfaces turn out to look a bit like this:

quasicrystalLook familiar? This is actually the surface of an Aluminium-Palladium-Manganese quasicrystal. Looks an awful lot like the Penrose Tiling, don’t you think?

But what about the 4th dimension? Well, it turns out that the best way to think about aperiodic tilings is to thing of them as a 3D slice through a regular 4D lattice. If the slice is parallel to one of the regular planes in 4D, it’s periodic. Any non-parallel slice is aperiodic.

What’s a 3D slice? Think of a line. A line is 1D, and you can cut it with a single point. Similarly, a square is 2D and can be cut with a line, which is 1D. A cube is 3D and can be cut with a plane, which is 2D. There’s a pattern here: an nD volume can be sliced with an (n-1)D object. So a 4D object (don’t worry about picturing it!) can be sliced with a 3D volume. Mathematically this a relatively easy thing to do, and so your 3D quasicrystal is a the 3D slice through a 4D object.

This idea of slicing through a 4D volume is something that crops up in my own work. I do a bit of work with 3D graphics. Here is turns out to be convenient to think of your 3D space containing your graphics as a slice through 4D. Why is that, you ask? Well, in 3D I want to move around and also rotate around a point. A rotation in 3D is written as a 3×3 matrix. In fact you can think of any rotation as the result of three other rotations: one about each of the x, y, and z axes. This is provable using – guess what – Group Theory.

The thing is, that once your 3×3 matrix is full of rotations (any other things like scale factors and sheers) there no room for the translations that move things around. The most natural thing to now is to make the matrix bigger to include extra elements for moving around, and this is equivalent to (guess what) using a fourth dimension. So, weirdly, 3D geometry turns out to be easier to think of in 4D. The 3D graphics in your favourite computer game are in fact a slice through a 4D space, just like the tilings and packings for the quasicrystals. The guy in the tiling shop seemed quite interested in this when I mentioned it, although I did steer clear of projective geometry and quasicrystals at the time. It was Sunday lunchtime.

One more thing before I go. I was talking to Mrs of-Science about this later on that day, and she immediately asked if I was talking about time as the 4th dimension. I’m not – all of this is about a 4th space-like dimension, and is very different from the geometry of relativity’s 4D Spacetime, which is fascinating in a completely different way and not at all like regular Euclidean space. Interestingly, though, you can add a 4th space dimension to general relativity. Kaluza ad Klein did it back in the 30s. It turns out that if you do this you get a unification of gravity and classical electromagnetism, which was the start of all the work on grand unification of the the four fundamental interactions that we hear so much about. So the fourth dimension reveal secrets here too.

Life in 4D is pretty cool, and I fully expect my bathroom to look pretty good too.