Logo Search packages:      
Sourcecode: csound version File versions  Download package

pvlock.c
#include "csdl.h"
#include "pstream.h"
#define MAXOUTS 2

00005 typedef struct dats{
  OPDS h;
  MYFLT *out[MAXOUTS], *time, *kamp, *kpitch, *knum, *klock, *iN, *idecim, *konset, *offset, *dbthresh;
  int cnt, hsize, curframe, N, decim,tscale,nchans;
  double pos;
  MYFLT accum;
  AUXCH outframe[MAXOUTS], win, bwin[MAXOUTS], fwin[MAXOUTS], nwin[MAXOUTS], prev[MAXOUTS], framecount[MAXOUTS];
} DATASPACE;


static int sinit(CSOUND *csound, DATASPACE *p)
{

  int N =  *p->iN, i,size, nchans;
    int decim = *p->idecim;

    if (N) {
      for (i=0; N; i++){
        N >>= 1;
      }
      N = (int) pow(2., i-1.);
    } else N = 2048;
    if (decim == 0) decim = 4;

    p->hsize = N/decim;
    p->cnt = p->hsize;
    p->curframe = 0;

    nchans = csound->GetOutputArgCnt(p);
  if (UNLIKELY(nchans < 1 || nchans > MAXOUTS))
      csound->Die(csound, Str("invalid number of output arguments"));
    p->nchans = nchans;

    for (i=0; i < nchans; i++){
    size = (N+2)*sizeof(MYFLT);
    if (p->fwin[i].auxp == NULL || p->fwin[i].size < size)
      csound->AuxAlloc(csound, size, &p->fwin[i]);
    if (p->bwin[i].auxp == NULL || p->bwin[i].size < size)
      csound->AuxAlloc(csound, size, &p->bwin[i]);
    if (p->prev[i].auxp == NULL || p->prev[i].size < size)
      csound->AuxAlloc(csound, size, &p->prev[i]);
    size = decim*sizeof(int);
    if (p->framecount[i].auxp == NULL || p->framecount[i].size < size)
      csound->AuxAlloc(csound, size, &p->framecount[i]);
    memset(p->framecount[i].auxp,0,size);
    size = decim*sizeof(MYFLT)*N;
    if (p->outframe[i].auxp == NULL || p->outframe[i].size < size)
      csound->AuxAlloc(csound, size, &p->outframe[i]);
    memset(p->outframe[i].auxp,0,size);
    }
   size = N*sizeof(MYFLT);
    if (p->win.auxp == NULL || p->win.size < size)
      csound->AuxAlloc(csound, size, &p->win);

    for (i=0; i < N; i++)
      ((MYFLT *)p->win.auxp)[i] = 0.5 - 0.5*cos(i*2*PI/N);

    p->N = N;
    p->decim = decim;
    return OK;
}


static int sprocess(CSOUND *csound, DATASPACE *p)
{
    MYFLT pitch = *p->kpitch, *time = p->time, lock = *p->klock,
      *out, amp =*p->kamp;
    MYFLT *tab, frac,scale;
    FUNC *ft;
    int N = p->N, hsize = p->hsize, cnt = p->cnt, nchans = p->nchans;
    int ksmps = csound->GetKsmps(csound), n;
    int sizefrs, size, post, i, j, spos = p->pos;
    double pos;
    MYFLT *fwin, *bwin, in,
      *prev, *win = (MYFLT *) p->win.auxp;
    MYFLT *outframe;
    MYFLT ph_real, ph_im, tmp_real, tmp_im, div;
    int *framecnt;
    int curframe = p->curframe, decim = p->decim;

    for (n=0; n < ksmps; n++) {

      if (cnt == hsize) {
        /* audio samples are stored in a function table */
        ft = csound->FTnp2Find(csound,p->knum);
        tab = ft->ftable;
        size = ft->flen;
       
        if (UNLIKELY((int) ft->nchanls != nchans))
          return csound->PerfError(csound, Str("number of output arguments "
                                               "inconsistent with number of "
                                               "sound file channels"));

       

        /* spos is the reading position in samples, hsize is hopsize,
           time[n] is current read position in secs
           esr is sampling rate
        */
        spos  = hsize*(int)((time[n])*csound->esr/hsize);
      sizefrs = size/nchans;
      while(spos > sizefrs - N) spos -= sizefrs;
      while(spos <= hsize)  spos += sizefrs;
        pos = spos;

      for (j = 0; j < nchans; j++) {

          bwin = (MYFLT *) p->bwin[j].auxp;
          fwin = (MYFLT *) p->fwin[j].auxp;
          prev = (MYFLT *)p->prev[j].auxp;
          framecnt  = (int *)p->framecount[j].auxp;
          outframe= (MYFLT *) p->outframe[j].auxp;
          /* this loop fills two frames/windows with samples from table,
             reading is linearly-interpolated,
             frames are separated by 1 hopsize
          */
          for (i=0; i < N; i++) {
            /* front window, fwin */
            post = (int) pos;
            frac = pos  - post;
            post *= nchans;
            post += j;
            if (post >= 0 && post < size)
              in = tab[post] + frac*(tab[post+nchans] - tab[post]);
            else in =  (MYFLT) 0;
            fwin[i] = in * win[i]; /* window it */
            /* back windo, bwin */
            post = (int) (pos - hsize*pitch);
            post *= nchans;
            post += j;
            if (post >= 0 && post < size)
              in =  tab[post] + frac*(tab[post+nchans] - tab[post]);
            else in =  (MYFLT) 0;
            bwin[i] = in * win[i];  /* window it */
            /* increment read pos according to pitch transposition */
            pos += pitch;
          }

          /* take the FFT of both frames
             re-order Nyquist bin from pos 1 to N
          */
          csound->RealFFT(csound, bwin, N);
          bwin[N] = bwin[1];
          bwin[N+1] = 0.0;
          csound->RealFFT(csound, fwin, N);
          fwin[N] = fwin[1];
          fwin[N+1] = 0.0;

          /* phase vocoder processing */

          for (i=0; i < N + 2; i+=2) {
            /* phases of previous output frame in exponential format,
               obtained by dividing by magnitude */
            div =  FL(1.0)/(HYPOT(prev[i], prev[i+1]) + 1e-20);
            ph_real  =    prev[i]*div;
            ph_im =       prev[i+1]*div;

            /* back window magnitudes, phase differences between
               prev and back windows */
            tmp_real =   bwin[i] * ph_real + bwin[i+1] * ph_im;
            tmp_im =   bwin[i] * ph_im - bwin[i+1] * ph_real;
            bwin[i] = tmp_real;
            bwin[i+1] = tmp_im;
          }

          for (i=0; i < N + 2; i+=2) {
            if (lock) {  /* phase-locking */
              if (i > 0) {
                if (i < N){
                  tmp_real = bwin[i] + bwin[i-2] + bwin[i+2];
                  tmp_im = bwin[i+1] + bwin[i-1] + bwin[i+3];
                }
                else { /* Nyquist */
                  tmp_real = bwin[i] + bwin[i-2];
                  tmp_im = FL(0.0);
                }
              }
              else { /* 0 Hz */
                tmp_real = bwin[i] + bwin[i+2];
                tmp_im = FL(0.0);
              }
            } else { /* no locking */
              tmp_real = bwin[i];
              tmp_im = bwin[i+1];
            }

            tmp_real += 1e-15;
            div =  FL(1.0)/(HYPOT(tmp_real, tmp_im));

            /* phases of tmp frame */
            ph_real = tmp_real*div;
            ph_im = tmp_im*div;

            /* front window mags, phase sum of
               tmp and front windows */
            tmp_real =   fwin[i] * ph_real - fwin[i+1] * ph_im;
            tmp_im =   fwin[i] * ph_im + fwin[i+1] * ph_real;

            /* phase vocoder output */
            prev[i] = fwin[i] = tmp_real;
            prev[i+1] = fwin[i+1] = tmp_im;
          }
          /* re-order bins and take inverse FFT */
          fwin[1] = fwin[N];
          csound->InverseRealFFT(csound, fwin, N);
          /* frame counter */
          framecnt[curframe] = curframe*N;
          /* write to overlapped output frames */
          for (i=0;i<N;i++) outframe[framecnt[curframe]+i] = win[i]*fwin[i];
        
      }

        cnt=0;
        curframe++;
        if (curframe == decim) curframe = 0;
      }
      
      for (j=0; j < nchans; j++) {
        framecnt  = (int *) p->framecount[j].auxp;
        outframe  = (MYFLT *) p->outframe[j].auxp;
        out = p->out[j];
        out[n] = (MYFLT)0;
        /* write output */
        for (i = 0; i < decim; i++) {
          out[n] += outframe[framecnt[i]];
          framecnt[i]++;
        }
        /* scale output */
        out[n] *= amp*(2./3.);
      }
      cnt++;
    }
 
    p->cnt = cnt;
    p->curframe = curframe;
    return OK;

}

static int sinit2(CSOUND *csound, DATASPACE *p)
{
    int size,i;
    sinit(csound, p);
    size = p->N*sizeof(MYFLT);
    for (i=0; i < p->nchans; i++)
    if (p->nwin[i].auxp == NULL || p->nwin[i].size < size)
      csound->AuxAlloc(csound, size, &p->nwin[i]);
    p->pos = *p->offset*csound->esr + p->hsize;
    p->tscale  = 0;
    p->accum = 0;
    return OK;
}

static int sprocess2(CSOUND *csound, DATASPACE *p)
{
    MYFLT pitch = *p->kpitch, time = *p->time, lock = *p->klock;
    MYFLT *out, amp =*p->kamp;
    MYFLT *tab,frac,scale,  dbtresh = *p->dbthresh;
    FUNC *ft;
    int N = p->N, hsize = p->hsize, cnt = p->cnt, sizefrs, nchans = p->nchans;
    int  ksmps = csound->GetKsmps(csound), n;
    int size, post, i, j;
    double pos, spos = p->pos;
    MYFLT *fwin, *bwin;
    MYFLT in, *nwin, *prev;
    MYFLT *win = (MYFLT *) p->win.auxp, *outframe;
    MYFLT powrat;
    MYFLT ph_real, ph_im, tmp_real, tmp_im, div;
    int *framecnt, curframe = p->curframe;
    int decim = p->decim;

    for (n=0; n < ksmps; n++) {
      
      if (cnt == hsize){
        ft = csound->FTnp2Find(csound,p->knum);
        tab = ft->ftable;
        size = ft->flen;

        if (time < 0 || time >= 1 || !*p->konset) {
          spos += hsize*time;
        }
        else if (p->tscale) {
          spos += hsize*(time/(1+p->accum));
          p->accum=0.0;
        }
        else  {
          spos += hsize;
          p->accum++;
          p->tscale = 1;
        }
        if (UNLIKELY((int) ft->nchanls != nchans))
          return csound->PerfError(csound, Str("number of output arguments "
                                               "inconsistent with number of "
                                               "sound file channels"));

        sizefrs = size/nchans;
        while(spos > sizefrs - N - hsize) spos -= sizefrs;
        while(spos <= hsize)  spos += sizefrs;
        pos = spos;

        for (j = 0; j < nchans; j++) {

          bwin = (MYFLT *) p->bwin[j].auxp;
          fwin = (MYFLT *) p->fwin[j].auxp;
          nwin = (MYFLT *) p->nwin[j].auxp;
          prev = (MYFLT *)p->prev[j].auxp;
          framecnt  = (int *)p->framecount[j].auxp;
          outframe= (MYFLT *) p->outframe[j].auxp;
        
          for (i=0; i < N; i++) {
            post = (int) pos;
            frac = pos  - post;
            post *= nchans;
            post += j;
            if (post >= 0 && post < size)
              in = tab[post] + frac*(tab[post+nchans] - tab[post]);
            else in =  (MYFLT) 0;
            fwin[i] = in * win[i]; 
         
            post = (int) (pos - hsize*pitch);
            post *= nchans;
            post += j;
            if (post >= 0 && post < size)
              in =  tab[post] + frac*(tab[post+nchans] - tab[post]);
            else in =  (MYFLT) 0;
            bwin[i] = in * win[i];  
            post = (int) pos + hsize;
            post *= nchans;
            post += j;
            if (post >= 0 && post < size) in =  tab[post];
            else in =  (MYFLT) 0;
            nwin[i] = in * win[i];
            pos += pitch;
          }
       
          csound->RealFFT(csound, bwin, N);
          bwin[N] = bwin[1];
          bwin[N+1] = FL(0.0);
          csound->RealFFT(csound, fwin, N);
          csound->RealFFT(csound, nwin, N);

          tmp_real = tmp_im = (MYFLT) 1e-20;
          for (i=2; i < N; i++) {
            tmp_real += nwin[i]*nwin[i] + nwin[i+1]*nwin[i+1];
            tmp_im += fwin[i]*fwin[i] + fwin[i+1]*fwin[i+1];
          }
          powrat = FL(20.0)*LOG10(tmp_real/tmp_im);
          if (powrat > dbtresh) p->tscale=0;
          /*else tscale=1;*/

          fwin[N] = fwin[1];
          fwin[N+1] = FL(0.0);

          for (i=0; i < N + 2; i+=2) {

            div =  FL(1.0)/(HYPOT(prev[i], prev[i+1]) + 1.0e-20);
            ph_real  =    prev[i]*div;
            ph_im =       prev[i+1]*div;

            tmp_real =   bwin[i] * ph_real + bwin[i+1] * ph_im;
            tmp_im =   bwin[i] * ph_im - bwin[i+1] * ph_real;
            bwin[i] = tmp_real;
            bwin[i+1] = tmp_im;
          }

          for (i=0; i < N + 2; i+=2) {
            if (lock) {
              if (i > 0) {
                if (i < N) {
                  tmp_real = bwin[i] + bwin[i-2] + bwin[i+2];
                  tmp_im = bwin[i+1] + bwin[i-1] + bwin[i+3];
                }
                else {
                  tmp_real = bwin[i] + bwin[i-2];
                  tmp_im = FL(0.0);
                }
              }
              else {
                tmp_real = bwin[i] + bwin[i+2];
                tmp_im = FL(0.0);
              }
            }
            else {
              tmp_real = bwin[i];
              tmp_im = bwin[i+1];
            }

            tmp_real += 1e-15;
            div =  FL(1.0)/(HYPOT(tmp_real, tmp_im));

            ph_real = tmp_real*div;
            ph_im = tmp_im*div;

            tmp_real =   fwin[i] * ph_real - fwin[i+1] * ph_im;
            tmp_im =   fwin[i] * ph_im + fwin[i+1] * ph_real;
          
            prev[i] = fwin[i] = tmp_real;
            prev[i+1] = fwin[i+1] = tmp_im;
          }

          fwin[1] = fwin[N];
          csound->InverseRealFFT(csound, fwin, N);

          framecnt[curframe] = curframe*N;

          for (i=0;i<N;i++) outframe[framecnt[curframe]+i] = win[i]*fwin[i];
        
        }
        cnt=0;
        curframe++;
        if (curframe == decim) curframe = 0;
      }
           
      for (j=0; j < nchans; j++) {
        out = p->out[j];
      framecnt  = (int *) p->framecount[j].auxp;
      outframe  = (MYFLT *) p->outframe[j].auxp;
 
        out[n] = (MYFLT) 0;
     
        for (i = 0; i < decim; i++) {
          out[n] += outframe[framecnt[i]];
          framecnt[i]++;      
        }
        out[n] *= amp*(2./3.);
      }
      cnt++;  
    }
    p->cnt = cnt;
    p->curframe = curframe;
    p->pos = spos;
    return OK;

}

00440 typedef struct {
  OPDS h;
  PVSDAT *fout;
  PVSDAT *fin;
  MYFLT *klock;
  MYFLT  *file;
  uint32 lastframe;
}PVSLOCK;

static int pvslockset(CSOUND *csound, PVSLOCK *p)
{
    int32    N = p->fin->N;

    if (UNLIKELY(p->fin == p->fout))
      csound->Warning(csound, Str("Unsafe to have same fsig as in and out"));
    p->fout->N = N;
    p->fout->overlap = p->fin->overlap;
    p->fout->winsize = p->fin->winsize;
    p->fout->wintype = p->fin->wintype;
    p->fout->format = p->fin->format;
    p->fout->framecount = 1;
    p->lastframe = 0;
   if (p->fout->frame.auxp == NULL ||
            p->fout->frame.size < sizeof(float) * (N + 2))
          csound->AuxAlloc(csound, (N + 2) * sizeof(float), &p->fout->frame);
   if (UNLIKELY(!(p->fout->format == PVS_AMP_FREQ) ))
          return csound->InitError(csound, Str("pvsfreeze: signal format "
                                               "must be amp-freq."));

   return OK;
}


static int pvslockproc(CSOUND *csound, PVSLOCK *p)
{
    int i,n,k,j, maxmagi;
    float mag,  maxmag, *fout = (float *) p->fout->frame.auxp, 
      *fin = (float *) p->fin->frame.auxp;
    int N = p->fin->N;
  
    if (p->lastframe < p->fin->framecount) {
      memcpy(fout,fin, sizeof(float)*(N+2));
  
      if (*p->klock) {
        for (i=2, n = 0; i < N-4; i+=2){
          float p2 = fin[i];
          float p3 = fin[i+2];
          float p1 = fin[i-2];
          float p4 = fin[i+4];
          float p5 = fin[i+6];
          if (p3 > p1 && p3 > p2 && p3 > p4 && p3 > p5) {       
            float freq = fin[i+3], d;
            d = 0.01*freq;
            if (FABS(fout[i-1] - freq) < d)fout[i-1]  = freq;
            if (FABS(fout[i+1] - freq) < d)fout[i+1]  = freq;
            if (FABS(fout[i+3] - freq) < d)fout[i+3]  = freq;
            if (FABS(fout[i+5] - freq) < d)fout[i+5]  = freq;
            if (FABS(fout[i+7] - freq) < d)fout[i+7]  = freq;
          }
        }
      }
      p->fout->framecount = p->lastframe = p->fin->framecount;
    }
    return OK;
}

static OENTRY localops[] = {
  {"mincer", sizeof(DATASPACE), 5, "mm", "akkkkoo",
                                               (SUBR)sinit, NULL,(SUBR)sprocess },
  {"temposcal", sizeof(DATASPACE), 5, "mm", "kkkkkooPOP",
                                               (SUBR)sinit2, NULL,(SUBR)sprocess2 },
  {"pvslock", sizeof(PVSLOCK), 3, "f", "fk", (SUBR) pvslockset,
         (SUBR) pvslockproc},
};

LINKAGE

Generated by  Doxygen 1.6.0   Back to index