Lagrange Interpolation Frequency Response Examples

Figure J.2: Lagrange Interpolator Frequency Responses: Order 4 -- $\Delta$ = 1.1 to 1.5 in steps of 0.1.

Figure J.3: Lagrange Interpolator Frequency Responses: Order 3 -- $\Delta$ = 1.1 to 1.5 in steps of 0.1.

Figure J.4: Lagrange Interpolator Frequency Responses: Order 2 -- $\Delta$ = 1.1 to 1.5 in steps of 0.1.

Figure J.5: Lagrange Interpolator Frequency Responses: Order 1 -- $\Delta$ = 1.1 to 1.5 in steps of 0.1.

Figure J.6: Lagrange Interpolator Frequency Responses: Orders 1, 2, and 3 -- $\Delta = 1.4$

Figure J.7: Faust program tlagrange.dsp used to generate Figures J.2 through J.5.

 

// tlagrange.dsp - test lagrange.lib

import("lagrange.lib");

N = 16;

// Change "processX" to "process" in one case below:

// Aggregate test:
processA = 1-1' <: (fdelay1(N,5.5),
                               fdelay2(N,5.5),
                               fdelay3(N,5.5),
                               fdelay4(N,5.5));
// To see results:
// [in a shell]:
//   make tlagrange.m
// [in Octave]:
//   plot(db(fft(faustout,1024)(1:512,:)));
// Conclusion: 4th order seems worth it, all considered.

// Individual tests:
process1 = 1-1' : fdelay1(N,5.5);
process2 = 1-1' : fdelay2(N,5.5);
process3 = 1-1' : fdelay3(N,5.5);
process4 = 1-1' : fdelay4(N,5.5);

// Test a range of 4th-order cases:
process = 1-1' <: (fdelay4(N,5.1),
                   fdelay4(N,5.2),
                   fdelay4(N,5.3),
                   fdelay4(N,5.4),
                   fdelay4(N,5.5));

// Test a range of 3rd-order cases:
process3R = 1-1' <: (fdelay3(N,5.1),
                   fdelay3(N,5.2),
                   fdelay3(N,5.3),
                   fdelay3(N,5.4),
                   fdelay3(N,5.5));

// Test a range of 2nd-order cases:
process2R = 1-1' <: (fdelay2(N,5.1),
                   fdelay2(N,5.2),
                   fdelay2(N,5.3),
                   fdelay2(N,5.4),
                   fdelay2(N,5.5));

// Test a range of 1st-order cases:
process1 = 1-1' <: (fdelay1(N,5.1),
                   fdelay1(N,5.2),
                   fdelay1(N,5.3),
                   fdelay1(N,5.4),
                   fdelay1(N,5.5));

Figure J.8: File lagrange.lib containing Lagrange interpolation utilities written in the Faust language.

 

// lagrange.lib - Lagrange interpolation utilities

declare name "Lagrange Interpolation Library";
declare author "Julius O. Smith (jos at ccrma.stanford.edu)";
declare copyright "Julius O. Smith III";
declare version "1.0";
declare license "STK-4.3"; // Synthesis Tool Kit 4.3 (MIT style license)
//declare reference 
//"http://ccrma.stanford.edu/~jos/pasp/Lagrange_Interpolation";

import("music.lib"); // define delay, frac

// NOTE: While the implementations below appear to use multiple delay lines,
//  they in fact use only one thanks to optimization by the Faust compiler.

// first-order case (linear), delay in [0,1]:
fdelay1(n,d,x)  = delay(n,id,x)*(1 - fd) + delay(n,id+1,x)*fd
with {
  id = int(d);
  fd = frac(d);
};

// second-order case (quadratic), delay in [0,1] ([1,2] same quality):
fdelay2(n,d,x) = delay(n,id,x)*(1-fd)*(2-fd)/2 
               + delay(n,id+1,x)*(2-fd)*fd
               + delay(n,id+2,x)*(fd-1)*fd
with {
  id = int(d);
  fd = frac(d);
};

// third-order (cubic) case, delay in [1,2]:
fdelay3(n,d,x) = delay(n,id,x) * (0-fdm1*fdm2*fdm3)/6
               + delay(n,id+1,x) * fd*fdm2*fdm3/2
               + delay(n,id+2,x) * (0-fd*fdm1*fdm3)/2
               + delay(n,id+3,x) * fd*fdm1*fdm2/6
with {
  id = int(d-1); 
  fd = 1+frac(d);
  fdm1 = fd-1;
  fdm2 = fd-2;
  fdm3 = fd-3;
};

// fourth-order (quartic) case, delay in [1,2] ([2,3] same quality):
fdelay4(n,d,x) = delay(n,id,x)   * fdm1*fdm2*fdm3*fdm4/24 
               + delay(n,id+1,x) * (0-fd*fdm2*fdm3*fdm4)/6
               + delay(n,id+2,x) * fd*fdm1*fdm3*fdm4/4
               + delay(n,id+3,x) * (0-fd*fdm1*fdm2*fdm4)/6
               + delay(n,id+4,x) * fd*fdm1*fdm2*fdm3/24
with {
  id = int(d-1);
  fd = 1+frac(d);
  fdm1 = fd-1;
  fdm2 = fd-2;
  fdm3 = fd-3;
  fdm4 = fd-4;
};

Figure J.9: Comparison of Lagrange and Optimal Chebyshev Fractional-Delay Filter Frequency Responses