1 /*
2     vpvoc.c:
3 
4     Copyright (C) 1992 Richard Karpen
5 
6     This file is part of Csound.
7 
8     The Csound Library is free software; you can redistribute it
9     and/or modify it under the terms of the GNU Lesser General Public
10     License as published by the Free Software Foundation; either
11     version 2.1 of the License, or (at your option) any later version.
12 
13     Csound is distributed in the hope that it will be useful,
14     but WITHOUT ANY WARRANTY; without even the implied warranty of
15     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16     GNU Lesser General Public License for more details.
17 
18     You should have received a copy of the GNU Lesser General Public
19     License along with Csound; if not, write to the Free Software
20     Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA
21     02110-1301 USA
22 */
23 
24 /**************************************************************/
25 /***********tableseg, tablexseg, voscili, vpvoc************/
26 /*** By Richard Karpen - July-October 1992************/
27 /************************************************************/
28 
29 #include "pvoc.h"
30 #include <math.h>
31 
tblesegset(CSOUND * csound,TABLESEG * p)32 int32_t tblesegset(CSOUND *csound, TABLESEG *p)
33 {
34     TSEG    *segp;
35     int32_t     nsegs;
36     MYFLT   **argp, dur;
37     FUNC    *nxtfunc, *curfunc;
38     int32    flength;
39 
40     if (UNLIKELY(!(p->INCOUNT & 1))) {
41       return csound->InitError(csound, "%s",
42                                Str("incomplete number of input arguments"));
43     }
44 
45     {
46       PVOC_GLOBALS  *p_ = PVOC_GetGlobals(csound);
47       p_->tbladr = p;
48     }
49 
50     nsegs = (p->INCOUNT >> 1);  /* count segs & alloc if nec */
51 
52     if ((segp = (TSEG *) p->auxch.auxp) == NULL ||
53         p->auxch.size<(nsegs+1)*sizeof(TSEG)) {
54       csound->AuxAlloc(csound, (size_t)(nsegs+1)*sizeof(TSEG), &p->auxch);
55       p->cursegp = segp = (TSEG *) p->auxch.auxp;
56       (segp+nsegs)->cnt = MAXPOS;
57     }
58     argp = p->argums;
59     if (UNLIKELY((nxtfunc = csound->FTnp2Finde(csound, *argp++)) == NULL))
60         return NOTOK;
61     flength = nxtfunc->flen;
62     p->outfunc =
63       (FUNC*) csound->Calloc(csound, sizeof(FUNC));
64     p->outfunc->ftable =
65       (MYFLT*)csound->Calloc(csound, (1 + flength) * sizeof(MYFLT));
66     p->outfunc->flen = nxtfunc->flen;
67     p->outfunc->lenmask = nxtfunc->lenmask;
68     p->outfunc->lobits = nxtfunc->lobits;
69     p->outfunc->lomask = nxtfunc->lomask;
70     p->outfunc->lodiv = nxtfunc->lodiv;
71     //memset(p->outfunc->ftable, 0, sizeof(MYFLT)*(flength+1)); not needed -- calloc
72     if (**argp <= 0.0)  return OK;         /* if idur1 <= 0, skip init  */
73     p->cursegp = segp;                      /* else proceed from 1st seg */
74     segp--;
75     do {
76         segp++;                 /* init each seg ..  */
77         curfunc = nxtfunc;
78         dur = **argp++;
79         if (UNLIKELY((nxtfunc = csound->FTnp2Finde(csound, *argp++)) == NULL))
80           return OK;
81         if (LIKELY(dur > FL(0.0))) {
82                 segp->d = dur * CS_EKR;
83                 segp->function =  curfunc;
84                 segp->nxtfunction = nxtfunc;
85                 segp->cnt = (int32) (segp->d + FL(0.5));
86         }
87         else break;             /*  .. til 0 dur or done */
88     } while (--nsegs);
89     segp++;
90     segp->d = FL(0.0);
91     segp->cnt = MAXPOS;         /* set last cntr to infin */
92     segp->function =  nxtfunc;
93     segp->nxtfunction = nxtfunc;
94     return OK;
95 }
96 
ktableseg(CSOUND * csound,TABLESEG * p)97 int32_t ktableseg(CSOUND *csound, TABLESEG *p)
98 {
99     TSEG        *segp;
100     MYFLT       *curtab, *nxttab,curval, nxtval, durovercnt=FL(0.0);
101     int32_t         i;
102     int32        flength, upcnt;
103 
104     /* RWD fix */
105     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1;
106     segp = p->cursegp;
107     curtab = segp->function->ftable;
108     nxttab = segp->nxtfunction->ftable;
109     upcnt = (int32)segp->d-segp->cnt;
110     if (upcnt > 0)
111       durovercnt = segp->d/upcnt;
112     while (--segp->cnt < 0)
113       p->cursegp = ++segp;
114     flength = segp->function->flen;
115     for (i=0; i<flength; i++) {
116       curval = curtab[i];
117       nxtval = nxttab[i];
118       if (durovercnt > FL(0.0))
119         p->outfunc->ftable[i] = (curval + ((nxtval - curval) / durovercnt));
120       else
121         p->outfunc->ftable[i] = curval;
122     }
123     return OK;
124  err1:
125     return csound->PerfError(csound, &(p->h), "%s",
126                              Str("tableseg: not initialised"));
127 }
128 
ktablexseg(CSOUND * csound,TABLESEG * p)129 int32_t ktablexseg(CSOUND *csound, TABLESEG *p)
130 {
131     TSEG        *segp;
132     MYFLT       *curtab, *nxttab,curval, nxtval, cntoverdur=FL(0.0);
133     int32_t         i;
134     int32        flength, upcnt;
135 
136     /* RWD fix */
137     if (UNLIKELY(p->auxch.auxp==NULL)) goto err1;
138     segp = p->cursegp;
139     curtab = segp->function->ftable;
140     nxttab = segp->nxtfunction->ftable;
141     upcnt = (int32)segp->d-segp->cnt;
142     if (upcnt > 0) cntoverdur = upcnt/ segp->d;
143     while(--segp->cnt < 0)
144       p->cursegp = ++segp;
145     flength = segp->function->flen;
146     for (i=0; i<flength; i++) {
147       curval = curtab[i];
148       nxtval = nxttab[i];
149       p->outfunc->ftable[i] =
150         (curval + ((nxtval - curval) * (cntoverdur*cntoverdur)));
151     }
152     return OK;
153  err1:
154     return csound->PerfError(csound, &(p->h), "%s",
155                              Str("tablexseg: not initialised"));
156 }
157 
158 /************************************************************/
159 /*****************VPVOC**************************************/
160 /************************************************************/
161 
162 #define WLN   1         /* time window is WLN*2*ksmps long */
163 #define OPWLEN (2*WLN*CS_KSMPS)    /* manifest used for final time wdw */
164 
vpvset_(CSOUND * csound,VPVOC * p,int32_t stringname)165 int32_t vpvset_(CSOUND *csound, VPVOC *p, int32_t stringname)
166 {
167     uint32_t i;
168     char     pvfilnam[MAXNAME];
169     PVOCEX_MEMFILE  pp;
170     int32_t  frInc, chans; /* THESE SHOULD BE SAVED IN PVOC STRUCT */
171 
172     p->pp = PVOC_GetGlobals(csound);
173     /* If optional table given, fake it up -- JPff  */
174     if (*p->isegtab == FL(0.0))
175       p->tableseg = p->pp->tbladr;
176     else {
177       csound->AuxAlloc(csound, sizeof(TABLESEG), &p->auxtab);
178       p->tableseg = (TABLESEG*) p->auxtab.auxp;
179       if (UNLIKELY((p->tableseg->outfunc =
180                     csound->FTnp2Find(csound, p->isegtab)) == NULL)) {
181         return csound->InitError(csound, "%s%f",
182                                  Str("vpvoc: Could not find ifnmagctrl table "),
183                                  *p->isegtab);
184       }
185     }
186     if (UNLIKELY(p->tableseg == NULL))
187       return csound->InitError(csound, "%s",
188                                Str("vpvoc: associated tableseg not found"));
189 
190     if (p->auxch.auxp == NULL) {              /* if no buffers yet, alloc now */
191       MYFLT *fltp;
192       csound->AuxAlloc(csound,
193                        (PVDATASIZE + PVFFTSIZE * 3 + PVWINLEN) * sizeof(MYFLT),
194                        &p->auxch);
195       fltp = (MYFLT *) p->auxch.auxp;
196       p->lastPhase = fltp;   fltp += PVDATASIZE;    /* and insert addresses */
197       p->fftBuf = fltp;      fltp += PVFFTSIZE;
198       p->dsBuf = fltp;       fltp += PVFFTSIZE;
199       p->outBuf = fltp;      fltp += PVFFTSIZE;
200       p->window = fltp;
201     }
202     if (stringname==0){
203       if (csound->ISSTRCOD(*p->ifilno))
204         strNcpy(pvfilnam,get_arg_string(csound, *p->ifilno), MAXNAME-1);
205       else csound->strarg2name(csound, pvfilnam, p->ifilno, "pvoc.",0);
206     }
207     else strNcpy(pvfilnam, ((STRINGDAT *)p->ifilno)->data, MAXNAME-1);
208 
209     if (UNLIKELY(csound->PVOCEX_LoadFile(csound, pvfilnam, &pp) != 0))
210       return csound->InitError(csound, Str("VPVOC cannot load %s"), pvfilnam);
211 
212     p->frSiz = pp.fftsize;
213     frInc    = pp.overlap;
214     chans    = pp.chans;
215     p->asr   = pp.srate;
216     if (UNLIKELY(p->asr != CS_ESR)) {                /* & chk the data */
217       csound->Warning(csound, Str("%s's srate = %8.0f, orch's srate = %8.0f"),
218                               pvfilnam, p->asr, CS_ESR);
219     }
220     if (UNLIKELY(p->frSiz > PVFRAMSIZE)) {
221       return csound->InitError(csound,
222                                Str("PVOC frame %ld bigger than %ld in %s"),
223                                (long) p->frSiz, (long) PVFRAMSIZE, pvfilnam);
224     }
225     if (UNLIKELY(p->frSiz < 128)) {
226       return csound->InitError(csound,
227                                Str("PVOC frame %ld seems too small in %s"),
228                                (long) p->frSiz, pvfilnam);
229     }
230     if (UNLIKELY(chans != 1)) {
231       return csound->InitError(csound, Str("%d chans (not 1) in PVOC file %s"),
232                                        (int32_t) chans, pvfilnam);
233     }
234     /* Check that pv->frSiz is a power of two too ? */
235     p->frPtr = (float*) pp.data;
236     p->baseFr = 0;  /* point to first data frame */
237     p->maxFr = pp.nframes - 1;
238     /* highest possible frame index */
239     p->frPktim = (MYFLT) CS_KSMPS / (MYFLT) frInc;
240     /* factor by which to mult expand phase diffs (ratio of samp spacings) */
241     p->frPrtim = CS_ESR / (MYFLT) frInc;
242     /* factor by which to mulitply 'real' time index to get frame index */
243     /* amplitude scale for PVOC */
244  /* p->scale = (MYFLT) pp.fftsize * ((MYFLT) pp.fftsize / (MYFLT) pp.winsize);
245   */
246     p->scale = (MYFLT) pp.fftsize * FL(0.5);
247     p->scale *= csound->GetInverseRealFFTScale(csound, pp.fftsize);
248     /* 2*incr/OPWLEN scales down for win ovlp, windo'd 1ce (but 2ce?) */
249     /* 1/frSiz is the required scale down before (i)FFT */
250     p->prFlg = 1;    /* true */
251     p->opBpos = 0;
252     p->lastPex = FL(1.0);   /* needs to know last pitchexp to update phase */
253     /* Set up time window */
254     memset(p->lastPhase, 0, sizeof(MYFLT)*pvdasiz(p));
255     /* for (i = 0; i < pvdasiz(p); ++i) {  /\* or maybe pvdasiz(p) *\/ */
256     /*   p->lastPhase[i] = FL(0.0); */
257     /* } */
258     if (UNLIKELY((OPWLEN / 2 + 1) > PVWINLEN)) {
259       return csound->InitError(csound, Str("ksmps of %d needs wdw of %d, "
260                                            "max is %d for pv %s"),
261                                        CS_KSMPS, (OPWLEN / 2 + 1),
262                                        PVWINLEN, pvfilnam);
263     }
264     for (i = 0; i < OPWLEN / 2 + 1; ++i)    /* time window is OPWLEN long */
265       p->window[i] = (FL(0.5) - FL(0.5) * COS(TWOPI_F*(MYFLT)i/(MYFLT)OPWLEN));
266     /* NB: HANNING */
267     memset(p->outBuf, 0, sizeof(MYFLT)*pvfrsiz(p));
268     /* for (i = 0; i < pvfrsiz(p); ++i) */
269     /*   p->outBuf[i] = FL(0.0); */
270     MakeSinc(p->pp);                    /* sinctab is same for all instances */
271     if (p->memenv.auxp == NULL || p->memenv.size < pvdasiz(p)*sizeof(MYFLT))
272         csound->AuxAlloc(csound, pvdasiz(p) * sizeof(MYFLT), &p->memenv);
273     return OK;
274 }
275 
vpvset(CSOUND * csound,VPVOC * p)276 int32_t vpvset(CSOUND *csound, VPVOC *p){
277   return vpvset_(csound,p,0);
278 }
279 
vpvset_S(CSOUND * csound,VPVOC * p)280 int32_t vpvset_S(CSOUND *csound, VPVOC *p){
281   return vpvset_(csound,p,1);
282 }
283 
vpvoc(CSOUND * csound,VPVOC * p)284 int32_t vpvoc(CSOUND *csound, VPVOC *p)
285 {
286     MYFLT     *ar = p->rslt;
287     MYFLT     frIndx;
288     MYFLT     *buf = p->fftBuf;
289     MYFLT     *buf2 = p->dsBuf;
290     int32_t       asize = pvdasiz(p); /* fix */
291     int32_t       size = pvfrsiz(p);
292     int32_t       buf2Size, outlen;
293     int32_t       circBufSize = PVFFTSIZE;
294     int32_t       specwp = (int32_t) *p->ispecwp;   /* spectral warping flag */
295     MYFLT     pex, scaleFac = p->scale;
296     TABLESEG  *q = p->tableseg;
297     int32     i, j;
298 
299     /* RWD fix */
300     if (UNLIKELY(p->auxch.auxp == NULL)) goto err1;
301 
302     pex = *p->kfmod;
303     outlen = (int32_t) (((MYFLT) size) / pex);
304     /* use outlen to check window/krate/transpose combinations */
305     if (UNLIKELY(outlen>PVFFTSIZE)) { /* Maximum transposition down is one octave */
306                             /* ..so we won't run into buf2Size problems */
307       goto err2;
308     }
309     if (UNLIKELY(outlen<(int32_t)(2*CS_KSMPS))) {
310       /* minimum post-squeeze windowlength */
311       goto err3;
312     }
313     buf2Size = OPWLEN;     /* always window to same length after DS */
314     if (UNLIKELY((frIndx = *p->ktimpnt * p->frPrtim) < 0)) {
315       goto err4;
316     }
317     if (frIndx > (MYFLT)p->maxFr) { /* not past last one */
318       frIndx = (MYFLT)p->maxFr;
319       if (UNLIKELY(p->prFlg)) {
320         p->prFlg = 0;   /* false */
321         csound->Warning(csound, Str("PVOC ktimpnt truncated to last frame"));
322       }
323     }
324 
325     FetchIn(p->frPtr, buf, size, frIndx);
326 
327 /**** Apply "spectral envelope" to magnitudes ********/
328     if (pex > FL(1.0))
329       scaleFac /= pex;
330     {
331       MYFLT *ftable = q->outfunc->ftable;
332       for (i = 0, j = 0; i <= size; i += 2, j++)
333         buf[i] *= ftable[j] * scaleFac;
334     }
335 /***************************************************/
336 
337     FrqToPhase(buf, asize, pex * (MYFLT) CS_KSMPS, p->asr,
338                (MYFLT) (0.5 * ((pex / p->lastPex) - 1)));
339     /* accumulate phase and wrap to range -PI to PI */
340     RewrapPhase(buf, asize, p->lastPhase);
341 
342     if (specwp == 0 || (p->prFlg)++ == -(int32_t)specwp) {
343       /* ?screws up when prFlg used */
344       /* specwp=0 => normal; specwp = -n => just nth frame */
345       if (UNLIKELY(specwp < 0))
346         csound->Warning(csound, Str("PVOC debug: one frame gets through\n"));
347       if (specwp > 0)
348         PreWarpSpec(buf, asize, pex, (MYFLT *)p->memenv.auxp);
349 
350       Polar2Real_PVOC(csound, buf, size);
351 
352       if (pex != FL(1.0))
353         UDSample(p->pp, buf,
354                  (FL(0.5) * ((MYFLT) size - pex * (MYFLT) buf2Size)),
355                  buf2, size, buf2Size, pex);
356       else
357         memcpy(buf2, buf + (int32_t) ((size - buf2Size) >> 1),
358                sizeof(MYFLT) * buf2Size);
359       if (specwp >= 0)
360         ApplyHalfWin(buf2, p->window, buf2Size);
361     }
362     else {
363       memset(buf2, 0, sizeof(MYFLT)*buf2Size);
364       /* for (n = 0; n < buf2Size; ++n) */
365       /*   buf2[n] = FL(0.0); */
366     }
367 
368     addToCircBuf(buf2, p->outBuf, p->opBpos, CS_KSMPS, circBufSize);
369     writeClrFromCircBuf(p->outBuf, ar, p->opBpos, CS_KSMPS, circBufSize);
370     p->opBpos += CS_KSMPS;
371     if (p->opBpos > circBufSize)
372       p->opBpos -= circBufSize;
373     addToCircBuf(buf2 + CS_KSMPS, p->outBuf, p->opBpos,
374                  buf2Size - CS_KSMPS, circBufSize);
375     p->lastPex = pex;        /* needs to know last pitchexp to update phase */
376 
377     return OK;
378  err1:
379     return csound->PerfError(csound, &(p->h),
380                              Str("vpvoc: not initialised"));
381  err2:
382     return csound->PerfError(csound, &(p->h),
383                              Str("PVOC transpose too low"));
384  err3:
385     return csound->PerfError(csound, &(p->h),
386                              Str("PVOC transpose too high"));
387  err4:
388     return csound->PerfError(csound, &(p->h),
389                              Str("PVOC timpnt < 0"));
390 }
391 
392