1 /*
2     ugens6.c:
3 
4     Copyright (C) 1991-2000 Barry Vercoe, John ffitch, Jens Groh,
5                             Hans Mikelson, Istvan Varga
6 
7     This file is part of Csound.
8 
9     The Csound Library is free software; you can redistribute it
10     and/or modify it under the terms of the GNU Lesser General Public
11     License as published by the Free Software Foundation; either
12     version 2.1 of the License, or (at your option) any later version.
13 
14     Csound is distributed in the hope that it will be useful,
15     but WITHOUT ANY WARRANTY; without even the implied warranty of
16     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17     GNU Lesser General Public License for more details.
18 
19     You should have received a copy of the GNU Lesser General Public
20     License along with Csound; if not, write to the Free Software
21     Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA
22     02110-1301 USA
23 */
24 
25 #include "csoundCore.h" /*                              UGENS6.C        */
26 #include "ugens6.h"
27 #include <math.h>
28 
29 #define log001 (-FL(6.9078))    /* log(.001) */
30 
downset(CSOUND * csound,DOWNSAMP * p)31 int32_t downset(CSOUND *csound, DOWNSAMP *p)
32 {
33     if (UNLIKELY((p->len = (uint32_t)*p->ilen) > CS_KSMPS))
34       return csound->InitError(csound, "ilen > ksmps");
35     return OK;
36 }
37 
downsamp(CSOUND * csound,DOWNSAMP * p)38 int32_t downsamp(CSOUND *csound, DOWNSAMP *p)
39 {
40     IGN(csound);
41     MYFLT       *asig, sum;
42     uint32_t offset = p->h.insdshead->ksmps_offset;
43     uint32_t early  = p->h.insdshead->ksmps_no_end;
44     int32_t len, n;
45 
46     if (p->len <= 1)
47       *p->kr = p->asig[offset];
48     else {
49       asig = p->asig;
50       sum = FL(0.0);
51       len = p->len;
52       if (len>(int32_t)(CS_KSMPS-early)) len = early;
53       for (n=offset; n<len; n++) {
54         sum += asig[n];
55       }
56       *p->kr = sum / p->len;
57     }
58     return OK;
59 }
60 
upsamp(CSOUND * csound,UPSAMP * p)61 int32_t upsamp(CSOUND *csound, UPSAMP *p)
62 {
63     IGN(csound);
64     MYFLT kval = *p->ksig;
65     MYFLT *ar = p->ar;
66     uint32_t offset = p->h.insdshead->ksmps_offset;
67     uint32_t early  = p->h.insdshead->ksmps_no_end;
68     uint32_t n, nsmps = CS_KSMPS;
69 
70     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
71     if (UNLIKELY(early)) {
72       nsmps -= early;
73       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
74     }
75     for (n = offset; n < nsmps; n++)
76       ar[n] = kval;
77     return OK;
78 }
79 
a_k_set(CSOUND * csound,INTERP * p)80 int32_t a_k_set(CSOUND *csound, INTERP *p)
81 {
82     IGN(csound);
83     p->prev = FL(0.0);
84     p->init_k = 0;              /* IV - Sep 5 2002 */
85     return OK;
86 }
87 
interpset(CSOUND * csound,INTERP * p)88 int32_t interpset(CSOUND *csound, INTERP *p)
89 {
90     IGN(csound);
91     if (*p->istor == FL(0.0)) {
92       p->prev = (*p->imode == FL(0.0) ? *p->istart : FL(0.0));
93       p->init_k = (*p->imode == FL(0.0) ? 0 : 1);       /* IV - Sep 5 2002 */
94     }
95 
96     return OK;
97 }
98 
interp(CSOUND * csound,INTERP * p)99 int32_t interp(CSOUND *csound, INTERP *p)
100 {
101     IGN(csound);
102     MYFLT *ar, val, incr;
103     uint32_t offset = p->h.insdshead->ksmps_offset;
104     uint32_t early  = p->h.insdshead->ksmps_no_end;
105     uint32_t n, nsmps = CS_KSMPS;
106 
107     ar = p->rslt;
108     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
109     if (UNLIKELY(early)) {
110       nsmps -= early;
111       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
112     }
113     if (p->init_k) {
114       p->init_k = 0;
115       p->prev = *p->xsig;
116     }
117     val = p->prev;
118     incr = (*p->xsig - val) / (nsmps-offset);
119     for (n=offset; n<nsmps; n++) {
120       ar[n] = val += incr;
121     }
122     p->prev = val;
123     return OK;
124 }
125 
indfset(CSOUND * csound,INDIFF * p)126 int32_t indfset(CSOUND *csound, INDIFF *p)
127 {
128     IGN(csound);
129     if (*p->istor == FL(0.0))   /* IV - Sep 5 2002 */
130       p->prev = FL(0.0);
131     return OK;
132 }
133 
kntegrate(CSOUND * csound,INDIFF * p)134 int32_t kntegrate(CSOUND *csound, INDIFF *p)
135 {
136     IGN(csound);
137     *p->rslt = p->prev += *p->xsig;
138     return OK;
139 }
140 
integrate(CSOUND * csound,INDIFF * p)141 int32_t integrate(CSOUND *csound, INDIFF *p)
142 {
143     IGN(csound);
144     MYFLT       *rslt, *asig, sum;
145     uint32_t offset = p->h.insdshead->ksmps_offset;
146     uint32_t early  = p->h.insdshead->ksmps_no_end;
147     uint32_t n, nsmps = CS_KSMPS;
148 
149     rslt = p->rslt;
150     if (UNLIKELY(offset)) memset(rslt, '\0', offset*sizeof(MYFLT));
151     if (UNLIKELY(early)) {
152       nsmps -= early;
153       memset(&rslt[nsmps], '\0', early*sizeof(MYFLT));
154     }
155     asig = p->xsig;
156     sum = p->prev;
157     for (n=offset; n<nsmps; n++) {
158       rslt[n] = sum += asig[n];
159     }
160     p->prev = sum;
161     return OK;
162 }
163 
kdiff(CSOUND * csound,INDIFF * p)164 int32_t kdiff(CSOUND *csound, INDIFF *p)
165 {
166     IGN(csound);
167     MYFLT       tmp;
168     tmp = *p->xsig;             /* IV - Sep 5 2002: fix to make */
169     *p->rslt = tmp - p->prev;   /* diff work when the input and */
170     p->prev = tmp;              /* output argument is the same  */
171     return OK;
172 }
173 
diff(CSOUND * csound,INDIFF * p)174 int32_t diff(CSOUND *csound, INDIFF *p)
175 {
176     IGN(csound);
177     MYFLT       *ar, *asig, prev, tmp;
178     uint32_t offset = p->h.insdshead->ksmps_offset;
179     uint32_t early  = p->h.insdshead->ksmps_no_end;
180     uint32_t n, nsmps = CS_KSMPS;
181 
182     ar = p->rslt;
183     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
184     if (UNLIKELY(early)) {
185       nsmps -= early;
186       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
187     }
188     asig = p->xsig;
189     prev = p->prev;
190     for (n=offset; n<nsmps; n++) {
191       tmp = asig[n];            /* IV - Sep 5 2002: fix to make */
192       ar[n] = tmp - prev;       /* diff work when the input and */
193       prev = tmp;               /* output argument is the same  */
194     }
195     p->prev = prev;
196     return OK;
197 }
198 
samphset(CSOUND * csound,SAMPHOLD * p)199 int32_t samphset(CSOUND *csound, SAMPHOLD *p)
200 {
201     IGN(csound);
202     if (!(*p->istor))
203       p->state = *p->ival;
204     p->audiogate = IS_ASIG_ARG(p->xgate) ? 1 : 0;
205     return OK;
206 }
207 
ksmphold(CSOUND * csound,SAMPHOLD * p)208 int32_t ksmphold(CSOUND *csound, SAMPHOLD *p)
209 {
210     IGN(csound);
211     if (*p->xgate > FL(0.0))
212       p->state = *p->xsig;
213     *p->xr = p->state;
214     return OK;
215 }
216 
samphold(CSOUND * csound,SAMPHOLD * p)217 int32_t samphold(CSOUND *csound, SAMPHOLD *p)
218 {
219     IGN(csound);
220     MYFLT       *ar, *asig, *agate, state;
221     uint32_t offset = p->h.insdshead->ksmps_offset;
222     uint32_t early  = p->h.insdshead->ksmps_no_end;
223     uint32_t n, nsmps = CS_KSMPS;
224 
225     ar = p->xr;
226     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
227     if (UNLIKELY(early)) {
228       nsmps -= early;
229       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
230     }
231     asig = p->xsig;
232     state = p->state;
233     if (p->audiogate) {
234       agate = p->xgate;
235       for (n=offset; n<nsmps; n++) {
236         if (agate[n] > FL(0.0))
237           state = asig[n];
238         ar[n] = state;
239       }
240     }
241     else {
242       if (*p->xgate > FL(0.0)) {
243         for (n=offset; n<nsmps; n++) {
244           ar[n] = state = asig[n];
245         }
246       }
247       else {
248         for (n=offset; n<nsmps; n++) {
249           ar[n] = state;
250         }
251       }
252     }
253     p->state = state;
254     return OK;
255 }
256 
delset(CSOUND * csound,DELAY * p)257 int32_t delset(CSOUND *csound, DELAY *p)
258 {
259     int32_t      npts;
260     char        *auxp;
261 
262     if (UNLIKELY(*p->istor && p->auxch.auxp != NULL))
263       return OK;
264     /* Round not truncate */
265     if (UNLIKELY((npts = MYFLT2LRND(*p->idlt * csound->esr)) <= 0)) {
266       return csound->InitError(csound, Str("illegal delay time"));
267     }
268 
269     if ((auxp = p->auxch.auxp) == NULL ||
270         npts != p->npts) { /* new space if reqd */
271       csound->AuxAlloc(csound, (int32_t)npts*sizeof(MYFLT), &p->auxch);
272       auxp = p->auxch.auxp;
273       p->npts = npts;
274     }
275     else if (!(*p->istor)) {                    /* else if requested */
276       memset(auxp, '\0', npts*sizeof(MYFLT));
277     }
278     p->curp = (MYFLT *) auxp;
279 
280     return OK;
281 }
282 
delrset(CSOUND * csound,DELAYR * p)283 int32_t delrset(CSOUND *csound, DELAYR *p)
284 {
285     uint32_t    npts;
286     MYFLT       *auxp;
287 
288     if (UNLIKELY(!IS_ASIG_ARG(p->ar)))
289       return csound->InitError(csound, Str("delayr: invalid outarg type"));
290     /* fifo for delayr pointers by Jens Groh: */
291     /* append structadr for delayw to fifo: */
292     if (csound->first_delayr != NULL)       /* fifo not empty */
293       ((DELAYR*) csound->last_delayr)->next_delayr = p;
294     else                                    /* fifo empty */
295       csound->first_delayr = (void*) p;
296     csound->last_delayr = (void*) p;
297     csound->delayr_stack_depth++;
298     p->next_delayr = NULL;
299     if (p->OUTOCOUNT > 1) {
300       /* set optional output arg if specified */
301       *(p->indx) = (MYFLT)-(csound->delayr_stack_depth);
302     }
303 
304     if (UNLIKELY(*p->istor != FL(0.0) && p->auxch.auxp != NULL))
305       return OK;
306     /* ksmps is min dely */
307     if (UNLIKELY((npts=(uint32_t)MYFLT2LRND(*p->idlt*csound->esr)) < CS_KSMPS)) {
308       return csound->InitError(csound, Str("illegal delay time"));
309     }
310     if ((auxp = (MYFLT*)p->auxch.auxp) == NULL ||       /* new space if reqd */
311         npts != p->npts) {
312       csound->AuxAlloc(csound, (int32_t)npts*sizeof(MYFLT), &p->auxch);
313       auxp = (MYFLT*)p->auxch.auxp;
314       p->npts = npts;
315     }
316     else if (*p->istor == FL(0.0)) {            /* else if requested */
317       memset(auxp, 0, npts*sizeof(MYFLT));
318     }
319     p->curp = auxp;
320     return OK;
321 }
322 
delwset(CSOUND * csound,DELAYW * p)323 int32_t delwset(CSOUND *csound, DELAYW *p)
324 {
325    /* fifo for delayr pointers by Jens Groh: */
326     if (UNLIKELY(csound->first_delayr == NULL)) {
327       return csound->InitError(csound,
328                                Str("delayw: associated delayr not found"));
329     }
330     p->delayr = (DELAYR*) csound->first_delayr;         /* adr delayr struct */
331     /* remove structadr from fifo */
332     if (csound->last_delayr == csound->first_delayr) {  /* fifo will be empty */
333       csound->first_delayr = NULL;
334     }
335     else    /* fifo will not be empty */
336       csound->first_delayr = ((DELAYR*) csound->first_delayr)->next_delayr;
337     csound->delayr_stack_depth--;
338     return OK;
339 }
340 
delayr_find(CSOUND * csound,MYFLT * ndx)341 static DELAYR *delayr_find(CSOUND *csound, MYFLT *ndx)
342 {
343     DELAYR  *d = (DELAYR*) csound->first_delayr;
344     int32_t     n = (int32_t)MYFLT2LRND(*ndx);
345 
346     if (UNLIKELY(d == NULL)) {
347       csound->InitError(csound, Str("deltap: associated delayr not found"));
348       return NULL;
349     }
350     if (!n)
351       return (DELAYR*) csound->last_delayr;     /* default: use last delayr */
352     else if (n > 0)
353       n = csound->delayr_stack_depth - n;       /* ndx > 0: LIFO index mode */
354     else
355       n = -n;                                   /* ndx < 0: FIFO index mode */
356     if (UNLIKELY(n < 1 || n > csound->delayr_stack_depth)) {
357       csound->InitError(csound,
358                         Str("deltap: delayr index %.0f is out of range"),
359                         (double)*ndx);
360       return NULL;
361     }
362     /* find delay line */
363     while (--n)
364       d = d->next_delayr;
365     return d;
366 }
367 
tapset(CSOUND * csound,DELTAP * p)368 int32_t tapset(CSOUND *csound, DELTAP *p)
369 {
370     p->delayr = delayr_find(csound, p->indx);
371     return (p->delayr != NULL ? OK : NOTOK);
372 }
373 
delay(CSOUND * csound,DELAY * p)374 int32_t delay(CSOUND *csound, DELAY *p)
375 {
376     MYFLT       *ar, *asig, *curp, *endp;
377     uint32_t offset = 0;
378     uint32_t n, nsmps = CS_KSMPS;
379 
380     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1;  /* RWD fix */
381     ar = p->ar;
382     if (csound->oparms->sampleAccurate) {
383       uint32_t early  = p->h.insdshead->ksmps_no_end;
384       offset = p->h.insdshead->ksmps_offset;
385 
386       if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
387       if (UNLIKELY(early)) {
388         nsmps -= early;
389         memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
390       }
391     }
392     asig = p->asig;
393     curp = p->curp;
394     endp = (MYFLT *) p->auxch.endp;
395     for (n=offset; n<nsmps; n++) {
396       MYFLT in = asig[n];       /* Allow overwriting form */
397       ar[n] = *curp;
398       *curp = in;
399       if (UNLIKELY(++curp >= endp))
400         curp = (MYFLT *) p->auxch.auxp;
401     }
402     p->curp = curp;             /* sav the new curp */
403 
404     return OK;
405  err1:
406     return csound->PerfError(csound, &(p->h),
407                              Str("delay: not initialised"));
408 }
409 
delayr(CSOUND * csound,DELAYR * p)410 int32_t delayr(CSOUND *csound, DELAYR *p)
411 {
412     MYFLT       *ar, *curp, *endp;
413     uint32_t offset = p->h.insdshead->ksmps_offset;
414     uint32_t early  = p->h.insdshead->ksmps_no_end;
415     uint32_t n, nsmps = CS_KSMPS;
416 
417     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1; /* RWD fix */
418     ar = p->ar;
419     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
420     if (UNLIKELY(early)) {
421       nsmps -= early;
422       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
423     }
424     curp = p->curp;
425     endp = (MYFLT *) p->auxch.endp;
426     for (n=offset; n<nsmps; n++) {
427       ar[n] = *curp++;
428       if (UNLIKELY(curp >= endp))
429         curp = (MYFLT *) p->auxch.auxp;
430     }
431     return OK;
432  err1:
433     return csound->PerfError(csound, &(p->h),
434                              Str("delayr: not initialised"));
435 }
436 
delayw(CSOUND * csound,DELAYW * p)437 int32_t delayw(CSOUND *csound, DELAYW *p)
438 {
439     DELAYR      *q = p->delayr;
440     MYFLT       *asig, *curp, *endp;
441     uint32_t offset = p->h.insdshead->ksmps_offset;
442     uint32_t early  = p->h.insdshead->ksmps_no_end;
443     uint32_t n, nsmps = CS_KSMPS;
444 
445     if (UNLIKELY(q->auxch.auxp==NULL)) goto err1; /* RWD fix */
446     asig = p->asig;
447     curp = q->curp;
448     endp = (MYFLT *) q->auxch.endp;
449     if (UNLIKELY(early)) nsmps -= early;
450     for (n=offset; n<nsmps; n++) {
451       *curp = asig[n];
452       if (UNLIKELY(++curp >= endp))
453         curp = (MYFLT *) q->auxch.auxp;
454     }
455     q->curp = curp;                                     /* now sav new curp */
456     return OK;
457  err1:
458     return csound->PerfError(csound, &(p->h),
459                              Str("delayw: not initialised"));
460 }
461 
deltap(CSOUND * csound,DELTAP * p)462 int32_t deltap(CSOUND *csound, DELTAP *p)
463 {
464     DELAYR      *q = p->delayr;
465     MYFLT       *ar, *tap, *endp;
466     uint32_t offset = p->h.insdshead->ksmps_offset;
467     uint32_t early  = p->h.insdshead->ksmps_no_end;
468     uint32_t n, nsmps = CS_KSMPS;
469 
470     if (UNLIKELY(q->auxch.auxp==NULL)) goto err1; /* RWD fix */
471     ar = p->ar;
472     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
473     if (UNLIKELY(early)) {
474       nsmps -= early;
475       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
476     }
477     tap = q->curp - MYFLT2LRND(*p->xdlt * csound->esr);
478     while (tap < (MYFLT *) q->auxch.auxp)
479       tap += q->npts;
480     endp = (MYFLT *) q->auxch.endp;
481     for (n=offset; n<nsmps; n++) {
482       if (UNLIKELY(tap >= endp))
483         tap -= q->npts;
484       ar[n] = *tap++;
485     }
486     return OK;
487  err1:
488     return csound->PerfError(csound, &(p->h),
489                              Str("deltap: not initialised"));
490 }
491 
deltapi(CSOUND * csound,DELTAP * p)492 int32_t deltapi(CSOUND *csound, DELTAP *p)
493 {
494     DELAYR      *q = p->delayr;
495     MYFLT       *ar, *tap, *prv, *begp, *endp;
496     uint32_t offset = p->h.insdshead->ksmps_offset;
497     uint32_t early  = p->h.insdshead->ksmps_no_end;
498     uint32_t n, nsmps = CS_KSMPS;
499     int32_t       idelsmps;
500     MYFLT       delsmps, delfrac;
501 
502     if (UNLIKELY(q->auxch.auxp==NULL)) goto err1;
503     ar = p->ar;
504     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
505     if (UNLIKELY(early)) {
506       nsmps -= early;
507       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
508     }
509     begp = (MYFLT *) q->auxch.auxp;
510     endp = (MYFLT *) q->auxch.endp;
511     if (!IS_ASIG_ARG(p->xdlt)) {
512       if (*p->xdlt == INFINITY) goto err2;
513       delsmps = *p->xdlt * csound->esr;
514       idelsmps = (int32_t)delsmps;
515       delfrac = delsmps - idelsmps;
516       tap = q->curp - idelsmps;
517       while (tap < begp) tap += q->npts;
518       for (n=offset; n<nsmps; n++) {
519         if (UNLIKELY(tap >= endp))
520           tap -= q->npts;
521         if (UNLIKELY((prv = tap - 1) < begp))
522           prv += q->npts;
523         ar[n] = *tap + (*prv - *tap) * delfrac;
524         tap++;
525       }
526     }
527     else {
528       MYFLT *timp = p->xdlt, *curq = q->curp;
529       for (n=offset; n<nsmps; n++) {
530         if (timp[n] == INFINITY) goto err2;
531         delsmps = timp[n] * csound->esr;
532         idelsmps = (int32_t)delsmps;
533         delfrac = delsmps - idelsmps;
534         tap = curq++ - idelsmps;
535         if (UNLIKELY(tap < begp)) tap += q->npts;
536         else if (UNLIKELY(tap >= endp))
537           tap -= q->npts;
538         if (UNLIKELY((prv = tap - 1) < begp))
539           prv += q->npts;
540         ar[n] = *tap + (*prv - *tap) * delfrac;
541       }
542     }
543     return OK;
544  err1:
545     return csound->PerfError(csound, &(p->h),
546                               Str("deltapi: not initialised"));
547   err2:
548     return csound->PerfError(csound, &(p->h),
549                               Str("deltapi: INF delaytime"));
550 }
551 
552 /* ***** From Hans Mikelson ************* */
553 /* Delay N samples */
deltapn(CSOUND * csound,DELTAP * p)554 int32_t deltapn(CSOUND *csound, DELTAP *p)
555 {
556     DELAYR *q = p->delayr;
557     MYFLT  *ar, *tap, *begp, *endp;
558     uint32_t offset = p->h.insdshead->ksmps_offset;
559     uint32_t early  = p->h.insdshead->ksmps_no_end;
560     uint32_t n, nsmps = CS_KSMPS;
561     int32_t  idelsmps;
562     MYFLT  delsmps;
563 
564     if (UNLIKELY(q->auxch.auxp==NULL)) goto err1;
565     ar = p->ar;
566     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
567     if (UNLIKELY(early)) {
568       nsmps -= early;
569       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
570     }
571     begp = (MYFLT *) q->auxch.auxp;
572     endp = (MYFLT *) q->auxch.endp;
573     if (!IS_ASIG_ARG(p->xdlt)) {
574       delsmps = *p->xdlt;
575       idelsmps = (int32_t)delsmps;
576       tap = q->curp - idelsmps;
577       while (tap < begp) tap += q->npts;
578       for (n=offset; n<nsmps; n++) {
579         while (UNLIKELY(tap >= endp ))
580           tap -= q->npts;
581         while (UNLIKELY(tap < begp))
582           tap += q->npts;
583         ar[n] = *tap;
584         tap++;
585       }
586     }
587     else {
588       MYFLT *timp = p->xdlt, *curq = q->curp;
589       for (n=offset; n<nsmps; n++) {
590         delsmps = timp[n];
591         idelsmps = (int32_t)delsmps;
592         if (UNLIKELY((tap = curq++ - idelsmps) < begp))
593           tap += q->npts;
594         else if (UNLIKELY(tap >= endp))
595           tap -= q->npts;
596         ar[n] = *tap;
597       }
598     }
599     return OK;
600  err1:
601     return csound->PerfError(csound, &(p->h),
602                              Str("deltapn: not initialised"));
603 }
604 
605 /* **** JPff **** */
deltap3(CSOUND * csound,DELTAP * p)606 int32_t deltap3(CSOUND *csound, DELTAP *p)
607 {
608     DELAYR      *q = p->delayr;
609     MYFLT       *ar, *tap, *prv, *begp, *endp;
610     uint32_t offset = p->h.insdshead->ksmps_offset;
611     uint32_t early  = p->h.insdshead->ksmps_no_end;
612     uint32_t n, nsmps = CS_KSMPS;
613     int32_t       idelsmps;
614     MYFLT       delsmps, delfrac;
615 
616     if (UNLIKELY(q->auxch.auxp==NULL)) goto err1;
617     ar = p->ar;
618     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
619     if (UNLIKELY(early)) {
620       nsmps -= early;
621       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
622     }
623     begp = (MYFLT *) q->auxch.auxp;
624     endp = (MYFLT *) q->auxch.endp;
625     if (!IS_ASIG_ARG(p->xdlt)) {
626       if (*p->xdlt == INFINITY) goto err2;
627       delsmps = *p->xdlt * csound->esr;
628       idelsmps = (int32_t)delsmps;
629       delfrac = delsmps - idelsmps;
630       tap = q->curp - idelsmps;
631       while (tap < begp) tap += q->npts;
632       for (n=offset; n<nsmps; n++) {
633         MYFLT ym1, y0, y1, y2;
634         if (UNLIKELY(tap >= endp))
635           tap -= q->npts;
636         if (UNLIKELY((prv = tap - 1) < begp))
637           prv += q->npts;
638         if (UNLIKELY(prv - 1 < begp))
639           y2 = *(prv-1+q->npts);
640         else
641           y2 = *(prv-1);
642         if (UNLIKELY(tap + 1 >= endp))
643           ym1 = *(tap+1-q->npts);
644         else
645           ym1 = *(tap+1);
646         y0 = *tap; y1 = *prv;
647         {
648           MYFLT w, x, y, z;
649           z = delfrac * delfrac; z--; z *= FL(0.16666666666667);
650           y = delfrac; y++; w = (y *= FL(0.5)); w--;
651           x = FL(3.0) * z; y -= x; w -= z; x -= delfrac;
652           ar[n] = (w*ym1 + x*y0 + y*y1 + z*y2) * delfrac + y0;
653         }
654         tap++;
655       }
656     }
657     else {
658       MYFLT *timp = p->xdlt, *curq = q->curp;
659       for (n=offset; n<nsmps; n++) {
660         MYFLT ym1, y0, y1, y2;
661         if (timp[n] == INFINITY) goto err2;
662         delsmps = *timp++ * csound->esr;
663         idelsmps = (int32_t)delsmps;
664         delfrac = delsmps - idelsmps;
665         if (UNLIKELY((tap = curq++ - idelsmps) < begp))
666           tap += q->npts;
667         else if (UNLIKELY(tap >= endp))
668           tap -= q->npts;
669         if (UNLIKELY((prv = tap - 1) < begp))
670           prv += q->npts;
671         if (UNLIKELY(prv - 1 < begp)) y2 = *(prv-1+q->npts);
672         else                          y2 = *(prv-1);
673         if (UNLIKELY(tap + 1 >= endp)) ym1 = *(tap+1-q->npts);
674         else                           ym1 = *(tap+1);
675         y0 = *tap; y1 = *prv;
676         {
677           MYFLT w, x, y, z;
678           z = delfrac * delfrac; z--; z *= FL(0.1666666667);
679           y = delfrac; y++; w = (y *= FL(0.5)); w--;
680           x = FL(3.0) * z; y -= x; w -= z; x -= delfrac;
681           ar[n] = (w*ym1 + x*y0 + y*y1 + z*y2) * delfrac + y0;
682         }
683       }
684     }
685     return OK;
686  err1:
687     return csound->PerfError(csound, &(p->h),
688                              Str("deltap3: not initialised"));
689   err2:
690     return csound->PerfError(csound, &(p->h),
691                               Str("deltapi: INF delaytime"));
692 
693 }
694 
695 
696 /* deltapx and deltapxw opcodes by Istvan Varga */
697 
tapxset(CSOUND * csound,DELTAPX * p)698 int32_t tapxset(CSOUND *csound, DELTAPX *p)
699 {
700     p->delayr = delayr_find(csound, p->indx);
701     if (UNLIKELY(p->delayr == NULL))
702       return NOTOK;
703     p->wsize = (int32_t)(*(p->iwsize) + FL(0.5));          /* window size */
704     p->wsize = ((p->wsize + 2) >> 2) << 2;
705     if (UNLIKELY(p->wsize < 4)) p->wsize = 4;
706     if (UNLIKELY(p->wsize > 1024)) p->wsize = 1024;
707     /* wsize = 4: d2x = 1 - 1/3, wsize = 64: d2x = 1 - 1/36 */
708     p->d2x = 1.0 - pow((double)p->wsize * 0.85172, -0.89624);
709     p->d2x /= (double)((p->wsize * p->wsize) >> 2);
710     return OK;
711 }
712 
deltapx(CSOUND * csound,DELTAPX * p)713 int32_t deltapx(CSOUND *csound, DELTAPX *p)                 /* deltapx opcode */
714 {
715     DELAYR  *q = p->delayr;
716     MYFLT   *out1, *del, *buf1, *bufp, *bufend;
717     uint32_t offset = p->h.insdshead->ksmps_offset;
718     uint32_t early  = p->h.insdshead->ksmps_no_end;
719     uint32_t n, nsmps = CS_KSMPS;
720     int32_t   indx, maxd, xpos;
721 
722     if (UNLIKELY(q->auxch.auxp == NULL)) goto err1; /* RWD fix */
723     out1 = p->ar; del = p->adlt;
724     if (UNLIKELY(offset)) memset(out1, '\0', offset*sizeof(MYFLT));
725     if (UNLIKELY(early)) {
726       nsmps -= early;
727       memset(&out1[nsmps], '\0', early*sizeof(MYFLT));
728     }
729     buf1 = (MYFLT *) q->auxch.auxp;
730     indx = (int32_t) (q->curp - buf1);
731     maxd = q->npts; bufend = buf1 + maxd;
732 
733     if (p->wsize != 4) {                /* window size >= 8 */
734       double  d, x1, n1, w, d2x;
735       int32_t     i2, i;
736       i2 = (p->wsize >> 1);
737       /* wsize = 4: d2x = 1 - 1/3, wsize = 64: d2x = 1 - 1/36 */
738       d2x = p->d2x;
739       for (n=offset; n<nsmps; n++) {
740         /* x1: fractional part of delay time */
741         /* x2: sine of x1 (for interpolation) */
742         /* xpos: integer part of delay time (buffer position to read from) */
743 
744         x1 = (double)indx - (double)*(del++) * (double)csound->esr;
745         while (x1 < 0.0) x1 += (double)maxd;
746         xpos = (int32_t)x1; x1 -= (double)xpos;
747         while (xpos >= maxd) xpos -= maxd;
748 
749         if (x1 > 0.00000001 && x1 < 0.99999999) {
750           xpos -= i2;
751           while (xpos < 0) xpos += maxd;
752           d = (double)(1 - i2) - x1;
753           bufp = buf1 + xpos;
754           i = i2;
755           n1 = 0.0;
756           do {
757             w = 1.0 - d * d * d2x;
758             if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
759             n1 += w * w * (double)*bufp / d; d++;
760             w = 1.0 - d * d * d2x;
761             if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
762             n1 -= w * w * (double)*bufp / d; d++;
763           } while (--i);
764           out1[n] = (MYFLT)(n1 * sin(PI * x1) / PI);
765         }
766         else {                                          /* integer sample */
767           xpos = MYFLT2LRND((double)xpos + x1);         /* position */
768           if (xpos >= maxd) xpos -= maxd;
769           out1[n] = buf1[xpos];
770         }
771         indx++;
772       }
773     }
774     else {                          /* window size = 4, cubic interpolation */
775       double  x, am1, a0, a1, a2;
776       for (n=offset; n<nsmps; n++) {
777         am1 = (double)indx - (double)*(del++) * (double)csound->esr;
778         while (am1 < 0.0) am1 += (double)maxd;
779         xpos = (int32_t) am1; am1 -= (double)xpos;
780 
781         a0  = am1 * am1; a2 = 0.16666667 * (am1 * a0 - am1);    /* sample +2 */
782         a1  = 0.5 * (a0 + am1) - 3.0 * a2;                      /* sample +1 */
783         am1 = 0.5 * (a0 - am1) - a2;                            /* sample -1 */
784         a0  = 3.0 * a2 - a0; a0++;                              /* sample 0  */
785 
786         bufp = (xpos ? (buf1 + (xpos - 1L)) : (bufend - 1));
787         while (bufp >= bufend) bufp -= maxd;
788         x = am1 * (double)*bufp;   if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
789         x += a0 * (double)*bufp;   if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
790         x += a1 * (double)*bufp;   if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
791         x += a2 * (double)*bufp;
792 
793         indx++; out1[n] = (MYFLT)x;
794       }
795     }
796     return OK;
797  err1:
798     return csound->PerfError(csound, &(p->h),
799                              Str("deltap: not initialised"));
800 }
801 
deltapxw(CSOUND * csound,DELTAPX * p)802 int32_t deltapxw(CSOUND *csound, DELTAPX *p)                /* deltapxw opcode */
803 {
804     DELAYR  *q = p->delayr;
805     MYFLT   *in1, *del, *buf1, *bufp, *bufend;
806     uint32_t offset = p->h.insdshead->ksmps_offset;
807     uint32_t early  = p->h.insdshead->ksmps_no_end;
808     uint32_t n, nsmps = CS_KSMPS;
809     int32_t   indx, maxd, xpos;
810 
811     if (UNLIKELY(q->auxch.auxp == NULL)) goto err1; /* RWD fix */
812     in1 = p->ar; del = p->adlt;
813     if (UNLIKELY(early)) nsmps -= early;
814     buf1 = (MYFLT *) q->auxch.auxp;
815     indx = (int32_t) (q->curp - buf1);
816     maxd = q->npts; bufend = buf1 + maxd;
817 
818     if (p->wsize != 4) {                /* window size >= 8 */
819       double  d, x1, n1, w, d2x;
820       int32_t     i2, i;
821       i2 = (p->wsize >> 1);
822       /* wsize = 4: d2x = 1 - 1/3, wsize = 64: d2x = 1 - 1/36 */
823       d2x = p->d2x;
824       for (n=offset; n<nsmps; n++) {
825         /* x1: fractional part of delay time */
826         /* x2: sine of x1 (for interpolation) */
827         /* xpos: integer part of delay time (buffer position to read from) */
828 
829         x1 = (double)indx - (double)*(del++) * (double)csound->esr;
830         while (x1 < 0.0) x1 += (double)maxd;
831         xpos = (int32_t) x1; x1 -= (double)xpos;
832         while (xpos >= maxd) xpos -= maxd;
833 
834         if (x1 > 0.00000001 && x1 < 0.99999999) {
835           n1 = (double)*in1 * (sin(PI * x1) / PI);
836           xpos -= i2;
837           while (xpos < 0) xpos += maxd;
838           d = (double)(1 - i2) - x1;
839           bufp = buf1 + xpos;
840           i = i2;
841           do {
842             w = 1.0 - d * d * d2x;
843             if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
844             *bufp = (MYFLT)((double)*bufp + w * w * n1 / d); d++;
845             w = 1.0 - d * d * d2x;
846             if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
847             *bufp = (MYFLT)((double)*bufp - w * w * n1 / d); d++;
848           } while (--i);
849         }
850         else {                                          /* integer sample */
851           xpos = MYFLT2LRND((double)xpos + x1);         /* position */
852           if (UNLIKELY(xpos >= maxd)) xpos -= maxd;
853           buf1[xpos] += in1[n];
854         }
855         indx++;
856       }
857     }
858     else {                          /* window size = 4, cubic interpolation */
859       double  x, am1, a0, a1, a2;
860       for (n=offset; n<nsmps; n++) {
861         am1 = (double)indx - (double)*(del++) * (double)csound->esr;
862         while (am1 < 0.0) am1 += (double)maxd;
863         xpos = (int32_t) am1; am1 -= (double)xpos;
864 
865         a0  = am1 * am1; a2 = 0.16666667 * (am1 * a0 - am1);    /* sample +2 */
866         a1  = 0.5 * (a0 + am1) - 3.0 * a2;                      /* sample +1 */
867         am1 = 0.5 * (a0 - am1) - a2;                            /* sample -1 */
868         a0  = 3.0 * a2 - a0; a0++;                              /* sample 0  */
869 
870         x = (double)in1[n];
871         bufp = (xpos ? (buf1 + (xpos - 1L)) : (bufend - 1));
872         while (bufp >= bufend) bufp -= maxd;
873         *bufp += (MYFLT)(am1 * x); if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
874         *bufp += (MYFLT)(a0 * x);  if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
875         *bufp += (MYFLT)(a1 * x);  if (UNLIKELY(++bufp >= bufend)) bufp = buf1;
876         *bufp += (MYFLT)(a2 * x);
877 
878         indx++;
879       }
880     }
881     return OK;
882  err1:
883     return csound->PerfError(csound, &(p->h),
884                              Str("deltap: not initialised"));
885 }
886 
del1set(CSOUND * csound,DELAY1 * p)887 int32_t del1set(CSOUND *csound, DELAY1 *p)
888 {
889     IGN(csound);
890     if (!(*p->istor))
891       p->sav1 = FL(0.0);
892     return OK;
893 }
894 
delay1(CSOUND * csound,DELAY1 * p)895 int32_t delay1(CSOUND *csound, DELAY1 *p)
896 {
897     IGN(csound);
898     MYFLT       *ar, *asig;
899     uint32_t offset = p->h.insdshead->ksmps_offset;
900     uint32_t early  = p->h.insdshead->ksmps_no_end;
901     uint32_t nsmps = CS_KSMPS;
902 
903     ar = p->ar;
904     /* asig = p->asig - 1; */
905     asig = p->asig;
906     ar[offset] = p->sav1;
907     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
908     if (UNLIKELY(early)) {
909       nsmps -= early;
910       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
911     }
912     memmove(&ar[offset+1], &asig[offset], sizeof(MYFLT)*(nsmps-1-offset));
913     p->sav1 = asig[nsmps-1];
914     return OK;
915 }
916 
cmbset(CSOUND * csound,COMB * p)917 int32_t cmbset(CSOUND *csound, COMB *p)
918 {
919     int32_t       lpsiz, nbytes;
920 
921     if (*p->insmps != 0) {
922       if (UNLIKELY((lpsiz = MYFLT2LRND(*p->ilpt))) <= 0) {
923         return csound->InitError(csound, Str("illegal loop time"));
924       }
925     }
926     else if (UNLIKELY((lpsiz = MYFLT2LRND(*p->ilpt * csound->esr)) <= 0)) {
927       return csound->InitError(csound, Str("illegal loop time"));
928     }
929     nbytes = lpsiz * sizeof(MYFLT);
930     if (p->auxch.auxp == NULL || (uint32_t)nbytes != p->auxch.size) {
931       csound->AuxAlloc(csound, (int32_t)nbytes, &p->auxch);
932       p->pntr = (MYFLT *) p->auxch.auxp;
933       p->prvt = FL(0.0);
934       p->coef = FL(0.0);
935     }
936     else if (!(*p->istor)) {
937       /* Does this assume sizeof(MYFLT)==sizeof(int32_t)?? */
938       p->pntr = (MYFLT *) p->auxch.auxp;
939       memset(p->auxch.auxp, '\0', nbytes);
940       p->prvt = FL(0.0);
941       p->coef = FL(0.0);
942     }
943     return OK;
944 }
945 
comb(CSOUND * csound,COMB * p)946 int32_t comb(CSOUND *csound, COMB *p)
947 {
948     uint32_t offset = p->h.insdshead->ksmps_offset;
949     uint32_t early  = p->h.insdshead->ksmps_no_end;
950     uint32_t n, nsmps = CS_KSMPS;
951     MYFLT       *ar, *asig, *xp, *endp;
952     MYFLT       coef = p->coef;
953 
954     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1; /* RWD fix */
955     if (p->prvt != *p->krvt) {
956       p->prvt = *p->krvt;
957       /*
958        * The argument to exp() in the following is sometimes a small
959        * enough negative number to result in a denormal (or worse)
960        * on Alpha. So if the result would be less than 1.0e-16, we
961        * just say it's zero and don't call exp().  heh 981101
962        */
963       double exp_arg = (double)(log001 * *p->ilpt / p->prvt);
964       if (UNLIKELY(exp_arg < -36.8413615))    /* ln(1.0e-16) */
965         coef = p->coef = FL(0.0);
966       else
967         coef = p->coef = (MYFLT)exp(exp_arg);
968     }
969     xp = p->pntr;
970     endp = (MYFLT *) p->auxch.endp;
971     ar = p->ar;
972     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
973     if (UNLIKELY(early)) {
974       nsmps -= early;
975       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
976     }
977     asig = p->asig;
978     for (n=offset; n<nsmps; n++) {
979       MYFLT out = *xp;
980       *xp *= coef;
981       *xp += asig[n];
982       ar[n] = out;
983       if (UNLIKELY(++xp >= endp))
984         xp = (MYFLT *) p->auxch.auxp;
985     }
986     p->pntr = xp;
987     return OK;
988  err1:
989     return csound->PerfError(csound, &(p->h),
990                              Str("comb: not initialised"));
991 }
992 
invcomb(CSOUND * csound,COMB * p)993 int32_t invcomb(CSOUND *csound, COMB *p)
994 {
995     int32_t n, nsmps = csound->ksmps;
996     MYFLT       *ar, *asig, *xp, *endp;
997     MYFLT       coef = p->coef;
998 
999     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1; /* RWD fix */
1000     if (p->prvt != *p->krvt) {
1001       p->prvt = *p->krvt;
1002       /*
1003        * The argument to exp() in the following is sometimes a small
1004        * enough negative number to result in a denormal (or worse)
1005        * on Alpha. So if the result would be less than 1.0e-16, we
1006        * just say it is zero and do not call exp().  heh 981101
1007        */
1008       double exp_arg = (double)(log001 * *p->ilpt / p->prvt);
1009       if (UNLIKELY(exp_arg < -36.8413615))    /* ln(1.0e-16) */
1010         coef = p->coef = FL(0.0);
1011       else
1012         coef = p->coef = (MYFLT)exp(exp_arg);
1013     }
1014     xp = p->pntr;
1015     endp = (MYFLT *) p->auxch.endp;
1016     ar = p->ar;
1017     asig = p->asig;
1018     MYFLT out;
1019     for (n=0; n<nsmps; n++) {
1020       out = *xp;
1021       ar[n] = (*xp = asig[n])-coef*out;
1022       if (UNLIKELY(++xp >= endp))
1023         xp = (MYFLT *) p->auxch.auxp;
1024     }
1025     p->pntr = xp;
1026     return OK;
1027  err1:
1028     return csound->PerfError(csound, &(p->h),
1029                              Str("combinv: not initialised"));
1030 }
1031 
alpass(CSOUND * csound,COMB * p)1032 int32_t alpass(CSOUND *csound, COMB *p)
1033 {
1034     uint32_t    offset = p->h.insdshead->ksmps_offset;
1035     uint32_t early  = p->h.insdshead->ksmps_no_end;
1036     uint32_t    n, nsmps = CS_KSMPS;
1037     MYFLT       *ar, *asig, *xp, *endp;
1038     MYFLT       y, z;
1039     MYFLT       coef = p->coef;
1040 
1041     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1; /* RWD fix */
1042     if (p->prvt != *p->krvt) {
1043       p->prvt = *p->krvt;
1044       coef = p->coef = EXP(log001 * *p->ilpt / p->prvt);
1045     }
1046     xp = p->pntr;
1047     endp = (MYFLT *) p->auxch.endp;
1048     ar = p->ar;
1049     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
1050     if (UNLIKELY(early)) {
1051       nsmps -= early;
1052       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
1053     }
1054     asig = p->asig;
1055     for (n=offset; n<nsmps; n++) {
1056       y = *xp;
1057       *xp++ = z = coef * y + asig[n];
1058       ar[n] = y - coef * z;
1059       if (UNLIKELY(xp >= endp))
1060         xp = (MYFLT *) p->auxch.auxp;
1061     }
1062     p->pntr = xp;
1063     return OK;
1064  err1:
1065     return csound->PerfError(csound, &(p->h),
1066                              Str("alpass: not initialised"));
1067 }
1068 
1069 static const MYFLT revlptimes[6] = {FL(0.0297), FL(0.0371), FL(0.0411),
1070                                     FL(0.0437), FL(0.0050), FL(0.0017)};
1071 
reverbinit(CSOUND * csound)1072 void reverbinit(CSOUND *csound)         /* called once by oload */
1073 {                                       /*  to init reverb data */
1074     const MYFLT *lptimp = revlptimes;
1075     int32_t     *lpsizp = csound->revlpsiz;
1076     int32_t n = 6;
1077 
1078     if (csound->revlpsum==0) {
1079       csound->revlpsum = 0;
1080       for (n=0; n<6; n++) {
1081         lpsizp[n] = MYFLT2LRND(lptimp[n] * csound->esr);
1082         csound->revlpsum += lpsizp[n];
1083       }
1084     }
1085 }
1086 
rvbset(CSOUND * csound,REVERB * p)1087 int32_t rvbset(CSOUND *csound, REVERB *p)
1088 {
1089     if (p->auxch.auxp == NULL) {                        /* if no space yet, */
1090       int32      *sizp = csound->revlpsiz;               /*    allocate it   */
1091       csound->AuxAlloc(csound, csound->revlpsum * sizeof(MYFLT), &p->auxch);
1092       p->adr1 = p->p1 = (MYFLT *) p->auxch.auxp;
1093       p->adr2 = p->p2 = p->adr1 + *sizp++;
1094       p->adr3 = p->p3 = p->adr2 + *sizp++;              /*    & init ptrs   */
1095       p->adr4 = p->p4 = p->adr3 + *sizp++;
1096       p->adr5 = p->p5 = p->adr4 + *sizp++;
1097       p->adr6 = p->p6 = p->adr5 + *sizp++;
1098       if (UNLIKELY(p->adr6 + *sizp != (MYFLT *) p->auxch.endp)) {
1099         return csound->InitError(csound, Str("revlpsiz inconsistent\n"));
1100       }
1101       p->prvt = FL(0.0);
1102     }
1103     else if (!(*p->istor)) {                    /* else if istor = 0 */
1104       memset(p->adr1, '\0', csound->revlpsum * sizeof(MYFLT));
1105       p->p1 = p->adr1;                          /*  and reset   */
1106       p->p2 = p->adr2;
1107       p->p3 = p->adr3;
1108       p->p4 = p->adr4;
1109       p->p5 = p->adr5;
1110       p->p6 = p->adr6;
1111       p->prvt = FL(0.0);
1112     }
1113     return OK;
1114 }
1115 
reverb(CSOUND * csound,REVERB * p)1116 int32_t reverb(CSOUND *csound, REVERB *p)
1117 {
1118     MYFLT       *asig, *p1, *p2, *p3, *p4, *p5, *p6, *ar, *endp;
1119     MYFLT       c1,c2,c3,c4,c5,c6;
1120     uint32_t offset = p->h.insdshead->ksmps_offset;
1121     uint32_t early  = p->h.insdshead->ksmps_no_end;
1122     uint32_t n, nsmps = CS_KSMPS;
1123 
1124     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1; /* RWD fix */
1125     if (UNLIKELY(p->prvt != *p->krvt)) {          /* People rarely change rvt */
1126       const MYFLT *lptimp = revlptimes;
1127       MYFLT       logdrvt = log001 / *p->krvt;
1128       c1=p->c1 = EXP(logdrvt * *lptimp++);
1129       c2=p->c2 = EXP(logdrvt * *lptimp++);
1130       c3=p->c3 = EXP(logdrvt * *lptimp++);
1131       c4=p->c4 = EXP(logdrvt * *lptimp++);
1132       c5=p->c5 = EXP(logdrvt * *lptimp++);
1133       c6=p->c6 = EXP(logdrvt * *lptimp++);
1134       p->prvt = *p->krvt;       /* JPff optimisation?? */
1135     }
1136     else {
1137       c1=p->c1;
1138       c2=p->c2;
1139       c3=p->c3;
1140       c4=p->c4;
1141       c5=p->c5;
1142       c6=p->c6;
1143    }
1144 
1145     p1 = p->p1;
1146     p2 = p->p2;
1147     p3 = p->p3;
1148     p4 = p->p4;
1149     p5 = p->p5;
1150     p6 = p->p6;
1151     endp = (MYFLT *) p->auxch.endp;
1152 
1153     ar = p->ar;
1154     if (UNLIKELY(offset)) memset(ar, '\0', offset*sizeof(MYFLT));
1155     if (UNLIKELY(early)) {
1156       nsmps -= early;
1157       memset(&ar[nsmps], '\0', early*sizeof(MYFLT));
1158     }
1159     asig = p->asig;
1160     for (n=offset; n<nsmps; n++) {
1161       MYFLT     cmbsum, y1, y2, z;
1162       MYFLT     sig = asig[n];
1163       cmbsum = *p1 + *p2 + *p3 + *p4;
1164       *p1 = c1 * *p1 + sig;
1165       *p2 = c2 * *p2 + sig;
1166       *p3 = c3 * *p3 + sig;
1167       *p4 = c4 * *p4 + sig;
1168       p1++; p2++; p3++; p4++;
1169       y1 = *p5;
1170       *p5++ = z = c5 * y1 + cmbsum;
1171       y1 -= c5 * z;
1172       y2 = *p6;
1173       *p6++ = z = c6 * y2 + y1;
1174       ar[n] = y2 - c6 * z;
1175       if (UNLIKELY(p1 >= p->adr2)) p1 = p->adr1;
1176       if (UNLIKELY(p2 >= p->adr3)) p2 = p->adr2;
1177       if (UNLIKELY(p3 >= p->adr4)) p3 = p->adr3;
1178       if (UNLIKELY(p4 >= p->adr5)) p4 = p->adr4;
1179       if (UNLIKELY(p5 >= p->adr6)) p5 = p->adr5;
1180       if (UNLIKELY(p6 >= endp))    p6 = p->adr6;
1181     }
1182     p->p1 = p1;
1183     p->p2 = p2;
1184     p->p3 = p3;
1185     p->p4 = p4;
1186     p->p5 = p5;
1187     p->p6 = p6;
1188     return OK;
1189  err1:
1190     return csound->PerfError(csound, &(p->h),
1191                              Str("reverb: not initialised"));
1192 }
1193 
panset(CSOUND * csound,PAN * p)1194 int32_t panset(CSOUND *csound, PAN *p)
1195 {
1196     FUNC  *ftp;
1197 
1198     if (UNLIKELY((ftp = csound->FTnp2Finde(csound, p->ifn)) == NULL))
1199       return NOTOK;
1200     p->ftp = ftp;
1201     p->xmul = (*p->imode == FL(0.0) ? FL(1.0) : (MYFLT)ftp->flen);
1202     p->xoff = (*p->ioffset == FL(0.0) ? (MYFLT)ftp->flen * FL(0.5) : FL(0.0));
1203 
1204     return OK;
1205 }
1206 
pan(CSOUND * csound,PAN * p)1207 int32_t pan(CSOUND *csound, PAN *p)
1208 {
1209     MYFLT   flend2, xndx_f, yndx_f, xt, yt, ch1, ch2, ch3, ch4;
1210     int32   xndx, yndx, flen;
1211     uint32_t offset = p->h.insdshead->ksmps_offset;
1212     uint32_t early  = p->h.insdshead->ksmps_no_end;
1213     uint32_t n, nsmps = CS_KSMPS;
1214     FUNC    *ftp;
1215 
1216     ftp = p->ftp;
1217     if (UNLIKELY(ftp == NULL)) goto err1;          /* RWD fix */
1218     xndx_f = (*p->kx * p->xmul) - p->xoff;
1219     yndx_f = (*p->ky * p->xmul) - p->xoff;
1220     flen = ftp->flen;
1221     flend2 = (MYFLT)flen * FL(0.5);
1222     xt = FABS(xndx_f);
1223     yt = FABS(yndx_f);
1224     if (xt > flend2 || yt > flend2) {
1225       if (xt > yt)
1226         yndx_f *= (flend2 / xt);
1227       else
1228         xndx_f *= (flend2 / yt);
1229     }
1230     xndx_f += flend2;
1231     yndx_f += flend2;
1232     xndx = MYFLT2LRND(xndx_f);
1233     yndx = MYFLT2LRND(yndx_f);
1234     xndx = (xndx >= 0L ? (xndx < flen ? xndx : flen) : 0L);
1235     yndx = (yndx >= 0L ? (yndx < flen ? yndx : flen) : 0L);
1236     ch1 = ftp->ftable[flen - xndx] * ftp->ftable[yndx];
1237     ch2 = ftp->ftable[xndx]        * ftp->ftable[yndx];
1238     ch3 = ftp->ftable[flen - xndx] * ftp->ftable[flen - yndx];
1239     ch4 = ftp->ftable[xndx]        * ftp->ftable[flen - yndx];
1240 
1241     if (UNLIKELY(offset)) {
1242       memset(p->r1, '\0', offset*sizeof(MYFLT));
1243       memset(p->r2, '\0', offset*sizeof(MYFLT));
1244       memset(p->r3, '\0', offset*sizeof(MYFLT));
1245       memset(p->r4, '\0', offset*sizeof(MYFLT));
1246     }
1247     if (UNLIKELY(early)) {
1248       nsmps -= early;
1249       memset(&p->r1[nsmps], '\0', early*sizeof(MYFLT));
1250       memset(&p->r2[nsmps], '\0', early*sizeof(MYFLT));
1251       memset(&p->r3[nsmps], '\0', early*sizeof(MYFLT));
1252       memset(&p->r4[nsmps], '\0', early*sizeof(MYFLT));
1253     }
1254     for (n=offset; n<nsmps; n++) {
1255       MYFLT sig = p->asig[n];
1256       p->r1[n] = sig * ch1;
1257       p->r2[n] = sig * ch2;
1258       p->r3[n] = sig * ch3;
1259       p->r4[n] = sig * ch4;
1260     }
1261 
1262     return OK;
1263  err1:
1264     return csound->PerfError(csound, &(p->h),
1265                              Str("pan: not initialised"));
1266 }
1267