Living Without the Bunnies

I'm happy to announce the release of Artifastring 1.0, a highly optimized physical simulation of a violin for sound synthesis. http://percival-music.ca/artifastring/

This is an implementation of modal bowed string physical modeling as described by Matthias Demoucron; references are given on the project's main page. The code is available under the GPL v3 or higher.

I had a blast optimising this stuff -- actually, I had far too much fun working on tiny speed improvements. Surprisingly, making the code multithreaded didn't help; the extra overhead of the mutex and condition variables nuked the benefit of modeling each violin string separately for any realistic length of discretized input.

I won't get into everything I did, but here's the main four points:

• As suggested in a conference paper, if a violin string wasn't vibrating and had no input, I didn't bother calculating new values for it -- I just returned 0.0. (this model doesn't include sympathetic vibrations, so there was absolutely no change in the output from doing this). I also detected if the string was barely vibrating (i.e. producing a force of less than 1e-6 on the bridge and each modal velocity was below 1e-4), and then "turned off" that string. This step did change the output, but only very slightly. I picked those values experimentally -- I plucked the string, then looked for the time when I really couldn't hear the sound any more. Admittedly, I was using a pair of £10 headphones... if I had a more expensive sound setup (whether headphones, speakers, a larger room, room with less background noise, etc), these values might be too low. If that were the case, it's easy enough to change those constants.

• Caching any values more complicated than addition or multiplication. In particular, the results of any trig functions were cached -- the amount of extra memory required is trivial (something like 20 kb ?).

• The code calculates eigenvectors of the string for each mode. This is a series of calculatations in the form

  y_i = sin(x*i)   , i = 0, 1, 2... N (where N=56 in my case)
  

  There's a well-known optimization for this: to calculate y_n = sin(nW+B), use y_n = 2*cos(W)*y_(n-1) - y(n-2).

  I could improve on that, though. In my case, there's no B, and I know that sin(0) = 0. So I only need to calculate sin(-W) and cos(W).

  But wait, it gets better! The quoted webpage includes this code:

  double y[3];
  double keep[N];
  for (int i=0; i<N; i++) {
      y[2] = p*y[1] - y[0];
      keep[i] = sqrt_two_div_l * y[2];
      y[0] = y[1];
      y[1] = y[2];
  }
  

  but this involves unnecessary copying. We can avoid that by using a ring buffer of size 2:

  double y[2];
  double keep[N];
  unsigned int latest = 0;
  for (int i=0; i<N; i++) {
      y[latest] = p*y[!latest] - y[latest];
      keep[i] = sqrt_two_div_l * y[latest];
      latest = !latest;
  }
  

  hmm... it just occurred to me that I could avoid using y if I postponed the multiplication by sqrt_two_div_l until after I'd finished calculating the entire array. I mean, go through the array once calculating the sin() values, then go through again multiplying everything by sqrt_two_div_l. I could also get rid of the "latest" variable, and just use i. Hmm... :)

• Fast ring buffer pointer calculation: the traditional way is to use modulo:

  i = (i++) % BUFFER_SIZE;
  

  But modulo is division -- it's slow. Another way is to use a comparison:

  i++;
  if (i == BUFFER_SIZE) {
      i = 0;
  }
  

  That's better, but still uses a bunch of CPU cycles. To improve it even more (courtesy of my brother), if you always use a buffer size which is a power of 2, you can just use an AND!

  i++;
  i &= BUFFER_SIZE_MINUS_ONE;
  

  Of course, we've defined BUFFER_SIZE and BUFFER_SIZE_MINUS_ONE at compile-time.

Sadly, none of this optimization is actually useful. Ok, the "inactive string" trick cut the processing time by more than half -- but the other tricks only gave 2-5% improvements in speed. One of the first lessons for people doing computer science is that running time doesn't matter. The real improvements are made at the level of algorithms... if you can do O(ln(n)) less operations than the next person, then it doesn't matter if you wrote your code in perl and he used hand-crafted assembly. Your code is guaranteed to be faster for any sufficiently large problem.

A long time ago in a city far, far away... I wrote a server in perl which transformed OSC messages into png's of sheet music.

This was back when I was doing my second Bachelor's degree; this time a B.Mus. It was my final project for the 4th-year computer music seminar at the University of Victoria, Canada. It wasn't a new idea, even back in 2006. But few people had done this kind of thing, and in any case, this was just one course in a bachelor degree. I got my 'A', finished the degree, and moved on.

However, I'd mentioned the program on the lilypond-user mailing list, and some people were interested. So interested that one of them recently emailed me to ask if my program was available -- I'd previously given him a "don't redistribute this" tarball a few years ago, with the understanding that I was still working on the infrastructure and was embarrassed about some of it.

Still, that's no reason to hold my earlier work hostage. So I'm announcing the release of Firelily 0.3.

This is old code which I wrote during my second degree (B.Mus). It does not represent my current ability level, with regards to programming ability, project organization, compositional techniques, or anything else.

My research interests have shifted; I now consider this project abandoned.

I'm happy to make it available under the BSD license. I'm not certain that it would be of any use to anybody, but I've had one such request, so here you go! If you would prefer it under a different license (GPL, Apache, artistic, whatever), then just ask and I'll (almost certainly) license it to you under that other license. http://percival-music.ca/software/firelily-0.3.tar.gz

Jamie Bullock has expressed an interest in maintaining and improving Firelily, so if anybody stumbles across this blog post in a few months/years, go look there. :)

Digital Signal Processing can be so awesome it's scary. That's true of any branch of mathematics, of course... but since I just reached this point with convolution to aid the violin physical model, that's what I'll talk about. With audio examples!

We've been working on a program that simulates a violin (called a "physical model" in this the jargon, but it's not physical, and only scientists use the word "model" in that fashion). We have a series of equations that describe how the violin reacts when you pluck it, bow it, put finger(s) on string(s), etc. We're using a lot of approximations and known-to-be-inaccurate physics in this program; we're not trying to push the boundaries of science in this area. We just want something that behaves (and thus sounds) vaguely like a violin, so that people can hear the difference that various bowings make.

In more precise terms, we give the computer a series of instructions about the physical actions of a violinist, and the computer does them. For example:

time  action  extra
0.0   finger  D string  from-nut 0.109
0.0   bow     D string  from-bridge 0.12  -0.1 m/s  0.4 N

This starts playing an E on the D string with a slow upbow (0.1 meters per second) with moderate pressure (0.4 Newtons).

After summing the forces that the four violin strings exert onto the bridge, we get the following sound:

twinkle-plain.wav.mp3

Why does this sound so bad? Well, the approximations we're making will hurt the sound quality to some extent, but the biggest reason is that I'm not doing anything fancy in terms of the playing (yet). I mean, even the most accurate physical simulation can produce horrible noises. Don't believe me? Give a violin to my father. A real violin is the most accurate physical "simulation" imagineable, but it won't produce nice sound unless it's played by a skilled violinist!

But I digress. Right now I'm excited about convolution of audio with an impulse response.

One neat trick I learned while doing my music degree at UVic was that you can make audio sound like it's played in a room by taking it's convolution with the impulse response of a room. Want to pretend that you're singing karaoke in the Sistine Chapel, or the Met opera hall? All you need to do is go to those places and record a clap. That's it -- a single clap ("an impulse") is enough to make any piece of audio sound like it was played in that environment.

When I heard that, I was like "woah!" [sic -- I'm trying to write like a UVic music undergraduate].

But it didn't really sink in, because I hadn't actually tried doing it myself. Like so many parts of mathematics, it doens't make sense until you do it yourself. (which makes this whole blog post rather pointless, of course)

The really cool thing, though, is that this applies to small things, not just big rooms. In particular, it applies to the violin body.

So my supervisor and I went down to the basement, wandered past all the liquid nitrogen (or whatever was in those tanks) and hopped in the anechoic chamber. After plugging in my supervisor's eee 901 netbook into the microphone that somebody left in the room, I hit the violin a few times with a dirty tea-spoon that we grabbed from our lounge.

Yeah. In the Centre for Music Technology, we're just that hardcore. A week ago previous, we ended up measuring the weight of the violin bow using a flat piece of metal (the leg of a speaker stand) as a lever, a ball-point pen as fulcrum, and a slightly-drank pint of milk as the counter-balance to the bow.

We went back upstairs, isolated one of those taps, and filtered it at 80 Hz to get rid of some background noise. (hey, there's noise everywhere, even in an anechoic chamber!)

Here's the resulting sound:

impulse-2048

If you think you missed something, nope -- that tiny click is all there is. It's 2048 samples long (slightly longer than 46 millisceonds), and that's actually way longer than it needs to be. Now we just need to convolve those two pieces of audio together.

The C code that does the magical convolution was simplicity in itself:

double ViolinInstrument::body_impulse(double sample)
{
    body_ringbuffer[body_write_index] = sample;
    double result = 0.0;
    unsigned int bi=body_read_index;
    for (unsigned int ki=0; ki < PC_KERNEL_SIZE; ki++) {
        result += body_ringbuffer[bi] * pc_kernel[ki];
        bi++;
        if (bi >= PC_KERNEL_SIZE) {
            bi -= PC_KERNEL_SIZE;
        }
    }
    update_body_pointers();
    return result;
}

... ok, it you're not a programmer, that probably doesn't look simple. The only important thing to note is that we're multiplying two numbers together, then repeating that process, and adding up all the results together. And also note that the two numbers we're multiply together are simply the raw samples.

(if you know the math and you're worried about the reverse+shift... don't be! That's done by the physics. The second function isn't actually the recorded "bonk" noise; the "bonk" noise itself is already the reverse+shifted version of the function)

The end result of the multiple+add is this:

twinkle-convolution.wav.mp3

Sounds much more like a violin, right? I think the mic placement wasn't ideal, and it was a really cheap violin, but we can make more records later with other violins. The important thing is that the basic system is working... and that multiplying a simple tea-spoon-bonk can magically change the sound in this way.

Oh, and I absolutely cannot resist adding a link to my favorite (and sadly no longer updated) web-comic, talking about convolution: a magical superpower. (the clown is the supervisor -- this will be obvious to anybody in academia)