Python: convert frequency response to impulse response

Question

I'm using an acoustic simulator that gives me a (complex) frequency response for any (positive real) frequency I feed in.

I would like to generate a 1024 point impulse response.

I believe I have the basic technique here: https://dsp.stackexchange.com/questions/13883/convert-frequency-response-to-impulse-response

However, I would like to implement it from the command line.

The simulator generates A file that looks like this (note I have shortened the numbers to make it more readable):

# Run name  = pro-ch0-pc1-gmr
# Run owner = umby
# Cfg file  = sphere_source.cfg
# Frequency = 0.000000 + i*1077.566280 (171.5000 Hz)
#
# 1    2  3  4  5       6       7        8        9        10        11      12      13    
# imic xm ym zm re(inc) im(inc) abs(inc) re(scat) im(scat) abs(scat) re(tot) im(tot) abs(tot)
0 +1.4E+00 +0.0E+00 +9.4E-02 -9.8E-04 -5.2E-02 +5.2E-02 -5.4E-03 +1.2E-02 +1.3E-02 -6.4E-03 -4.0E-02 +4.0E-02
1 +1.4E+00 +0.0E+00 +1.8E-01 -3.8E-03 -5.2E-02 +5.2E-02 -5.1E-03 +1.3E-02 +1.3E-02 -9.0E-03 -3.9E-02 +4.0E-02
:
:
etc

• The Second column is the angular frequency in radians/sec

• The third column is the frequency in Hz

• Columns 11 and 12 are Re(z) and Im(z), i.e. Real and imaginary part of the frequency response for this frequency

• Column 13 is the magnitude of the frequency response, I'm guessing this can be discarded along with one of the first two columns

So my question is: how can I process this file and get out an impulse response?

Answer 1

Here is the necessary code for transforming the frequency response to the impulse response (it contains a test signal; the signal is transformed to the frequency domain and then recovered, demonstrating that the algorithm does work):

import numpy

#http://stackoverflow.com/questions/9062387/ifft-of-symmetric-spectrum/9062837#9062837

# test
waveform = [1,2,3,4,5,6]
print( waveform )

fullSpectrum = list( numpy.fft.fft( foo ) )
for w in fullSpectrum:
    print w

# In my practical application, I would be receiving frequencies DC through Nyquist,
#  i.e. the half of the spectrum containing the positive frequencies, and I would 
#  need to reconstitute the negative frequencies using symmetry:
#
DC_to_Nyquist = list( numpy.fft.fft( foo ) )[0:4]


# RESTORE FULL SPECTRUM BY SYMMETRY
DC      = DC_to_Nyquist[0].real + 0j
nyquist = DC_to_Nyquist[-1].real + 0j

# discard DC and Nyquist from list
Z = DC_to_Nyquist[1:-1]

Z_reverse = Z[::-1]
Z_conj_reverse = [ w.conjugate() for w in Z_reverse ]

spectrum = [DC] + Z + [nyquist] + Z_conj_reverse

print( "\nSpectrum" )
for w in spectrum:
    print( w )

# INVERSE FFT
timedomain = list( numpy.fft.ifft( spectrum ) )

#ir = [s.real for s in timedomain]

print( "\nTimedomain" )
for w in timedomain:
    print w

And here is the boilerplate code that extracts the complex numbers from the input file:

#!/usr/bin/python

import sys
import csv

infile = sys.argv[1]

Z = [ ]

with open( infile ) as f_in:
    lines = csv.reader( f_in, delimiter=' ', skipinitialspace=True )
    for n, line in enumerate( lines ):
        if n > 6:
            re = float( line[11] )
            im = float( line[12] )

            Z.append( complex( re, im ) )

print( Z )

EDIT: Better is:

Z_conj_reverse = [ w.conjugate() for w in reversed( Z ) ]