The possibilities offered by the algorithm can be extended using allpass filters. The allpass filters destroy the inharmonicity of the string.
t_int *kastro_perform(t_int *w) { t_float *fb = (t_float *)(w[1]); t_float *freq = (t_float *)(w[2]); t_float *damp= (t_float *)(w[3]); t_float *atten= (t_float *)(w[4]); t_float *inharm = (t_float *)(w[5]); t_float *out = (t_float *)(w[6]); t_kastro *x = (t_kastro *)(w[7]); int n_tick = (int)(w[8]); float yn,xn,ym1 = x->b_ym1,ym2 = x->b_ym2; float ym1all=x->b_ym1all,ym2all=x-> b_ym2all,xm2=x->b_xm2, xm1=x->b_xm1; float *vinut; float *vinbridge; int posr=x->v_posright; int posl=x->v_posleft; long i; double wn, sn, cs, alpha, a0, a1, a2, b0, b1, b2; double sr =x->srate; double gain = 1.0; // x->f_gain; float a0l,a1l,a2l,b1l,b2l; float alphar; float ynout; float ynn; float vin, vh; float srate=x->srate; float del_tot; int samp_rperiod; int samp_lperiod; vinut=x->vinut; if (x->v_obj.z_disabled) // if the user disactivates the object from the dsp chain goto out; for(i=0; i<n_tick; i++) { del_tot=srate/(freq[i]); //total length of the delay line samp_rperiod=(int)del_tot; //integer part if(samp_rperiod<0) samp_rperiod = 0; if(samp_rperiod>BUFSIZE-1) samp_rperiod = BUFSIZE-1; alphar=(float)(del_tot-samp_rperiod); //fractional part vin = vinut[(posr+i + BUFSIZE - samp_rperiod)% BUFSIZE] *(1-alphar)+(vinut[(posr+i+ BUFSIZE - samp_rperiod-1)% BUFSIZE])*(alphar) ; //check if damping is >= 1, to avoid to distroy the filter if (damp[i]>0.99999) damp[i]=0.999; if (damp[i]<0.000001) damp[i]=0.000001; //check for attenuation too if (atten[i]>1) atten[i]=1; if (atten[i]<0.00001) atten[i]=0.00001; //low-pass filter for damping ynn = atten[i]*(damp[i]*vin +((1-damp[i]) *ym1)); //update ym1 = vin; vh=fb[i]+ynn; //all-pass filter for inharmonicity if (inharm[i]>0.99999) inharm[i]=0.999; if (inharm[i]<0) inharm[i]=0; ynout = -inharm[i] * vh + xm1+ inharm[i] * ym1all; xm1 = vh; ym1all = ynout; vinut[(posr+i+BUFSIZE)%BUFSIZE] =ynout; out[i]=ynout; } posr=(posr+n_tick)%BUFSIZE; x->b_ym1 = ym1; x->b_ym2 = ym2; x->b_ym2all = ym2all; x->b_xm1 = xm1; x->b_ym1all = ym1all; x->v_posright=posr; x->v_posleft=posl; out: return (w+9); }