Sound waves as equations

Author: Mario Román
GitHub: https://github.com/M42/sum-of-waves

Introduction

This notebook contains code to a sound file (specifically, a 16bit PCM WAV file) into a pdf file containing the mathematical formula of the sound wave. For an example of input and output, you can check the files "test2.wav" and "output.pdf" in the GitHub repository.

We use the Fast Fourier Transform algorithm to express the wave as a sum of sinusoidal functions, and the complete sum is translated to a LaTeX file, which can be used to produce the final pdf.

Plotting wave files

In this first step, we are going to open the sound file and make a simple plot of the signal, using matplotlib. The sound file:

Must be a 16bit PCM WAV file.
Must be mono; not stereo.

To convert multiple sound files to this format, you can use sound editing programs, such as Audacity.



In [1]:

    
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import wave
import sys



In [2]:

    
# Only opens 16bit PCM WAV files.
# The format can be changed in Audacity.
spf = wave.open('test2.wav','r')



In [3]:

    
# Only opens Mono files
if spf.getnchannels() == 2:
    print 'Just mono files'
    sys.exit(0)



In [4]:

    
# Extracts the signal
signal = spf.readframes(-1)
signal = np.fromstring(signal, 'Int16')



In [5]:

    
# Prints the signal
plt.figure(1)
plt.title('Signal Wave...')
plt.plot(signal)









    Out[5]:





[<matplotlib.lines.Line2D at 0x7f5366af7c10>]

Fast Fourier Transform

We use the Fast Fourier Transform algorithm from the numpy library.



In [6]:

    
Y = np.fft.fft(signal)
plt.figure(2)
plt.title('Fourier transform')
plt.plot(Y)









    



/usr/lib/python2.7/dist-packages/numpy/core/numeric.py:460: ComplexWarning: Casting complex values to real discards the imaginary part
  return array(a, dtype, copy=False, order=order)






    Out[6]:





[<matplotlib.lines.Line2D at 0x7f5366621e90>]



In [7]:

    
Y









    Out[7]:





array([  3.79386760e+07    +0.j        ,  -1.11978943e+04-25406.69537752j,
        -1.67714729e+04 +5782.91817767j, ...,
        -1.09229381e+04+21810.02214398j,  -1.67714729e+04 -5782.91817766j,
        -1.11978943e+04+25406.69537752j])

Writing the Fast Fourier Transform in Latex

We output the formula to a document, using the Latex syntax. The final pdf file will be too big to be handled by one latex file. Therefore, we split it in multiple files, included in a main latex file.



In [14]:

    
lfile = open("latex.tex","w")
N = str(len(Y))

lfile.write("\\documentclass{scrartcl}\n\\usepackage{amsmath}\n\\begin{document}\n")
lfile.write("\\allowdisplaybreaks[1]\n")
lfile.write("\\title{Hungarian March}\n\\subtitle{First 8 seconds!}\n")


# Huge files cause a memory error in pdflatex
# We split in small files

j = 0
ffile = open("latex"+str(0)+".tex","w")
ffile.write("\\begin{align*}\n")
for k in range(len(Y)):
    j = j+1
    if j==28:
        j=0
        ffile.write("\\end{align*}\n")
        lfile.write("\\input{latex"+str((k-1)/28)+".tex}\n")
        ffile.close()
        ffile = open("latex"+str(k/28)+".tex","w")
        ffile.write("\\begin{align*}\n")
    ffile.write(str(Y[k]) + "e^{ \\frac{j \\tau}{" + N +  "} " + str(k) +  "t } && + \\\\ \n")
ffile.write("\\end{align*}\n")

lfile.write("\\end{document}")
lfile.close()

Inverse Fast Fourier Transform

We will check the formula is correct performing the inverse transformation on the formula to obtain the original signal and outputting it into a new .wav file, which will sound exactly as the original one.



In [8]:

    
# Inverses the fast fourier transform
yinv = np.fft.ifft(Y)
newsignal = yinv.real.astype(np.int16)
newsignal = newsignal.copy(order='C')



In [9]:

    
# Recreates the wave file
spfw = wave.open('test3.wav','w')
spfw.setnchannels(spf.getnchannels())
spfw.setsampwidth(spf.getsampwidth())
spfw.setframerate(spf.getframerate())
spfw.writeframes(newsignal)