Sister blog of Physicists of the Caribbean. Shorter, more focused posts specialising in astronomy and data visualisation.
Showing posts with label Visualisation. Show all posts
Showing posts with label Visualisation. Show all posts

Monday, 17 February 2025

Sports Stars Can Save Humanity

I know, I know, I get far less than my proverbial five-a-day so far as reading papers goes. Let me try and make some small amends.

Today, a brief overview of a couple of visualisation papers I read while I was finishing off my own on FRELLED, plus a third which is somewhat tangentially related.


The first is a really comprehensive review of the state of astronomical visualisation tools in 2021. Okay, they say it isn't comprehensive, which is strictly speaking true, but that would be an outright impossible task. In terms of things at a product-level state, with useable interfaces, few bugs and plenty of documentation, this is probably as close as anyone can realistically get.

Why is a review needed ? Mainly because with the "digital tsunami" of data flooding our way, we need to know which tools already exist before we go about reinventing the wheel. As they say, there are data-rich but technique-poor astronomers and data-poor but technique-rich visualisation experts, so giving these groups a common frame of reference is a big help. And as they say, "science not communicated is science not done". The same is true for science ignored as well, of which I'm extremely guilty... you can see from the appallingly-low frequency of posts here how little time I manage to find for reading papers. 

So yeah, having everything all together in one place makes things very much easier. They suggest a dedicated keyword in papers "astrovis" to make everything easier to find. As far as I know this hasn't been adopted anywhere, but it's a good idea all the same.

Most of the paper is given to summarising the capabilities of assorted pieces of software, some of which I still need to check out properly (and yes, they include mine, so big brownie points to them for that !). But they've also thought very carefully about how to organise all this into a coherent whole. For them there are five basic categories for their selected tools : data wrangling (turning data into something suitable for general visualisation), exploration, feature identification, object reconstruction, and outreach. They also cover the lower-level capabilities (e.g. graph plotting, uncertainty visualisation, 2D/3D, interactivity) without getting bogged-down in unproductively pigeon-holing everything. 

Perhaps the best bit of pigeon-unholing is something they quote from another paper : the concept of explornation, an ugly but useful word meaning the combination of exploration and explanation. This, I think, has value. It's possible to do both independently, to go out looking at stuff and never getting any understanding of it at all, or conversely to try and interpret raw numerical data without ever actually looking at it. But how much more powerful is the combination ! Seeing can indeed be believing. The need for good visualisation tools is not only about making pretty pictures (although that is a perfectly worthwhile end in itself) but also in helping us understand and interpret data in different ways, every bit as much as developing new techniques for raw quantification. 

I also like the way they arrange things here because we too often tend to ignore tools developed for different purposes other than our own field of interest. And they're extraordinarily non-judgemental, both about individual tools and different techniques. From personal experience it's often difficult to remain so aloof, to avoid saying, "and we should all do it this way because it's just better". Occasionally this is true, but usually what's good for one person or research topic just isn't useful at all for others.

On the "person" front I also have to mention that people really do have radically different preferences for what they want out of their software. Some, inexplicably, genuinely want everything to do be done via text and code and nothing else, with only the end result being shown graphically. Far more, I suspect, don't like this. We want to do everything interactively, only using code when we need to do something unusual that has to be carefully customised. And for a long time astronomy tools have been dominated too much by the interface-free variety. The more that's done to invert the situation, the better, so far as I'm concerned.


The second paper presents a very unusual overlap between the world of astronomy and... professional athletes. I must admit this one languished in my reading list for quite a while because I didn't really understand what it was about from a quick glance at the abstract or text, mostly because of my own preconceptions : I was expecting it to be about evaluating the relative performance of different people at source-finding. Actually this is (almost) only tangential to the main thrust of the paper, though it's my own fault for misreading what they wrote.

Anyway, professional sports people train themselves and others by reviewing their behaviour using dedicated software tools. One of the relatively simple features that one of these (imaginatively named "SPORTSCODE") has is the ability to annotate videos. This means that those in training can go back over past events and see relevant features, e.g. an expert can point out exactly what and where something of interest happened – and thereby, one hopes, improve their own performance.

What the authors investigate is whether astronomers can use this same technique, even using the same code, to accomplish the same thing. If an expert marks on the position of a faint source in a data cube, can a non-expert go back and gain insight into how they made that identification ? Or indeed if they mark something they think is spurious, will that help train new observers ? The need for this, they say, is that ever-larger data volumes threaten to make training more difficult, so having some clear strategy for how to proceed would be nice. They also note that medical data, where the stakes are much, much higher, relies on visual extraction, while astronomical algorithms have traditionally been... not great. "Running different source finders on the same data set rarely generates the same set of candidates... at present, humans have pattern recognition and feature identification skills that exceed those of any automated approach."

Indeed. This is a sentiment I fully endorse, and I would advocate using as much visual extraction as possible. Nevertheless, my own tests have found that more modern software can approach visual performance in some limited cases, but a full write-up on that is awaiting the referee's verdict.

While this paper asks all the right questions, it presents only limited answers. I agree that it's an interesting question as to whether source finding is a largely inherent or learned (teachable) skill, but most of the paper is about the modifications they made to SPORTSCODE and its setup to make this useful. The actual result is a bit obvious : yes indeed, annotating features is useful for training, and subjectively this feels like a helpful thing to do. I mean... well yeah, but why would you expect it to be otherwise ? 

I was hoping for some actual quantification of how users perform before and after training – to my knowledge nobody has ever done this for astronomy. We muddle through training users as best we can, but we don't quantify which technique works best. That I would have found a lot more interesting. As it is, it's an interesting proof of concept, and it asks all the right questions, but the potential follow-up is obvious and likely much more interesting and productive. I also have to point out that FRELLED comes with all the tools they use for their training methods, without having to hack any professional athletes (or their code) to get them to impart their pedagogical secrets.


The final paper ties back into the question of whether humans can really outperform algorithms. I suppose I should note that these algorithms are indeed truly algorithms in the traditional, linear, procedural sense, and nothing at all to do with LLMs and the like (which are simply no good at source finding). What they try to do here is use the popular SoFiA extractor in combination with a convolutional neural network. SoFiA is a traditional algorithm, which for bright sources can give extremely reliable and complete catalogues, but it doesn't do so well for fainter sources. So to go deeper, the usual approach is to use a human to vet its initial catalogues to reject all the likely-spurious identifications.

The authors don't try to replace SoFiA with a neural network. Instead they use the network to replace this human vetting stage. Don't ask me how neural networks work but apparently they do. I have to say that while I think this is a clever and worthwhile idea, the paper itself leaves me with several key questions. Their definition of signal to noise appears contradictory, making it hard to know exactly how well they've done : it isn't clear to me if they're really used the integrated S/N (as they claim) or the peak S/N (as per their definition). The two numbers mean very different things. It doesn't help that the text is replete with superlatives, which did annoy me quite a bit.

The end result is clear enough though, at least at a qualitative level : this method definitely helps, but not as much as visual inspection. It's interesting to me that they say this can fundamentally only approach but not surpass humans. I would expect that a neural network could be trained on data containing (artificial) sources so faint a human wouldn't spot them, but knowing they were there, the program could be told when it found them and thereby learn their key features. If this isn't the case, then it's possible we've already hit a fundamental limit, that when humans start to dig into the noise, they're doing about as well as it's ever possible to do by any method. When you get to the faintest features we can find, there simply aren't any clear traits that distinguish signal from noise. Actually improving, in any significant way, on human vision, might be a matter of a radically different approach... but it might even be an  altogether hopeless challenge.

And that's nice, isn't it ? Cometh the robot uprising, we shall make ourselves useful by doing astronomical source-finding under the gentle tutelage of elite footballers. 

Or not, because that algorithms can be thousands of times faster can more than offset their lower reliability levels, but that's another story.

Phew ! Three papers down, several hundred more to go.

Thursday, 11 July 2024

ChatGPT Is Not A Source Extractor

When ChatGPT-4o came along I was pretty keen to try out its shiny new features, especially since some of the shine has rubbed off the chatbots of late. Oh, ignore the hype trains completely : those who are saying it's going to cause the apocalypse or usher in the Utopian end of history are equally deluded. I'm talking about actual use cases for LLMs. This situation remains pretty much as it has been since they were first unleashed. That is...
  • Decent enough if you want free-form discussions (especially if you need new ideas and don't care too much about factual accuracy)
  • Genuinely actually very useful indeed for coding (brilliant at doing boiler-plate work, a serious time-saver !)
  • Largely crap if you need facts, and even worse if you need those facts to be reliably accurate
Pretending that those first two are unimportant is in my view quite silly, legitimate concerns about energy expenditure notwithstanding. But that third one... nothing much seems to have shifted on that at all. They're still plagued with frequent hallucinations, since they're not grounded in anything so that they have no internal distinction between a verifiable, observable truth and the CPU-equivalent of a random brain-fart firing of the neurons.

Unfortunately GPT-4o just seems to extend this into a multi-modal world, giving results basically consistent with my earlier tests of chatbots. But I was intrigued by its apparent accuracy when supplying image files. It seemed to be, albeit from limited testing, noticeably more accurate when asked questions about image files than, say, PDFs. So I had a passing thought : could I use ChatGPT-4o to find sources in my data ?

Spoiler : no. It doesn't work.

It's not possible to share the chat itself because it contains images, but basically what I did was this. I uploaded an image of a typical data set I would customarily trawl look looking for galaxies. The very short version is that the HI detections of galaxies typically look like elongated blobs, sometimes appearing saturated and sometimes as mere enhancements in the noise. You can find a much more thorough explanation on my website, but that's the absolute basics. For example, in the image below, there are seven very obvious detections and one which is a bit fainter. 

I began by giving ChatGPT a detailed description of the image and the task at hand. This is the kind of thing that takes a few minutes to explain to a new observer; the actual training of data inspection can take a few days, but the explanations need be only very short indeed. And finding the bright sources is trivial : almost anyone can do that almost immediately. The bright galaxies are inherently obvious in the data when presented like this. Even if you have no idea what the axes labels refer to, it's clear that some parts of the image are very different to the others.

I asked ChatGPT to mark the location of the sources or otherwise describe their position. It didn't mark them but instead gave descriptions. Its world coordinates weren't precise enough to verify what it had identified, however, being limited to only values directly readable in the image and not doing any interpolation. I also gave it a broad alpha-numeric grid (A-J along the x-axis and 1-6 along the y-axis), but this was too coarse to properly confirm what it thought it had found. 

Its results were ambiguous at best. Even with this coarse grid it was clear some of its results were simply wrong. So I did what I'd do with new observers. I marked the sources with red outlines and numbers, uploaded the new image and described what I'd done, so it would have some kind of reference image. I also described the sources in more detail, e.g. which ones were bright and which were faint, and whether they extended into adjacent cells.

Next I gave it a new image with a finer grid (A-O and 1-14). This time, two sources (out of the ten or so visible) were reported correctly while the rest were wrong.  By mistake, I missed out the "D" cell in the coordinate labels, but ChatGPT reported a source at D4 ! Its revised claims were still wrong though, with once again getting only two correct.

This wasn't going well. I decided to dial it back and try something simpler. Maybe ChatGPT was able to "see" the features but not was accurately reading the coordinates, or perhaps hallucinating its answers and so mangling its results. So now I uploaded an image devoid of any coordinates and asked it for a simple count of the number of bright blobs. It got the answer right ! Okay, better... I asked it if it could mark the locations directly on the image, but it said it couldn't edit images. Instead it suggested giving the coordinates of the sources as a percentage of the axis length from the top left. Fair enough, but when comparing its reported coordinates it had again two near-misses and got all the rest simply wrong.

Finally I decided to check if at least the reported number count wasn't just a fluke. I uploaded three images in one file (thus circumventing OpenAI's painfully-limited restrictions on the free plan), each labelled with a number, and asked for the number of sources in each. It got one right and the rest wrong. It also gave descriptions of where it thought the sources were (i.e. upper left, middle, that sort of thing) and these were all wrong. Then, rather surprisingly and quite unprompted, it decided that it actually could edit images to mark the positions after all. The result came back :


Well... it's less than stellar. 

The upshot is that nothing much has changed about chatbot use cases at all. Good for discussions,  useless for facts. Whether it is "seeing" the images in some sense I don't know : possibly at some level it does recognise the sources but hallucinates both when trying to mark them and describe their positions, or possibly it's just making stuff up and nothing else. The latter seems rather unlikely though. Too often in other tests it was capable of giving results from figures in PDFs and image files which could not have been obtained from reading any of the text, that required actually "looking" at the images. 

Regardless of what it's actually doing, in terms of using ChatGPT as a source extractor, it's a non-starter. It doesn't matter why it gets things wrong, for practical application it only matters that it does. Maybe there's something capable under there, maybe there isn't. For now it's just an energy-intensive way of getting the wrong answers. Well, I could have done that anyway !

Friday, 11 December 2020

Data visualisation : the next level

This is going to be a very strange post in which I describe pretty pictures but don't actually show you what I'm looking at. Why ? Because I'm using a VR headset, and I can't yet make the final media into a shareable format.

I now have a computer capable of VR. This is almost as big a jump as getting the standalone Quest headset itself, because the graphical capabilities of the PC far exceed the high-end smartphone level of the headset alone. 

One of the first things I tried was to examine the tiny handful of models I've uploaded to Sketchfab. On the Quest by itself, these are barely functional. The framerate and/or tracking are lousy, and though you can get the general idea, the experience is unpleasant. Not so with the PC, which easily handles much more complex models than these. So I can walk around my model as though it was actually there in my living room. I can even interactively rescale it with the thumbsticks. But of course, using Sketchfab isn't very convenient, especially given the pathetic limitations imposed on uploads.

That's where FRELLED comes in*. I used Blender 2.79 for this, partly out of ignorance. When I first looked at Eevee, back during the early test builds for 2.80, it wasn't up to much. Textures loaded slowly and even simple tests of FRELLED weren't at all successful, with the view essentially re-rendering whenever anything changed at all. That completely breaks the main benefit of FRELLED, which is that the view should update instantaneously (that is, at > 25 fps) and in real time. Since Blender 2.8+ lacks any other realtime capability, I stuck with 2.79.

* I've now managed to confirm that this works both on Windows and Mac. There's still quite a bit to do, but it's getting closer and closer to being released into the wild.

But nowadays Eevee is massively more powerful. For scenes as simple as those FRELLED constructs, Eevee renders in true real time, not the pseudo-realtime of before. And it doesn't have the problem of transparent materials needing to be ordered correctly, and - maybe best of all - you can adjust the brightness and contrast of materials in realtime as well. That makes it a dramatic and wholesale improvement over Blender 2.79's capabilities, not a poor substitute with a few fringe benefits.

Unfortunately, the switch from the Python in Blender 2.78 to 2.8+ is significant enough that I can't just directly convert everything. So FRELLED version 5 is already looking obsolete compared to a planned version 6, though, mercifully, the conversion will be nowhere near what I've had to do in the upgrade from Blender 2.49 (essentially a complete re-write of all 11,000 lines of code). But there's one capability of Eevee which I simply had to try out : virtual reality.

Getting this to work was remarkably easy. I started with importing isosurfaces and adding a few lights. Then I plugged in the headset, enabled the Link, and hit "start VR server" in Blender. And it just worked. I had a greyscale surface of M33 floating in front of me that I could walk around. Or at least partway, due to the limited length of the cable (I've got the wireless version of Virtual Desktop running, but so far that only works with Steam and I haven't figured out how to run it with Oculus software yet).

From there it was a simple matter of playing with Eevee's materials to get something more shiny. For my purposes, I barely need any lights - I can do everything with the material preview. In an hour or so I had this, floating in front of me. I could even stick my head inside it, though Blender gets sluggish if I do :

Now of course that needs some more effort to come up with nicer materials, but the proof of concept is solid. I was so impressed by how well this worked I began to wonder if it might even be possible to do the full volumetric display of FRELLED. And soon I found that yes, yes it is. I wrote a couple of short Python scripts to automate most of the process. So now I get to see M33 in its full volumetric glory, rendered as a cube about half a metre across that I can walk right round.

The main limitation appears to be proximity. From around 0.5m away the frame rate is very good. Get much closer, though, and it drops sharply. You can stick your head inside, but it's not much fun. I'm not sure why this is but I guess it's a limitation of Blender. Still, even this is more than sufficient for outreach.

And in some ways this is even simpler than the old process. Blender 2.79 had problems ordering materials, so that transparent surfaces weren't rendered correctly from behind. This meant an elaborate series of forward and reverse images with a background script deciding which ones to show based on the viewing angle. This isn't necessary in 2.91, where I can just show everything at once. I might eventually use a simpler version of the script to deal with larger data sets, but for small ones it isn't needed at all.

Is there any practical, scientific benefit to this though ? Honestly, I dunno. Personally I think the more visualisation techniques we have access to, the better. When you actually see it, I don't think there's any question that this is inherently better. The greater immersion helps focus on features you might never have noticed (but I won't really know this until it's developed enough to use in anger). Granted, VR could still benefit from lighter, cheaper headsets at higher resolution, but this was once true of television as well. If there's one piece of tech that's come closest to the sci-fi predictions of the last few decades, then VR is surely a leading contender.

Of course, at this stage it's nice for outreach but useless for science. Still, the current technology appears adequate to the point that it's the software which is now the chief bottleneck. Greater native integration of VR hardware in Blender would be nice, though I'd prefer some format which could be easily shared online without having to give away the actual .blend file*. But it's entirely feasible to conceive of sticking on a headset and doing all the standard analysis in VR, with negligible additional effort compared to using an ordinary monitor. In fact it's already possible, in that this could be developed on a timescale of weeks or months - certainly not years.

The only major practical issue, though, may not be the hardware so much as the space requirement. We're not going to be giving up 2D screens any time soon. Whether anyone will feel that the capabilities of VR are so beneficial as to give it dedicated areas (at least in astronomy) is something we're just going to have to find out by experiment. Personally I think it's something well worth exploring.

Wednesday, 14 October 2020

Coming soon : FRELLED version 5

I had two lockdown projects. One was to develop an interactive model of Arecibo, which I more-or-less have working but just have to find the time and inclination to get into a useable format (which is tedious and boring). The second was to recode FRELLED, my Python script that imports 3D FITS files into Blender. This was originally written for Blender 2.49, released back in 2009 (!), but Blender 2.5 has a completely different Python syntax - more like using another language than making minor modifications. So it took a global pandemic to force me to re-write the bloody thing in a modern version of Blender.

After several months, I'm pleased to announce that this is done. Well, sort of.

The new version uses Blender 2.79. This isn't the very latest version, but for a very good reason. Blender 2.8+ doesn't support the OpenGL realtime shaders that FRELLED relies on, and unfortunately neither Cycles nor Eevee are suitable replacements. Apparently it will get a modern OpenGL equivalent at some point though, and the Python syntax is almost identical to that used in 2.79. This means the next update won't be anything like burdensome as recoding the entire thing again.

FRELLED version 5 looks like this :


This is the main display section with an example cube loaded. Blender's GUI now enables adjustable panels, so the user won't be overwhelmed with information. Presets are now such that loading a cube should be a matter of about five mouse clicks. This, I hope, will be easy enough to persuade people that it's worth installing and using.

Incidentally, installation should now be MUCH simpler. Blender 2.79 comes with its own internal Python and PIP kept completely separate from system Python (you can download it in a zip file, no other installation needed). So installing the modules FRELLED needs is now trivial... at least it was for me on Windows. It should even work on Linux networks.

Loading cubes is now much faster, hence there's generally not so much need to worry about which projections are being imported, so this is all hidden by default. But all that is still there for enthusiasts and those using very large or weird data sets. In particular, the "sparse sampling" option now lets you import only every nth slice of the data, adaptive to the size of the cube in different directions, so in principle even arbitrarily large data sets should be no problem. Contrary to expectations, loading in less of the data often makes the appearance better rather than worse.

Not every feature in the GUI is currently functional - the major one being multi-component rendering, but also the quick import and preview buttons. Both of these are actually fairly simple -  they just requires me to work out the most efficient way to do it (for multi-compment/volume rendering, the GUI buttons will greatly simplify what used to be a rather ugly, hacky process that worked but was unpleasant to use). 2D mode, though, is fully functional.

The Analysis menu does look a bit scary, but you can hide any panels you're not using and most of them should be simple enough. A big advance is that you can change the spectral axis units and it's no longer hardcoded to assume the data is HI, so the velocity of any molecular line just needs the rest frequency (a drop-down menu provides a few preset values and also access to Spatalogue). It's also possible to hide the axes with a single button, which used to be a much sillier process.

Region analysis tools remain much the same as in the original FRELLED but with improvements. Contours are now much faster and true isosurfaces are supported (and are fast enough that you could even show these for an entire cube, as an alternative to volume renders). You can also show velocity maps as well as much nicer-looking flux maps that use the requested colour scheme rather than only greyscale. There's a simple toggle for using a geometrical progression for contours (or logarithmic display for maps), with built-in safeguards to stop the user trying to display unfeasibly large and complex contours. And SDSS maps are now opaque, making them very much easier to see.

Isosurfaces are functional, though currently with only limited display capabilities. The mbspect section also has limited (but significantly improved) capabilities : it can only produce the input files and not run interactively, but it does allow all the options to be set directly in the GUI. The interactive version will be restored once I have access to a Linux system to test it on. Finally, users can also set some options for NED queries instead of just having it return absolutely everything.

I've tried as much as possible to test everything and test again. But as you can imagine, it's just not possible to test everything to destruction. So before making an official release, it's time for some beta testing. Volunteers are welcome ! Preferably those who aren't scared of working with FITS files. I can provide instructions and example cubes to try, but I'm especially keen to see what happens with data sets I've never tried, and with using it in anger : doing things in odd sequences and using features in unexpected ways. So if anyone out there wants to help, do get in touch. Leave a comment on social media or this blog, or contact me directly at feedback @ rhysy . net, and I'll add you to a beta-testing email list for next week.

Monday, 4 May 2020

Bringing wide open spaces indoors

My lockdown side-project is reaching the stage where I feel confident enough of finishing it that I can show what I've got so far.

Seven years ago (!) I made a model of Arecibo observatory, which became a pretty glass cube. Some time later this got heavily updated to the standards of having enough detail for a passable rendering (even if the materials were never that great) for a pre-rendered VR video. But wouldn't this be much cooler if it was an interactive game-like thing you could walk around in ? Answer : yes, yes it would.

After playing around with a few realtime data-based experiments of the amazing Blend4Web plugin, I decided it was time to take on the more laborious task of converting Arecibo to something people can explore for themselves. This is a lot more work than the previous experiments, which are largely a matter of tweaking Python scripts - essentially they're automatic except for some manual window-dressing. With Arecibo, everything has to be done manually.

A quick test proved that the plugin was easily capable of handling the mesh at an extremely high frame rate. But the materials were designed for the world of raytracing, and for various reasons they mostly look pretty horrendous in the realtime view. The easiest solution would be to remove them all and replace them with very simple plain colours, but this would look meh. So instead I learned about texture baking, which essentially does the rendering and stores the image on the mesh so it looks just like the rendered view but in real time.

Unfortunately this can't easily be scripted. Meshes have to have clean geometry and unwrapped in a reasonably decent way. Large parts of the mesh are, for reasons best known to my younger self, not clean at all. Large parts have had to be remodelled completely just so texture baking will work. Even then it isn't perfect, but I've decided that version 1.0 of Half Life Arecibo (working title) will be a quick(ish) learning experiment. For version 2 I'll learn the nice new materials available in Blender 2.8 and make everything look way better.

Anyway, screenshot time ! Here's the whole site. The landscape will be extended using a Blender plugin that automatically extracts textured terrain meshes from Google Earth (I actually already did that but I seem to have misplaced the landscape somewhere, as one does). The grey sky is what you see in Blender internally, although Blend4Web can replace this with something sensible.


The ground screen (that protects the beam from receiving the hot ground) was a particular challenge to re-texture. The baking isn't perfect, but I think it's probably good enough. The supporting struts do not yet have correct materials so they don't display nicely yet. 


The towers are all fully textured, though one of them has wrong material settings on the ladders and suchlike. Again, the cable materials are not yet correct.



The Gregorian dome was surprisingly easy to convert. For now, all textures are 4k resolution, even the really small objects, but this will probably change to save memory.


Having found a texture baking solution that generally works well, I was rather annoyed to find this method didn't work at all for the triangle. Then I found out the mesh has - lord knows why - far too many faces, so I completely remodelled it and it worked (there are a lot of minor defects, though I think they don't notice much). The same needs to be done for the azimuth arm. 


You can see a few rendering artifacts in a few places. These can probably be fixed, but for this first test I won't bother. In general they don't cause any serious problems.


Getting this building to work took up a good few hours and it still has some nasties. And the signs (not yet shown) are being particularly uncooperative. I don't know why this little shack is being such a pain, but I'll get there in the end.


That's for now. Most parts of the mesh are relatively straightforward. It's the few bits that don't work so easily that slow everything down. Still, if I can do a little bit every day, it shouldn't be too long before everyone can recreate Sean Bean's infamous death-by-telescope scene... well, maybe.

Monday, 23 March 2020

A field guide to mapping the Milky Way

How do you go about mapping the galaxy you happen to live inside of ? There's a hell of a lot of information on GalaxyMap.org, but it's not quite what I'm after. So in this post I'll do a step-by-step guide as to how to construct a map using 21 cm neutral atomic hydrogen data. If you actually want to try this for yourself, I'm going to be assuming some familiarity with Python (especially numpy and astropy/pyfits). Otherwise this post should describe the theoretical aspect well enough to be of interest by itself.


Working in Galactic coordinates

HI data has two main advantages : first, it's not at all subject to extinction by intervening stars and dust, and second, it gives us an easy way to measure velocities. Using trigonometry and a few reasonable assumptions, we can convert velocity into distance, with some limitations.

But before that we need to define a coordinate system. The convention for all-sky HI data is to use Galactic coordinates. In this system, the centre of the Galaxy is defined to be at longitude and latitude of 0. Galactic latitude l and longitude b are defined as follows :

Optical image of the Milky Way overlaid with all-sky data from LAB, with the Magellanic Stream highlighted in orange.
There are two main all-sky Galactic HI surveys : the Leiden-Argentine-Bonn survey and the HI4PI survey*. Both have been gridded in a nice friendly way, such that the pixel size is fixed in latitude, longitude, and velocity. Thus once you know the pixel size in each dimension and the world coordinates of any given pixel, you can very easily calculate the exact coordinates of every other pixel.

* It's 4Ï€ steradians, but I always read it to mean "HI For Principal Investigator".

For reference, both surveys have the origin at the bottom left  (l = +180, b = -90). For the LAB survey, the maximum x-pixel (l) range is 720 and the y-pixel (b) range is 360. The pixel size is 0.5 degrees for both axes. For HI4PI, the x and y ranges are 4320 and 2160 respectively, while the pixel size is 5 arcminutes.

And velocity ? For LAB data this spans the velocity range -458.6 km/s (z = 0) to +458.6 km/s (z = 890), with a channel size of 1.03 km/s. For HI4PI the velocity range is -600.0 km/s  (z = 0) to +600 km/s  (z = 945*), with a channel size of 1.288 km/s.

* This value may be slightly off. At the time of writing, I can't access the files I need to check.

Note that these values refer explicitly to the gridded pixel values. These are slightly different from the true resolution values, which are more often quoted in the papers. For the researcher, the real resolution is what matters, but for the data visualiser, it's all about the pixels.


Converting to distance

Once we've found the world coordinates of a pixel, and its flux value, we can then convert this to true 3D position, as follows. First we'll need some assumptions. The Sun is reckoned to be rotating around the centre of the Galaxy with speed V= 220.0 km/s, at a distance R= 8.5 km/s. The rotation curve of the Milky Way we can approximate to be totally flat, so that the velocity at any point Vpnt is also always 220.0 km/s. Given the velocity (vel) of any point, we can then calculate is distance R from the Galactic centre :
Note that this further assumes that this is independent of galactic latitude.This is reasonable because the disc is quite thin, but causes problems for structures which are outside the disc completely.

Next we can calculate the distance of the point from the observer :

Where d1,2 refers to the fact that the equation has two solutions. However, it turns out that this is only really true within side the solar circle, so we don't need to do the calculations twice. Rather we should accept that this region of distance ambiguity is inaccessible to us, so we should only do this calculation if R > R (we'll see what happens if we disregard this sage advice later on).

If that's so, we can proceed to calculate Cartesian coordinates of our pixel in PPP (position-position-position) space :
These will be in kpc since those are the units we've been working with. We can now iterate over every pixel in our data set and create a full PPP map from our original PPV (position-position-velocity) cube. This is relatively easy to do, and you can find the Python code to do so for LAB data here (note that some simple extra transform is applied to these final equations, just to ensure the data appears in a sensible position in our PPP cube). We just have to specify the size of the cube we want to make and hardly have to worry about anything else at all. Sounds great ! We'll be using the full information from the original data, so we should get a nice, clean, super detailed map at the end, without even having to specify the pixel resolution or anything even slightly complicated, right ?

Wrong. The problem is that PPV maps to PPP in a very strange way, which is not at all intuitive from the equations (unless you're some sort of trigonometric super freak, I guess). There's no guarantee that every pixel in our PPP cube even corresponds to one in our PPV cube. And not all our PPV pixels will contribute anything, since many of them will lie well outside the Galactic disc where our equations are invalid.

Here's what we get from the LAB data if we do this :

Slice through a PPP cube created from LAB data.
In some regions things are relatively good and we can see some nice astrophysical structures, but other parts are hugely undersampled while others have downright weird artifacts. We can do quite a bit better with HI4PI, which has higher resolution and so more fully samples PPP space, but it's still far from perfect.


Why mapping from PPV to PPP is a bad idea

Let's start with the artifacts. Our observations give us velocity along our line of sight, that is, how fast the gas is moving towards or away from us. In reality, the gas is also moving across the sky, but we can't measure that. We can get these "proper motions" for stars with considerable effort, but we just can't get it for gas.

The first problem is not so much that we'd like to know the proper motion (although that'd be nice), but that the equations assume our line of sight velocity measurements are accurate. But because we're inside the disc of the Galaxy, this is not always true. When we look towards the Galactic centre, or in the opposite direction, the only motion of the gas is across the sky - except for a little bit of random motion (~10 km/s). This means that in those regions of low measured velocity, our data is just too inaccurate for our equations to properly convert line of sight to true velocity. Better instrumentation won't help, it's a fundamental limitation of the structure of the Galaxy and our method.

Velocity vectors relative to the centre in green. Blue and red show the components
towards and away from us, respectively.
The second problem is that the equations have that annoying distance ambiguity within the solar radius, where the equation gives two solutions. Although we might be able to break this degeneracy using other data (e.g. by associating the gas with stars of known distances), by itself there's nothing much we can do to save the HI data. So this region, like the low velocity regions, has to be thrown away.

It might help to visualise how velocity maps to distance. One way of doing this is to plot isovelocity lines : lines of constant velocity.


Being inside the disc has weird consequences for what we can detect and where. Since everything's so darn close, sensitivity is extraordinarily high : we can detect essentially all Galactic gas. Looking through our original data cube, we see tonnes of stuff at very low velocities across the entire sky, because gas at high latitudes is only found when it's close to us and, therefore, moving slowly relative to us, due to the thin nature of the disc. But at the same velocities we can also be detecting material on the far side of the Galaxy !

The bottom line for visualisation is that if we start with the PPV map, we don't necessarily fully sample the PPP cube. This explains the other even more serious problem of the image - all those ugly black lines. How can we fix this ? One answer would be to interpolate extra velocity channels and/or spatial pixels in the PPV cube, so that we'd have more points that map to the PPP data. This does help, but it's inefficient and far from perfect. Even using the enormous HI4PI data set, which has vastly better spatial resolution (though similar velocity resolution) gives only a modest improvement in the sampling.


Alternatively, go directly from PPP to PPV

A much better approach is to work backwards. Beginning with a blank PPP cube, we can calculate the corresponding pixel in the PPV cube and use that to fill in the flux values. This essentially knocks all the problems on the head. By iterating every pixel in the PPP cube, we guarantee that we'll sample the whole thing. Although our calculated pixel positions in the PPV cube won't be integer values, all we have to do is simple rounding and we effectively interpolate the missing data (there are more sophisticated ways to do this, but they can wait for another time).

How exactly do we go about this ? We define the coordinate system of the PPP cube arbitrarily. Then, knowing our Galactic coordinate system, we can use some basic trig to calculate the longitude and latitude of any given pixel. We need to get R first, but this is easy because we know the position relative to the galactic centre gc :

I work in degrees, hence the +90 for convention. The pixel positions xyz must be in physical units (kpc). The tan2 function is a wonderful programmatic convention that simplifies things enormously. Using the usual arctan returns values ±90 degrees, since there's a degeneracy in the tan function. Atan2 gets around this by providing two values, returning values ±180 degrees, which is exactly in accordance with Galactic data gridding conventions (we could convery this easily enough to the range 0-360 if we wanted to, but there's absolutely no need).

All we need now is the line of sight velocity. We can get that by rearranging the earlier equation to calculate R :
Since the original PPV data is gridded in a nice friendly way, the hard part's over. Now that we know the longitude, latitude, and velocity of a pixel in the PPP cube, it's trivial to convert this to the pixel in the PPV cube - remember, the original data has pixel size of constant latitude/longitude/velocity.

Voila. We can now extract the corresponding flux value and create a fully sampled PPP cube.


What to do if your data set is feckin' enormous

But wait ! There's one extra complication. If we want to use the LAB data, we can go right ahead and use the final code. Of course, it's always going to be better to use the HI4PI data, but this is difficult to work with because of its gargantuan 35 GB size. The astropy "pyfits" module is not at all good at dealing with large data sets, so we'll need to convert it into a format it can handle. We can do this using the much older miriad software, which was written from an era when a 100 MB file was considered hefty. Consequently it's massively superior in terms of memory management and can process 35 GB files without breaking a sweat, even on a system with 16 GB RAM.

This bit is trivial. First, we convert the FITS file to miriad's own format using the FITS task. Next, we use the same task to extract individual channels, by setting the region parameter to give single-channel slices (unfortunately, this only works when converting from miriad->FITS, not the other way around, which is why we had to convert the file). Annoyingly, miriad insists on producing FITS files of not one but two channels. We can either accept this and have the script only work with the first velocity channel of each cube (it's super simple to slice the data for this), or first process the files and save ourselves from an extra 35 GB of data we don't actually need. Both steps are easy anyway.

Right. We're pretty much there. We've converted our 35 GB single data file into 944 smaller, more manageable files. All we need to do now is modify our PPP code to work with multiple files - and here it is. As a result of all this, we get the following :


Ta-da ! Lovely. Doesn't quite eliminate all of the artifacts, but we can get rid of those in the visualisation stage.

(For the enthusiast : we could in principle go through every pixel in the PPP cube and open the necessary FITS file every time, but this is hugely inefficient - it means opening the files hundreds of millions of times. I estimated it would take the code 4 months to complete, which made me sad. So instead the code precalculates which pixels it needs to extract for each file. It then orders the list, opens each file, extracts all the pixels it needs from that file, and moves on to the next one. This means it only needs a maximum number of 945 file-opening operations and runs in an hour or so. An extra complication is that a list of hundreds of millions of entries starts to cause memory issues, so the code can work in chunks. The user then has to specify the pixel range of the PPP cube they want to search.)


Making things look pretty

For visuals, we have two options. We can either make a volume render using FRELLED, or we can try and pick out the major features using isosurfaces. Volume renders look prettier and use all of the data, but isosurfaces can make it easier to reveal the important structures and are small enough to display as interactive on-line sections of a web page.

Volume rendering via FRELLED need not be explained in any detail here. Let's skip straight to the render :

Sun position in yellow and the galactic centre in green.

What about isosurfaces ? These are something I've struggled with for a while. Eventually I stumbled on a couple of different Python modules that specialise in this. The one I'm using is scikit-image, which can produce meshes in a format that Blender can recognise. It's also fast and deals with large, complicated meshes pretty darn well. It's not 100% foolproof -  sometimes meshes have their normal vectors pointing the wrong way, but generally this is easy to fix manually.

(I tried other solutions - extensively - like using metaballs and point cloud skinning scripts, and I'd strongly advise everyone else not to try this. They just don't work very well, at best giving ugly results and at worst being useless and inaccurate. Use a dedicated module and save your sanity !)

The code to generate isosurfaces is a bit of a hack at this stage - incorporating it into the next generation of FRELLED is definitely happening, but that's still a ways off. But for now, it works. You'll need this .blend file (containing an internal script and some pre-set materials) and this external Python script, plus this readme file. Eventually I'll make something less hacky, but not today.

Last but not least - exporting to the web. I wanted to use Sketchfab, which I've been very impressed with, but then I discovered it has a stupid 50 MB file size limit. So I spent a weekend investigating Blend4Web, which is free and totally awesome. It's one of those nice things that just works. So here's an interactive 3D model of the HI content of the Milky Way. It'll work on mobiles (it even has a VR mode option !) but it's better on PC - the labels tend to get cut off on a phone. Click on the buttons in the top left panel to toggle different components, and on the annotations for a bit more information.

It's a 28 MB file, may take a few minutes to load, and Blogger won't let me
embed it correctly, so click here for the interactive version.
It's far from perfect yet  - the bloom effect is a bit strong and the anti-aliasing is lousy. But this is my first attempt, so I'm pretty pleased with it. Expect updates on this and lots more interactive content.

So that's it : you now know absolutely everything about how to map the hydrogen content of the Milky Way disc. In a future post I'll look at how to do something similar for the surrounding clouds, which are a lot more fun in 3D because they're found across the entire sky.

Friday, 14 February 2020

Accidental optical illusions

I have a fun little side-project to make a 3D model of the Milky Way using all-sky HI data. By measuring how fast the gas is moving and doing some trigonometry, it's possible to convert velocity into distance. The equations are a bit awkward, and if you get things a bit off, the end result looks very strange. They also have a limitation that they give a meaningless double solution for any point closer to the centre of the Galaxy than the Sun.

This meant there was quite a bit of trial and error involved until I got the correct result (more on that in a future post). To check I where things were going wrong, I had the code output the calculated galactic coordinates (latitude and longitude across the sky, measured from the galactic centre, as well as velocity along the line of sight), with the data set to zero inside the solar circle where the solutions would be garbage. Actually I'm pretty sure I got the position of the solar circle wrong, so this is just a complete mistake.

But it did produce a couple of fun little optical illusions :

Raw image here.

Both are quite similar. The colour in the grey circle looks like its varies, but it doesn't : it's completely uniform.

The effect is strongest with the left figure, which is velocity. The right side of the circle appears significantly brighter and the left significantly darker, a bit like looking at a crater in partial shadow (or, if you take the inverse perspective, a dome). Cover everything except the circle with your hands and you'll see this is entirely the result of your brain inventing stuff.

The figure on the right (galactic longitude) can be subtle at first, but once you see it, it's very hard indeed to make it stop. This time the right side of the circle appears darker and the left brighter, especially when you focus on the edges. I find that I can more-or-less control how strong this appears by concentrating on different parts of the circle, but sometimes it becomes so strong that I can barely make it stop even by covering the edges.

If we apply an animated mask to the regions outside the circles then things get even more fun (apologies for the small radial artifacts caused by gif compression) :

Raw image here.
It really is quite hard to graphically prove that the circles are always of constant colour. The only way to show it for sure is to download the images and an examine them in extreme close-up for yourself.

Monday, 18 November 2019

Wibbly-wobbly spacey-wacey

A few weeks ago I posted some first tests of displaying isosurfaces in Blender. This got interrupted by the receipt of a referee report, so now back to the wibbly-wobbly spacey-wacey stuff.

I used to poo-pooh isosurfaces as being inherently inferior to volumetric renders because they lose information. They also feel somewhat like cheating, because rendering a surface is easier and uses a lot less memory. While all this is true, I'm somewhat reconciled to their uses : isosurfaces depend far less strongly on viewing angle than volumetric renders. That makes them much easier to highlight features in an objective way, and contrary to my intuition, this can make them better for finding faint structures rather than worse. The eye tends to get very confused when you have really bright and faint sources together, but slap on a fixed level surface and BAM you can see if what you're looking at is significant or not. You might still miss the very faintest stuff, but there can be a surprising amount to see at relatively bright levels.

The other nice thing about isosurfaces is that because they're very fast to generate (typically seconds per surface in these examples, if not less), they're easy to animate. Last time I showed some fixed-level turntable animations. Here's M33 from the AGES HI survey shown at a variable flux level, starting with the brightest gas and ending with the faintest. Each frame decreases the flux to 97% of the previous value.


There's a whole bunch of clouds around M33, most of which show up quite well in the animation. It's not as good as manually tweaking the level of each region, but it basically works. I probably should have frozen the final level and rendered a full rotation, but never mind.

Here's the same sort of animation but using NGC 4361 from the WSRT HALOGAS survey, which also has a bunch of weird stuff going on around it :


And finally, here's a test of a true time series from a simulation. This is an old one from when I was trying to learn the FLASH hydrocode and my galaxy got all unstable because of a bug in the boundary conditions. Four fixed levels all animated, showing how the nice stable galaxy gradually goes fully wibbly-wobbly and eventually gives up. Does it help analyse what went wrong ? Not really, but it looks nice.


Thursday, 24 October 2019

Isosurfaces for fun and profit

The best way to look at at 3D volumetric data is, in my opinion, in its original 3D volumetric form. To this end I've spent several years years developing the 11,000 lines of Python code that is FRELLED (albeit all of which is just a script to load astronomical FITS files in Blender, which is what does the hard work). I've always been a bit skeptical of other ways of visualising the data... 2D slices are fine, and often necessary, but things like isosurfaces seem to me to be throwing away a lot of really pretty* information.

* I don't much care if it's meaningful or not.

To be fair, FRELLED does already include the ability to display renzograms, which are essentially contours of each slice of the data. While it's true that viewing the data at a fixed level, reducing it from a full volume to a thin surface, does remove a lot of information, I've been realising that by cutting away a lot of the noise it can actually become a lot easier to see interesting features. With the full volume, sometimes the noise just gets in the way. Using the renzogram facility of FRELLED, we've found a bunch of hydrogen streams in the Virgo cluster we'd just never noticed before (paper submitted).

So renzograms are super useful. But while 3D renzograms are sort-of isosurfaces, they're not proper 3D fits to the data. That's harder to do in Blender - I tried to get this a while back, and it sort of worked but it was very, very hacky. That method used someone's old Python script that skins a point cloud of vertices. It works well in some situations but not in others - a lot of manual cleaning is needed on complex data sets, and that's not much fun. The experience is a bit like using a half-broken toaster : you're never quite sure if you're going to have a nice breakfast or burn your house down.

But now I've found that the Python scikit-image module includes a "marching cubes" algorithm that generates proper isosurfaces. Fast, effective, and no mucking about with cleaning up artifacts at all. I've quickly hacked this into Blender, using another module to convert the vertex data generated into a Blender-readable format. The basic code is just a few lines long and it works without complaint.

So, time for some examples ! This first one is a bog-standard Virgo cluster galaxy with no interesting features whatsoever - it's just a long, cigar-like blob, with the long axis being velocity. Colours indicate brightness of the emisssion (purple, blue, green, yellow and red going from bright to faint).


For a second example, here's another Virgo galaxy which does seem to have a distinct protuberence on one side. It's probably losing gas as it moves through the cluster.


And then there are oddballs like this one, which seem to have a distinctly noisy appearance and a "tail" in velocity space :


Just to prove how incredibly easy this is, here's the whole data set of 102 galaxies. Even in this zoomed-out view, you can see that most galaxies are quite smooth and symmetrical, but some have pretty clear extensions and other weirdness (after a laborious statistical analysis we're highly confident these are real and not just due to variations in the noise).




You may be thinking that that's all very nice, but what about some nice resolved high resolution data ? No problem, here's one of my favourites - the M33 galaxy and its many associated clouds :



The M33 system is so complicated that I cheated a bit with the contours on that one. In all other cases, the same colour is used for identical brightness levels, but in the case of M33 I set the levels manually for each cloud - otherwise you start being totally dominated by noise in some cases, while not seeing anything at all in others.

Finally, here's an isosurface of Medusa a simulated galaxy undergoing strong ram pressure stripping. No noise to worry about at all for simulations.


All this is part of a larger effort to recode FRELLED in modern Blender. FRELLED currently relies on Blender 2.49, which is more than 10 years old. Blender 2.8 has a lot more features and comes with its own Python and Pip install, making it waaay easier to install the necessary modules. Perhaps that will help catapult FRELLED from obscurity to total domination of the astronomical community having more than a dozen users. That'd be nice.

Friday, 11 October 2019

Pretty things are pretty

Time for some more pretty pictures of hydrogen....

It's proving surprisingly difficult to convince people that some streams I've found in the Virgo cluster radio data are real. So to settle the matter once and for all, I've resorted to creating synthetic galaxies and adding fake streams and noise extracted from real data. Then I blindly search the data, labelling what I think looks like a stream and what doesn't. Since I don't know in advance which galaxies have streams or not, this should be a good way to quantify very rigorously how many false positives can occur just due to the noise, as well as measuring how many of the known fake streams would actually be detected by the search technique.

Each of these 100 pillars is a synthetic galaxy, with the vertical axis being velocity. Each of the "segments" is a contour at a different velocity channel, extended into 3D. Real galaxies would look a bit more complicated than this, but these are good enough to search for features as extended as the ones in the real data. Arranging them into a neat grid makes it easy to search the whole data set very quickly and isn't just for the sake of making something minimalist.


You can't really see the fake extensions from this angle - they're more visible from underneath. This particular data set has a stream in every galaxy, always pointing in the same direction. For the real search I vary the length, angle, brightness, presence, and number of velocity channels of the streams.


These pretty pictures aren't going in the paper - for that, I'm showing a boring but more easily comprehensible 2D plot. They look nice though.


Sunday, 30 June 2019

Unleash the render from hell

Ever since I got a VR headset I've been meaning to make more content for it but never quite manage to get started. Now I'm trying to make amends, but the process has been considerably... less smooth than I would have liked. So come, CGI enthusiasts of the internet, and let me regale you with tales of daring renders, dashing rogues, beautiful princesses and murder on the high seas ! Or, well, some of those, anyway...

(This is also going to be, if not a fully-fledged tutorial, then at least a highly practical guide with many hints and tips for those interested in doing VR content in Blender, complete with Python scripts.)

One of the first things I wanted to do was to recreate the ALFALFA Sky videos I did in the glory of 360 3D VR. Here's one of the originals for reference :


Clearly very well-suited indeed to VR. A bunch of galaxies, scientific authenticity, and an incredibly simple setup. Great ! Let's just update the file and turn on VR, right ? Wrong. In practise, this became something of an albatross, but unlike the Ancient Mariner I would happily see it dead.


Galaxies in 3D that are really in 2D

Making the 2D version was trivial. I even used the venerable Blender 2.49 for the first one, since I was much more comfortable with its Python operations and materials settings at the time. It just needs a few things :
1) The ALFALFA catalogue itself, found on the survey website.
2) A query to the SDSS to get the optical size of each galaxy. I had to split the galaxy table into seven because the SDSS doesn't large "large" queries, but honestly in this era 30,000 positions and radii shouldn't be large by anyone's standards, let alone the SDSS.
3) Using the catalogue file from step two I obtained images of the optical counterparts of each detection via this script. The optical size isn't terribly accurate but it's good enough.
3) A Blender script to produce a plane mesh for each galaxy, with a shadeless material and the correct image texture, scaled according to the optical size of the galaxy. The old 2.49 script can be found here and the version for Cycles in 2.79 (may need minor changes to run in 2.8) is here.
4) A constraint applied to each galaxy object to face the camera.

Okay, there are quite a few steps to follow, but none of them are very difficult or take very long. And that's really it. To do this in VR would, you'd think, simply be a matter of using the correct camera settings. Strictly speaking this is true, but there are a number of whacking great ugly-as-being-hit-in-the-face-with-a-sock-full-of-Ebola complications.


Galaxies in proper 3D

The first is that Blender only supports spherical stereo rendering using Cycles. And believe you me, using Cycles is necessary. Yes, you can render conventional side-by-side 3D in Blender Internal, but spherical stereo is different. The thing is that when you look around, the position of your eyeballs changes as you rotate. The Cycles Spherical Stereo mode accounts for this, but doing it in BI requires ancient eldritch knowledge the like of which has long since gone out of the world. Maybe if you sacrifice enough chickens to B'aalzathak, Demon God of Renders, you can make it work, but don't. No really, don't. Just bite the damn bullet and accept that spherical stereo rendering requires Cycles.

Fortunately it was quite easy to produce the script to make Cycles-compatible galaxy meshes and materials. You can find them in the above list.

You would think that's the hard part. You would think that learning how to do the Cycles materials and the new internal Python syntax in Blender was difficult enough, and that surely by now we should have reached the green fields and sunny uplands of VR utopia. You'd be wrong, poor fool.

Now normally to render VR content all you have to do is is enable "spherical stereo" for the camera and its associated options, and set which "view" you want to render. For some reason Blender still* can't render left and right views and automatically composite them, so you have to render each one separately and composite them together later. This is very easy - here's an example of how to do this using Blender's sequencer.

* Caveat - we'll get back to this later.

Where it gets complicated is not with VR rendering, but with a particularly horrid combination of local circumstances. Although my gaming laptop is powerful enough, I tend to use it quite a lot. So I don't like to use it for much in the way of system-hogging rendering unless it's not going to take too long (I'm also a little worried about leaving it on continuously). Fortunately, my desktop in work is even more powerful, so I use my laptop at home to create files and my work machine to render them. It can chug away for days on end without breaking a a sweat, if need be.

But I can't do that here. See, there was a stupid limit in Blender <= 2.78 where Cycles couldn't handle more than 1,024 image textures. That's been removed in 2.79, but the glibc version my work machine uses is too old to support 2.79. Can I fix it ? I guess so... but as far as I can tell this is a big job and could end up with me needing to reinstall Linux from scratch. And I really don't want to do that while I'm at the revising stage of a paper. Afterwards, fine, but not during the process. That would be silly.

(Quite honestly I haven't got a soddin' clue what a "glibc" is. When anyone tells me about libraries and suchlike, my eyes glaze over and I come out in a nasty rash).

So if a straightforward render on my laptop is out, and my work PC is unusable, what are the options ?


1) Get it rendering faster on the laptop
This is an incredibly simple scene - just a bunch of transparent textured planes. So absolutely all settings can be set to the bare minimum, with even Samples set right down to 1. And once rendering begins, it's blazing fast. The problem is not the objects themselves but simply their vast number, which Blender does not like. Not one bit.

One thing I notice about Cycles is that it spends friggin' ages on the "synchronising objects" stage before even starting to render anything. Typically the total render time is about 5.5 minutes or so (per half-frame, that is, for a Left or Right image, so the actual time is double that), of which almost all consists of synchronising objects, updating shaders and loading images. Only a few seconds are needed for the render itself. Frustratingly, in Preview mode there's no need to do the preliminary stages, but you can't screen capture that like you can with OpenGL view for BI. So you have to render.

Somehow I found that there's a way to significantly decrease the synchronising speed : create a blank file with linked (not appended !) copies of all objects. I have no idea how or why, but it works - it gets the total rendering speed down to maybe 3.5 minutes per half-frame. Rendering the linked file in 2.8 gets things down to 2 min 50 seconds per half-frame, with the synchronisation stage almost eliminated but the "updating shaders" phase taking somewhat longer.

The penalty for this is that it makes the file even slower to open. It was already slow, maybe 10 minutes or so, but now that approximately doubles. Not nice, but worth it for the render times. And in 2.8 the files open very much more quickly, so that's not much of a factor.

Other tricks are less successful. GPU rendering isn't an option due to memory limitations. Command line rendering is more annoying, because Blender insists on printing out the status of the object synchronisation. And because there are 30,000 objects, the print statements cause a significant slowdown, making it actually slower than rendering with the GUI. Googling how to suppress the command line output for Blender didn't help because all that comes back was for Linux. Eventually I realised I should be less specific and look for how to suppress terminal output in general under Windows. That turns out to be easy : just add >$null to the end of a command sequence in PowerShell (I don't think this works in the regular command prompt).

Does this actually help with render times ? Annoyingly, no it does not. Command line rendering is still slower, inexplicably.

Another trick I found online somewhere was to use PNG images instead of jpegs. This does reduce the synchronisation stage (for some reason), but significantly increases the updating shader time, and the two approximately balance out so it wasn't worth doing. But I found it quite frustratingly difficult to figure out how to access Cycles nodes (e.g. image texture nodes) via Python, so in case anyone wants it, here's a script that will go through a file and replace Cycles image textures with PNG versions (if you need to convert the files in the first place, this script can be easily modified - but watch out as it'll overwrite existing files).


2) Make it render in 2.78 and run it in work
Option 1 having only limited success, is there any hope at all for rendering the file using a version of Blender my beefy work desktop can handle ?

The answer to that one is "yeah but don't". The only way to overcome the texture limit is to render in passes, about 30 per frame. So I wrote a script that hides (or better yet changes the layer) of the most distant 1,000, 2,000, 3000 etc. objects, progressively rendering each set of galaxies such that compositing is a simple matter of alpha-overing every different pass. This does work. The problem is there's a huge speed penalty, such that I estimate the total rendering time to be about six weeks. And if anything went wrong, that'd be bad.

While removing objects from view (this works best on layers rather than hiding them) does cause a big decrease in the synchronisation time, there are two reasons this doesn't help. The first is the increased number of passes. The second is that each rendered image must be saved with transparency, and to get a decent-looking result (in 2.79 - this is not so much an issue in 2.8) necessitates using at least 15-20 samples, slowing down the rendering stage substantially. So this does work, but it's really not worth doing.

A somewhat related, less serious issue is that rendering an image with transparency and then overlaying on a black background gives different, substantially paler results than rendering directly with a black background. I don't know why this is, but I guess it relates to the emission node and how it computes colour values based on what's behind it.


3) What about a render farm ?
Well now I mean this should be a project ideally suited to such a thing. It's not a complicated file, it just has a bit of a long loading time but thereafter renders quite quickly. All it needs it a bunch of a computers thrown at it.

A very nice gentleman offered me the use of his SheepIt renderfarm account. SheepIt is an absolutely amazing service and once my work desktop is properly updated I plan to leave the client running in the background in perpetuity, unless I need to run something intensive. Unfortunately, my file is more than double the 500 MB upload limit. And they're not fooled by zip files easier, frustratingly (as the .blend file compresses down to just 200 MB). Sheepit doesn't allow Python scripts either, but nor would I want for anyone else to suffer the several hours it takes to recreate the original file from the catalogue (Blender doesn't handle large numbers of objects well, no matter how good your system specs are).

I briefly considered a truly insane plan : for every frame, create a new Blender file with sets of objects at different distances. That would get the file size down and possibly increase the rendering time nonlinearly due to the reduced file size needing less memory. While deleting objects didn't take too long on my work machine, on my laptop this could run for hours and still not complete (at least in 2.79). So I wrote a script to delete objects one by one and print its progress to the terminal. This didn't work in 2.79 (or rather it took about ten minutes to delete one object) but did work in 2.8 (deleting thousands of objects in about half an hour). But I realised this would be far too slow and far too labour intensive, so I abandoned this hair-brained scheme on the grounds of Nope.

I briefly tried a commerical renderfarm, Blendergrid. This looks very nice, doesn't have a file upload limit, and the staff were very helpful and pro-active in dealing with issues. I deem them to be first rate when it comes to customer service. It has an especially nice feature that you can upload a file to get a price quote, so I tried that. Alas ! After about an hour it was still trying to run the file, but giving generic error messages. But no biggie : it can email you the result when it's done... but not, it seems, if the file just keeps crashing.

Ultimately we might have been able to solve that one and do it commercially. How much it would have cost, though, we'll never know.


4) Work the problem
Alack poor Rhys ! At this stage we find him distraught, having tried every trick under the Sun and still stuck with an absurdly long render time. The only solutions appeared to be six weeks running on a work computer or 150 hours on a much-needed laptop. Most frustrating of all was that the render time was absolutely dominated by unnecessary processes, not complex calculations : once the preliminaries have been done, the render completes in about five seconds or less. Never have I felt so near and yet so far.

Then I had two crucial breakthroughs. First, because the camera is moving continuously through a three-dimensional cloud of galaxies, there are no obvious ways to cut out a single set of objects and render multiple groups separately for compositing later. That's why I had to use a script for my distance-based rendering approach earlier. But then I realised that actually, thanks to the nature of the survey, there is a natural break point where the render can be split :

This early version looks very dull nowadays, sorry about that.

When the camera is on the left side, galaxies in the left-hand wedge will always be in front of those on the right (and vice-versa), even given the spherical stereo camera. That means the render can be split after all. Now the split isn't even, unfortunately, but it's still substantial, with about 20,000 galaxies on one side and 10,000 on the other. That means that the smaller file is easily small enough to render on SheepIt, and the reduced size of the file for the bigger half makes it substantially faster to render on the laptop : about 1.5 minutes per half-frame. So my object-deletion script came in useful after all.

And second, this business of rendering left and right separately is (somewhat) nonsense. It's true that for some godforsaken reason Blender can't display top/down left/right renders in the rendered image window... but if you save the image to a file then it's correct. That's totally counter-intuitive but it works. Similarly when animated it works as well. The settings needed for this are simple :
Stereo mode : Stereo 3D
Left and right both enabled
Views format : stereo 3D
Stereo mode : top/bottom
Squeezed frame : enabled
The "squeezed frame" is important. If you don't do this it will render both images at the size specified and stick them together, so the final image will be twice the vertical size requested. That gives odd results, so the squeezed option should be enabled to keep everything in correct proportion. The image size should be square, i.e. accounting for the fact that you'll have two rectangular 2:1 images stuck together.

This discovery might not have been so fortuitous for most other scenes. But this one is massively dominated by the preliminary stages (synching objects and updating shaders), which, it turns out, do not have to be recalculated for each viewpoint. This means that the render time is very nearly halved, with the full image still only taking 1.5 minutes to render. And that gets me exactly where I want to be : the render isn't as fast as it could be in principle, but it's certainly fast enough for my purposes.

Later, I found that some simple scenes (like the opening text) render very much slower at 4k than at 2k for no obvious reason. So it turned out to be much faster to render the left and right images separately and composite them in the sequencer, so that sequencer setup file is still quite useful. Why this doesn't apply to the main rendering sequence I have absolutely no clue.

So that was the solution : splitting the file, rending part on a render farm and part locally, rendering both left and right together, and using 2.8. Whew.


Thirty Thousand Galaxies At Five Hundred Trillion Times The Speed Of Light

Without further ado then, here's the final result that took an insane amount of time :


The final rendering was bestraught with all the usual difficulties. I rendered it once, taking about 37 hours, and it looked very nice except that the camera froze for the last ten seconds so I had to render the whole thing again. Annoyingly that increased the rendering time non-linearly to about 48 hours. But that version worked. And the short sequence where the galaxies first appear made use of a modified version of the the multi-pass rendering script, so I got some use out of that after all.

The final sting in the tale was that my VR headeset (an obscure but affordable device by the name of Magicsee M1) decided, as I was about the view the movie for the very first time, that now would be a good time for a firmware update. It's had only one other update since I got it, and since it's obscenely obscure I presumed it was now dead. Apparently not. Foolishly I decided to accept the update since the render wasn't quite done yet... and that broke everything. All my apps, gone. The button to enable VR mode no longer working.

I was not happy. Not at all.

Fortunately it wasn't nearly as bad as it appeared. The built-in movie player simply had a different option to enable VR mode, and my apps were still there but took longer to show up. The VR mode button still enables VR in apps that need it. Christmas was saved, joy was unconfined, there was dancing in the streets, I decided the movie looked quite nice, and that was the end of the whole sorry business.

Why Bother ?

It's rare that I manage to read any longer pieces on arXiv that aren't strictly about galaxy evolution, but today I indulge myself. ...