diff options
Diffstat (limited to '2.3-1/src/c/signalProcessing')
-rw-r--r-- | 2.3-1/src/c/signalProcessing/faurre/dfaurrea.c | 112 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/hank/dhanka.c | 4 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/includes/faurre.h | 29 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/includes/kalm.h | 30 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/includes/wiener.h | 30 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/interfaces/int_faurre.h | 20 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/interfaces/int_kalm.h | 20 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/interfaces/int_weiner.h | 20 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/interfaces/int_wiener.h | 20 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/kalm/dkalma.c | 118 | ||||
-rw-r--r-- | 2.3-1/src/c/signalProcessing/wiener/dwienera.c | 166 |
11 files changed, 567 insertions, 2 deletions
diff --git a/2.3-1/src/c/signalProcessing/faurre/dfaurrea.c b/2.3-1/src/c/signalProcessing/faurre/dfaurrea.c new file mode 100644 index 00000000..04d5a697 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/faurre/dfaurrea.c @@ -0,0 +1,112 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Funtion to find faurre */ + +#include "lapack.h" +#include "stdlib.h" +#include "string.h" +#include "faurre.h" +#include "matrixTranspose.h" +#include "matrixMultiplication.h" +#include "matrixInversion.h" +#include "addition.h" +#include "subtraction.h" +#include "eye.h" +#include "matrixDivision.h" + +void dfaurrea(int n, double* H, int H_row, int H_col, double* F, int F_row, int F_col, double* G, int G_row, int G_col, double * R0, int R0_row, int R0_col, double* P, double *R , double* T ) +{ + +/*Variable Declaration Start*/ + double R0_inv[R0_row* R0_col]; + double temp1[G_row* R0_col]; + double Pn[G_row* G_row]; + double Pn_temp[G_row*G_row]; + double G_trans[G_col*G_row]; + double H_trans[H_col*H_row]; + double F_trans[F_col*F_row]; + double A_trans[G_col*G_row]; + double temp2[F_row*G_row]; + double temp3[F_row*H_row]; + double A[G_row*G_col]; + double temp4[H_row*G_row]; + double temp5[H_row*H_row]; + double temp6[R0_row*R0_col]; + double temp7[G_row*R0_row]; + double temp8[G_row*G_row]; + double temp9[F_row*F_col]; + double temp10[H_row*G_row]; + double temp11[H_row*H_row]; + double temp12[F_row*G_row]; + double temp13[F_row*H_row]; + double temp14[G_row*G_col]; +/*Variable Declaration End*/ + + dinverma(R0, R0_inv, R0_row); //R0_inv= inv(R0) + dmulma(G, G_row, G_col, R0_inv, R0_row, R0_col, temp1); //temp1= G* inv(R0) + dtransposea(G, G_row, G_col, G_trans); //G_trans= G' + dmulma(temp1, G_row, R0_col, G_trans, G_col, G_row, Pn ); // Pn= G* inv(R0)* G' + dtransposea(H, H_row, H_col, H_trans); //H_trans= H' + dtransposea(F, F_row, F_col, F_trans); //F_trans= F' + +/*//Debug Only*/ +/* for(int i=0; i< G_row*G_row; i++)*/ +/* printf("Before loop Pn[%d]= %lf \t", i, Pn[i]);*/ +/* printf("\n");*/ + + for(int i=1;i <= n ; i++) + { + dmulma(F, F_row, F_col, Pn, G_row, G_row, temp2); //temp2= F*Pn; + dmulma(temp2, F_row, G_row, H_trans, H_col, H_row, temp3); //temp3= F*Pn*H' + ddiffa( G, G_row*G_col, temp3, F_row*H_row, A); //A = G- F*Pn*H' + dtransposea(A, G_row, G_col, A_trans); //A_trans= A' + dmulma(H, H_row, H_col, Pn, G_row, G_row, temp4); //temp4= H*Pn + dmulma(temp4, H_row, G_row, H_trans, H_col, H_row, temp5 ); //temp5= H*Pn*H' + ddiffa(R0, R0_row*R0_col, temp5, H_row*H_row,temp6 ); //temp6= R0-H*Pn*H' + drdivma(A,G_row, G_col, temp6, R0_row, R0_col, temp7); //temp7= A / (R0-H*Pn*H') + dmulma(temp7, G_row, R0_row, A_trans, G_col, G_row, temp8); //temp8 = A / (R0-H*Pn*H')*A' + dmulma(temp2, F_row, G_row, F_trans, F_col, F_row, temp9); //temp9= F*Pn*F' + dadda(temp9, F_row*F_col,temp8,G_row*G_row, Pn_temp); //Pn_temp = F*Pn*F'+ A / (R0-H*Pn*H')*A' + for(int j=0; j< G_row*G_row; j++) + Pn[j]= Pn_temp[j]; + + } + +/* //Debug Only*/ +/* for(int i=0; i< G_row*G_row; i++)*/ +/* printf("Afterloop Pn[%d]= %lf \t", i, Pn[i]);*/ +/* printf("\n");*/ + + for(int j=0; j< G_row*G_row; j++) + P[j]= Pn[j]; + + dmulma(H, H_row, H_col, P, G_row, G_row, temp10); //temp10= H*P + dmulma(temp10, H_row, G_row, H_trans, H_col, H_row, temp11 ); //temp11= H*P*H' + ddiffa(R0, R0_row*R0_col, temp11, H_row*H_row, R ); //R= R0-H*P*H' + +/* //Debug Only*/ +/* for(int i=0; i< H_row*H_row; i++)*/ +/* printf("R[%d]= %lf \t", i, R[i]);*/ +/* printf("\n");*/ + + dmulma(F, F_row, F_col, P, G_row, G_row, temp12); //temp12= F*P; + dmulma(temp12, F_row, G_row, H_trans, H_col, H_row, temp13); //temp13= F*P*H' + ddiffa( G, G_row*G_col, temp13, F_row*H_row, temp14); //temp14 = G- F*P*H' + drdivma(temp14,G_row, G_col, R, R0_row, R0_col, T); //T= (G-F*P*H')/ R + +/* //Debug Only*/ +/* for(int i=0; i< G_row*R0_row; i++)*/ +/* printf("T[%d]= %lf \t", i, T[i]);*/ +/* printf("\n");*/ + +} diff --git a/2.3-1/src/c/signalProcessing/hank/dhanka.c b/2.3-1/src/c/signalProcessing/hank/dhanka.c index 37008fd6..5cef7f3b 100644 --- a/2.3-1/src/c/signalProcessing/hank/dhanka.c +++ b/2.3-1/src/c/signalProcessing/hank/dhanka.c @@ -43,7 +43,7 @@ for(int i=1; i< mr+1; i= i+row) { temp1[j]=i; j++; - printf("%d", i); + //printf("%d", i); } for(int j=0;j< nr; j++) temp2[j]=j; @@ -65,7 +65,7 @@ for(int i=0; i< m*nr; i++) } } - +//Debug Only /* printf("ones1\n"); for(int k=0; k<1*nr; k++) diff --git a/2.3-1/src/c/signalProcessing/includes/faurre.h b/2.3-1/src/c/signalProcessing/includes/faurre.h new file mode 100644 index 00000000..e4a2de3c --- /dev/null +++ b/2.3-1/src/c/signalProcessing/includes/faurre.h @@ -0,0 +1,29 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __FAURRE_H__ +#define __FAURRE_H__ +#include "types.h" + +#ifdef __cplusplus +extern "C" { +#endif + +void dfaurrea(int n, double* H, int H_row, int H_col, double* F, int F_row, int F_col, double* G, int G_row, int G_col, double * R0, int R0_row, int R0_col, double * P, double* R, double* T); + + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /* __FAURRE_H__ */ + diff --git a/2.3-1/src/c/signalProcessing/includes/kalm.h b/2.3-1/src/c/signalProcessing/includes/kalm.h new file mode 100644 index 00000000..6f73ac15 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/includes/kalm.h @@ -0,0 +1,30 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __KALM_H__ +#define __KALM_H__ +#include "types.h" + +#ifdef __cplusplus +extern "C" { +#endif + +void dkalma(double* y, int y_row, int y_col, double* x0, int x0_row, int x0_col, double* p0, int p0_row, int p0_col, double* f, int f_row, int f_col, double* g, int g_row, int g_col, double* h, int h_row, int h_col, double* q, int q_row, int q_col, double* r, int r_row, int r_col, double* x1, double* p1, double* x, double* p); + + + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /* __KALM_H__ */ + diff --git a/2.3-1/src/c/signalProcessing/includes/wiener.h b/2.3-1/src/c/signalProcessing/includes/wiener.h new file mode 100644 index 00000000..ba775266 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/includes/wiener.h @@ -0,0 +1,30 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __WEINER_H__ +#define __WEINER_H__ +#include "types.h" + +#ifdef __cplusplus +extern "C" { +#endif + +void dwienera(double* y, int y_row, int y_col, double* x0, int x0_row, int x0_col, double* p0, int p0_row, int p0_col, double* f, int f_row, int f_col, double* g, int g_row, int g_col, double* h, int h_row, int h_col, double* q, int q_row, int q_col, double* r, int r_row, int r_col, double* xs, double* ps, double* xf, double* pf); + + + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /* __WEINER_H__ */ + diff --git a/2.3-1/src/c/signalProcessing/interfaces/int_faurre.h b/2.3-1/src/c/signalProcessing/interfaces/int_faurre.h new file mode 100644 index 00000000..e527a639 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/interfaces/int_faurre.h @@ -0,0 +1,20 @@ +/* Copyright (C) 2016 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __INT_FAURRE_H__ +#define __INT_FAURRE_H__ + +#define d0d2d2d2d2faurred2d2d2(in1,in2, size2, in3, size3, in4, size4, in5, size5,out1, out2, out3) \ + dfaurrea(in1,in2, size2[0], size2[1],in3, size3[0], size3[1],in4, size4[0], size4[1],in5, size5[0], size5[1], out1, out2, out3) + + +#endif /* !INT_FAURRE_H__! */ diff --git a/2.3-1/src/c/signalProcessing/interfaces/int_kalm.h b/2.3-1/src/c/signalProcessing/interfaces/int_kalm.h new file mode 100644 index 00000000..b4369c90 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/interfaces/int_kalm.h @@ -0,0 +1,20 @@ +/* Copyright (C) 2016 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __INT_KALM_H__ +#define __INT_KALM_H__ + +#define d2d2d2d2d2d2d2d2kalmd2d2d2d2(in1, size1, in2, size2, in3, size3, in4, size4, in5, size5, in6, size6, in7, size7, in8, size8, out1, out2, out3, out4) \ + dkalma(in1, size1[0], size1[1],in2, size2[0], size2[1],in3, size3[0], size3[1],in4, size4[0], size4[1],in5, size5[0], size5[1] ,\ + in6, size6[0], size6[1],in7, size7[0], size7[1],in8, size8[0], size8[1], out1, out2, out3, out4) + +#endif /* !INT_KALM_H__! */ diff --git a/2.3-1/src/c/signalProcessing/interfaces/int_weiner.h b/2.3-1/src/c/signalProcessing/interfaces/int_weiner.h new file mode 100644 index 00000000..daefcf03 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/interfaces/int_weiner.h @@ -0,0 +1,20 @@ +/* Copyright (C) 2016 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __INT_WIENER_H__ +#define __INT_WIENER_H__ + +#define d2d2d2d2d2d2d2d2wienerd2d2d2d2(in1, size1, in2, size2, in3, size3, in4, size4, in5, size5, in6, size6, in7, size7, in8, size8, out1, out2, out3, out4) \ + dwienera(in1, size1[0], size1[1],in2, size2[0], size2[1],in3, size3[0], size3[1],in4, size4[0], size4[1],in5, size5[0], size5[1] ,\ + in6, size6[0], size6[1],in7, size7[0], size7[1],in8, size8[0], size8[1], out1, out2, out3, out4) + +#endif /* !INT_WEINER_H__! */ diff --git a/2.3-1/src/c/signalProcessing/interfaces/int_wiener.h b/2.3-1/src/c/signalProcessing/interfaces/int_wiener.h new file mode 100644 index 00000000..daefcf03 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/interfaces/int_wiener.h @@ -0,0 +1,20 @@ +/* Copyright (C) 2016 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +#ifndef __INT_WIENER_H__ +#define __INT_WIENER_H__ + +#define d2d2d2d2d2d2d2d2wienerd2d2d2d2(in1, size1, in2, size2, in3, size3, in4, size4, in5, size5, in6, size6, in7, size7, in8, size8, out1, out2, out3, out4) \ + dwienera(in1, size1[0], size1[1],in2, size2[0], size2[1],in3, size3[0], size3[1],in4, size4[0], size4[1],in5, size5[0], size5[1] ,\ + in6, size6[0], size6[1],in7, size7[0], size7[1],in8, size8[0], size8[1], out1, out2, out3, out4) + +#endif /* !INT_WEINER_H__! */ diff --git a/2.3-1/src/c/signalProcessing/kalm/dkalma.c b/2.3-1/src/c/signalProcessing/kalm/dkalma.c new file mode 100644 index 00000000..d71a86a9 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/kalm/dkalma.c @@ -0,0 +1,118 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Function to find kalm */ + +#include "lapack.h" +#include "stdlib.h" +#include "stdio.h" +#include "string.h" +#include "kalm.h" +#include "matrixTranspose.h" +#include "matrixMultiplication.h" +#include "matrixInversion.h" +#include "addition.h" +#include "subtraction.h" +#include "eye.h" + +void dkalma(double* y, int y_row, int y_col, double* x0, int x0_row, int x0_col, double* p0, int p0_row, int p0_col, double* f, int f_row, int f_col, double* g, int g_row, int g_col, double* h, int h_row, int h_col, double* q, int q_row, int q_col, double* r, int r_row, int r_col, double* x1, double* p1, double* x, double* p) + +{ + double k[p0_row*h_row]; + double h_trans[h_col*h_row]; + double temp1[p0_row*h_row]; + double temp2[h_row*p0_col]; + double temp3[h_row*h_row]; + double temp4[h_row*h_row]; + double temp5[h_row*h_row]; + double temp6[p0_row*p0_col]; + double temp7[p0_row*h_col]; + double temp8[p0_row*p0_col]; + double f_trans[f_col*f_row]; + double g_trans[g_col*g_row]; + double temp9 [f_row*p0_col]; + double temp10[f_row*f_row]; + double temp11 [g_row*q_col]; + double temp12[g_row*g_row]; + double temp13[h_row*x0_col]; + double temp14[h_row*x0_col]; + double temp15[p0_row*x0_col]; + + + + dtransposea(h, h_row, h_col, h_trans); + dmulma(p0, p0_row, p0_col, h_trans, h_col, h_row, temp1); //temp1= p0*h' + dmulma(h, h_row, h_col, p0, p0_row, p0_col, temp2); //temp2= h*p0 + dmulma(temp2, h_row, p0_col, h_trans, h_col, h_row, temp3); //temp3= h*p0*h' + dadda(temp3, h_row*h_row, r, r_row*r_col, temp4); //temp4= h*p0*h'+r + +/*//Debug Only*/ +/* for(int i=0; i< h_row*h_row; i++)*/ +/* printf("temp4[%d] = %lf \t",i, temp4[i]);*/ +/* printf("\n");*/ + + dinverma(temp4, temp5, h_row); //temp5= (h*p0*h'+r)**(-1) + + dmulma(temp1, p0_row, h_row, temp5, h_row, h_row, k); //k=p0*h'*(h*p0*h'+r)**(-1) + +/*//Debug Only*/ +/* for(int i=0; i< p0_row*h_row; i++)*/ +/* printf("k[%d] = %lf \t",i, k[i]);*/ +/* printf("\n");*/ + + + deyea(temp6, p0_row, p0_col); //temp6 = eye(p0) + dmulma(k, p0_row, h_row, h, h_row, h_col, temp7); //temp7 = k*h + ddiffa(temp6, p0_row*p0_col, temp7, p0_row*h_col, temp8); //temp8= eye(p0)- k*h + dmulma(temp8, p0_row, p0_col, p0, p0_row, p0_col, p); //p=(eye(p0)-k*h)*p0 + +/*//Debug Only*/ +/* for(int i=0; i< p0_row*p0_col; i++)*/ +/* printf("p[%d] = %lf \t",i, p[i]);*/ +/* printf("\n");*/ + + + dtransposea(g, g_row, g_col, g_trans); + dtransposea(f, f_row, f_col, f_trans); + dmulma(f, f_row, f_col, p, p0_row, p0_col, temp9); //temp9= f*p + dmulma(temp9, f_row, p0_col, f_trans, f_col, f_row, temp10); //temp10= f*p*f' + dmulma(g, g_row, g_col, q, q_row, q_col, temp11); //temp11= g*q + dmulma(temp11, g_row, q_col, g_trans, g_col, g_row, temp12); //temp12= g*q*g' + dadda(temp10, f_row*f_row, temp12, g_row*g_row, p1); //p1=f*p*f'+g*q*g' + +/* //Debug Only*/ +/* for(int i=0; i< f_row*f_row; i++)*/ +/* printf("p1[%d] = %lf \t",i, p1[i]);*/ +/* printf("\n");*/ + + + dmulma(h, h_row, h_col, x0, x0_row, x0_col, temp13); //temp13= h*x0 + ddiffa(y, y_row*y_col, temp13, h_row*x0_col, temp14); //temp14= y-h*x0 + dmulma(k, p0_row, h_row, temp14, h_row, x0_col, temp15); //temp15= k*(y-h*x0) + dadda(x0, x0_row*x0_col, temp15, p0_row*x0_col, x); //x=x0+k*(y-h*x0) + +/* //Debug Only*/ +/* for(int i=0; i< x0_row*x0_col; i++)*/ +/* printf("x[%d] = %lf \t",i, x[i]);*/ +/* printf("\n");*/ + + + dmulma(f, f_row, f_col, x, x0_row, x0_col, x1) ; //x1= f*x + +/* //Debug Only*/ +/* for(int i=0; i< f_row*x0_col; i++)*/ +/* printf("x1[%d] = %lf \t",i, x1[i]);*/ +/* printf("\n");*/ + + + +} diff --git a/2.3-1/src/c/signalProcessing/wiener/dwienera.c b/2.3-1/src/c/signalProcessing/wiener/dwienera.c new file mode 100644 index 00000000..18f49bc9 --- /dev/null +++ b/2.3-1/src/c/signalProcessing/wiener/dwienera.c @@ -0,0 +1,166 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Abhinav Dronamraju + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Function to find kalm */ + +#include "lapack.h" +#include "stdlib.h" +#include "stdio.h" +#include "string.h" +#include "wiener.h" +#include "kalm.h" +#include "matrixTranspose.h" +#include "matrixMultiplication.h" +#include "matrixInversion.h" +#include "addition.h" +#include "subtraction.h" +#include "eye.h" + +void dwienera(double* y, int y_row, int y_col, double* x0, int x0_row, int x0_col, double* p0, int p0_row, int p0_col, double* f, int f_row, int f_col, double* g, int g_row, int g_col, double* h, int h_row, int h_col, double* q, int q_row, int q_col, double* r, int r_row, int r_col, double* xs, double* ps, double* xf, double* pf) + +{ + +/* dtransposea(h, h_row, h_col, h_trans);*/ +/* dmulma(p0, p0_row, p0_col, h_trans, h_col, h_row, temp1); //temp1= p0*h'*/ +/* dmulma(h, h_row, h_col, p0, p0_row, p0_col, temp2); //temp2= h*p0*/ +/* dmulma(temp2, h_row, p0_col, h_trans, h_col, h_row, temp3); //temp3= h*p0*h'*/ +/* dadda(temp3, h_row*h_row, r, r_row*r_col, temp4); //temp4= h*p0*h'+r*/ + + int n= x0_row, x0j= x0_col, m= y_row, tf= y_col, to=1, k; + double ind_nk[n]; + double ind_mk[m]; + double xf1[x0_row*x0_col]; + double pf1[p0_row*p0_col]; + double yk[y_row*1]; + double fk[f_row*n]; + double gk[g_row*n]; + double hk[h_row*n]; + double qk[q_row*n]; + double rk[r_row*m]; + double x1[f_row*x0_col]; + double p1[f_row*f_row]; + double x[x0_row*x0_col]; + double p[p0_row*p0_col]; + + +for(int i=0; i< x0_row*x0_col; i++) + xf1[i]= x0[i]; + +for(int i=0; i< p0_row*p0_col; i++) + pf1[i]= p0[i]; + +for(k=t0; k<=tf; k++) +{ + int j=0; + for(int i=1+(k-1)*n; i<= k*n; i++) + { + ind_nk[j]= i; + j++; + } + + int j=0; + for(int i=1+(k-1)*m; i<= k*m; i++) + { + ind_mk[j]= i; + j++; + } + + for(int i=0; i< y_row; y++) + { + int j=k-1; + + yk[i]= y[i+j*y_row]; + + } +int l=0; + for(int i=0; i<n; i++) + { + for(int j=0; j< f_row; j++) + { + int k=ind_nk[i]-1; + + fk[l]= f[j+k*y_row]; + l++; + } + + + } + +int l=0; + for(int i=0; i<n; i++) + { + for(int j=0; j< g_row; j++) + { + int k=ind_nk[i]-1; + + gk[l]= g[j+k*y_row]; + l++; + } + + + } + +int l=0; + for(int i=0; i<n; i++) + { + for(int j=0; j< h_row; j++) + { + int k=ind_nk[i]-1; + + hk[l]= h[j+k*y_row]; + l++; + } + + + } + +int l=0; + for(int i=0; i<n; i++) + { + for(int j=0; j< q_row; j++) + { + int k=ind_nk[i]-1; + + qk[l]= q[j+k*y_row]; + l++; + } + + + } + +int l=0; + for(int i=0; i<m; i++) + { + for(int j=0; j< r_row; j++) + { + int k=ind_mk[i]-1; + + rk[l]= r[j+k*y_row]; + l++; + } + + + } + + +dkalma(yk, y_row, 1, x0, x0_row, x0_col, p0, p0_row, p0_col, fk, f_row, n, gk, g_row, n, hk, h_row, n, qk, q_row, n, rk, r_row, m, x1, p1, x, p); + + + + +} + + + + + +} |