Thanks Bassik for pointing to the Farina article. Particularly interesting is the section on Kirkeby inverse filters for compensating spectral deficiencies in measurement equipment.

Meanwhile, I have made a patch to inspect characteristics of [chirp~] output in detail, with horizontal and vertical zoom options on the IR plot, and spectrum magnitude and phase plots. Real-signal chirps can never have a completely white spectrum. It turns out that logarithmic chirps are quite sensitive to parameter settings, while linear chirps are less so. Therefore, a log chirp may have considerable ripple near start and stop frequency, but I have also found settings which induce low ripple (+/- 0.3 or 0.4 db). It may have been coincidence that I stumbled upon such favourable settings. There is probably no simple rule to derive them numerically.

Because the log chirp is so sensitive to parameter settings, it needs special attention. For an IR measurement tool, we could provide 'approved' chirps in stereo wave files (32 bit floating point) together with their inverse. Chirp generator / tester and IR measurement tool would then be two different Pd root patches which can be developed independantly.

So, the chirp generator / tester is attached to this post. Again, there is only an OSX binary for [chirp~] (same one as before). Sorry. Building on other platforms is easy pie if you have developer tools installed and know how to operate them, but otherwise, don't spoil your time on it. I'll do it later on a friend's computer.

Katja

http://www.pdpatchrepo.info/hurleur/chirp~test.zip

Hello Katjav,

hope you are well.

I have performed some test on a real case and I have to say that despite the first 2 measurements worked very well, after that the deconvolution process stopped to work properly.

Attached there is one of the wrong results I have obtained from the soft.

you can clearly hear that instead of the IR you get the inverse filter at the same time where the IR should be.

Is it a user problem?
Is it a multi measurement issue with the software?

I will investigate further before any conclusion.

thanks

Bassik, there may be some confusion. This last patch chirp~test.pd is only a patch to test the qualities of the generated chirp itself.

I am kind of swallowed by details of the method. Instead of quickly producing a working IR measurement tool, I am now trying to refine the method, starting with refining the chirp itself. Last time I mentioned the sensitivity for parameter settings. Now I have found a systematic method to optimize a chirp. The formula is different from the regular ones as found in texts. This has lead to yet another chirp object [expochirp~] (yes it's exponential, not logarithmic, also see Farina's 2007 article).

The regular exponential chirp formula's make a chirp start at sine phase to avoid a step, but they do not care about the phase at the end of the function, therefore it may well end with a step and produce excessive ripple at the high spectrum end. This is traditionally relieved by windowing, or by manually cutting the chirp at the last zero-crossing as Farina suggests. My new formula computes the chirp in such a way that it will start and stop in sine phase. No need for windowing or cutting, and the ripples are minimized.

Still there will always remain ripples in a chirp's spectrum, and I am now trying to produce inversion-filters to clean up a system's measured IR. For the high spectrum end, this is fairly well possible. But in the low spectrum end, where it is more important, the ripple character of exponential chirps is such that conventional methods fail. I did find some solution but it's not without compromise. Therefore I am investigating alternative methods.

I am sorry that proceedings are slow now. But my hope is to contribute to a good quality tool in the end. I'll post details on the new chirp formula soon.

Katja

Hello Katjav,

I haven't used the latest chirp~test patch but the previous version of the soft with teh old chirp computation.

So it was supposed to work at least in terms of recording the chirp and deconvolve it; but as i said it needs more testing as it was just a one off mesurements set the one i performed.

I am sorry as well I haven't contributed much in teh last time but I am full of stuff at work.

Anyway I will post a comparison between the PD measurement and a professional soft measurement in teh same positions soon.

Also I am not able to use your latest chirp~test patch as I am on windows.

thank you for keeping the study alive!!!

HEllo Katjav,

As promised I have performed some comaprisons between the IR measured by EASERA software and the one done with the PD patch we are developing.

As said, I have used the old chirp computation patch and have obtained useful IRs only for the first 2 measurements.

comparison in RT results are in the attached file showing:

• RT results
• Waterfall plots (where possible to compute)
• Magnitude plots

Interesting is to see that EASERA is not able to compute the dB SPL of the wav file as it doesn't have any reference of the microphone calibration.
I need to investigate more if I can insert this info for the wav file so for the moment the waterfall plots are not comparable (EASERA measurements are in dB SPL while PD measurements in dB FS).

Anyway this is just a start and you can see that the RT measurements are quite similar apart for the 125 Hz results.

Hope this helps

More to come

Cheers

http://www.pdpatchrepo.info/hurleur/RT.xls

Hey Bassik, thanks so much for your test results.

Apart from the 125 Hz the found RT are comparable indeed. But the magnitude plots of the second case? Should M0002_S01_R02.etm be comparable to Test 2b proc.wav? Even without a reference, the figure should show a similar shape because it is in deciBels, no?

But well, at least it shows that the Pd patch is not complete nonsense. We're on the right track.

I wonder how the magnitude of a full IR is computed in EASERA. For a 20 second chirp the IR it is at least 40 seconds, anticausal and causal part taken together. Or is the analysis limited to the causal part? No, I would think that full IR means the whole thing with the left side included. Pd can not do such long FFT's, 131072 points is the maximum. Probably the magnitude could be computed as an average of several overlapping frames. But the phase information is lost then.

I consider writing an FFT object for Pd that operates on a named buffer, and works with 64 bit floats internally. Maybe it is then possible to do larger FFT's. This is also important for calculating inverse filters (correction filters), since the full IR spectrum is sometimes needed for that. There is also something wrong with [fft~] and [rfft~] phase interpretation: it does not do a centered FFT. So with an IR centered you get these alternating phases instead of a zero-phase representation. Therefore, complex multiplication of spectra will give wrong results. Quite annoying. A new FFT object should have an option to select centered and non-centered FFT so you could work well with centered IR, and with causal part only. Hmm, plans....

Katja

Here is my new formula / method for an exponential test chirp.

sin((((pi/2^p)*L)/ln(2^p)) * e^((n/N)*ln(2^p)))

condition: ((pi/2^p)*L/ln(2^p)) = M*2*pi

M = a positive integer
L = a decimal number (floating point number) expressing chirpsize in unity samples
N = ceil(L), that is L rounded up to an integer, the real arraysize
p = the integer number of octaves in the chirp (user argument)

With this method, you get a chirp with evenly-spaced octaves. The chirp starts in sine-phase and also ends in (approximate) sine-phase, thus minimizing the ripple height. It is even the case that each octave starts and ends in approximate sine-phase, whatever use that may have.

The crucial element in the method is the exponential curve at index 0 and N.

e^((0/N)*ln(2^p)) = e^0 = 1
e^((N/N)*ln(2^p)) = e^ln(2^p) = 2^p

Notice that ((pi/2^p)*L/ln(2^p)) is a constant. If this constant is an integer multiple of 2*pi, the exponential curve will start with 2*pi and end with an integer multiple of 2*pi, since 2^p is also an integer.

Using this method, the start frequency depends on the sampling rate and number of octaves. The start frequency is (sr/2)/2^p. For example: (44100/2)/(2^10) = 21.53 Hz. The stop frequency is always near Nyquist, with a small deviance resulting from the rounding of L to N. For large N, the deviance is negligible for all practical purposes.

If y[n] is the chirp, the inverse chirp is formulated:

y[N-n] * -ln(1/2^p) * (2^(p/N))^-n

This calculates an inverse chirp with average power equal to the forward chirp. That means, more than unity peak level in the high frequencies, and less than unity in the low frequencies.

The user supplies a maximum length argument in nr of samples, and the required number of octaves. From this, a new object [expochirp~] calculates the first ideal length L below the maximum length. It returns the integer chirp size N as a message (so you can resize arrays), and starts feeding chirp and inverse chirp to the signal outputs.

Chirp and inverse chirp can be stored in a 32 bit float stereo .wav file, but the risk is that someone will play the stereo file in a regular soundfile player, which for the inverse chirp could lead to severe clipping. The reason for normalising to amplitude instead of unity peak level is that no further information is now needed about the chirp for normalisation after deconvolution. That normalisation will be simply N/2, like it would be the case for a linear chirp. For safety, chirp&inverse could be stored as .txt files instead of .wav. I don't know yet what is wise.

Although my expochirp method so far guarantuees a straightforward chirp with minimized ripple, the ripple is still considerable, and at the low frequency side quite problematic in my view. Some speaker systems may already have their own resonance frequencies in the region of the start frequency. Therefore this is exactly the region where a flat chirp spectrum is most important, to not push eventual resonances to further height. I have experimented with correction filters on the IR, but it would be better to relieve the problem at it's origin again, the shape of the chirp itself. I am now trying to find a fade-in window for that purpose. Windowing means amplitude modulation, causing sum- and difference frequencies. A proper window may cancel ripples, but a bad window could easily amplify them, or not? For me this is a new research topic. Only when results are acceptable I will finish [expochirp~] and share it.

Katja

@katjav said:

Hey Bassik, thanks so much for your test results.

Apart from the 125 Hz the found RT are comparable indeed. But the magnitude plots of the second case? Should M0002_S01_R02.etm be comparable to Test 2b proc.wav? Even without a reference, the figure should show a similar shape because it is in deciBels, no?

But well, at least it shows that the Pd patch is not complete nonsense. We're on the right track.

I wonder how the magnitude of a full IR is computed in EASERA. For a 20 second chirp the IR it is at least 40 seconds, anticausal and causal part taken together. Or is the analysis limited to the causal part? No, I would think that full IR means the whole thing with the left side included. Pd can not do such long FFT's, 131072 points is the maximum. Probably the magnitude could be computed as an average of several overlapping frames. But the phase information is lost then.

I consider writing an FFT object for Pd that operates on a named buffer, and works with 64 bit floats internally. Maybe it is then possible to do larger FFT's. This is also important for calculating inverse filters (correction filters), since the full IR spectrum is sometimes needed for that. There is also something wrong with [fft~] and [rfft~] phase interpretation: it does not do a centered FFT. So with an IR centered you get these alternating phases instead of a zero-phase representation. Therefore, complex multiplication of spectra will give wrong results. Quite annoying. A new FFT object should have an option to select centered and non-centered FFT so you could work well with centered IR, and with causal part only. Hmm, plans....

Katja

Hello Katja,

Magnitude is expressed in dBFS and for definition our recorded signal should not go above 0 dBFS hence clipping in the analog domain.
So there is something wrong in EASERA when computing our IRs; I haven't figure out yet as the manual does not contain any reference to wav file loading and analysis.

Regarding which part is computed in Easera:
For acoustic measurements (RT, STI etc) only the casual part should be used.
The purpose of the swept sine technique is to get all the reproduction system non-linearities in the anti casual part (leftmost part) of the deconvoluted IR in order to then eliminate them.
So all the analysis should be done on a IR that has the time 0 where the peak of the IR is; also it will be then windowed to eliminate the late part that is not useful to be investigated as it might contain other artifacts of the deconvolution process due to e.g. the noise in the environment during measurements.

I hope I have understood your question.

Glad this keeps you thinking about new objects for PD.

I will review your post on the new chirp formula.

Thank you

Bassik

Hello Katja,

please have a look at the attached file!
Magnitude plots and waterfall plots with our wav file processed (arrival to IR peak moved and windowing) made in Easera; the magnitude is now about right!!!

Probably Audacity applies some processing to the file when exported after processing (just cutting) that lead EASERA to misleading the IR characteristics.

You can see that our IR has its peak normalized to 0 dBFS, did you apply already all our thoughts on IR scaling??

Still need to find out how to import mic calibration in EASERA....

Bassik

http://www.pdpatchrepo.info/hurleur/RT.xls

Hi Bassik,

These comparisons look quite promising! And this is done with the old wiggly chirp. The low frequency ripples over which I'm cracking my brains hardly play a role because the system's response starts at around 100 Hz. I'm curious, what is actually the system under test?

The promising results are of course no reason to stop the quest for the ideal chirp. I have done some fade-in window tests, and indeed they can completely remove the low frequency ripple. But the fade-in window has to be crazy long before it has good result, stretching fully over the first octave. The ripples actually result from the chirp's sudden start. A regular smooth fade-in like the half Hann window I am using now does not effectively counteract this effect. So now I'm thinking of a critically damped oscillator, applied time-reversed as a window. Or maybe a windowed sinc.

About normalization I had some more thoughts. It occurred to me that the direct response will (in most cases) be spreaded over two neighbouring indexes. Ideally there should be just one index representing the indentity diagonal in a Toeplitz convolution matrix. With the two largest peaks in the center of an FFT frame, the group delay can be estimated over a region which has linear phase. Then, the complete IR phase spectrum should be rotated, to center the IR exactly on one index. After this, peak normalization of the IR should be done. Then, the maximum spectrum magnitude must be found. The multiplicative inverse of this maximum provides the normalization factor to prevent clipping. To do all these things, we really need that FFT-XXL object yet to write.

For acoustic measurements (RT, STI etc) only the casual part should be used. The purpose of the swept sine technique is to get all the reproduction system non-linearities in the anti casual part (leftmost part) of the deconvoluted IR in order to then eliminate them.

This is clear. So 'full IR' actually means the causal part for this purpose. I can't remember having heard about STI. Found this DIY article on the topic:

http://www.synaudcon.com/site/articles/a-do-it-yourselfers-guide-to-computing-the-speech-transmission-index/

Seems rather complicated. Should we build this in Pd too?

Katja

More on chirp windowing, good news.

Yesterday I reported about a Hann-fade-in window eating so much of the chirp, before producing good result. Actually that need not be a problem. The chirp can be extended with an extra octave, and this octave is then sacrificed to the fade-in. With 11 octaves in the chirp, and the first octave fully fade-in-windowed, a beautiful spectrum results, with a fairly steep transition band climbing up to precisely the start of the second octave (around 21 Hz at sr 44K1) and with only the faintest remainders of a ripple.

1/11 of the chirp is now spent on a fade-in. For example, almost 2 seconds of a 20 seconds chirp. So what? The spectral result is excellent, and the chirp length may be extended with a few seconds to compensate the lost time. I see no disadvantages so far.

Think I'll implement this long Hann-fade-in by default in my chirp generator.

Katja

So I have this object [expochirp~] finished and done some simple IR measurement on my MacBook internal soundcard (out >> in loopback). Earlier tests had already shown a substantial dip around Nyquist but that seems reasonable as the signal must be filtered to not generate aliases in the discrete domain.

But now the surprise, at least for me. A couple of days ago I uttered the assumption that the direct impulse response will be spreaded over two neighbouring samples. After all, when a pulse is received from acoustic or analog source, it would be very coincidental if it falls exactly on a single discrete index. From my first experiments however, I learned that the pulse is not spreaded over two neighbouring samples, but... over a whole bunch of samples. It is a sin(x)/x figure! (edit: sin(pi*x)/(pi*x) actually). Like spectral smearing in frequency domain. Should have known that....

Ok then, time zero lies somewhere inbetween two indexes. And because of the sin(x)/x smearing, there is no way of exactly taking the causal part of the IR. Everywhere you try to cut it near the highest coefficient, it gives a completely different spectrum. This can not be resolved until the impulse response is time-shifted by a fraction of a sample-interval, to position the direct response exactly on a single index. It's funny that Farina did not mention this issue in his texts on the ESS method.

Katja

p.s. I'll post [expochirp~] and patches soon, whenever I find opportunity to compile on different systems.

@katjav said:

Hi Bassik,

These comparisons look quite promising! And this is done with the old wiggly chirp. The low frequency ripples over which I'm cracking my brains hardly play a role because the system's response starts at around 100 Hz. I'm curious, what is actually the system under test?

Katja

Hello Katja,

the system under test is a big exhibition hall in london (http://www.vam.ac.uk/collections/fashion/galleries/40/index.html) measured using a genelec 8030A monitor and an omni directional mic only.

Keeping it simple to understand teh results quicker.

Regarding your exp chirp computation, I am afraid this goes above my current knowledge so I need to study a bit before giving you an answer or even a comment.

But we should bear in mind that exp Sweep sine technique is done for measuring acoustic spaces (hence the storing of non linearities of the reproduction and recording system in the anti casual part of the IR) and to obtain high quality IR for convolution and auralisation use.

In this case I don't really get the point why you are researching the exact sample value of the peak value.
In my experience, analysis of simple IR don't give substantial different results if the time zero of the IR is not exactly on the peak but probably I am missing something here expecially regarding the signal processing.

base of the exponentianl or logaritmic sweep is that the instantaneous frequency vary exponentially with time (see http://www.acoustics.net/objects/pdf/article_farina02.pdf)
Hope this can be of help.

Also Farina described a method for fast convolution for ambiophonics that can maybe be of inspiration for the deconvolution engine (see http://www.acoustics.net/objects/pdf/article_farina04.pdf)

Regarding STI, it is a measure of the speech intelligibility in a given space and takes into account many factors of which th emost important are the S/N ratio and the frequency masking.

Voice is a very complex mathematical beast as it has a characteristics frequency respone that vary with spl and also a frequency modulation mechanism that we use to articulate words. this is obvioulsy different for different languags but has been standardised in a certain way.

I will send to you the ISO standard that clearly explains th emethod of calculation impplemented in any commercial software in teh next days.

Also I would like to discuss possible development after the chirp issues are resolved.

I believe that it can be useful to implent the following in the patch:

• Sound editor to delete the anti casual part of the IR and move the arrival time of the peak to zero

• Room acoustics parameters calculation (I will send you the relative ISO and I will write more info later)

• Multichannel IR recording and processing including an Ambisonics module; IEM in Gratz (http://iem.at/Members/noisternig/bin_ambi) has done a great work on this. Our implementation would be much simpler and would include A-format to B-Format encoding and additional analysis.

Hope everything is well in life in the meantime

speak to you soon

Bassik

ExpoChirpToolbox.

This will be the name of an IR measurement toolchain, based on the [expochirp~] object. The first two patches are built now:

• [expochirp-generator], generates an exponential chirp of excelllent quality, together with inverse chirp. Deconvolution-test and inspection of the chirp's spectral characteristics within the patch. Chirp & inverse chirp can be stored with metadata in a textfile.

• [IRrecorder], loads a chirpfile as created with [expochirp-generator], applies the chirp to the system under test. Deconvolution, and inspection of spectral characteristics of the impulse response. The raw IR can be stored as 32 bit floating point .wav file.

Next component on my program will be [IReditor], where you can load a raw IR as produced with [IRrecorder] and edit it: trim the useful portion, and eventually window the tail if necessary. The editor should apply the IR as a filter in a soundfile player, so you can check the resulting sound. When these patches work satisfactory, we can think of additional tools: [IRinverter] to produce correction filters, [IRanalyser] for extracting acoustics parameters, [IRambisonic] for multichannel IR...

Attached is ExpoChirpToolbox.zip containing [expochirp~] binaries for OSX and Windows (I'll do Linux later) plus patches [expochirp-generator] and [IRrecorder]. Please post comments or suggestions.

edit: can't attach a file at the moment, download from here;

http://www.katjaas.nl/temp/ExpoChirpToolbox.zip

Katja

Thank you very much for this Katjav.

Much appreciated.

Hope that my post yesterday have been of some help.
I was just reading it today and it contains some incongruence.

Anyway, I will start some testing on the patch asap and let you know.

I will also prepare some more sketches for B-format encoding and more info for the analyzer.

Cheers

Bassik

@bassik said:

In this case I don't really get the point why you are researching the exact sample value of the peak value. In my experience, analysis of simple IR don't give substantial different results if the time zero of the IR is not exactly on the peak

Some more experimenting has learned me that it depends on the system under test. For an IR of the MacBook internal soundcard loopback, the problem is manifest. For an IR of my room, that is not the case. The difference is probably this: the soundcard's spectrum is near to perfect and there is only a direct response. The room test response missed a lot of the low spectrum end (also due to lousy 'test equipment'), and the direct response is buried in reflections because of the small room size. It's funny that the room response is much easier to trim than the loopback response.

Searching the internet I found the name for the phenomenon: fractional delay. And of course implementations for fractional delay, which can also be used to get rid of fractional delay. As far as I can see now, we'll only need it for cases where the direct response has a prominent role. For example, when we want to produce inverse filters for correction of spectral deficiencies in test equipment. We're not at that point yet.

Anyway, I did a patch with a model of fractional delay, which clearly illustrates the problem. (See attachment). And if you do a loopback test with [IRrecorder], you'll see what I mean as well.

edit: bah, again impossible to attach file. (Did I exceed my upload quotum?). It's here:

http://www.katjaas.nl/temp/fractional-delay.pd.zip

Katja

Patches IReditor and IRanalyser are added to the ExpoChirpToolbox.

The editor is a basic sample-editor for sample-precise trimming of the IR, windowing the tail(s), and normalising the maximum spectrum magnitude to 0 dBFS (full scale). The analyser presents spectral info of an edited IR, and you can test the IR as FIR filter with test signal, sound file or soundcard input. Later, the analyser will be extended to do parameter extraction, but since this will again require a lot of research, I want to share the patches now. Together with [expochirp-generator] and [IR-recorder], they already form a complete toolchain for generating your own FIR filters from acoustic environments.

Katja

ps .zip file attached this time (many thanks to the people who repaired this forum feature)

http://www.pdpatchrepo.info/hurleur/ExpoChirpToolbox.zip

Hi Guys,
I was referred to this thread by Serafino Di Rosario, and I will test this PD patch for performing ESS measurements. Everything is very interesting for me, and it seems that Katja did a very good job!

Regarding the problems encountered, I give you these infos:

1. sine-phase-matched sweep. This method is very useful when performing distortion measurements, or computing multiple-order IRS to be used in not-linear convolution processor (for emulating the nonlinearities of a device). For the method to work, it is mandatory that the sine sweep is sine-phased not only at the beginning, but also at the end of each octave. This way, each harmonic-order IR will be phase matched with the linear IR. The provided formulation solves this problem, and it is very good to see it explained here so simply.
   The importance of using a phase-synced exponential sweep was first discovered by Antonin Novak, a Ph.D. student of the universities of Prague and Le Mans.

2. The ripple at low frequencies can be controlled by proper fade-in. The choice of the "optimal" fade law is still a big subject under scientific discussion. Hann windowing is just a very initial, suboptimal approach. I plan to investigate further the choice of the optimal fade law, and publish something on this topic, soon.

3. The concept of cutting away everything before the arrival of the direct sound is wrong, in my opinion. The "silence" before the arrival of the direct sound has a very important physical meaning, it is the "time-of-flight" of the sound, and provides an accurate measurement of the distance between the source and the receiver. Furthermore, it contains "background noise", which is a very important quantity to know, for example when deriving STI from the IR measurement.
   So PLEASE, do not cut away this initial silence! If the IR has to be used as a filter for a convolution-based reverb plugin, the plugin must be intelligent enough for analyzing the IR, and giving the user the possibility to keep this initial silence or cut it away. For example, IR-1 from Waves gives these possibilities. In any case, a measured IR of a room should always contain the time-of-flight... Publishing "pre-cutted" IRs is wrong, and in the long run will cause a lot of troubles...

4. The "Fractional delay", for the same reasons, should NOT be corrected! If the time-of-flight is fractional, good, let's stay with this fact. As pointed out, cutting (time-shifting) improperly the measured IR can alter its spectrum. So, please, keep every measured IR as it comes out from the convolution with the inverse sweep... If the higher-order distortion products are not needed, it make sense to only keep the linear part only, but always starting from the true "zero time". Let's make an example: I generate a 20s-long IR at 48 kHz, that is 960,000 samples.
   The inverse sweep will also be 960,000 samples long.
   I play the sweep, and record the room response for, say, 1,200,000 samples, for being sure of capturing the complete reverberant tail even at the higher frequencies.
   Now I convolve the recorded signal (1,200,000 samples) with the inverse sweep (960,000 samples), and I get a convolved signal which is long 2,159,999 samples.
   If I want to keep a 4s long IR, containing only the linear response, you should throw away the first 959,999 samples, and keep the following 192,000 samples.
   As this signal starts from the true "zero time", the main peak will not be at the very beginning, but delayed of an amount corresponding to the source-receiver distance. If it was 10m, it will be 10/340=0.0294 seconds...

5. For performing efficiently convolution of very long filters (in the example above, the Inverse Sweep was nearly 1 million points) it is advisable to employ a partitioned convolution scheme. That is, the filter is splitted in a number of blocks, so that instead of performing a single, very long FFT , a number of shorter FFTs is performed instead. On my web site you find a couple of papers explaining the partitioned convolution algorithm. This is the same algorithm employed in the well-known BruteFIR open-source program by Anders Torger.

Bye!

Angelo Farina

Hello Angelo,

Welcome to this discussion. It is a joy to hear the opinion and advice of the foremost expert in this field. Our explorations are based on your work. I want to thank you for all the publications which you have so generously published online.

With the info in your post, a major point has now become clear to me: time zero in an IR refers to the start of the measurement, not to the arrival of the (first) response peak. I have been confused about this all the time, and happy to learn that it is actually so simple and logical.

As you may have guessed from the discussion, we learn from trial and error, just as much as from literature and education. I was 'initiated' in the IR measurement topic some years ago by the dsp teacher of the Sonology Department of the Conservatory of The Hague, Peter Pabon. Some of the mathematics were revealed, and I was fascinated. But now that we are trying to build our own Pd implementation in this forum cooperation, we run into the practical details. One of these was the undetermined phase of the chirp end in the regular mathematical formulation, causing excessive frequency ripple in many cases. Unaware of the work of Antonin Novak, I developed the earlier described solution, which now turns out to be a re-invention, like it happens so often.

Regarding the partioned convolution method: fortunately this is implemented in Pd, by Benjamin Saylor. We use this for the deconvolution of the long chirp responses (~ 20 seconds).

Hopefully you'll find time to try ExpoChirpToolbox and post your comments. It is a work in progress and we are working on further extensions.

Katja

Hello Katjav,

I have been chasing the swept sine technique with Farina's Aurora tools on the Cool Edit/Adobe Audition platform. I too found that octave extension of low frequency side only to use as fade in made possible much better result. Also small fade out at the end for non zero end state and all was greatly improved. But applied to getting IR of sound card loopback produces wrong answer. AC coupling of output and input stages produces min phase Butterworth type answer and swept sine with fades always contaminates with its sinc like IR. MLS gives correct answer, but of course brings its own baggage.

Your webpage: http://www.katjaas.nl/expochirp/expochirp.html shows IR result that I've seen on many of my screens while working on the fade issue. Very nicely done indeed.

Farina released plugins for Audacity for convolution and sweep generation. Audacity display doesn't provide zoom tools for display of fine detail. The plugins do provide ready reference.

The solution to the fade issue is to not use them in producing the sweep pair. The answer is in generating at least one of the sweeps in the frequency domain.

Here is link to a great paper covering some of this: "Transfer Function Measurement with Sweeps," SWEN MÜLLER and PAULO MASSARANI I found it here: http://www.melaudia.net/zdoc/comparisonMesure.PDF

Chapter 5 has got the goods.

Turns out that Kirkeby inverse is perfect tool for the job, and many more! To bad it is trapped in Cool Edit/Audition domain. Perhaps Farina will liberate it soon. The perfect place to apply it is before all the trouble begins. An arbitrary broadband signal of 2^n length may be inverted with Kirkeby, and convolution of the pair returns extremely good result. Samples directly next to pulse of 1e-6 easy to obtain. Starting with white noise, results emulating MLS are possible. Starting with exponential sweep a sweep is returned by Kirkeby with circular convolution properties: To get right convolution result, two copies of forward sweep need to be used and cropped to central result. Alternately Kirkeby may be directed to return result 4x longer for direct convolution. I found that a 2x long Kirkeby result produced artifacts. One may also use 2 copies of a 1x Kirkeby to recover IR with single instance of the seed wave.

I communicated this with Farina, and he honored me by placing his confirmations to: http://pcfarina.eng.unipr.it/Public/Wolkoff/

There are six wave files: 2^20 sample length forward sweep and inverse sweep generated without fades; and 1x,2x,4x Kirkeby inverses of the forward sweep. Additionally a time shift Kirkeby is posted. This is a hack I proposed as a partial work around for straight convolution. It is simply cutting first half of 1x Kirkeby and pasting it back at the end. Fading the ends controls bandwidth, and provides kernel for auralization applications.

Cheers!