<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><div>I think now would be a good place to pause the implementation and discuss </div><div>if there is a need for this kind of algorithm in the PD community. </div><div><br class=""></div><div>We saw some use cases for the Whittaker-Shannon interpolation where</div><div>we gain in quality and/or speed. Namely waveguides and pitch shifters. </div><div><br class=""></div><div>Is there anything else, where we would like to use this interpolation kernel?</div><div>Like in general resampling?</div><div>In the tabread~ for sample playback?</div><div>Are there really quality or speed issues that we could solve there?</div><div><br class=""></div><div><div>For our use case the code I submitted is good enough. But I would be happy </div><div>to spend more time optimizing if there is a need / a broad use of that algorithm.</div><div>In that case we need another restructuring and some help from somebody</div><div>who is very proficient in writing pd source code. :)</div></div><div><br class=""></div><div>Chuck, I commented your last message below!</div><div><br class=""></div><div><blockquote type="cite" class=""><div class="">Am 2021-05-03 um 03:52 schrieb Charles Z Henry <<a href="mailto:czhenry@gmail.com" class="">czhenry@gmail.com</a>>:</div><br class="Apple-interchange-newline"><div class=""><div class="">On Sun, May 2, 2021 at 4:28 PM Clemens <<a href="mailto:clemens@chair.audio" class="">clemens@chair.audio</a>> wrote:<br class=""><blockquote type="cite" class=""><br class="">When there is no aliasing to worry about, i might set the cutoff to pi<br class="">again...<br class="">On low sample rates (22kHz), the lower cutoff is definitely noticeable.<br class=""></blockquote><br class="">I like this thesis you posted earlier<br class=""><a href="https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf" class="">https://www2.spsc.tugraz.at/www-archive/downloads/Mueller11_DopplerSRC_0.pdf</a><br class=""><br class="">f(at) <-> 1/|a|*F(w/a)<br class="">This (and two formulas that follow) is listed as Smith's algorithm.  I<br class="">actually got to speak with Julius Smith at the 2012 LAC about this<br class="">formula.  I asked, can't we do any better in terms than O(a*n) number<br class="">of computations for a>1?  He said nope!, but I still have some<br class="">questions there.  Playback with speeds less than 1 always use O(n),<br class="">rather than O(a*n).  I wrote an anti-aliasing external tabread4a~ that<br class="">implements this formula per sample and works pretty well, except it<br class="">becomes expensive when you transpose a few octaves up.<br class=""></div></div></blockquote><div><br class=""></div><div>I guess that’s the code you are referring to (?):</div><div><a href="https://lists.puredata.info/pipermail/pd-list/2007-03/048397.html" class="">https://lists.puredata.info/pipermail/pd-list/2007-03/048397.html</a></div><div><br class=""></div><div>Is it using a lookup of the sync table at certain fixed points?</div><div>If so, would it be „compatible“ in the sense that these points are part </div><div>of the lookup table in delreadsinc~?</div><br class=""><blockquote type="cite" class=""><div class=""><div class=""><br class="">For a variable delay line (like vd~), the paper's contents are more<br class="">relevant and maybe you should consider writing a vdsinc~ object next<br class="">(once you've optimized this one).  delread is more like a fixed,<br class="">stationary delay line, and I think it's better to default as a literal<br class="">interpolator (LP_SCALE = 1).<br class=""></div></div></blockquote><div><br class=""></div><div>Ok. That’s what I thought as well. Still I would need to read a bit more  </div><div>about the subject to understand how the variable delay case relates to </div><div>down-/upsampling. In tabread4a~ you derive the sampling factor from </div><div>the difference in delay time, right?</div><br class=""><blockquote type="cite" class=""><div class=""><div class=""><br class=""><blockquote type="cite" class="">Btw. I just implemented sharing of the interpolation table of the<br class="">delreadsinc~ object according to your suggestions.<br class="">It counts the number of references and frees the pointer when no object<br class="">uses it anymore.<br class=""></blockquote><br class="">    I read through the changes.  For those reading along see the rest<br class="">at line 348 of d_delay.c at<br class=""><a href="https://github.com/chairaudio/pure-data/blob/feature/delreadsinc/src/d_delay.c" class="">https://github.com/chairaudio/pure-data/blob/feature/delreadsinc/src/d_delay.c</a><br class=""> (in progress, current code re: the tables pasted below).<br class="">    1. I wonder is this global table properly scoped for Pd.  I<br class="">seriously don't know.  But it could be used for multiple objects,<br class="">vdsinc~, tabreadsinc~, etc.... So you ought to think about how it<br class="">could be re-used between classes and choose a scalable approach now.<br class=""></div></div></blockquote><div><br class=""></div><div>Good question. :) I’m not sure if Miller has suggestions or examples how </div><div>global variables are handled in PD.</div><br class=""><blockquote type="cite" class=""><div class=""><div class="">    2. I think the declared LP_SCALE variable is a bad approach as<br class="">currently implemented.  You can't change it at runtime.  It can be<br class="">used in the perform routine instead, with an additional argument to<br class="">set the cutoff frequency per object.  This is also good for your<br class="">testing, as you'll be able to put objects with different LP_SCALE<br class="">values side-by-side for comparison.<br class=""></div></div></blockquote><div><br class=""></div><div>That’s right, LP_SCALE should be tunable in real-time. </div><div>I think that would need only little change to the code. </div><br class=""><blockquote type="cite" class=""><div class=""><div class="">    3. I'm convinced this is the right approach for making an optimal<br class="">interpolator.  Adding more length pays off over upsampling, because<br class="">you evaluate the convolution at a fewer number of points.  Your<br class="">implementation needs a lot of tuning---you can squeeze out more<br class="">performance.  Focus on the inner most loops.  Also, you'll have to<br class="">compare some arm vs intel/amd platforms at some point.  I think you're<br class="">closer to the beginning than the end.  This is maybe not what you<br class="">wanted to research in the first place<br class=""></div></div></blockquote><div><br class=""></div><div>When it gets down to processor specific optimization, I’m not sure how </div><div>much the compiler does already and if it’s really worth the effort. </div><div>I would investigate this if we find that we want to put this code into pd.</div><div><br class=""></div><div>Best wishes, </div><div>Clemens</div><br class=""><blockquote type="cite" class=""><div class=""><div class=""><br class=""><br class="">...<br class="">typedef struct _sigvdsinc<br class="">{<br class="">t_object x_obj;<br class="">t_symbol *x_sym;<br class="">t_float x_sr; /* samples per msec */<br class="">int x_zerodel; /* 0 or vecsize depending on read/write order */<br class="">t_float x_f;<br class="">} t_sigvdsinc;<br class=""><br class="">typedef struct _sigvdsinc_sharedtable<br class="">{<br class="">// sinc table is held in global variable<br class="">t_sample sinc_array[SINC_LEN];<br class="">// derivative of sinc funtion for interpolation of sinc function table<br class="">t_sample sinc_diff_array[SINC_LEN];<br class="">// reference count for shared pointer<br class="">int ref_count;<br class="">} t_sigvdsinc_sharedtable;<br class=""><br class="">t_sigvdsinc_sharedtable *gSharedTable = NULL;<br class=""><br class="">void sigvdsinc_initialize_sinc_table(t_sample* sinc_array, t_sample*<br class="">sinc_diff_array) {<br class="">sinc_array[0] = LP_SCALE;<br class="">sinc_diff_array[SINC_LEN] = 0;<br class="">float a0 = 0.35875;<br class="">float a1 = 0.48829;<br class="">float a2 = 0.14128;<br class="">float a3 = 0.01168;<br class="">for (int i=1; i<SINC_LEN; i++) {<br class="">float idx = 0.5*(M_PI*i/SINC_LEN - M_PI);<br class="">// four term blackmanharris after<br class=""><a href="https://www.mathworks.com/help/signal/ref/blackmanharris.html" class="">https://www.mathworks.com/help/signal/ref/blackmanharris.html</a><br class="">float blackmanharris= a0 - a1*cos(2*idx) + a2*cos(4*idx) - a3*cos(6*idx);<br class="">// blackmanharris windowed sinc function<br class="">sinc_array[i] =<br class="">sin(LP_SCALE*M_PI*(float)i/STEPS_ZC)/(LP_SCALE*M_PI*(float)i/STEPS_ZC)<br class="">* blackmanharris * LP_SCALE;<br class="">// calculate derivative for table interpolation<br class="">sinc_diff_array[i-1] = sinc_array[i]-sinc_array[i-1];<br class="">}<br class="">}<br class=""><br class="">...<br class="">void sigvdsinc_free(t_sigvdsinc *x)<br class="">{<br class="">// delete shared interpolation table if need be<br class="">if(gSharedTable != NULL){<br class="">if(gSharedTable->ref_count==1){<br class="">free((void*)(gSharedTable));<br class="">gSharedTable = NULL;;<br class="">}<br class="">else{<br class="">gSharedTable->ref_count--;<br class="">}<br class="">}<br class="">}<br class=""></div></div></blockquote></div><br class=""></body></html>