Interesting. I've included that in the test patch. It exposes the same samplerate dependent timbre change. So far only the Sinc function solves the issue.

On 25.04.21 23:15, Sebastian Shader wrote:

I have a [vdhs~] object in my library that does hermite spline interpolation for a delay line (like tabread4c~). https://github.com/sebshader/shadylib (also on deken)

I had to re-make delread~ as delread4c~ because delread~s delay lines are not exposed in the .h files

-seb

-----Original Message----- From: Max abonnements@revolwear.com To: pd-list@lists.iem.at pd-list@lists.iem.at Sent: Sun, Apr 25, 2021 5:45 am Subject: [PD] 4-point interpolation changes timbre depending on sample rate

Hi list,

the 4-point interpolation in tabread4~ has been a popular topic in the past, going back to at least 2008. [1]

A similar issue is in delread4~. In fact a simple resonator changes its timbre quite drastically by just changing the sample rate of the audio interface. Attached is a test patch.

The issue becomes audible when choosing an odd delay time and compare the result between an odd and an even sample rate (e.g. 44.1k vs. 48k).

This is not good. In fact this is a serious defect. Imagine you want to market a product like a synth plugin (based on libpd) which sounds different depending on if the daw is opened in 44.1 or 48 kHz.

Cyrille Henry has coded tabread4c~ with a four-point cubic interpolation in the nusmuk library, but there is no delread4c~ equivalent in nusmuk.

Clemens Wegener has coded delreadsinc~ which implements a Whittaker–Shannon interpolation (Sinc function). This implementation does sounds consistent in any sample rate. It also sounds much better at very slow speeds inside a pitch shifter where delread4~ produces serious artefacts. The Sinc function requires a larger padding for the interpolation.

There are a couple of questions regarding on how to contribute this to Pd. Vanilla objects are currently:

delwrite~ = the sink delread~ = control rate tap delread4~ (vd~) = audio rate tap with a four-point interpolation

1. the code in delwrite~ isn't agnostic towards the interpolation since

it already provides the padding for the 4 point interpolation

2. if we add more interpolation methods to delread / tabread, the

cleanest way would be to just have one tap object and the interpolation can be changed by a message and argument to it. currently there are implementations for the miller 4pt, cubic 4pt and Sinc.

Unfortunately delread4 already carries the 4 from 4-point interpolation in its name, so probably it would be the best to deprecate that and find a new name like varidelay~ or so.

[1] Review of tabread4~ threads in the archive https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html

[2] https://github.com/chairaudio/pure-data/tree/feature/delreadsinc https://github.com/chairaudio/pure-data/tree/feature/delreadsinc

---

Pd-list@lists.iem.at mailto:Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list https://lists.puredata.info/listinfo/pd-list

I read through the code in your d_delay.c for the delread_sinc implementation. It's a 22-point interpolator--very high quality and it looks like a good implementation

The Hermite spline is the best 4-point cubic interpolator, but the point of those polynomial interpolators is just to be cheap, and good enough for regular usage. It doesn't require a table, so it also doesn't need as much cache space

I don't think it's *the* optimal 4-point interpolator. It's just the best cubic polynomial 4-point interpolator. You could re-use your sinc interpolator code for other shorter interpolation kernels, ones that are better than cubics but more expensive to calculate.

At least for me, it's still an open mathematical problem that I've come back to, from time to time. I started to wonder if I could expand the error terms in the truncated sinc spectrum in the form of 1/w, 1/w^2, 1/w^3, .... and then cancel them out by polynomial correctors. An alternative to windowing the sinc function, this could be a more precise way to set the rate of attenuation in the stop-band, while still retaining all the smoothness of the sinc shape. I get a break from teaching in about 3 weeks. I'm looking at the math right now, and I don't think I can solve it quickly (or maybe at all).

How short does your sinc interpolator with Blackman-Harris window need to be when you start to get the timbre change effect?

Best, Chuck

On Sun, Apr 25, 2021 at 7:18 PM Max abonnements@revolwear.com wrote:

Interesting. I've included that in the test patch. It exposes the same samplerate dependent timbre change. So far only the Sinc function solves the issue.

On 25.04.21 23:15, Sebastian Shader wrote:

...

I have a [vdhs~] object in my library that does hermite spline interpolation for a delay line (like tabread4c~). https://github.com/sebshader/shadylib (also on deken)

I had to re-make delread~ as delread4c~ because delread~s delay lines are not exposed in the .h files

-seb

-----Original Message----- From: Max abonnements@revolwear.com To: pd-list@lists.iem.at pd-list@lists.iem.at Sent: Sun, Apr 25, 2021 5:45 am Subject: [PD] 4-point interpolation changes timbre depending on sample rate

Hi list,

the 4-point interpolation in tabread4~ has been a popular topic in the past, going back to at least 2008. [1]

A similar issue is in delread4~. In fact a simple resonator changes its timbre quite drastically by just changing the sample rate of the audio interface. Attached is a test patch.

The issue becomes audible when choosing an odd delay time and compare the result between an odd and an even sample rate (e.g. 44.1k vs. 48k).

This is not good. In fact this is a serious defect. Imagine you want to market a product like a synth plugin (based on libpd) which sounds different depending on if the daw is opened in 44.1 or 48 kHz.

Cyrille Henry has coded tabread4c~ with a four-point cubic interpolation in the nusmuk library, but there is no delread4c~ equivalent in nusmuk.

Clemens Wegener has coded delreadsinc~ which implements a Whittaker–Shannon interpolation (Sinc function). This implementation does sounds consistent in any sample rate. It also sounds much better at very slow speeds inside a pitch shifter where delread4~ produces serious artefacts. The Sinc function requires a larger padding for the interpolation.

There are a couple of questions regarding on how to contribute this to Pd. Vanilla objects are currently:

delwrite~ = the sink delread~ = control rate tap delread4~ (vd~) = audio rate tap with a four-point interpolation

1. the code in delwrite~ isn't agnostic towards the interpolation since

it already provides the padding for the 4 point interpolation

2. if we add more interpolation methods to delread / tabread, the

cleanest way would be to just have one tap object and the interpolation can be changed by a message and argument to it. currently there are implementations for the miller 4pt, cubic 4pt and Sinc.

Unfortunately delread4 already carries the 4 from 4-point interpolation in its name, so probably it would be the best to deprecate that and find a new name like varidelay~ or so.

[1] Review of tabread4~ threads in the archive https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html

[2] https://github.com/chairaudio/pure-data/tree/feature/delreadsinc https://github.com/chairaudio/pure-data/tree/feature/delreadsinc

---

Pd-list@lists.iem.at mailto:Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

On Wed, Apr 28, 2021 at 11:14 AM Miller Puckette msp@ucsd.edu wrote:

On Wed, Apr 28, 2021 at 10:56:58AM -0500, Charles Z Henry wrote: My 2 cents...

The 4-point interpolation scheme gets radically better if the signal it's used on is oversampled (error goes down asymptotically by 24 dB for each doubling of sample rate) - so my own strategy is simply to 4x upsample everything I send through tabread4~ or delread4~. This moves the "problem" to that od designing an upsampling filter, which is much easier than a general interpolator.

I found a related problem this fall, teaching a University Physics I class that got thrown my way at the last minute.... and I decided to make a couple lessons about numerically integrating systems of differential equations

The first one (Euler's method) worked just fine, but later on, I wanted to show we could move past the 1st order derivative approximation and get better results than just upsampling. Then, my simulations started blowing up during class and I realized I don't understand implicit methods as well as I thought I did. Whoops

So, the related problem is what is the best truncated differentiation kernel on [-a, 1]? Once sampled, you'd get the coefficients of a numerical derivative scheme that can be re-arranged into an implicit method.

Either case, it goes back to the spectrum of the kernel. The derivative approximation used in Euler's method is only good over a small range of frequencies near zero, so the signal has to be upsampled in order to produce good results. By making that kernel longer, you can get higher frequencies and a spectrum that more closely matches 2*pi*i*f (although I haven't found a best scheme for it either). I think this question is different from the one Runge-Kutta methods answer. I think there's something here to find that's relevant

I think this is a bit too application-specific to put in Pd vanilla... but anyway, I don't understand why you can't just run the entire waveguide at 4x sample rate - am I missing something?

cheers M

On Thu, Apr 29, 2021 at 02:06:11PM +0200, Max wrote:

Hi Miller,

I think the oversampling workaround is great for cases like slow playback of a sound file, but I'm not sure how this would work in a waveguide, ans particularly our situation where we have different timbre depending on if the number of samples in the delay line is odd or even. Since we vary the delay time in ms it's kind of a gamble and inconsistent depending on samplerate.

Do you think delreadsinc~ could have a chance to be merged into vanilla and if so how?

I've outlined some ideas here: https://urldefense.com/v3/__https://github.com/pure-data/pure-data/issues/13...

Max

On 28.04.21 18:14, Miller Puckette via Pd-list wrote:

...

On Wed, Apr 28, 2021 at 10:56:58AM -0500, Charles Z Henry wrote: My 2 cents...

The 4-point interpolation scheme gets radically better if the signal it's used on is oversampled (error goes down asymptotically by 24 dB for each doubling of sample rate) - so my own strategy is simply to 4x upsample everything I send through tabread4~ or delread4~. This moves the "problem" to that of designing an upsampling filter, which is much easier than a general interpolator.

cheers Miller

...

I read through the code in your d_delay.c for the delread_sinc implementation. It's a 22-point interpolator--very high quality and it looks like a good implementation

The Hermite spline is the best 4-point cubic interpolator, but the point of those polynomial interpolators is just to be cheap, and good enough for regular usage. It doesn't require a table, so it also doesn't need as much cache space

I don't think it's *the* optimal 4-point interpolator. It's just the best cubic polynomial 4-point interpolator. You could re-use your sinc interpolator code for other shorter interpolation kernels, ones that are better than cubics but more expensive to calculate.

At least for me, it's still an open mathematical problem that I've come back to, from time to time. I started to wonder if I could expand the error terms in the truncated sinc spectrum in the form of 1/w, 1/w^2, 1/w^3, .... and then cancel them out by polynomial correctors. An alternative to windowing the sinc function, this could be a more precise way to set the rate of attenuation in the stop-band, while still retaining all the smoothness of the sinc shape. I get a break from teaching in about 3 weeks. I'm looking at the math right now, and I don't think I can solve it quickly (or maybe at all).

How short does your sinc interpolator with Blackman-Harris window need to be when you start to get the timbre change effect?

Best, Chuck

On Sun, Apr 25, 2021 at 7:18 PM Max abonnements@revolwear.com wrote:

...

Interesting. I've included that in the test patch. It exposes the same samplerate dependent timbre change. So far only the Sinc function solves the issue.

On 25.04.21 23:15, Sebastian Shader wrote:

...

I have a [vdhs~] object in my library that does hermite spline interpolation for a delay line (like tabread4c~). https://urldefense.com/v3/__https://github.com/sebshader/shadylib__;!!Mih3wA... (also on deken)

I had to re-make delread~ as delread4c~ because delread~s delay lines are not exposed in the .h files

-seb

-----Original Message----- From: Max abonnements@revolwear.com To: pd-list@lists.iem.at pd-list@lists.iem.at Sent: Sun, Apr 25, 2021 5:45 am Subject: [PD] 4-point interpolation changes timbre depending on sample rate

Hi list,

the 4-point interpolation in tabread4~ has been a popular topic in the past, going back to at least 2008. [1]

A similar issue is in delread4~. In fact a simple resonator changes its timbre quite drastically by just changing the sample rate of the audio interface. Attached is a test patch.

The issue becomes audible when choosing an odd delay time and compare the result between an odd and an even sample rate (e.g. 44.1k vs. 48k).

This is not good. In fact this is a serious defect. Imagine you want to market a product like a synth plugin (based on libpd) which sounds different depending on if the daw is opened in 44.1 or 48 kHz.

Cyrille Henry has coded tabread4c~ with a four-point cubic interpolation in the nusmuk library, but there is no delread4c~ equivalent in nusmuk.

Clemens Wegener has coded delreadsinc~ which implements a Whittaker???Shannon interpolation (Sinc function). This implementation does sounds consistent in any sample rate. It also sounds much better at very slow speeds inside a pitch shifter where delread4~ produces serious artefacts. The Sinc function requires a larger padding for the interpolation.

There are a couple of questions regarding on how to contribute this to Pd. Vanilla objects are currently:

delwrite~ = the sink delread~ = control rate tap delread4~ (vd~) = audio rate tap with a four-point interpolation

1. the code in delwrite~ isn't agnostic towards the interpolation since

it already provides the padding for the 4 point interpolation

2. if we add more interpolation methods to delread / tabread, the

cleanest way would be to just have one tap object and the interpolation can be changed by a message and argument to it. currently there are implementations for the miller 4pt, cubic 4pt and Sinc.

Unfortunately delread4 already carries the 4 from 4-point interpolation in its name, so probably it would be the best to deprecate that and find a new name like varidelay~ or so.

[1] Review of tabread4~ threads in the archive https://urldefense.com/v3/__https://lists.puredata.info/pipermail/pd-list/20... <https://urldefense.com/v3/__https://lists.puredata.info/pipermail/pd-list/20... >

[2] https://urldefense.com/v3/__https://github.com/chairaudio/pure-data/tree/fea... <https://urldefense.com/v3/__https://github.com/chairaudio/pure-data/tree/fea... >

---

Pd-list@lists.iem.at mailto:Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!... <https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!... >

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!...

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!...

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!...

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://urldefense.com/v3/__https://lists.puredata.info/listinfo/pd-list__;!...

Hi Chuck,

Clemens isn't on the list (yet), he sent me some measurements he made in response to your post.

All at 48kHz.

https://cloud.yuca.biz/s/WZseAmBQm9z6NdL

he says he can't hear a difference with 10 samples cutoff@pi and 14 samples cutoff@0.85*pi It sounds better than than delread4~ from 6 samples cutoff@pi upwards.

greetings, Max & Clemens

On 28.04.21 17:56, Charles Z Henry wrote:

I read through the code in your d_delay.c for the delread_sinc implementation. It's a 22-point interpolator--very high quality and it looks like a good implementation

The Hermite spline is the best 4-point cubic interpolator, but the point of those polynomial interpolators is just to be cheap, and good enough for regular usage. It doesn't require a table, so it also doesn't need as much cache space

I don't think it's *the* optimal 4-point interpolator. It's just the best cubic polynomial 4-point interpolator. You could re-use your sinc interpolator code for other shorter interpolation kernels, ones that are better than cubics but more expensive to calculate.

At least for me, it's still an open mathematical problem that I've come back to, from time to time. I started to wonder if I could expand the error terms in the truncated sinc spectrum in the form of 1/w, 1/w^2, 1/w^3, .... and then cancel them out by polynomial correctors. An alternative to windowing the sinc function, this could be a more precise way to set the rate of attenuation in the stop-band, while still retaining all the smoothness of the sinc shape. I get a break from teaching in about 3 weeks. I'm looking at the math right now, and I don't think I can solve it quickly (or maybe at all).

How short does your sinc interpolator with Blackman-Harris window need to be when you start to get the timbre change effect?

Best, Chuck

On Sun, Apr 25, 2021 at 7:18 PM Max abonnements@revolwear.com wrote:

...

Interesting. I've included that in the test patch. It exposes the same samplerate dependent timbre change. So far only the Sinc function solves the issue.

On 25.04.21 23:15, Sebastian Shader wrote:

...

I have a [vdhs~] object in my library that does hermite spline interpolation for a delay line (like tabread4c~). https://github.com/sebshader/shadylib (also on deken)

I had to re-make delread~ as delread4c~ because delread~s delay lines are not exposed in the .h files

-seb

-----Original Message----- From: Max abonnements@revolwear.com To: pd-list@lists.iem.at pd-list@lists.iem.at Sent: Sun, Apr 25, 2021 5:45 am Subject: [PD] 4-point interpolation changes timbre depending on sample rate

Hi list,

the 4-point interpolation in tabread4~ has been a popular topic in the past, going back to at least 2008. [1]

A similar issue is in delread4~. In fact a simple resonator changes its timbre quite drastically by just changing the sample rate of the audio interface. Attached is a test patch.

The issue becomes audible when choosing an odd delay time and compare the result between an odd and an even sample rate (e.g. 44.1k vs. 48k).

This is not good. In fact this is a serious defect. Imagine you want to market a product like a synth plugin (based on libpd) which sounds different depending on if the daw is opened in 44.1 or 48 kHz.

Cyrille Henry has coded tabread4c~ with a four-point cubic interpolation in the nusmuk library, but there is no delread4c~ equivalent in nusmuk.

Clemens Wegener has coded delreadsinc~ which implements a Whittaker–Shannon interpolation (Sinc function). This implementation does sounds consistent in any sample rate. It also sounds much better at very slow speeds inside a pitch shifter where delread4~ produces serious artefacts. The Sinc function requires a larger padding for the interpolation.

There are a couple of questions regarding on how to contribute this to Pd. Vanilla objects are currently:

delwrite~ = the sink delread~ = control rate tap delread4~ (vd~) = audio rate tap with a four-point interpolation

1. the code in delwrite~ isn't agnostic towards the interpolation since

it already provides the padding for the 4 point interpolation

2. if we add more interpolation methods to delread / tabread, the

cleanest way would be to just have one tap object and the interpolation can be changed by a message and argument to it. currently there are implementations for the miller 4pt, cubic 4pt and Sinc.

Unfortunately delread4 already carries the 4 from 4-point interpolation in its name, so probably it would be the best to deprecate that and find a new name like varidelay~ or so.

[1] Review of tabread4~ threads in the archive https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html

[2] https://github.com/chairaudio/pure-data/tree/feature/delreadsinc https://github.com/chairaudio/pure-data/tree/feature/delreadsinc

---

Pd-list@lists.iem.at mailto:Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

@Miller, I just saw your message - we will try that asap!

Hello Chuck,

Good point! Thanks for these thoughts, there is no good reason to duplicate the table with each object, apart from initially being too lazy to code it correctly.

I will have a look at the code again and see which solution I can arrive at! Method 3. sounds a bit like the unique pointer idea in modern C++…

There is another thing that I struggle to understand: (I already posed this question on Github)

Should the sinc interpolation change its filter characteristics when you change the delay time? If that is the case, the code would need a little bit of adaption.

From what I understand, increasing the delay time is like downsampling (is it?). Now, I wonder how exactly delay time and sample rate are related mathematically. E.g. if I took the derivative of the delay time, would I get something that is proportional to the change in sample rate?

I read about the filter curve adaption during downsampling here: https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf ...but how would that apply to delay time interpolation?

Also: Is it really mandatory to choose a filter cutoff below pi, when the input and output Sample rates match? For a simple static delay that should be the case.

Best wishes, Clemens

Am 2021-04-29 um 17:55 schrieb Charles Z Henry czhenry@gmail.com:

It's sort of disappointing that it has to be so expensive, but you do what you have to to get the quality and consistency of the sound. The question of whether there is an optimal interpolating kernel (which I'd like to work on) isn't really going to help you in the near term. It's enough just to find a solution that works well and performs well.

You've got a scalability issue when using multiple instances of delreadsinc~: Each new instance recalculates the table when created. The table is stored in the struct _sigvdsinc so it is also duplicated for each instance beyond the first. So, I propose a few other choices:

1. Calculate the table during setup and add a class method that

returns a pointer to that table, so that a new object doesn't calculate the table. It just points at the existing table. This option would use memory, even if the object is never used, and would not free the table under any circumstances (not great, but easy) 2. Add a class method that returns NULL if the table has not been created or a pointer to the existing table and a class method to set the pointer. The first instance of the object used would run the initialize_sinc_table function, but each additional instance would not. The table hangs around for future use, even if the instances are deleted. 3. Add the class methods in 2. plus methods to keep a count of the number of instances. Then, when the last instance is deleted, the table is freed and the pointer re-assigned to NULL.

By my approximation, that table uses about (4*4096*10) 160 kB (or 320 kB in double precision). I'm not real clear on the differences in stack, static array allocation, and dynamic allocation in behavior/performance. I know it's there, just can't seem to digest it when I look at it. You might have had some reasons in mind to allocate the array the way you did--I'd like to learn more about this

Best, Chuck

On Wed, Apr 28, 2021 at 2:27 PM Max abonnements@revolwear.com wrote:

...

Hi Chuck,

Clemens isn't on the list (yet), he sent me some measurements he made in response to your post.

All at 48kHz.

https://cloud.yuca.biz/s/WZseAmBQm9z6NdL

he says he can't hear a difference with 10 samples cutoff@pi and 14 samples cutoff@0.85*pi It sounds better than than delread4~ from 6 samples cutoff@pi upwards.

greetings, Max & Clemens

On 28.04.21 17:56, Charles Z Henry wrote:

...

I read through the code in your d_delay.c for the delread_sinc implementation. It's a 22-point interpolator--very high quality and it looks like a good implementation

The Hermite spline is the best 4-point cubic interpolator, but the point of those polynomial interpolators is just to be cheap, and good enough for regular usage. It doesn't require a table, so it also doesn't need as much cache space

I don't think it's *the* optimal 4-point interpolator. It's just the best cubic polynomial 4-point interpolator. You could re-use your sinc interpolator code for other shorter interpolation kernels, ones that are better than cubics but more expensive to calculate.

At least for me, it's still an open mathematical problem that I've come back to, from time to time. I started to wonder if I could expand the error terms in the truncated sinc spectrum in the form of 1/w, 1/w^2, 1/w^3, .... and then cancel them out by polynomial correctors. An alternative to windowing the sinc function, this could be a more precise way to set the rate of attenuation in the stop-band, while still retaining all the smoothness of the sinc shape. I get a break from teaching in about 3 weeks. I'm looking at the math right now, and I don't think I can solve it quickly (or maybe at all).

How short does your sinc interpolator with Blackman-Harris window need to be when you start to get the timbre change effect?

Best, Chuck

On Sun, Apr 25, 2021 at 7:18 PM Max abonnements@revolwear.com wrote:

...

Interesting. I've included that in the test patch. It exposes the same samplerate dependent timbre change. So far only the Sinc function solves the issue.

On 25.04.21 23:15, Sebastian Shader wrote:

...

I have a [vdhs~] object in my library that does hermite spline interpolation for a delay line (like tabread4c~). https://github.com/sebshader/shadylib (also on deken)

I had to re-make delread~ as delread4c~ because delread~s delay lines are not exposed in the .h files

-seb

-----Original Message----- From: Max abonnements@revolwear.com To: pd-list@lists.iem.at pd-list@lists.iem.at Sent: Sun, Apr 25, 2021 5:45 am Subject: [PD] 4-point interpolation changes timbre depending on sample rate

Hi list,

the 4-point interpolation in tabread4~ has been a popular topic in the past, going back to at least 2008. [1]

A similar issue is in delread4~. In fact a simple resonator changes its timbre quite drastically by just changing the sample rate of the audio interface. Attached is a test patch.

The issue becomes audible when choosing an odd delay time and compare the result between an odd and an even sample rate (e.g. 44.1k vs. 48k).

This is not good. In fact this is a serious defect. Imagine you want to market a product like a synth plugin (based on libpd) which sounds different depending on if the daw is opened in 44.1 or 48 kHz.

Cyrille Henry has coded tabread4c~ with a four-point cubic interpolation in the nusmuk library, but there is no delread4c~ equivalent in nusmuk.

Clemens Wegener has coded delreadsinc~ which implements a Whittaker–Shannon interpolation (Sinc function). This implementation does sounds consistent in any sample rate. It also sounds much better at very slow speeds inside a pitch shifter where delread4~ produces serious artefacts. The Sinc function requires a larger padding for the interpolation.

There are a couple of questions regarding on how to contribute this to Pd. Vanilla objects are currently:

delwrite~ = the sink delread~ = control rate tap delread4~ (vd~) = audio rate tap with a four-point interpolation

1. the code in delwrite~ isn't agnostic towards the interpolation since

it already provides the padding for the 4 point interpolation

2. if we add more interpolation methods to delread / tabread, the

cleanest way would be to just have one tap object and the interpolation can be changed by a message and argument to it. currently there are implementations for the miller 4pt, cubic 4pt and Sinc.

Unfortunately delread4 already carries the 4 from 4-point interpolation in its name, so probably it would be the best to deprecate that and find a new name like varidelay~ or so.

[1] Review of tabread4~ threads in the archive https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html https://lists.puredata.info/pipermail/pd-list/2019-06/125391.html

[2] https://github.com/chairaudio/pure-data/tree/feature/delreadsinc https://github.com/chairaudio/pure-data/tree/feature/delreadsinc

---

Pd-list@lists.iem.at mailto:Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

---

Pd-list@lists.iem.at mailing list UNSUBSCRIBE and account-management -> https://lists.puredata.info/listinfo/pd-list

When there is no aliasing to worry about, i might set the cutoff to pi again... On low sample rates (22kHz), the lower cutoff is definitely noticeable.

Btw. I just implemented sharing of the interpolation table of the delreadsinc~ object according to your suggestions. It counts the number of references and frees the pointer when no object uses it anymore.

Best wishes, Clemens

Am 29.04.21 um 21:03 schrieb Charles Z Henry:

On Thu, Apr 29, 2021 at 12:08 PM Clemens Wegener clemens@chair.audio wrote:

...

Also: Is it really mandatory to choose a filter cutoff below pi, when the input and output Sample rates match? For a simple static delay that should be the case.

There is an exact requirement when what you're making is an interpolator: it has to pass through all the original points! For that to happen, your interpolating kernel i(t) must have i(t)=1, t=0 i(t)=0, for non-zero integer values of t

An interpolator will have a cutoff freq at exactly pi rad/sec

So, if you're changing the cutoff frequency or amplitudes, it won't be an interpolator. It will just be some kind of table reader. And that's good for performance, I'd say, but your users will have to know that (and you'll have to justify the choice if it's for an academic purpose)

I think now would be a good place to pause the implementation and discuss if there is a need for this kind of algorithm in the PD community.

We saw some use cases for the Whittaker-Shannon interpolation where we gain in quality and/or speed. Namely waveguides and pitch shifters.

Is there anything else, where we would like to use this interpolation kernel? Like in general resampling? In the tabread~ for sample playback? Are there really quality or speed issues that we could solve there?

For our use case the code I submitted is good enough. But I would be happy to spend more time optimizing if there is a need / a broad use of that algorithm. In that case we need another restructuring and some help from somebody who is very proficient in writing pd source code. :)

Chuck, I commented your last message below!

Am 2021-05-03 um 03:52 schrieb Charles Z Henry czhenry@gmail.com:

On Sun, May 2, 2021 at 4:28 PM Clemens clemens@chair.audio wrote:

...

When there is no aliasing to worry about, i might set the cutoff to pi again... On low sample rates (22kHz), the lower cutoff is definitely noticeable.

I like this thesis you posted earlier https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf

f(at) <-> 1/|a|*F(w/a) This (and two formulas that follow) is listed as Smith's algorithm. I actually got to speak with Julius Smith at the 2012 LAC about this formula. I asked, can't we do any better in terms than O(a*n) number of computations for a>1? He said nope!, but I still have some questions there. Playback with speeds less than 1 always use O(n), rather than O(a*n). I wrote an anti-aliasing external tabread4a~ that implements this formula per sample and works pretty well, except it becomes expensive when you transpose a few octaves up.

I guess that’s the code you are referring to (?): https://lists.puredata.info/pipermail/pd-list/2007-03/048397.html https://lists.puredata.info/pipermail/pd-list/2007-03/048397.html

Is it using a lookup of the sync table at certain fixed points? If so, would it be „compatible“ in the sense that these points are part of the lookup table in delreadsinc~?

For a variable delay line (like vd~), the paper's contents are more relevant and maybe you should consider writing a vdsinc~ object next (once you've optimized this one). delread is more like a fixed, stationary delay line, and I think it's better to default as a literal interpolator (LP_SCALE = 1).

Ok. That’s what I thought as well. Still I would need to read a bit more
about the subject to understand how the variable delay case relates to down-/upsampling. In tabread4a~ you derive the sampling factor from the difference in delay time, right?

...

Btw. I just implemented sharing of the interpolation table of the delreadsinc~ object according to your suggestions. It counts the number of references and frees the pointer when no object uses it anymore.

I read through the changes. For those reading along see the rest at line 348 of d_delay.c at https://github.com/chairaudio/pure-data/blob/feature/delreadsinc/src/d_delay... (in progress, current code re: the tables pasted below).

1. I wonder is this global table properly scoped for Pd. I

seriously don't know. But it could be used for multiple objects, vdsinc~, tabreadsinc~, etc.... So you ought to think about how it could be re-used between classes and choose a scalable approach now.

Good question. :) I’m not sure if Miller has suggestions or examples how global variables are handled in PD.

2. I think the declared LP_SCALE variable is a bad approach as

currently implemented. You can't change it at runtime. It can be used in the perform routine instead, with an additional argument to set the cutoff frequency per object. This is also good for your testing, as you'll be able to put objects with different LP_SCALE values side-by-side for comparison.

That’s right, LP_SCALE should be tunable in real-time. I think that would need only little change to the code.

3. I'm convinced this is the right approach for making an optimal

interpolator. Adding more length pays off over upsampling, because you evaluate the convolution at a fewer number of points. Your implementation needs a lot of tuning---you can squeeze out more performance. Focus on the inner most loops. Also, you'll have to compare some arm vs intel/amd platforms at some point. I think you're closer to the beginning than the end. This is maybe not what you wanted to research in the first place

When it gets down to processor specific optimization, I’m not sure how much the compiler does already and if it’s really worth the effort. I would investigate this if we find that we want to put this code into pd.

Best wishes, Clemens

... typedef struct _sigvdsinc { t_object x_obj; t_symbol *x_sym; t_float x_sr; /* samples per msec */ int x_zerodel; /* 0 or vecsize depending on read/write order */ t_float x_f; } t_sigvdsinc;

typedef struct _sigvdsinc_sharedtable { // sinc table is held in global variable t_sample sinc_array[SINC_LEN]; // derivative of sinc funtion for interpolation of sinc function table t_sample sinc_diff_array[SINC_LEN]; // reference count for shared pointer int ref_count; } t_sigvdsinc_sharedtable;

t_sigvdsinc_sharedtable *gSharedTable = NULL;

void sigvdsinc_initialize_sinc_table(t_sample* sinc_array, t_sample* sinc_diff_array) { sinc_array[0] = LP_SCALE; sinc_diff_array[SINC_LEN] = 0; float a0 = 0.35875; float a1 = 0.48829; float a2 = 0.14128; float a3 = 0.01168; for (int i=1; i<SINC_LEN; i++) { float idx = 0.5*(M_PI*i/SINC_LEN - M_PI); // four term blackmanharris after https://www.mathworks.com/help/signal/ref/blackmanharris.html float blackmanharris= a0 - a1*cos(2*idx) + a2*cos(4*idx) - a3*cos(6*idx); // blackmanharris windowed sinc function sinc_array[i] = sin(LP_SCALE*M_PI*(float)i/STEPS_ZC)/(LP_SCALE*M_PI*(float)i/STEPS_ZC)

• blackmanharris * LP_SCALE;

// calculate derivative for table interpolation sinc_diff_array[i-1] = sinc_array[i]-sinc_array[i-1]; } }

... void sigvdsinc_free(t_sigvdsinc *x) { // delete shared interpolation table if need be if(gSharedTable != NULL){ if(gSharedTable->ref_count==1){ free((void*)(gSharedTable)); gSharedTable = NULL;; } else{ gSharedTable->ref_count--; } } }

Am 2021-05-03 um 03:52 schrieb Charles Z Henry czhenry@gmail.com:

On Sun, May 2, 2021 at 4:28 PM Clemens clemens@chair.audio wrote:

When there is no aliasing to worry about, i might set the cutoff to pi again... On low sample rates (22kHz), the lower cutoff is definitely noticeable.

I like this thesis you posted earlier https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf

f(at) <-> 1/|a|*F(w/a) This (and two formulas that follow) is listed as Smith's algorithm. I actually got to speak with Julius Smith at the 2012 LAC about this formula. I asked, can't we do any better in terms than O(a*n) number of computations for a>1? He said nope!, but I still have some questions there. Playback with speeds less than 1 always use O(n), rather than O(a*n). I wrote an anti-aliasing external tabread4a~ that implements this formula per sample and works pretty well, except it becomes expensive when you transpose a few octaves up.

I guess that’s the code you are referring to (?): https://lists.puredata.info/pipermail/pd-list/2007-03/048397.html

Is it using a lookup of the sync table at certain fixed points? If so, would it be „compatible“ in the sense that these points are part of the lookup table in delreadsinc~?

It was more a test of this basic formula, rather than anything else. It's the tabread4~ polynomial, not a sinc. So I could do a direct comparison to tabread4~, and show how smoothly it makes the transition from low speeds to high speeds. Any choice of the interp(x) function yields the same basic trade-off: when you make your interpolating function longer and flatter by factor a, proportional to playback speeds greater than 1, the cutoff frequency goes to +pi/a rad/s. Then, resampling at speed a is an evaluation of the convolution with (1/a*interp(t/a)) with points spaced a apart and it stretches out that spectrum back to +pi rad/s

The test patch was an obvious case, sweep playback of a 256 point table with loud partials at f, 2f, 4f, 8f, and you'd see the partials in an fft graph march right up to the Nyquist frequency and fall off the edge. At 8 octaves transposition, this object was doing a 1024-point interpolation (table was wrapped). At that point, it was completely silent and using 50% cpu doing very little else. Meanwhile, you'd hear the aliasing effect from tabread4~ and all it's partials still bouncing around in the passband. So, the formula scales well. With a little bit of headroom from your LP_SCALE parameter, there would be some improvement also.

float interp(float x) { float absx=fabsf(x); return ((absx<2.0f)*((absx<1.0f)?(1-absx*(0.5f+absx*(1-0.5f*absx))):(1-absx*(1.833333f-absx*(1.0f-0.1666666f*absx))))); } ...... and the innermost loop that uses interp(x) just uses an extra term "diff" for difference between last two indexes ... while((a++)<b) sum_left+=(*(fp++))*interp((findex-a)/diff); ... (there's also a sum_right and they are added together to form the output per sample) ... *out++=(sum_left+sum_right)/diff;

Your implementation is better as a starting point anyhow. I've suspected that tabulating the interpolation function would yield better performance.