/* $Log$ Revision 1.15 2004/06/26 03:50:14 markster Merge source cleanups (bug #1911) Revision 1.14 2003/02/12 13:59:15 matteo mer feb 12 14:56:57 CET 2003 Revision 1.1.1.1 2003/02/12 13:59:15 matteo mer feb 12 14:56:57 CET 2003 Revision 1.2 2000/01/05 08:20:39 markster Some OSS fixes and a few lpc changes to make it actually work * Revision 1.1 1996/08/19 22:32:00 jaf * Initial revision * */ /* -- translated by f2c (version 19951025). You must link the resulting object file with the libraries: -lf2c -lm (in that order) */ #include "f2c.h" #ifdef P_R_O_T_O_T_Y_P_E_S extern int invert_(integer *order, real *phi, real *psi, real *rc); #endif /* **************************************************************** */ /* INVERT Version 45G */ /* $Log$ * Revision 1.15 2004/06/26 03:50:14 markster * Merge source cleanups (bug #1911) * * Revision 1.14 2003/02/12 13:59:15 matteo * mer feb 12 14:56:57 CET 2003 * * Revision 1.1.1.1 2003/02/12 13:59:15 matteo * mer feb 12 14:56:57 CET 2003 * * Revision 1.2 2000/01/05 08:20:39 markster * Some OSS fixes and a few lpc changes to make it actually work * * Revision 1.1 1996/08/19 22:32:00 jaf * Initial revision * */ /* Revision 1.3 1996/03/18 20:52:47 jaf */ /* Just added a few comments about which array indices of the arguments */ /* are used, and mentioning that this subroutine has no local state. */ /* Revision 1.2 1996/03/13 16:51:32 jaf */ /* Comments added explaining that none of the local variables of this */ /* subroutine need to be saved from one invocation to the next. */ /* Eliminated a comment from the original, describing a local array X */ /* that appeared nowhere in the code. */ /* Revision 1.1 1996/02/07 14:47:20 jaf */ /* Initial revision */ /* **************************************************************** */ /* Invert a covariance matrix using Choleski decomposition method. */ /* Input: */ /* ORDER - Analysis order */ /* PHI(ORDER,ORDER) - Covariance matrix */ /* Indices (I,J) read, where ORDER .GE. I .GE. J .GE. 1.*/ /* All other indices untouched. */ /* PSI(ORDER) - Column vector to be predicted */ /* Indices 1 through ORDER read. */ /* Output: */ /* RC(ORDER) - Pseudo reflection coefficients */ /* Indices 1 through ORDER written, and then possibly read. */ /* Internal: */ /* V(ORDER,ORDER) - Temporary matrix */ /* Same indices written as read from PHI. */ /* Many indices may be read and written again after */ /* initially being copied from PHI, but all indices */ /* are written before being read. */ /* NOTE: Temporary matrix V is not needed and may be replaced */ /* by PHI if the original PHI values do not need to be preserved. */ /* Subroutine */ int invert_(integer *order, real *phi, real *psi, real *rc) { /* System generated locals */ unsigned i__2; integer phi_dim1, phi_offset, i__1, i__3; real r__1, r__2; /* Local variables */ real save; integer i__, j, k; real v[100] /* was [10][10] */; /* Arguments */ /* $Log$ * Revision 1.15 2004/06/26 03:50:14 markster * Merge source cleanups (bug #1911) * * Revision 1.14 2003/02/12 13:59:15 matteo * mer feb 12 14:56:57 CET 2003 * * Revision 1.1.1.1 2003/02/12 13:59:15 matteo * mer feb 12 14:56:57 CET 2003 * * Revision 1.2 2000/01/05 08:20:39 markster * Some OSS fixes and a few lpc changes to make it actually work * * Revision 1.1 1996/08/19 22:32:00 jaf * Initial revision * */ /* Revision 1.3 1996/03/29 22:03:47 jaf */ /* Removed definitions for any constants that were no longer used. */ /* Revision 1.2 1996/03/26 19:34:33 jaf */ /* Added comments indicating which constants are not needed in an */ /* application that uses the LPC-10 coder. */ /* Revision 1.1 1996/02/07 14:43:51 jaf */ /* Initial revision */ /* LPC Configuration parameters: */ /* Frame size, Prediction order, Pitch period */ /* Parameters/constants */ /* Local variables that need not be saved */ /* Decompose PHI into V * D * V' where V is a triangular matrix whose */ /* main diagonal elements are all 1, V' is the transpose of V, and */ /* D is a vector. Here D(n) is stored in location V(n,n). */ /* Parameter adjustments */ --rc; --psi; phi_dim1 = *order; phi_offset = phi_dim1 + 1; phi -= phi_offset; /* Function Body */ i__1 = *order; for (j = 1; j <= i__1; ++j) { i__2 = *order; for (i__ = j; i__ <= i__2; ++i__) { v[i__ + j * 10 - 11] = phi[i__ + j * phi_dim1]; } i__2 = j - 1; for (k = 1; k <= i__2; ++k) { save = v[j + k * 10 - 11] * v[k + k * 10 - 11]; i__3 = *order; for (i__ = j; i__ <= i__3; ++i__) { v[i__ + j * 10 - 11] -= v[i__ + k * 10 - 11] * save; } } /* Compute intermediate results, which are similar to RC's */ if ((r__1 = v[j + j * 10 - 11], abs(r__1)) < 1e-10f) { goto L100; } rc[j] = psi[j]; i__2 = j - 1; for (k = 1; k <= i__2; ++k) { rc[j] -= rc[k] * v[j + k * 10 - 11]; } v[j + j * 10 - 11] = 1.f / v[j + j * 10 - 11]; rc[j] *= v[j + j * 10 - 11]; /* Computing MAX */ /* Computing MIN */ r__2 = rc[j]; r__1 = min(r__2,.999f); rc[j] = max(r__1,-.999f); } return 0; /* Zero out higher order RC's if algorithm terminated early */ L100: i__1 = *order; for (i__ = j; i__ <= i__1; ++i__) { rc[i__] = 0.f; } /* Back substitute for PC's (if needed) */ /* 110 DO J = ORDER,1,-1 */ /* PC(J) = RC(J) */ /* DO I = 1,J-1 */ /* PC(J) = PC(J) - PC(I)*V(J,I) */ /* END DO */ /* END DO */ return 0; } /* invert_ */