Merge branch 'uk_branch'

11 files changed, 627 insertions, 81 deletions

diff --git a/2.3-1/src/c/linearAlgebra/includes/sva.h b/2.3-1/src/c/linearAlgebra/includes/sva.h
new file mode 100644
index 00000000..ea628a32
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/includes/sva.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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+ */
+
+#ifndef __SVA_H__
+#define __SVA_H__
+#include "types.h"
+
+#ifdef  __cplusplus
+extern "C" {
+#endif
+
+void dsvaa(int ninp,double *in1,int row,int col,double in2,double *out1, \
+	double *out2,double *out3);	
+
+#ifdef  __cplusplus 
+} /* extern "C" */
+#endif
+
+#endif /*__SVA_H__*/
+
diff --git a/2.3-1/src/c/linearAlgebra/includes/svd.h b/2.3-1/src/c/linearAlgebra/includes/svd.h
new file mode 100644
index 00000000..dea681fc
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/includes/svd.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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+ */
+#ifndef __SVD_H__
+#define __SVD_H__
+#include "types.h"
+#include "doubleComplex.h"
+
+#ifdef  __cplusplus
+extern "C" {
+#endif
+
+double dsvda(double tol,double *in1,int row,int col,double in2,double nout,double *out1, \
+	double *out2,double *out3);
+void zsvda(doubleComplex *in1,int row,int col,int in2,int nout, doubleComplex *out1,\
+		doubleComplex *out2,doubleComplex *out3);	
+
+#ifdef  __cplusplus 
+} /* extern "C" */
+#endif
+
+#endif /*__SVD_H__*/
diff --git a/2.3-1/src/c/linearAlgebra/interfaces/int_sva.h b/2.3-1/src/c/linearAlgebra/interfaces/int_sva.h
new file mode 100644
index 00000000..f1f8260a
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/interfaces/int_sva.h
@@ -0,0 +1,28 @@
+ /* 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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+ */
+
+#ifndef __INT_SVA_H__
+#define __INT_SVA_H__
+
+#ifdef  __cplusplus
+extern "C" {
+#endif
+
+#define d2svad2d2d2(in1,size,out1,out2,out3) dsvaa(1,in1,size[0],size[1],0,out1,out2,out3);
+#define d2d0svad2d2d2(in1,size1,in2,out1,out2,out3) dsvaa(2,in1,size1[0],size1[1],in2,out1,out2,out3);
+
+#ifdef  __cplusplus
+} /* extern "C" */
+#endif
+
+#endif /*__INT_SVA_H__*/
+
diff --git a/2.3-1/src/c/linearAlgebra/interfaces/int_svd.h b/2.3-1/src/c/linearAlgebra/interfaces/int_svd.h
new file mode 100644
index 00000000..8f40bffe
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/interfaces/int_svd.h
@@ -0,0 +1,35 @@
+ /* 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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+ */
+
+#ifndef __INT_SVD_H__
+#define __INT_SVD_H__
+
+#ifdef  __cplusplus
+extern "C" {
+#endif
+
+#define d2svdd2(in1,size1,out1) dsvda(0,in1,size1[0],size1[1],0,1,out1,NULL,NULL)
+#define d2g2svdd2d2d2(in1,size1,in2,size2,out1,out2,out3) dsvda(0,in1,size1[0],size1[1],1,3,out1,out2,out3);
+#define d2svdd2d2d2(in1,size1,out1,out2,out3) dsvda(0,in1,size1[0],size1[1],0,3,out1,out2,out3);
+
+#define d2svdd2d2d2d0(in1,size1,out1,out2,out3) dsvda(0,in1,size1[0],size1[1],0,4,out1,out2,out3);
+#define d2d0svdd2d2d2d0(in1,size1,tol,out1,out2,out3) dsvda(tol,in1,size1[0],size1[1],0,4,out1,out2,out3);
+
+#define z2svdz2(in1,size1,out2) zsvda(in1,size1[0],size1[1],0,1,NULL,out2,NULL);
+#define z2g2svdz2z2z2(in1,size1,in2,size2,out1,out2,out3) zsvda(in1,size1[0],size1[1],1,3,out1,out2,out3);
+#define z2svdz2z2z2(in1,size1,out1,out2,out3) zsvda(in1,size1[0],size1[1],0,3,out1,out2,out3);
+
+#ifdef  __cplusplus
+} /* extern "C" */
+#endif
+
+#endif /*__INT_SVD_H__*/
diff --git a/2.3-1/src/c/linearAlgebra/sva/dsvaa.c b/2.3-1/src/c/linearAlgebra/sva/dsvaa.c
new file mode 100644
index 00000000..b7d07d8c
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/sva/dsvaa.c
@@ -0,0 +1,92 @@
+/* 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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+
+ */
+#include "sva.h"
+#include "svd.h"
+#include "lapack.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include "string.h"
+#include "matrixTranspose.h"
+
+#define eps 2.22044604925e-16
+
+void dsvaa(int ninp,double *in1,int row,int col,double in2,double *out1, \
+	double *out2,double *out3){
+
+	double tol;
+	double rk=0;
+	int N=col,M=row;
+	if(row == 0 && col == 0) return;
+	int i,j;
+	int arow;		/*Actual row of given matrix*/
+	int acol;		/*Actual col of given matrix*/
+
+	/* Calculation of svd of a given matrix */
+	double *U,*S,*V;
+	U = (double *)malloc((double)row*min(row,col)*sizeof(double));
+	S = (double *)malloc((double)min(row,col)*min(row,col)*sizeof(double));
+	V = (double *)malloc((double)col*min(row,col)*sizeof(double));
+		
+	dsvda(0,in1,M,N,1,3,U,S,V);
+
+	if (ninp == 1){ 		/* [u,s,v] = sva(A) when input is only matrix */
+		tol = max(row,col)*S[0]*eps;
+		rk = 0;		
+		for(i=0;i<col;i++){
+			if(S[i+i*row] > tol){
+				rk+=1;			
+			}		
+		}
+	}
+	else{				/*[u,s,v] = sva(A,tol) when two input's are there */
+		tol = in2;
+		if(tol > 1){
+			rk = tol;
+			if(rk > min(row,col)){
+				printf("ERROR:  Wrong value for input argument !");	
+				out1 = NULL;
+				out2 = NULL;
+				out3 = NULL;
+				return;	
+			}		
+		}
+		else{
+			rk = 0;
+			for(i=0;i<col;i++){
+				if(S[i+i*row] > tol){
+					rk+=1;
+				}
+			}	
+		}	
+	}
+	arow = M;
+	acol = min(M,N);	
+	for(i=0;i<arow;i++){
+		for(j=0;j<rk;j++){
+			out1[i+j*row]=U[i+j*arow];		
+		}	
+	}
+	arow = min(M,N);
+	for(i=0;i<rk;i++){
+		for(j=0;j<rk;j++){
+			out2[i+j*(int)rk] = S[i+j*arow];	
+		}			
+	}
+	arow = N;
+	acol = min(M,N);
+	for(i=0;i<arow;i++){
+		for(j=0;j<rk;j++){
+			out3[i+j*arow] = V[i+j*arow];		
+		}
+	}	
+}
diff --git a/2.3-1/src/c/linearAlgebra/svd/dsvda.c b/2.3-1/src/c/linearAlgebra/svd/dsvda.c
new file mode 100644
index 00000000..c3bcfc29
--- /dev/null
+++ b/2.3-1/src/c/linearAlgebra/svd/dsvda.c
@@ -0,0 +1,184 @@
+/* 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: Sandeep Gupta
+ Organization: FOSSEE, IIT Bombay
+ Email: toolbox@scilab.in
+
+ */
+
+/*Funtion to find singular value decomposition of given matrix */
+
+#include "lapack.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include "string.h"
+#include <math.h>
+#include "svd.h"
+#include "matrixTranspose.h"
+
+int min(int a,int b);
+int max(int a,int b);
+
+extern double dgesvd_(char*,char*,int*,int*,double*,int*,double*,double*,int*,\
+	double*,int*,double *,int*,int*);
+
+#define eps 2.22044604925e-16 /* pow(2,-52) */
+
+/*  DGESVD computes the singular value decomposition (SVD) of a real
+ M-by-N matrix A, optionally computing the left and/or right singular
+ vectors. The SVD is written
+
+      A = U * SIGMA * transpose(V) */
+
+/*Function support -
+
+s=svd(X)
+[U,S,V]=svd(X)
+[U,S,V]=svd(X,0) (obsolete)
+[U,S,V]=svd(X,"e")
+[U,S,V,rk]=svd(X [,tol])
+
+*/
+
+double dsvda(double tol,double *in1,int row,int col,double in2,double nout,double *out1, \
+	double *out2,double *out3){
+	
+	char JOBU,JOBVT;
+	int i,j,k;
+	int LDU=1;			/*Leading Dimension of U */
+	int LDVT=1;			/*Leading Dimension of VT */
+	int M = row;
+	int N = col;
+	double *buf;
+	double *S,*U,*VT;
+	double *WORK;	
+	
+	int rk; /*Fourth output if needed */
+
+	/*if((nout > 1 && in2 == 1) && (M != N)){   // [U,S,VT] = svd(x,'e')
+		if(M > N){		
+			JOBU = 'S';
+			JOBVT = 'A';
+			LDVT = N;
+		}
+		else{
+			JOBU = 'A';
+			JOBVT = 'S';		
+			LDVT = min(M,N);	
+		}
+		LDU = M; 		
+		U = (double*) malloc((double) (LDU)*min(M,N)*sizeof(double));
+		VT = (double*) malloc((double) (LDVT)*N*sizeof(double));
+	}
+	else */if(nout > 1){		/* [U,S,VT = svd(x)] */
+		JOBU = 'A'; 		/*If JOBU = 'A', U contains the M-by-M orthogonal matrix U */
+		JOBVT = 'A';		/*JOBVT = 'A':  all N rows of V**T are returned in the array VT;*/
+		LDU = M;
+		LDVT = N;
+		U = (double*) malloc((double) M*M*sizeof(double));
+		VT = (double*) malloc((double) N*N*sizeof(double));
+	}	
+	else{				/* ans = svd(x) */
+		JOBU = 'N';	
+		JOBVT = 'N';
+	}
+	int LDA = max(1,M);
+	
+	/* Making a copy of input matrix */
+	buf = (double*) malloc((double)M*N*sizeof(double));
+	memcpy(buf,in1,M*N*sizeof(double));
+	
+	S = (double*)malloc((double)min(col,row)*sizeof(double));
+
+	int LWORK = 5*min(M,N);
+	WORK = (double*)malloc((double)LWORK*sizeof(double));
+	int INFO = 0;	/*For successful exit */
+
+	dgesvd_(&JOBU,&JOBVT,&M,&N,buf,&LDA,S,U,&LDU,VT,&LDVT,WORK,&LWORK,&INFO);
+	/*Subroutine DGESVD from Lapack lib. */
+
+	if (nout  == 1){	/* ans = svd(x)*/
+		memcpy(out1,S,min(row,col)*sizeof(double));
+		//printf("%lf %lf %lf",*(S),*(S+1),*(S+2));
+	}			/* [U,S,VT] = svd(x) */
+	else if(in2 == 0 && nout > 1){
+		memcpy(out1,U,LDU*M*sizeof(double));
+		//memcpy(out3,VT,LDVT*min(row,col)*sizeof(double));
+		for(j=0;j<M;j++){
+			for(k=0;k<N;k++){
+				if(j == k) *((out2+j*(min(M,N)))+k) = *(S+j);
+				else *((out2+j*(min(M,N)))+k) = 0;					
+			}
+		}
+							
+		//dtransposea(VT,LDVT,N,out3);
+		/*As there is some patch of error in SVD, these lines are added */
+
+		for(j=1;j<=N;j++){
+			for(i=j;i<=N;i++){
+				*(out3+i+(j-1)*N-1) = VT[j+(i-1)*N-1];
+				*(out3+j+(i-1)*N-1) = VT[i+(j-1)*N-1];
+			}		
+		}
+		/*for(i=0;i<N;i++){
+			for(j=0;j<N;j++){
+				printf("%lf ",VT[i*row+j]);		
+			}	
+			printf("\n");
+		}*/		
+	}
+	else{
+		memcpy(out1,U,M*min(M,N)*sizeof(double));
+		for(j=0;j<min(M,N);j++){
+			for(k=0;k<min(M,N);k++){
+				if(j == k) *((out2+j*(min(M,N)))+k) = *(S+j);
+				else *((out2+j*(min(M,N)))+k) = 0;					
+			}
+		}
+		//dtransposea(VT,LDVT,N,out3);
+		/*As there is some patch of error in DGESVD, these lines are added */
+		/* out3 first taken in some array then will be copied from it. */
+		double *outV;
+		outV = (double *)malloc(N*N*sizeof(double));
+		for(j=1;j<=N;j++){
+			for(i=j;i<=N;i++){
+				*(outV+i+(j-1)*N-1) = VT[j+(i-1)*N-1];
+				*(outV+j+(i-1)*N-1) = VT[i+(j-1)*N-1];
+			}		
+		}
+		
+		for(j=0;j<min(M,N)*N;j++){
+			*(out3+j) = *(outV+j);		
+		}		
+	}	
+	
+	/* From the fortran file of scilab code - if(tol.eq.0.0d0) tol=dble(max(M,N))*eps*stk(lSV) */
+	if(tol == 0){
+		tol = (double)max(M,N)*eps*S[0];
+	}
+	if(nout == 4){	/*[U,S,VT,rk] = svd(X,tol) where tol - tolerance*/
+		rk = 0;
+		for(i=0;i<min(M,N);i++){
+			if(S[i] > tol){
+				rk = i+1;	
+			}	
+		}
+		return rk;	
+	}
+	return 0;   
+}
+
+int min(int a,int b){
+	if(a > b) return b;
+	return a;
+}
+
+int max(int a,int b){
+	if(a > b) return a;
+	return b;
+}
diff --git a/2.3-1/src/c/signalProcessing/transforms/dct/cdcta.c b/2.3-1/src/c/signalProcessing/transforms/dct/cdcta.c
index 5bc27929..7ff83645 100644
--- a/2.3-1/src/c/signalProcessing/transforms/dct/cdcta.c
+++ b/2.3-1/src/c/signalProcessing/transforms/dct/cdcta.c
@@ -15,6 +15,7 @@
 #include "addition.h"
 #include "types.h"
 #include "floatComplex.h" 
+#include "multiplication.h"
 /*#include "matrixMultiplication"*/
 /*#include <fftw3.h>*/
 #include <math.h>
@@ -24,10 +25,10 @@ void cdcta(floatComplex *in,int row,int col,int sign,floatComplex *out)
 	int i,j,k,u,v;
 	int n;
 	int x,y;
-	float res,ress;
+	float res,ress,vv,ff;
 	float re,z,q,m;
-	floatComplex accu = DoubleComplex(0, 0);
-	floatComplex temp,mm;
+	floatComplex accu = FloatComplex(0, 0);
+	floatComplex temp,mm,aa,bb,cc;
 	if(sign==-1)
 	{
 		if(row==1)
@@ -43,17 +44,25 @@ void cdcta(floatComplex *in,int row,int col,int sign,floatComplex *out)
 					{
 						for(j=0;j<col;j++)
 						{
-							y=row*j+i;		
-							temp=in[y]*(cos(((M_PI)*(y+1-1./2.)*(x))/n));
+							y=row*j+i;
+							vv = cos(((M_PI)*(y+1-1./2.)*(x))/n);
+							aa = FloatComplex(vv,0);		
+							temp=cmuls(in[y],aa);
 							out[x]=cadds(out[x],temp);
 						}
 					}
 					if(x==0)
-					out[x]*=1./(sqrt(n));
+					{
+						vv = 1./(sqrt(n));
+						aa = FloatComplex(vv,0);
+						out[x]=cmuls(out[x],aa);
+					}
 					else
 					{
 						float res=2./n;
-						out[x]*=sqrt(res);
+						res = sqrt(res);
+						aa = FloatComplex(res,0);
+						out[x]=cmuls(out[x],aa);
 					}
 				}
 			}
@@ -76,29 +85,53 @@ void cdcta(floatComplex *in,int row,int col,int sign,floatComplex *out)
 							y=j*row+i;
 							z=(float)(((float)j+1.0/2.0)*(float)v);
 							q=(float)(M_PI/(float)col);
-							mm=in[y]*(cos(q*z));	
+							vv = cos(q*z);
+							aa = FloatComplex(vv,0);
+							mm=cmuls(in[y],aa);
 							temp=cadds(temp,mm);
 						}
 						z=(float)(((float)i+1.0/2.0)*(float)u);
 						q=(float)(M_PI/(float)row);
-						temp*=cos(q*z);
+						ff = cos(q*z);
+						bb = FloatComplex(ff,0);
+						temp=cmuls(temp,bb);
 						out[x]=cadds(out[x],temp);
 						}
 					if(u==0)
 					{
-						out[x]*=1./sqrt((float)row);
+						vv = 1./sqrt((float)row);
+						aa = FloatComplex(vv,0);
+						out[x]=cmuls(out[x],aa);
 						if(v==0)
-						out[x]*=1./sqrt((float)col);
+						{
+							vv = 1./sqrt((float)col);
+							aa = FloatComplex(vv,0);
+							out[x]=cmuls(out[x],aa);
+						}
 						else
-						out[x]*=sqrt(2./col);
+						{
+							vv = sqrt(2./col);
+							aa = FloatComplex(vv,0);
+							out[x]=cmuls(out[x],aa);
+						}
 					}
 					else
 					{
-						out[x]*=sqrt(2./row);
+						vv = sqrt(2./row);
+						aa = FloatComplex(vv,0);
+						out[x]=cmuls(out[x],aa);
 						if(v==0)
-						out[x]*=1./sqrt((float)col);
+						{
+							vv = 1./sqrt((float)col);
+							aa = FloatComplex(vv,0);
+							out[x]=cmuls(out[x],aa);
+						}
 						else
-						out[x]*=sqrt(2./col);
+						{
+							vv = sqrt(2./col);
+							aa = FloatComplex(vv,0);
+							out[x]=cmuls(out[x],aa);
+						}
 					}
 				}
 			}
@@ -125,12 +158,14 @@ void cdcta(floatComplex *in,int row,int col,int sign,floatComplex *out)
 							if(y==0)
 							{
 								q=res*(cos(((M_PI)*(j)*(v+1./2.))/n));
-								out[x]=cadds(out[x],in[y]*q);
+								aa = FloatComplex(q,0);
+								out[x]=cadds(out[x],cmuls(in[y],aa));
 							}
 							else
 							{
 								q=ress*(cos(((M_PI)*(j)*(v+1./2.))/n));
-								out[x]=cadds(out[x],in[y]*q);
+								aa = FloatComplex(q,0);
+								out[x]=cadds(out[x],cmuls(in[y],aa));
 							}
 						}
 					}
@@ -156,19 +191,37 @@ void cdcta(floatComplex *in,int row,int col,int sign,floatComplex *out)
 							y=row*j+i;
 							mm=in[j*row+i];
 							z=(float)(((float)v+1.0/2.0)*(float)j);
-							q=(float)(M_PI/(float)col);	
-							mm=mm*(cos(q*z));
+							q=(float)(M_PI/(float)col);
+							vv = cos(q*z);
+							aa = FloatComplex(vv,0);		
+							mm=cmuls(mm,aa);
 							if(j==0)
-							temp=cadds(temp,mm*(1./sqrt((float)col)));
+							{
+								vv = 1./sqrt((float)col);
+								aa = FloatComplex(vv,0);
+								temp=cadds(temp,cmuls(mm,aa));
+							}		
 							else
-							temp=cadds(temp,mm*sqrt(2./col));
+							{
+								vv = sqrt(2./col);
+								aa = FloatComplex(vv,0);
+								temp=cadds(temp,cmuls(mm,aa));
+							}
 						}
 						z=(float)(((float)u+1.0/2.0)*(float)i);
 						q=(float)(M_PI/(float)row);	
 						if(i==0)
-						out[x]=cadds(out[x],temp*((cos(z*q))*(1./sqrt(row))));
+						{	
+							vv = (cos(z*q))*(1./sqrt(row));
+							aa = FloatComplex(vv,0);
+							out[x]=cadds(out[x],cmuls(temp,aa));
+						}
 						else
-						out[x]=cadds(out[x],temp*((cos(z*q))*sqrt(2./row)));
+						{
+							vv = (cos(z*q))*sqrt(2./row);
+							aa = FloatComplex(vv,0);
+							out[x]=cadds(out[x],cmuls(temp,aa));
+						}
 					}
 				}
 			}
diff --git a/2.3-1/src/c/signalProcessing/transforms/dct/zdcta.c b/2.3-1/src/c/signalProcessing/transforms/dct/zdcta.c
index 0204d682..3ae2e333 100644
--- a/2.3-1/src/c/signalProcessing/transforms/dct/zdcta.c
+++ b/2.3-1/src/c/signalProcessing/transforms/dct/zdcta.c
@@ -15,6 +15,7 @@
 #include "addition.h"
 #include "types.h"
 #include "doubleComplex.h" 
+#include "multiplication.h"
 /*#include "matrixMultiplication"*/
 /*#include <fftw3.h>*/
 #include <math.h>
@@ -24,10 +25,10 @@ void zdcta(doubleComplex *in,int row,int col,int sign,doubleComplex *out)
 	int i,j,k,u,v;
 	int n;
 	int x,y;
-	double res,ress;
-	double re,z,q,m;
+	double  res,ress,vv,ff;
+	double  re,z,q,m;
 	doubleComplex accu = DoubleComplex(0, 0);
-	doubleComplex temp,mm;
+	doubleComplex temp,mm,aa,bb,cc;
 	if(sign==-1)
 	{
 		if(row==1)
@@ -43,22 +44,30 @@ void zdcta(doubleComplex *in,int row,int col,int sign,doubleComplex *out)
 					{
 						for(j=0;j<col;j++)
 						{
-							y=row*j+i;		
-							temp=in[y]*(cos(((M_PI)*(y+1-1./2.)*(x))/n));
+							y=row*j+i;
+							vv = cos(((M_PI)*(y+1-1./2.)*(x))/n);
+							aa = DoubleComplex(vv,0);		
+							temp=zmuls(in[y],aa);
 							out[x]=zadds(out[x],temp);
 						}
 					}
 					if(x==0)
-					out[x]*=1./(sqrt(n));
+					{
+						vv = 1./(sqrt(n));
+						aa = DoubleComplex(vv,0);
+						out[x]=zmuls(out[x],aa);
+					}
 					else
 					{
-						double res=2./n;
-						out[x]*=sqrt(res);
+						double  res=2./n;
+						res = sqrt(res);
+						aa = DoubleComplex(res,0);
+						out[x]=zmuls(out[x],aa);
 					}
 				}
 			}
-		}	
-		else	
+		}
+		else
 		{	
 			n=col*row;
 			for(u=0;u<row;u++)
@@ -74,31 +83,55 @@ void zdcta(doubleComplex *in,int row,int col,int sign,doubleComplex *out)
 						for(j=0;j<col;j++)
 						{
 							y=j*row+i;
-							z=(double)(((double)j+1.0/2.0)*(double)v);
-							q=(double)(M_PI/(double)col);
-							mm=in[y]*(cos(q*z));	
+							z=(double )(((double )j+1.0/2.0)*(double )v);
+							q=(double )(M_PI/(double )col);
+							vv = cos(q*z);
+							aa = DoubleComplex(vv,0);
+							mm=zmuls(in[y],aa);
 							temp=zadds(temp,mm);
 						}
-						z=(double)(((double)i+1.0/2.0)*(double)u);
-						q=(double)(M_PI/(double)row);
-						temp*=cos(q*z);
+						z=(double )(((double )i+1.0/2.0)*(double )u);
+						q=(double )(M_PI/(double )row);
+						ff = cos(q*z);
+						bb = DoubleComplex(ff,0);
+						temp=zmuls(temp,bb);
 						out[x]=zadds(out[x],temp);
-					}
+						}
 					if(u==0)
 					{
-						out[x]*=1./sqrt((double)row);
+						vv = 1./sqrt((double )row);
+						aa = DoubleComplex(vv,0);
+						out[x]=zmuls(out[x],aa);
 						if(v==0)
-						out[x]*=1./sqrt((double)col);
+						{
+							vv = 1./sqrt((double )col);
+							aa = DoubleComplex(vv,0);
+							out[x]=zmuls(out[x],aa);
+						}
 						else
-						out[x]*=sqrt(2./col);
+						{
+							vv = sqrt(2./col);
+							aa = DoubleComplex(vv,0);
+							out[x]=zmuls(out[x],aa);
+						}
 					}
 					else
 					{
-						out[x]*=sqrt(2./row);
+						vv = sqrt(2./row);
+						aa = DoubleComplex(vv,0);
+						out[x]=zmuls(out[x],aa);
 						if(v==0)
-						out[x]*=1./sqrt((double)col);
+						{
+							vv = 1./sqrt((double )col);
+							aa = DoubleComplex(vv,0);
+							out[x]=zmuls(out[x],aa);
+						}
 						else
-						out[x]*=sqrt(2./col);
+						{
+							vv = sqrt(2./col);
+							aa = DoubleComplex(vv,0);
+							out[x]=zmuls(out[x],aa);
+						}
 					}
 				}
 			}
@@ -125,12 +158,14 @@ void zdcta(doubleComplex *in,int row,int col,int sign,doubleComplex *out)
 							if(y==0)
 							{
 								q=res*(cos(((M_PI)*(j)*(v+1./2.))/n));
-								out[x]=zadds(out[x],in[y]*q);
+								aa = DoubleComplex(q,0);
+								out[x]=zadds(out[x],zmuls(in[y],aa));
 							}
 							else
 							{
 								q=ress*(cos(((M_PI)*(j)*(v+1./2.))/n));
-								out[x]=zadds(out[x],in[y]*q);
+								aa = DoubleComplex(q,0);
+								out[x]=zadds(out[x],zmuls(in[y],aa));
 							}
 						}
 					}
@@ -155,20 +190,38 @@ void zdcta(doubleComplex *in,int row,int col,int sign,doubleComplex *out)
 						{
 							y=row*j+i;
 							mm=in[j*row+i];
-							z=(double)(((double)v+1.0/2.0)*(double)j);
-							q=(double)(M_PI/(double)col);	
-							mm=mm*(cos(q*z));
+							z=(double )(((double )v+1.0/2.0)*(double )j);
+							q=(double )(M_PI/(double )col);
+							vv = cos(q*z);
+							aa = DoubleComplex(vv,0);		
+							mm=zmuls(mm,aa);
 							if(j==0)
-							temp=zadds(temp,mm*(1./sqrt((double)col)));
+							{
+								vv = 1./sqrt((double )col);
+								aa = DoubleComplex(vv,0);
+								temp=zadds(temp,zmuls(mm,aa));
+							}		
 							else
-							temp=zadds(temp,mm*sqrt(2./col));
+							{
+								vv = sqrt(2./col);
+								aa = DoubleComplex(vv,0);
+								temp=zadds(temp,zmuls(mm,aa));
+							}
 						}
-						z=(double)(((double)u+1.0/2.0)*(double)i);
-						q=(double)(M_PI/(double)row);	
+						z=(double )(((double )u+1.0/2.0)*(double )i);
+						q=(double )(M_PI/(double )row);	
 						if(i==0)
-						out[x]=zadds(out[x],temp*((cos(z*q))*(1./sqrt(row))));
+						{	
+							vv = (cos(z*q))*(1./sqrt(row));
+							aa = DoubleComplex(vv,0);
+							out[x]=zadds(out[x],zmuls(temp,aa));
+						}
 						else
-						out[x]=zadds(out[x],temp*((cos(z*q))*sqrt(2./row)));
+						{
+							vv = (cos(z*q))*sqrt(2./row);
+							aa = DoubleComplex(vv,0);
+							out[x]=zadds(out[x],zmuls(temp,aa));
+						}
 					}
 				}
 			}
diff --git a/2.3-1/src/c/signalProcessing/transforms/idct/cidcta.c b/2.3-1/src/c/signalProcessing/transforms/idct/cidcta.c
index ec0df0c7..ae98ba19 100644
--- a/2.3-1/src/c/signalProcessing/transforms/idct/cidcta.c
+++ b/2.3-1/src/c/signalProcessing/transforms/idct/cidcta.c
@@ -14,7 +14,8 @@
 #include "idct.h"
 #include "addition.h"
 #include "types.h"
-#include "floatComplex.h" 
+#include "floatComplex.h"
+#include "multiplication.h" 
 /*#include "matrixMultiplication"*/
 /*#include <fftw3.h>*/
 #include <math.h>
@@ -24,10 +25,10 @@ void cidcta(floatComplex *in,int row,int col,floatComplex *out)
 	int i,j,k,u,v;
 	int n=col;
 	int x,y;
-	float res,ress;
+	float res,ress,vv,ff;
 	float re,z,q,m;
-	floatComplex accu = DoubleComplex(0, 0);
-	floatComplex temp,mm;
+	floatComplex accu = FloatComplex(0, 0);
+	floatComplex temp,mm,aa,bb;
 	if(row==1)
 	{
 		res=1./sqrt(n);
@@ -46,12 +47,14 @@ void cidcta(floatComplex *in,int row,int col,floatComplex *out)
 						if(y==0)
 						{
 							q=res*(cos(((M_PI)*(j)*(v+1./2.))/n));
-							out[x]=cadds(out[x],in[y]*q);
+							aa=FloatComplex(q,0);
+							out[x]=cadds(out[x],cmuls(in[y],aa));
 						}
 						else
 						{
 							q=ress*(cos(((M_PI)*(j)*(v+1./2.))/n));
-							out[x]=cadds(out[x],in[y]*q);
+							aa=FloatComplex(q,0);
+							out[x]=cadds(out[x],cmuls(in[y],aa));
 						}
 					}
 				}
@@ -77,19 +80,37 @@ void cidcta(floatComplex *in,int row,int col,floatComplex *out)
 						y=row*j+i;
 						mm=in[j*row+i];
 						z=(float)(((float)v+1.0/2.0)*(float)j);
-						q=(float)(M_PI/(float)col);	
-						mm=mm*(cos(q*z));
+						q=(float)(M_PI/(float)col);
+						vv=cos(q*z);
+						aa=FloatComplex(vv,0);	
+						mm=cmuls(mm,aa);
 						if(j==0)
-						temp=cadds(temp,mm*(1./sqrt((float)col)));
+						{
+							vv=1./sqrt((float)col);
+							aa=FloatComplex(vv,0);
+							temp=cadds(temp,cmuls(mm,aa));
+						}
 						else
-						temp=cadds(temp,mm*sqrt(2./col));
+						{
+							vv=sqrt(2./col);
+							aa=FloatComplex(vv,0);
+							temp=cadds(temp,cmuls(mm,aa));
+						}
 					}
 					z=(float)(((float)u+1.0/2.0)*(float)i);
 					q=(float)(M_PI/(float)row);	
 					if(i==0)
-					out[x]=cadds(out[x],temp*((cos(z*q))*(1./sqrt(row))));
+					{
+						vv=(cos(z*q))*(1./sqrt(row));
+						aa=FloatComplex(vv,0);
+						out[x]=cadds(out[x],cmuls(temp,aa));
+					}
 					else
-					out[x]=cadds(out[x],temp*((cos(z*q))*sqrt(2./row)));
+					{
+						vv=(cos(z*q))*sqrt(2./row);
+						aa=FloatComplex(vv,0);
+						out[x]=cadds(out[x],cmuls(temp,aa));
+					}
 				}
 			}
 		}
diff --git a/2.3-1/src/c/signalProcessing/transforms/idct/zidcta.c b/2.3-1/src/c/signalProcessing/transforms/idct/zidcta.c
index 2177b18c..cc01c966 100644
--- a/2.3-1/src/c/signalProcessing/transforms/idct/zidcta.c
+++ b/2.3-1/src/c/signalProcessing/transforms/idct/zidcta.c
@@ -15,6 +15,7 @@
 #include "addition.h"
 #include "types.h"
 #include "doubleComplex.h" 
+#include "multiplication.h"
 /*#include "matrixMultiplication"*/
 /*#include <fftw3.h>*/
 #include <math.h>
@@ -24,10 +25,10 @@ void zidcta(doubleComplex *in,int row,int col,doubleComplex *out)
 	int i,j,k,u,v;
 	int n=col;
 	int x,y;
-	double res,ress;
+	double res,ress,vv,ff;
 	double re,z,q,m;
 	doubleComplex accu = DoubleComplex(0, 0);
-	doubleComplex temp,mm;
+	doubleComplex temp,mm,aa,bb;
 	if(row==1)
 	{
 		res=1./sqrt(n);
@@ -46,17 +47,19 @@ void zidcta(doubleComplex *in,int row,int col,doubleComplex *out)
 						if(y==0)
 						{
 							q=res*(cos(((M_PI)*(j)*(v+1./2.))/n));
-							out[x]=zadds(out[x],in[y]*q);
+							aa=DoubleComplex(q,0);
+							out[x]=zadds(out[x],zmuls(in[y],aa));
 						}
 						else
 						{
 							q=ress*(cos(((M_PI)*(j)*(v+1./2.))/n));
-							out[x]=zadds(out[x],in[y]*q);
+							aa=DoubleComplex(q,0);
+							out[x]=zadds(out[x],zmuls(in[y],aa));
 						}
 					}
 				}
 			}
-				
+			
 		}
 	}	
 	else	
@@ -77,19 +80,37 @@ void zidcta(doubleComplex *in,int row,int col,doubleComplex *out)
 						y=row*j+i;
 						mm=in[j*row+i];
 						z=(double)(((double)v+1.0/2.0)*(double)j);
-						q=(double)(M_PI/(double)col);	
-						mm=mm*(cos(q*z));
+						q=(double)(M_PI/(double)col);
+						vv=cos(q*z);
+						aa=DoubleComplex(vv,0);	
+						mm=zmuls(mm,aa);
 						if(j==0)
-						temp=zadds(temp,mm*(1./sqrt((double)col)));
+						{
+							vv=1./sqrt((double)col);
+							aa=DoubleComplex(vv,0);
+							temp=zadds(temp,zmuls(mm,aa));
+						}
 						else
-						temp=zadds(temp,mm*sqrt(2./col));
+						{
+							vv=sqrt(2./col);
+							aa=DoubleComplex(vv,0);
+							temp=zadds(temp,zmuls(mm,aa));
+						}
 					}
 					z=(double)(((double)u+1.0/2.0)*(double)i);
 					q=(double)(M_PI/(double)row);	
 					if(i==0)
-					out[x]=zadds(out[x],temp*((cos(z*q))*(1./sqrt(row))));
+					{
+						vv=(cos(z*q))*(1./sqrt(row));
+						aa=DoubleComplex(vv,0);
+						out[x]=zadds(out[x],zmuls(temp,aa));
+					}
 					else
-					out[x]=zadds(out[x],temp*((cos(z*q))*sqrt(2./row)));
+					{
+						vv=(cos(z*q))*sqrt(2./row);
+						aa=DoubleComplex(vv,0);
+						out[x]=zadds(out[x],zmuls(temp,aa));
+					}
 				}
 			}
 		}
diff --git a/2.3-1/src/c/string/disp/ddispa.c b/2.3-1/src/c/string/disp/ddispa.c
index 5e6bb84e..29aaceeb 100644
--- a/2.3-1/src/c/string/disp/ddispa.c
+++ b/2.3-1/src/c/string/disp/ddispa.c
@@ -16,7 +16,7 @@ double ddispa (double* in, int rows, int columns){
 	int i = 0,j = 0;
 	
 	for (i = 0; i < rows; ++i) {
-		for (j=0;j<columns;j++) printf ("  %1.20f  ", in[i+j*rows]);
+		for (j=0;j<columns;j++) printf ("  %e  ", in[i+j*rows]);
 		printf("\n");
 	}
 	return 0;