Sonification - sound of sand - 7

An interesting error - In this program I tried to map the shape to the frequency domain while preserving its shape in the amplitude domain. For each different frequency I resampled the shape. This yields an interesting result, both auditively:

mp3 - test shape

mp3 - star shape

And visually:

The error is in the quantization of frequency. I use this formula in the integer domain. The rounding from float to integer maps many different frequencies (especially the higher frequencies) to the same number of samples. This causes the "blockiness" of the spectrum: 

The sampling itself works nicely. For lower frequencies the sampling distortion is invisible. For higher frequencies the shape gets distorted but is still recognizable

This is how it looks in Audacity:

from Nsound import *
import math

def convert_chaincode_to_x(c,x):
    if c == 1 or c == 0 or c == 7:
        x = x + 1
        return (x)
    elif ( c == 2 or c == 6):
        x = x
        return (x)
    else:
        x = x - 1
        return (x)

def convert_chaincode_to_y(c,y):
    if c == 1 or c == 2 or c == 3:
        y = y + 1
        return(y)
    elif c == 4 or c == 0:
        y = y
        return(y)
    else:
        y = y - 1
        return(y)

def convert_xy_frequency(xy, minxy, maxxy, minf, maxf):
    r = float(maxf - minf)/float(maxxy - minxy)
    f = (xy - minxy)*r + minf
    return (f)

def resample_list(nr_samples, list):
    resampled = Buffer()
    for i in range(nr_samples):
        resampled << list[int(i*(len(list)-1)/(nr_samples -1))]
    return(resampled)

# ==============================

debug1  = False
debug2  = False
debug3  = True
# read a chaincode .chc file that has been generated by SHAPE
infile = open("C:\\Users\\user\\Documents\\shape\\shape\\tiny_test.chc")
instr = infile.read()
infile.close()
if debug1:
    print instr
# parse the input file - split it into words
inwords = instr.split(' ')
if debug1:
    print inwords
# delete anything except the chain code
i = 0
for str in inwords:
    if str.find('0E+0') > -1 :
        break
    i = i + 1
inwords = inwords[i+2:len(inwords)-1]
if debug1:
    print inwords

# fill the x and y buffers with the chaincode values
b_x_chaincode = Buffer()
b_y_chaincode = Buffer()
x = 0
y = 0
for str in inwords:
    c = int(str)
    x = convert_chaincode_to_x(c,x)
    b_x_chaincode << x
    y = convert_chaincode_to_y(c,y)
    b_y_chaincode << y
b_x_chaincode = b_x_chaincode - b_x_chaincode.getMean()
b_y_chaincode = b_y_chaincode - b_y_chaincode.getMean()
if debug2:
    b_x_chaincode.plot("x value from chaincode")
    Plotter.show()
    b_y_chaincode.plot("y value from chaincode")
    Plotter.show()

# copy buffer to list so we can point to it by index
list_x = b_x_chaincode.toList()
list_y = b_y_chaincode.toList()

# generate a frequency modulated x and y signal
b_x_long = Buffer()
b_y_long = Buffer()
min_f = 40.0
max_f = 2000.0
sampling_rate = 44100.0
sound_pixel_length = sampling_rate * 0.01
min_x = b_x_chaincode.getMin()
max_x = b_x_chaincode.getMax()
if debug3:
    i=0
for x in b_x_chaincode:
    frequency = convert_xy_frequency(x, min_x, max_x, min_f, max_f)
    nr_samples = int(round(sampling_rate / frequency))
    temp = Buffer()
    temp << resample_list(nr_samples, list_x)
    for j in range(int(math.ceil(sound_pixel_length/nr_samples))):
        b_x_long << temp
    if debug3:
        i = i+1
        p = int(len(b_x_chaincode)/5)
        if i%p == 0:
            b = Buffer()
            b = resample_list(nr_samples, list_x)
            b.plot(frequency)
            Plotter.show()
b_x_long.normalize()

min_y = b_y_chaincode.getMin()
max_y = b_y_chaincode.getMax()
for y in b_y_chaincode:
    frequency = convert_xy_frequency(y, min_y, max_y, min_f, max_f)
    nr_samples = int(round(sampling_rate / frequency))
    temp = Buffer()
    temp << resample_list(nr_samples, list_y)
    for i in range(int(math.ceil(sound_pixel_length/nr_samples))):
        b_y_long << temp
b_y_long.normalize()

# code the x and y signal into the left and right channel of an audio stream
# write the audio stream into a .wav file
a = AudioStream(44100.0, 2)
a[0] = b_x_long
a[1] = b_y_long
a.writeWavefile("C:\\Users\\user\\Documents\\shape\\shape\\tiny_test_xy_freq_shape.wav")

Sonification - sound of sand - 6

Frequency domain - What we did with volume in the previous post we now do with frequencies. The python program below converts the X- and Y-components of the shape into varying frequencies. The input shapes are the same as in the previous example. You can hear the sound samples here:

mp3 of small test shape

mp3 of big star shape

Again the sounds are very technical and abstract. I'll leave it like that for the moment. I'm still learning about what I can do. You can see the spectrum of the signal as it appears in Audacity:

 

It is easy to see that the spectrum now has the shape of the X- and Y- components of the shapes.

Next time I'll try to map the shape to the frequency domain while using the original waveform and not a pure sine signal. To do this I'll have to resample the waveform. This is a bit tricky to program but should give a more interesting signal.

In the meantime I've found a lot of interesting articles about shape sonification. In the future I'll try to do things with the curvature of the shape. And I haven't done anything with the spectral components yet.

from Nsound import *
import math

def convert_chaincode_to_x(c,x):
    if c == 1 or c == 0 or c == 7:
        x = x + 1
        return (x)
    elif ( c == 2 or c == 6):
        x = x
        return (x)
    else:
        x = x - 1
        return (x)

def convert_chaincode_to_y(c,y):
    if c == 1 or c == 2 or c == 3:
        y = y + 1
        return(y)
    elif c == 4 or c == 0:
        y = y
        return(y)
    else:
        y = y - 1
        return(y)

def sine_duration_frequency(duration, frequency):
    g = Generator(44100.0)
    length = math.ceil(float(duration) * float(frequency))/float(frequency)
    return g.drawSine(length, frequency)

def convert_xy_frequency(xy, minxy, maxxy, minf, maxf):
    r = float(maxf - minf)/float(maxxy - minxy)
    f = (xy - minxy)*r + minf
    return (f)

# ==============================

debug1  = False
debug2  = False
debug2a = False
debug3  = True

# read a chaincode .chc file that has been generated by SHAPE
infile = open("C:\\Users\\user\\Documents\\shape\\shape\\tiny_test.chc")
instr = infile.read()
infile.close()
if debug1:
    print instr

# parse the input file - split it into words
inwords = instr.split(' ')
if debug1:    print inwords

# delete anything except the chain code
i = 0
for str in inwords:
    if str.find('0E+0') > -1 :
        break
    i = i + 1
inwords = inwords[i+2:len(inwords)-1]
if debug1:
    print inwords

# fill the x and y buffer with the chaincode values
b_x_chaincode = Buffer()
b_y_chaincode = Buffer()

x = 0
y = 0
for str in inwords:
    c = int(str)
    x = convert_chaincode_to_x(c,x)
    b_x_chaincode << x
    y = convert_chaincode_to_y(c,y)
    b_y_chaincode << y

    if debug2a:
        print c
        print x

if debug2:
    b_x_chaincode.plot("x value from chaincode")
    Plotter.show()
    b_y_chaincode.plot("y value from chaincode")
    Plotter.show()

b_x_chaincode = b_x_chaincode - b_x_chaincode.getMean()
b_y_chaincode = b_y_chaincode - b_y_chaincode.getMean()

# generate a frequency modulated x and y signal
temp = Buffer()
b_x_long = Buffer()
b_y_long = Buffer()

soundpixel_length = 0.02
min_f = 40.0
max_f = 5000.0

min_x = b_x_chaincode.getMin()
max_x = b_x_chaincode.getMax()

if debug3:
    i=0

for x in b_x_chaincode:
    frequency = convert_xy_frequency(x, min_x, max_x, min_f, max_f)
    b_x_long << sine_duration_frequency(soundpixel_length, frequency)

    if debug3:
        i = i+1
        p = int(len(b_x_chaincode)/4)
        if i%p == 0:
            b = Buffer()
            b = sine_duration_frequency(soundpixel_length, frequency)
            b.plot(frequency)
            Plotter.show()

b_x_long.normalize()

min_y = b_y_chaincode.getMin()
max_y = b_y_chaincode.getMax()
for y in b_y_chaincode:
    frequency = convert_xy_frequency(y, min_y, max_y, min_f, max_f)
    b_y_long << sine_duration_frequency(soundpixel_length, frequency)
b_y_long.normalize()

# code the x and y signal into the left and right channel of an audio stream
# write the audio stream into a .wav file
a = AudioStream(44100.0, 2)
a[0] = b_x_long
a[1] = b_y_long
a.writeWavefile("C:\\Users\\user\\Documents\\shape\\shape\\tiny_test xy_freq_chaincode.wav")



Sonification - sound of sand - 5

First sonification results - I've managed to sonify my first experimental shapes. The results are encouraging. They are not very musical yet but with some optimization they might get interesting.

Small random test shape - We use the small experimental test shape for our first sonification. It has been used in my previous blog post so we're quite familiar with its properties:

The python software (see below) first extracts the X and Y values that we get while traversing the outline in a counterclockwise direction. It transforms these values into two sound waves:

Then it concatenates these sound waves (they are extremely short at a sampling rate of 44100 Hz) into a longer sound sample and it modulates the amplitude of this signal using the same sound shape.

This is a nice fractal twist and it feels like a very natural thing to do with the signal. This way the signal is made self-similar on two levels. (I don't think it would be feasible to add more than two levels of self-similarity, the signal would get too long.) 

Then we put the X and Y signal into the left and right channels of an audio stream. Again this feels like a very natural thing to do with the signal.

Big star - In the same way we generate a sound sample for the star shape of one of the previous experiments. This gives comparable results.

Notice how recognizable the X and Y values are in the signal. In the X-direction the star has two points. In the Y-direction it has only one point.

Discussion of the results - You can listen to the resulting sound files here:

mp3 file of test shape

mp3 file of star shape

The outline of a shape has been mapped directly to a 44100 Hz sampling rate. This means that a smaller shape will generate a higher note and a shorter sample. For the moment I will leave it like that because it's the most natural mapping. Later I will explore other possibilities. This means that the test shape produces a high mosquito like drone. And the star shape produces a low atmospheric, almost inaudible soundscape. Both sounds are quite abstract and unmusical and this is how it should be for the moment.

This is how the two sound files look in Audacity:

test shape

star

And here you see how the signal-in-signal looks in audacity if you zoom into the details.

Now we've used the amplitude domain to map shapes into sound. I'll also try to use the frequency domain for this mapping.

Python program - I'm not sure this will run correctly if you copy it directly into your Python environment. I'm using Python-XY and this has all the necessary modules pre-installed. And Blogger may destroy some of the whitespace. So this could explain some unexpected bugs.
If I could do things in a more Pythonesque way then I'm open for comments.

from Nsound import *
debug = True

# read a chaincode .chc file that has been generated by SHAPE
infile = open("C:\\Users\\user\\Documents\\shape\\shape\\04 star.chc")
instr = infile.read()
infile.close()
if debug:
    print instr

# parse the input file - split it into words
inwords = instr.split(' ')
if debug:
    print inwords

# delete anything except the chain code
i = 0
for str in inwords:
    if str.find('0E+0') > -1 :
        break
    i = i + 1
inwords = inwords[i+2:len(inwords)-1]
if debug:
    print inwords

# fill the x and y buffer with the chaincode values
b_x_chaincode = Buffer()
b_y_chaincode = Buffer()
x = 0
y = 0
for str in inwords:
    c = int(str)

    # convert a chaincode into a plot of the x value against time
    if ((c == 1 or c == 0) or c == 7):
        x = x + 1
    elif ( c == 2 or c == 6):
        x = x
    else:
        x = x - 1
    b_x_chaincode << x

    # convert a chaincode into a plot of the x value against time
    if ((c == 1 or c == 2) or c == 3):
        y = y + 1
    elif ( c == 4 or c == 0):
        y = y
    else:
        y = y - 1
    b_y_chaincode << y

b_x_chaincode = b_x_chaincode - b_x_chaincode.getMean()
b_y_chaincode = b_y_chaincode - b_y_chaincode.getMean()
if debug:
    b_x_chaincode.plot("x plot from .chc file")
    Plotter.show()
    b_y_chaincode.plot("y plot from .chc file")
    Plotter.show()

# generate an amplitude modulated x and y signal
b_x_long = Buffer()
b_y_long = Buffer()
for level in b_x_chaincode:
    b_x_long << b_x_chaincode * level
for level in b_y_chaincode:   
    b_y_long << b_y_chaincode * level

# normalize to prevent clipping of the output signal
b_x_long.normalize()
b_y_long.normalize()
if debug:
    b_x_long.plot("x plot from .chc file")
    Plotter.show()
    b_y_long.plot("y plot from .chc file")
    Plotter.show()

# make sure that the sound sample is long enough to hear anything
while len(b_x_long) < 200000:
    b_x_long << b_x_long
    b_y_long << b_y_long

# code the x and y signal into the left and right channel of an audio stream
# write the audio stream into a .wav file
a = AudioStream(44100.0, 2)
a[0] = b_x_long
a[1] = b_y_long
a.writeWavefile("C:\\Users\\user\\Documents\\shape\\shape\\04 star xy_chaincode.wav")

Sonification - sound of sand - 4

Chain coding - What does the SHAPE software do with a shape? It is interesting to get a clear understanding of what's happening. The transformation of a shape into a fourier spectrum goes in two steps. The first one is chain coding. The software determines the pixellated outline of a shape and transforms it into a chain-coded string of numbers. Each step along the contour is translated into a number in this way:

 For example the leftmost edge of the pixellated shape above is translated into this number string:

...... 4 4 4 4 4 4   5 5 5 5 5 5 5   6 6 6 6 6 6 6 6 6 6 6 6 6 6 6   0 0 0 0 0 0 0 0 .......

We can use this a the starting point for a very direct, very primitive kind of sonification. Any shape can be decomposed into pixel-sized increments in the X and Y directions. Using ChcViewer.exe it is possible to visualize the X and Y components of a 2-D figure. For example if you start going counterclockwise from the position of the arrow:

Then you get this plot of the X-values of the shape as you go anticlockwise around its contour:

And this plot of the Y-values of the shape as you go around its contour:

You can see immediately that these shapes can be transformed into soundwaves easily. And we have a lot of degrees of freedom while combining the X and Y waveforms: we can add them together in different ratios and we can time-shift them with respect to each other.

Next time we'll see if we can get some sounds using python and nsound.

Note
You can make an X-plot by substituting these values in the chain code:
1,0,7 => 1
2,6   => 0
3,4,5 => 7

And an Y-plot by substituting these values:
1,2,3 => 1
4,0   => 0
5,6,7 => 7

References
Shape software - http://lbm.ab.a.u-tokyo.ac.jp/~iwata/shape/index.html
Python XY - http://code.google.com/p/pythonxy/
Nsound (included in Python XY) - http://nsound.sourceforge.net/users_guide/index.html

Sonification - sound of sand - 3

Testing - To get a feeling for elliptic fourier analysis with the SHAPE software I made a set of test shapes:

Then I determined their outline with ChainCoder.exe and  calculated the elliptic spectra with CHC2NEF.exe. The I made a plot with OpenOffice Calc. These are the results with some analysis:

Note: Below I've replaced the primary component of 1.0 with 0 otherwise the "harmonics" would be totally invisible. That's why you see nothing in place 1.

Circle - Theoretically a perfect circle should not have a spectrum. It should only have the first component. But my hand-drawn digtized circle made in MS-paint is not perfectly symmetrical and it has rough edges. So there are still some "harmonics" but these are much fainter than the harmonics of the other shapes (a factor of 100: 10^-3 instead of 10^-1). I assume that I'm just seeing "random noise" and "quantization noise" in this spectrum.

Triangle - One would expect the order-3 harmonics to be dominant for a triangular shape but that is not the case at all. There is no immediately visible correlation between a figure and its spectrum. This is even more obvious if one looks at the shape as more harmonics are added. The second harmonic is already sufficient for a nice triangular shape. I edited the .nef output files by hand and plotted them with NefViewer.exe.

Note: Below I've replaced the primary component of 1.0 with 0 otherwise the "harmonics" would be totally invisible. That's why you see nothing in place 1.

But if we look at the spectral components then we see that 10 harmonics is not really sufficient for a nice sharp triangle. We need at least 20 harmonics.

Star - Here you would expect that the order 2, 3 and 6 harmonics would be dominant. And you would expect that the spectrum of the triangle looks limilar to the spectrum of the star. But things are not that intuitive.

Square - Here there's a surprising similarity between normal fourier spectra and 2-D fourier spectra. A square wave spectrum has only odd harmonics. And this 2-D spectrum has also only odd harmonics!

Rectangle - There is a lot of similarity between the square and the rectangle. Also only odd harmonics, but the signs and ratios of the harmonic components are different from the square. It is interesting to see how the rectangle is constructed from the different harmonics.

Sonification - sound of sand - 2

Progress - Last week I did some more research on the Internet and I found new and better keywords. This lead me to better articles and promising software. The techniques I'm looking for are popular in paleontology and they're used for comparing species.
Norman MacLeod of the The Natural History Museum in London has put a series of articles online about quantitative analysis of paleontological data. These describe how shapes can be quantified, compared and classified mathematically. Link here.

The most interesting technique is "Elliptic Fourier Analysis". Using this technique any shape, even the most complex one, can be described by a series of 2-dimensional harmonics. This is demonstrated by this very cool picture:

The Centre Cannot Hold II: Elliptic Fourier Analysis by Norman MacLeod

The mathematics is easy to understand if you've got a technical background. And even better, the author has a list of relevant software packages that can be downloaded from the Internet.

More background literature - The morphometric techniques are explained in other documents as well:

• A very short introduction to circular harmonics.
• The morphometrics website of Stony Brook. It also has a long list of morphometric software packages.
• A personal site with morphometric research (used in archaeology) and some software.

Most promising image software - The most promising software is SHAPE, written by Hiroyoshi Iwata. It is a free software for quantitative evaluation of biological shapes based on elliptic fourier analysis. Link here. If I can get this working it will solve the spectrum analysis part of my project. Then I will only have to find software for the sonification.

Other relevant image software - I've found more interesting software that could be useful for my project but I'll only investigate that further if SHAPE does not work as expected:

• Potrace and autotrace : these software packages translate a bitmap picture into a scalable vector drawing. I don't know how easy it would be to do fourier analysis on their output. (E: Thanks anyway for the tip!)
• PyNGL : a Python language module used to visualize scientific data. A very interesting software package but it has no image processing functionality. But it might be useful for GPS mapping plots.
• SHTOOLS : an archive of software that can be used for spectral analyses on the sphere. This may do the job of spherical harmonic analysis. But it looks too powerful and too complex for my purposes.
• Fiji and ImageJ : image processing software written in Java. Fiji has a GUI with menus.

Less progress with sonification - I haven't found a software package that does simple spectral synthesis just like I need it. But I have found some interesting websites:

• SonEnvir : a research project that investigates sonification in a number of scientific disciplines. It uses SuperCollider.
• MAX/MSP was mentioned as useful software. In the meantime I realize that PureData is a branch of this development and that it has fourier transforms. So maybe I'll have to dive into that too.

Next time maybe real results!

Sonification project - sounds of sand

Sand spectrum - I want to listen to the sound of the shapes of sand grains. This seems a weird idea but it is a standard geological technique. The outline of a sand grain is first extracted by image processing:

And then the roughness of the grain outline is translated into a series of harmonics by Fourier analysis. This is a more formal description of sphericity, roundness or roughness:

Different grains with different provenance have different spectra:

The technique is described in these research papers:

• On the shapes of natural sand grains - David R. Barclay and Michael J. Buckingham - Published 21 February 2009 - JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 114 - Link
• Particle Shape Characterisation using Fourier Analysis - Elisabeth T. Bowman, Kenichi Saga & Tom W. Drummond - CUED/D-SoWTR315 (2000) - Link

Software tools discarded - I've spent a lot of time searching for tools. I need tools for image processing and sound synthesis. These don't often go together. I've discarded the following options:

• Matlab - Has image and sound and mathematics. But is prohibitively expensive and probably a steep learning curve.
• Scilab and SIP (Scilab Image Processing Toolbox) - Free and functionality looks sufficient. But I'm worried about stability, reliability and the steep learning curve.
• Processing - Not many examples of image processing, I'm not sure it's powerful enough.
• Supercollider and PureData - I'm afraid of the learning curve.

Software tools chosen - Finally I've settled on tools written in Python. I already know the language, it is portable to both Windows and Linux and it has both image processing and sound synthesis tools. These are the tools I'm going forward with:

• Python Vision - This seems a workable set of image processing tools. I'll just try it for a first start. It consists of a set of image processing tools, described here.
• Python XY - Most of the tool set can be installed in Windows 7 in one go from here. You have to deinstall all other Python software first.
• Two additional modules:
• Pymorph - Downloadable from here. Install by: python setup.py install. No problems here.
• Mahotas - Downloadable from here. This one is more tricky to get working. You get a well known error that is described here. But the solution is different from what is written here. You don't need to install MinGW and you don't need to set the path. I think Python-XY has done that already. You just need to install the module with a different command: python setup.py install build --compiler=mingw32. This worked for me after some googling and experimentation.

I've not yet chosen which tools I'll use for sonification of the images. But there seem to be sufficient tools available. At first sight this looks sufficient:

• Nsound - Seems sufficient. I'm still a bit dubious about the 'buffer' format. Can I write random data into a buffer or am I limited by the functions in the nsound library? The usage of the library looks a bit weird to me.

We'll see if these choices work. I think it would be difficult to find a better match between learning curve and functional power. But I'm open for advice!