lpc10/pitsyn.c

   1 /*
   2
   3  * Revision 1.2  1996/08/20  20:40:12  jaf
   4  * Removed all static local variables that were SAVE'd in the Fortran
   5  * code, and put them in struct lpc10_decoder_state that is passed as an
   6  * argument.
   7  *
   8  * Removed init function, since all initialization is now done in
   9  * init_lpc10_decoder_state().
  10  *
  11  * Revision 1.1  1996/08/19  22:31:12  jaf
  12  * Initial revision
  13  *
  14
  15 */
  16
  17 /*  -- translated by f2c (version 19951025).
  18    You must link the resulting object file with the libraries:
  19         -lf2c -lm   (in that order)
  20 */
  21
  22 #include "f2c.h"
  23
  24 extern int pitsyn_(integer *order, integer *voice, integer *pitch, real *rms, real *rc, integer *lframe, integer *ivuv, integer *ipiti, real *rmsi, real *rci, integer *nout, real *ratio, struct lpc10_decoder_state *st);
  25
  26 /* ***************************************************************** */
  27
  28 /*      PITSYN Version 53 */
  29
  30 /*
  31  * Revision 1.2  1996/08/20  20:40:12  jaf
  32  * Removed all static local variables that were SAVE'd in the Fortran
  33  * code, and put them in struct lpc10_decoder_state that is passed as an
  34  * argument.
  35  *
  36  * Removed init function, since all initialization is now done in
  37  * init_lpc10_decoder_state().
  38  *
  39  * Revision 1.1  1996/08/19  22:31:12  jaf
  40  * Initial revision
  41  * */
  42 /* Revision 1.2  1996/03/25  18:49:07  jaf */
  43 /* Added commments about which indices of array arguments are read or */
  44 /* written. */
  45
  46 /* Rearranged local variable declarations to indicate which need to be */
  47 /* saved from one invocation to the next.  Added entry INITPITSYN to */
  48 /* reinitialize local state variables, if desired. */
  49
  50 /* Added lots of comments about proving that the maximum number of pitch */
  51 /* periods (NOUT) that can be returned is 16.  The call to STOP that */
  52 /* could happen if NOUT got too large was removed as a result. */
  53
  54 /* Also proved that the total number of samples returned from N calls, */
  55 /* each with identical values of LFRAME, will always be in the range */
  56 /* N*LFRAME-MAXPIT+1 to N*LFRAME. */
  57
  58 /* Revision 1.1  1996/02/07 14:48:18  jaf */
  59 /* Initial revision */
  60
  61
  62 /* ***************************************************************** */
  63
  64 /*   Synthesize a single pitch epoch */
  65
  66 /* Input: */
  67 /*  ORDER  - Synthesis order (number of RC's) */
  68 /*  VOICE  - Half frame voicing decisions */
  69 /*           Indices 1 through 2 read. */
  70 /*  LFRAME - Length of speech buffer */
  71 /* Input/Output: */
  72 /*  PITCH  - Pitch */
  73 /*           This value should be in the range MINPIT (20) to MAXPIT */
  74 /*           (156), inclusive. */
  75 /*           PITCH can be modified under some conditions. */
  76 /*  RMS    - Energy  (can be modified) */
  77 /*           RMS is changed to 1 if the value passed in is less than 1. */
  78 /*  RC     - Reflection coefficients */
  79 /*           Indices 1 through ORDER can be temporarily overwritten with */
  80 /*           RCO, and then replaced with original values, under some */
  81 /*           conditions. */
  82 /* Output: */
  83 /*  IVUV   - Pitch epoch voicing decisions */
  84 /*           Indices (I) of IVUV, IPITI, and RMSI are written, */
  85 /*           and indices (J,I) of RCI are written, */
  86 /*           where I ranges from 1 to NOUT, and J ranges from 1 to ORDER. */
  87 /*  IPITI  - Pitch epoch length */
  88 /*  RMSI   - Pitch epoch energy */
  89 /*  RCI    - Pitch epoch RC's */
  90 /*  NOUT   - Number of pitch periods in this frame */
  91 /*           This is at least 0, at least 1 if MAXPIT .LT. LFRAME (this */
  92 /*           is currently true on every call), and can never be more than */
  93 /*           (LFRAME+MAXPIT-1)/PITCH, which is currently 16 with */
  94 /*           LFRAME=180, MAXPIT=156, and PITCH .GE. 20, as SYNTHS */
  95 /*           guarantees when it calls this subroutine. */
  96 /*  RATIO  - Previous to present energy ratio */
  97 /*           Always assigned a value. */
  98
  99 /* Subroutine */ int pitsyn_(integer *order, integer *voice,
 100         integer *pitch, real *rms, real *rc, integer *lframe, integer *ivuv,
 101         integer *ipiti, real *rmsi, real *rci, integer *nout, real *ratio,
 102                                struct lpc10_decoder_state *st)
 103 {
 104     /* Initialized data */
 105
 106     real *rmso;
 107     logical *first;
 108
 109     /* System generated locals */
 110     integer rci_dim1 = 0, rci_offset, i__1, i__2;
 111     real r__1;
 112
 113     /* Builtin functions */
 114     double log(doublereal), exp(doublereal);
 115
 116     /* Local variables */
 117     real alrn, alro, yarc[10], prop;
 118     integer i__, j, vflag, jused, lsamp;
 119     integer *jsamp;
 120     real slope;
 121     integer *ipito;
 122     real uvpit;
 123     integer ip, nl, ivoice;
 124     integer *ivoico;
 125     integer istart;
 126     real *rco;
 127     real xxy;
 128
 129 /*       Arguments */
 130
 131 /*   LPC Configuration parameters: */
 132 /* Frame size, Prediction order, Pitch period */
 133 /*       Local variables that need not be saved */
 134 /*       LSAMP is initialized in the IF (FIRST) THEN clause, but it is */
 135 /*       not used the first time through, and it is given a value before
 136 */
 137 /*       use whenever FIRST is .FALSE., so it appears unnecessary to */
 138 /*       assign it a value when FIRST is .TRUE. */
 139 /*       Local state */
 140 /* FIRST  - .TRUE. only on first call to PITSYN. */
 141 /* IVOICO - Previous VOICE(2) value. */
 142 /* IPITO  - Previous PITCH value. */
 143 /* RMSO   - Previous RMS value. */
 144 /* RCO    - Previous RC values. */
 145
 146 /* JSAMP  - If this routine is called N times with identical values of */
 147 /*          LFRAME, then the total length of all pitch periods returned */
 148 /*          is always N*LFRAME-JSAMP, and JSAMP is always in the range 0
 149 */
 150 /*          to MAXPIT-1 (see below for why this is so).  Thus JSAMP is */
 151 /*          the number of samples "left over" from the previous call to */
 152 /*          PITSYN, that haven't been "used" in a pitch period returned */
 153 /*          from this subroutine.  Every time this subroutine is called,
 154 */
 155 /*          it returns pitch periods with a total length of at most */
 156 /*          LFRAME+JSAMP. */
 157
 158 /* IVOICO, IPITO, RCO, and JSAMP need not be assigned an initial value */
 159 /* with a DATA statement, because they are always initialized on the */
 160 /* first call to PITSYN. */
 161
 162 /* FIRST and RMSO should be initialized with DATA statements, because */
 163 /* even on the first call, they are used before being initialized. */
 164     /* Parameter adjustments */
 165     if (rc) {
 166         --rc;
 167         }
 168     if (rci) {
 169         rci_dim1 = *order;
 170         rci_offset = rci_dim1 + 1;
 171         rci -= rci_offset;
 172         }
 173     if (voice) {
 174         --voice;
 175         }
 176     if (ivuv) {
 177         --ivuv;
 178         }
 179     if (ipiti) {
 180         --ipiti;
 181         }
 182     if (rmsi) {
 183         --rmsi;
 184         }
 185
 186     /* Function Body */
 187     ivoico = &(st->ivoico);
 188     ipito = &(st->ipito);
 189     rmso = &(st->rmso);
 190     rco = &(st->rco[0]);
 191     jsamp = &(st->jsamp);
 192     first = &(st->first_pitsyn);
 193
 194     if (*rms < 1.f) {
 195         *rms = 1.f;
 196     }
 197     if (*rmso < 1.f) {
 198         *rmso = 1.f;
 199     }
 200     uvpit = 0.f;
 201     *ratio = *rms / (*rmso + 8.f);
 202     if (*first) {
 203         lsamp = 0;
 204         ivoice = voice[2];
 205         if (ivoice == 0) {
 206             *pitch = *lframe / 4;
 207         }
 208         *nout = *lframe / *pitch;
 209         *jsamp = *lframe - *nout * *pitch;
 210
 211 /*          SYNTHS only calls this subroutine with PITCH in the range
 212 20 */
 213 /*          to 156.  LFRAME = MAXFRM = 180, so NOUT is somewhere in th
 214 e */
 215 /*          range 1 to 9. */
 216
 217 /*          JSAMP is "LFRAME mod PITCH", so it is in the range 0 to */
 218 /*          (PITCH-1), or 0 to MAXPIT-1=155, after the first call. */
 219
 220         i__1 = *nout;
 221         for (i__ = 1; i__ <= i__1; ++i__) {
 222             i__2 = *order;
 223             for (j = 1; j <= i__2; ++j) {
 224                 rci[j + i__ * rci_dim1] = rc[j];
 225             }
 226             ivuv[i__] = ivoice;
 227             ipiti[i__] = *pitch;
 228             rmsi[i__] = *rms;
 229         }
 230         *first = FALSE_;
 231     } else {
 232         vflag = 0;
 233         lsamp = *lframe + *jsamp;
 234         slope = (*pitch - *ipito) / (real) lsamp;
 235         *nout = 0;
 236         jused = 0;
 237         istart = 1;
 238         if (voice[1] == *ivoico && voice[2] == voice[1]) {
 239             if (voice[2] == 0) {
 240 /* SSUV - -   0  ,  0  ,  0 */
 241                 *pitch = *lframe / 4;
 242                 *ipito = *pitch;
 243                 if (*ratio > 8.f) {
 244                     *rmso = *rms;
 245                 }
 246             }
 247 /* SSVC - -   1  ,  1  ,  1 */
 248             slope = (*pitch - *ipito) / (real) lsamp;
 249             ivoice = voice[2];
 250         } else {
 251             if (*ivoico != 1) {
 252                 if (*ivoico == voice[1]) {
 253 /* UV2VC2 - -  0  ,  0  ,  1 */
 254                     nl = lsamp - *lframe / 4;
 255                 } else {
 256 /* UV2VC1 - -  0  ,  1  ,  1 */
 257                     nl = lsamp - *lframe * 3 / 4;
 258                 }
 259                 ipiti[1] = nl / 2;
 260                 ipiti[2] = nl - ipiti[1];
 261                 ivuv[1] = 0;
 262                 ivuv[2] = 0;
 263                 rmsi[1] = *rmso;
 264                 rmsi[2] = *rmso;
 265                 i__1 = *order;
 266                 for (i__ = 1; i__ <= i__1; ++i__) {
 267                     rci[i__ + rci_dim1] = rco[i__ - 1];
 268                     rci[i__ + (rci_dim1 << 1)] = rco[i__ - 1];
 269                     rco[i__ - 1] = rc[i__];
 270                 }
 271                 slope = 0.f;
 272                 *nout = 2;
 273                 *ipito = *pitch;
 274                 jused = nl;
 275                 istart = nl + 1;
 276                 ivoice = 1;
 277             } else {
 278                 if (*ivoico != voice[1]) {
 279 /* VC2UV1 - -   1  ,  0  ,  0 */
 280                     lsamp = *lframe / 4 + *jsamp;
 281                 } else {
 282 /* VC2UV2 - -   1  ,  1  ,  0 */
 283                     lsamp = *lframe * 3 / 4 + *jsamp;
 284                 }
 285                 i__1 = *order;
 286                 for (i__ = 1; i__ <= i__1; ++i__) {
 287                     yarc[i__ - 1] = rc[i__];
 288                     rc[i__] = rco[i__ - 1];
 289                 }
 290                 ivoice = 1;
 291                 slope = 0.f;
 292                 vflag = 1;
 293             }
 294         }
 295 /* Here is the value of most variables that are used below, depending
 296 on */
 297 /* the values of IVOICO, VOICE(1), and VOICE(2).  VOICE(1) and VOICE(2
 298 ) */
 299 /* are input arguments, and IVOICO is the value of VOICE(2) on the */
 300 /* previous call (see notes for the IF (NOUT .NE. 0) statement near th
 301 e */
 302 /* end).  Each of these three values is either 0 or 1.  These three */
 303 /* values below are given as 3-bit long strings, in the order IVOICO,
 304 */
 305 /* VOICE(1), and VOICE(2).  It appears that the code above assumes tha
 306 t */
 307 /* the bit sequences 010 and 101 never occur, but I wonder whether a
 308 */
 309 /* large enough number of bit errors in the channel could cause such a
 310  */
 311 /* thing to happen, and if so, could that cause NOUT to ever go over 1
 312 1? */
 313
 314 /* Note that all of the 180 values in the table are really LFRAME, but
 315  */
 316 /* 180 has fewer characters, and it makes the table a little more */
 317 /* concrete.  If LFRAME is ever changed, keep this in mind.  Similarly
 318 , */
 319 /* 135's are 3*LFRAME/4, and 45's are LFRAME/4.  If LFRAME is not a */
 320 /* multiple of 4, then the 135 for NL-JSAMP is actually LFRAME-LFRAME/
 321 4, */
 322 /* and the 45 for NL-JSAMP is actually LFRAME-3*LFRAME/4. */
 323
 324 /* Note that LSAMP-JSAMP is given as the variable.  This was just for
 325 */
 326 /* brevity, to avoid adding "+JSAMP" to all of the column entries. */
 327 /* Similarly for NL-JSAMP. */
 328
 329 /* Variable    | 000  001    011,010  111       110       100,101 */
 330 /* ------------+-------------------------------------------------- */
 331 /* ISTART      | 1    NL+1   NL+1     1         1         1 */
 332 /* LSAMP-JSAMP | 180  180    180      180       135       45 */
 333 /* IPITO       | 45   PITCH  PITCH    oldPITCH  oldPITCH  oldPITCH */
 334 /* SLOPE       | 0    0      0        seebelow  0         0 */
 335 /* JUSED       | 0    NL     NL       0         0         0 */
 336 /* PITCH       | 45   PITCH  PITCH    PITCH     PITCH     PITCH */
 337 /* NL-JSAMP    | --   135    45       --        --        -- */
 338 /* VFLAG       | 0    0      0        0         1         1 */
 339 /* NOUT        | 0    2      2        0         0         0 */
 340 /* IVOICE      | 0    1      1        1         1         1 */
 341
 342 /* while_loop  | once once   once     once      twice     twice */
 343
 344 /* ISTART      | --   --     --       --        JUSED+1   JUSED+1 */
 345 /* LSAMP-JSAMP | --   --     --       --        180       180 */
 346 /* IPITO       | --   --     --       --        oldPITCH  oldPITCH */
 347 /* SLOPE       | --   --     --       --        0         0 */
 348 /* JUSED       | --   --     --       --        ??        ?? */
 349 /* PITCH       | --   --     --       --        PITCH     PITCH */
 350 /* NL-JSAMP    | --   --     --       --        --        -- */
 351 /* VFLAG       | --   --     --       --        0         0 */
 352 /* NOUT        | --   --     --       --        ??        ?? */
 353 /* IVOICE      | --   --     --       --        0         0 */
 354
 355
 356 /* UVPIT is always 0.0 on the first pass through the DO WHILE (.TRUE.)
 357  */
 358 /* loop below. */
 359
 360 /* The only possible non-0 value of SLOPE (in column 111) is */
 361 /* (PITCH-IPITO)/FLOAT(LSAMP) */
 362
 363 /* Column 101 is identical to 100.  Any good properties we can prove
 364 */
 365 /* for 100 will also hold for 101.  Similarly for 010 and 011. */
 366
 367 /* SYNTHS calls this subroutine with PITCH restricted to the range 20
 368 to */
 369 /* 156.  IPITO is similarly restricted to this range, after the first
 370 */
 371 /* call.  IP below is also restricted to this range, given the */
 372 /* definitions of IPITO, SLOPE, UVPIT, and that I is in the range ISTA
 373 RT */
 374 /* to LSAMP. */
 375
 376         while(TRUE_) {
 377
 378 /*             JUSED is the total length of all pitch periods curr
 379 ently */
 380 /*             in the output arrays, in samples. */
 381
 382 /*             An invariant of the DO I = ISTART,LSAMP loop below,
 383  under */
 384 /*             the condition that IP is always in the range 1 thro
 385 ugh */
 386 /*             MAXPIT, is: */
 387
 388 /*             (I - MAXPIT) .LE. JUSED .LE. (I-1) */
 389
 390 /*             Note that the final value of I is LSAMP+1, so that
 391 after */
 392 /*             the DO loop is complete, we know: */
 393
 394 /*             (LSAMP - MAXPIT + 1) .LE. JUSED .LE. LSAMP */
 395
 396             i__1 = lsamp;
 397             for (i__ = istart; i__ <= i__1; ++i__) {
 398                 r__1 = *ipito + slope * i__;
 399                 ip = r__1 + .5f;
 400                 if (uvpit != 0.f) {
 401                     ip = uvpit;
 402                 }
 403                 if (ip <= i__ - jused) {
 404                     ++(*nout);
 405
 406 /*                   The following check is no longer nece
 407 ssary, now that */
 408 /*                   we can prove that NOUT will never go
 409 over 16. */
 410
 411 /*                  IF (NOUT .GT. 16) STOP 'PITSYN: too many epochs'
 412 */
 413
 414                     ipiti[*nout] = ip;
 415                     *pitch = ip;
 416                     ivuv[*nout] = ivoice;
 417                     jused += ip;
 418                     prop = (jused - ip / 2) / (real) lsamp;
 419                     i__2 = *order;
 420                     for (j = 1; j <= i__2; ++j) {
 421                         alro = log((rco[j - 1] + 1) / (1 - rco[j - 1]));
 422                         alrn = log((rc[j] + 1) / (1 - rc[j]));
 423                         xxy = alro + prop * (alrn - alro);
 424                         xxy = exp(xxy);
 425                         rci[j + *nout * rci_dim1] = (xxy - 1) / (xxy + 1);
 426                     }
 427                     rmsi[*nout] = log(*rmso) + prop * (log(*rms) - log(*rmso));
 428                     rmsi[*nout] = exp(rmsi[*nout]);
 429                 }
 430             }
 431             if (vflag != 1) {
 432                 goto L100;
 433             }
 434
 435 /*             I want to prove what range UVPIT must lie in after
 436 the */
 437 /*             assignments to it below.  To do this, I must determ
 438 ine */
 439 /*             what range (LSAMP-ISTART) must lie in, after the */
 440 /*             assignments to ISTART and LSAMP below. */
 441
 442 /*             Let oldLSAMP be the value of LSAMP at this point in
 443  the */
 444 /*             execution.  This is 135+JSAMP in state 110, or 45+J
 445 SAMP in */
 446 /*             states 100 or 101. */
 447
 448 /*             Given the loop invariant on JUSED above, we know th
 449 at: */
 450
 451 /*             (oldLSAMP - MAXPIT + 1) .LE. JUSED .LE. oldLSAMP */
 452
 453 /*             ISTART is one more than this. */
 454
 455 /*             Let newLSAMP be the value assigned to LSAMP below.
 456  This */
 457 /*             is 180+JSAMP.  Thus (newLSAMP-oldLSAMP) is either 4
 458 5 or */
 459 /*             135, depending on the state. */
 460
 461 /*             Thus, the range of newLSAMP-ISTART is: */
 462
 463 /*             (newLSAMP-(oldLSAMP+1)) .LE. newLSAMP-ISTART */
 464 /*             .LE. (newLSAMP-(oldLSAMP - MAXPIT + 2)) */
 465
 466 /*             or: */
 467
 468 /*             46 .LE. newLSAMP-ISTART .LE. 133+MAXPIT .EQ. 289 */
 469
 470 /*             Therefore, UVPIT is in the range 23 to 144 after th
 471 e first */
 472 /*             assignment to UVPIT below, and after the conditiona
 473 l */
 474 /*             assignment, it is in the range 23 to 90. */
 475
 476 /*             The important thing is that it is in the range 20 t
 477 o 156, */
 478 /*             so that in the loop above, IP is always in this ran
 479 ge. */
 480
 481             vflag = 0;
 482             istart = jused + 1;
 483             lsamp = *lframe + *jsamp;
 484             slope = 0.f;
 485             ivoice = 0;
 486             uvpit = (real) ((lsamp - istart) / 2);
 487             if (uvpit > 90.f) {
 488                 uvpit /= 2;
 489             }
 490             *rmso = *rms;
 491             i__1 = *order;
 492             for (i__ = 1; i__ <= i__1; ++i__) {
 493                 rc[i__] = yarc[i__ - 1];
 494                 rco[i__ - 1] = yarc[i__ - 1];
 495             }
 496         }
 497 L100:
 498         *jsamp = lsamp - jused;
 499     }
 500 /*       Given that the maximum pitch period MAXPIT .LT. LFRAME (this is
 501 */
 502 /*       currently true on every call, since SYNTHS always sets */
 503 /*       LFRAME=180), NOUT will always be .GE. 1 at this point. */
 504     if (*nout != 0) {
 505         *ivoico = voice[2];
 506         *ipito = *pitch;
 507         *rmso = *rms;
 508         i__1 = *order;
 509         for (i__ = 1; i__ <= i__1; ++i__) {
 510             rco[i__ - 1] = rc[i__];
 511         }
 512     }
 513     return 0;
 514 } /* pitsyn_ */