diff options
author | jofret | 2009-04-28 06:57:25 +0000 |
---|---|---|
committer | jofret | 2009-04-28 06:57:25 +0000 |
commit | 5bb72c97b02e1cea0dbf008ec16029491956f971 (patch) | |
tree | a3b38e2ad33af5514091e67f57e8a60e0fe0cef3 /src/c/signalProcessing/lev/dleva2.c | |
parent | 5868a7f86c42335cdad7252ea55117acf7bafe83 (diff) | |
download | scilab2c-5bb72c97b02e1cea0dbf008ec16029491956f971.tar.gz scilab2c-5bb72c97b02e1cea0dbf008ec16029491956f971.tar.bz2 scilab2c-5bb72c97b02e1cea0dbf008ec16029491956f971.zip |
Moving source code
Diffstat (limited to 'src/c/signalProcessing/lev/dleva2.c')
-rw-r--r-- | src/c/signalProcessing/lev/dleva2.c | 68 |
1 files changed, 68 insertions, 0 deletions
diff --git a/src/c/signalProcessing/lev/dleva2.c b/src/c/signalProcessing/lev/dleva2.c new file mode 100644 index 00000000..8b2ab8ae --- /dev/null +++ b/src/c/signalProcessing/lev/dleva2.c @@ -0,0 +1,68 @@ +/* + * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab + * Copyright (C) 2008 - INRIA - Arnaud TORSET + * + * This file must be used under the terms of the CeCILL. + * This source file is licensed as described in the file COPYING, which + * you should have received as part of this distribution. The terms + * are also available at + * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt + * + */ + + +/*Resolve the Yule-Walker equations: + + |r(0) r(1) ... r(N-1)|| a(1) | |sigma2| + |r(1) r(0) ... r(n-1)|| a(2) | | 0 | + | : : ... : || : |=| 0 | + | : : ... : || : | | 0 | + |r(N-1) r(N-2) ... r(0) ||a(N-1)| | 0 | + + using Levinson's algorithm. + r :Correlation coefficients + ar :Auto-Regressive model parameters + sigma2 :Scale constant + rc :Reflection coefficients +*/ + + +#include "lev.h" +#include "malloc.h" +#include "stdlib.h" + + +double dleva2(double* in,int size, double* ar){ + int i=0, j=0; + double accu=0; + double* ak1; + double sigma2; + + /* FIXME : malloc here */ + ak1=(double*)malloc((unsigned int)size*sizeof(double)); + + /* initialize levinson's algorithm */ + ar[0]=-in[1]/in[0]; + sigma2=(1-ar[0]*ar[0])*in[0]; + ak1[0]=0; + + /* iterative solution to yule-walker equations */ + for (i=1;i<size-1;i++){ + accu=0; + for (j=0;j<i;j++){ + accu+=ar[j]*in[i-j]; + } + ak1[i]=-(in[i+1]+accu)/(sigma2); + for (j=0;j<i;j++){ + ak1[j] = ar[j]+ak1[i]*ar[i-1-j]; + } + sigma2=(1-ak1[i]*ak1[i])*(sigma2); + for (j=0;j<=i;j++){ + ar[j]=ak1[j]; + } + } + free(ak1); + return sigma2; +} + + |