Starry, starry night

The holidays are nearing, and it always makes me ponder the universe. Maybe it's because the stars come out in the cold and clear dark winter evenings, or maybe it's because of the unusually bright star in the Christmas story. Anyway, it gave me the idea of creating something like a classic 3D star field. What if you used actual star positions for this star field, and would you see something like you do on Star Trek if you flew through it at warp speed?

I started out by searching for star catalogs. These are plain text databases with lots and lots of scientific data about the stars.
Star databases are a bit of an annoyance; each has its own format, and they all contain different fields. Some contain parallax data (which is needed for calculating the distance of stars), but most don't. Some contain the visible color of stars, others don't. But most annoying of all, they don't use the same unique identifiers; star Procyon A is named HIP-37279 in the Hipparcos catalog, HD 61421 in the Henry Draper catalog, Gliese 280 in the Gliese catalog, while it's Bayer Flamsteed designation is 10 Alpha Canis Minoris.
The Tycho-2 star catalog contains over 2.5 million stars, but unfortunately it is worthless for making a 3D model because it lacks parallax data. So I took the HYG database, which is a modified version of the Hipparcos database, and combined it with usable data from the SKY2000 star catalog, plus a short list of well-known brightest stars. This involved some awk and Python scripting to cross reference the databases and calculate the 3D positions.

Star positions are usually given in right ascension (RA) and declination (Dec). These are angles, like latitude and longitude, but are often described in hours, minutes, and seconds. Given the angles and the parallax, you can calculate an x,y,z position in space, and this is what we need for a 3D star field. These coordinates have the earth at its origin. You might want to use galactic coordinates instead, which puts the Sun at the origin.
Now, because all of this stuff is turning, turning in outer space, the positions of stars change (unlike what I was told as a child). Star databases contain the position at a point in time (referred to as the epoch) and you can adjust the position using the proper motion of stars. Frankly, I didn't bother, but if you want to know the correct position of a star as it is today, you are going to have to do some extra math.

The constellation Orion

At a distance of 4 kpc

So, I put all coordinates into OpenGL and the result is pretty neat. In the left image you can clearly see Orion. The red star in the top-left of Orion is Betelgeuse. All the way to the left from Betelgeuse is Procyon. The image on the right shows what you see from a distance of 4 kiloparsecs. It's like a bubble with a disc through the center. This is just what the data set looks like though; in reality, it probably does not have this shape as we are part of a much, much larger Milky Way.

This is what almost 220 thousand stars look like on a computer screen, and you can look around and point at stars in real-time.

Fun facts:

• The furthest bright star, just outside the bubble, is Deneb at 4754.8 light years.
• Most stars around us are red dwarfs. Hence the pinkish color of the bubble.
• When you travel to another star, the night sky looks different.
• OpenGL vertex buffer objects (VBO) kick ass.

I was amazed that my (now three years old) computer had no problems whatsoever in dealing with this data. Once again, I find that the modern PC is an amazingly powerful computer. Although I'm no scientist and no astronomer, it enabled me to put together a pretty realistic 3D model of our stellar neighborhood in my attic.

A Pool Allocator or The (un)importance of Code Optimization

In school I learned to use dynamic memory allocation, and to use linked lists to keep track of dynamically allocated objects. The idea was that arrays are static and fixed-size, while linked lists can virtually limitless continue to grow dynamically. The linked list also doesn't use memory that it doesn't need, and there are some nice performance characteristics when inserting and removing items from a linked list.

The malloc() library function takes memory from the heap, which is a region of memory. Whenever malloc() runs out of heap space, it will call sbrk() to grow the heap by asking the operating system kernel for extra memory pages. Modern operating systems actually use mmap() instead, which works by virtue of virtual memory and demand paging. It's all pretty sophisticated.

How different were things for classic arcade game machines, tiny computers with very small system memory. There was practically no operating system at all, it was just the game code itself running on the bare metal. Arcade games like Space Invaders and Galaga did not have the luxury of a malloc() library function. They were not even written in C, but in the machine's assembly language. The way to allocate monsters and bullets was to use an array in the bss segment—an area of system RAM memory sitting just outside where the program code was loaded, or rather, copied from eeprom to RAM when the machine was switched on.

Things were so much simpler then, but even today I like allocating monsters and bullets from a ‘free’ list backed by a statically allocated array. Such a pool allocator gives a very clear picture of what is going on and how things work.

Even though the pool allocator has a hard maximum on how many objects may be allocated, it will never deny you memory if you stay within that limit. No other task in a multitasking operating system can steal this memory away from you, making it a very robust memory allocator.
In terms of performance, the pool allocator always completes in a fixed time, which is important for games that need constant performance.

Discussions about memory usage and performance are becoming pointless in this day and age, when we have massive amounts of memory and CPU power at our disposal even on mobile devices. As computers before were frustratingly slow machines, you had to push the limits to get the most out of it. Writing efficient code was the challenge, and every time you managed to squeeze out a few cycles more it was an achievement well earned.

Nowadays this hardly applies anymore. Although I have to say, lately I had some scripting code for doing large data processing. I rewrote it in C++ with multi-threading, and boy, did that make a difference in run time.

File wrapper class in C++ (2)

Last time I showed how to implement a File class. In the solution Brian (pictured on the right) and me used a second class as a wrapper around a standard FILE pointer. The two classes worked together by means of a shared_ptr. The shared_ptr ensured that the file would not be closed until there were no more references to the file wrapper object. Reread the last post to see what I'm talking about.

A reader of this blog pointed out that there is another solution, one that turns out far more idiomatic. The constructor of shared_ptr accepts an instance of a functor, named a Deleter. This functor will be called when the reference count drops to zero, and the object will be destructed. The resulting code looks like this:

#include <cstdio>
#include <tr1/memory>

class File {
public:

    // conversion constructor
    File(FILE *f) : stream_(std::tr1::shared_ptr<FILE>(f,
        FileDeleter() )) { }

    // copy constructor copies the shared_ptr
    File(const File& f) : stream_(f.stream_) { }

private:
    std::tr1::shared_ptr<FILE> stream_;

    // local functor closes the FILE upon destruction   
    class FileDeleter {
    public:
        void operator()(FILE *f) const {
            if (f != NULL)
                std::fclose(f);
        }
    };
};

Right here I inlined the FileDeleter class. You may also choose to put it outside the enclosing class, or it can even be implemented as a struct with operator()().

All in all, this solution is more elegant and more idiomatic than the previous one.

File wrapper class in C++

Last week my good friend Brian (see picture on the right) and I encountered an interesting problem: wrapping up a file descriptor in a File class. Making class wrappers in C++ generally isn't all that hard, but this particular one was a lot more tricky than we had expected.

What we wanted to have was a File class that would act the following way:

• based on standard C's FILE*
• factory function open() creates a File, and opens it
• destructor of File closes the open file

So that you could write code like this:

    File f = open("filename");
    f.write("hello\n");

    // destructor closes the file

You might code open() like this:

File open(const char *filename) {
    File f;
    f.open(filename);
    return f;
}

Unfortunately, this code doesn't work, as it has a bug. Do you see it already?
The thing is, when you return f, the value f gets copied [by means of the copy constructor] and then the destructor gets called in the original f, because it goes out of scope. Since our destructor closes the open file, the result is that the file is closed once open() returns. Read that last sentence again, argh.

There are two issues at play here:
1. The copy constructor doesn't work because you can't properly copy an open FILE* object.
2. Exiting the function causes the destructor to run, closing the file — despite the fact that we are returning the instance.

Can I have the attention of the class
How to solve this problem, and no ugly hacks, please. A good solution lies in using shared_ptr. You can copy a shared_ptr and safely return it, and its internal reference count will cause the File instance to stick around. The File object will not destructed until it really needs to be.
[More specifically, the shared_ptr keeps a pointer to dynamically allocated memory, which never auto-destructs like instances on the stack do. The copy constructor and operator=() of shared_ptr happily copy the pointer, but while doing so it increments a reference counter that keeps track of how many shared_ptr instances are currently referencing the dynamically allocated object. The destructor of the shared_ptr decrements the reference count, but does not delete the pointer until the refcount drops to zero. So, a shared_ptr can be used to 'keep objects around' for as long as they're needed].

Beware, you can't simply dump a FILE* into a shared_ptr, since the FILE* may be allocated with malloc() rather than with new. And besides, a FILE* opened (or allocated) with fopen() must be deallocated using fclose(). Therefore we wrap the FILE* into another new class named FileStream and use that class with the shared_ptr.

#include <cstdio>
#include <tr1/memory>

class FileStream {
public:
    FileStream(FILE *f) : f_(f) { }

    ~FileStream() {
        if (f_ != NULL)
            std::fclose(f_);
    }

    ...

private:
    FILE *f_;
};

class File {
public:

    // copy constructor copies the shared_ptr
    File(const File& f) : stream_(f.stream_) { }

    ...

private:
    std::tr1::shared_ptr<Filestream> stream_;
};

File open(const char *);

Wrapping up
This is one of the reasons why C++ is so hard. In any other language, you'd just return the basic File object and be done with it. Even simple things like returning a value is hard in C++. You need a degree in Computer Science to even attempt it.
[It makes sense when you realize that upon exiting the subroutine, the stack is cleaned and therefore the destructor is called. But we're returning a value that is on that stack ... so the value needs to be copied to another place, first].

Of course, I should've used the std::fstream class to begin with ... but std::fstream doesn't give me fdopen() nor popen(). It feels so unfinished. The documentation doesn't say whether its destructor closes the stream. I didn't feel like fixing fstream, so I went with FILE*.

C++ is a nice language but sometimes trivial things can get quite complex. Situations like this remind you that coding in C++ is really hard.