ocn_f.h

图像处理的压缩算法

	
	Matrix Q
	
	-0.5000 -0.5000  0.0000 -0.5000 -0.5000
	-0.5000 -0.5000  0.0000  0.5000  0.5000
	-0.4000  0.4000 -0.6000 -0.4000  0.4000
	-0.5000  0.5000  0.0000  0.5000 -0.5000
	-0.3000  0.3000  0.8000 -0.3000  0.3000



Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	70: On entry, parameter wheret had an illegal value.
	17: On entry, m = _value_ while n = _value_. These parameters must satisfy m = n.  On entry, tda = _value_ while max(n,ncolq) = _value_. These parameters must satisfy tda = max(n,ncolq).
	19: On entry, ncolq = _value_ while m = _value_. These parameters must satisfy ncolq = m.
	11: On entry, n must not be less than 0: n = _value_.  On entry, ncolq must not be less than 0: ncolq = _value_.
	73: Memory allocation failed.

	successfully call of the nag_real_form_q function.


*/
	int  nag_real_form_q(
	Nag_WhereElements wheret, 
	int m, // the number of rows of A.
	int n, // the number of columns of A.
	int ncolq, // the required number of columns of Q.  
	double a[], // Input: the required number of columns of Q.  When ncolq = 0 then an immediate return is effected.  Output: the first ncolq columns of the array a are overwritten by the first ncolq columns of the m by m orthogonal matrix Q. 
	int tda, // the second dimension of the array a as declared in the function from which nag real form q is called.
	double zeta[] // contain the elements of zeta.
	);


/**	f01rcc
		finds the QR factorization of the complex m by n matrix A, where m = n.

Example:
	To obtain the QR factorization of the 5 by 3 matrix
	A =
	.
	....
	0.5i -0.5 + 1.5i -1.0+ 1.0i
	0.4 + 0.3i 0.9 + 1.3i 0.2+ 1.4i
	0.4 -0.4 + 0.4i 1.8
	0.3 - 0.4i 0.1 + 0.7i 0.0
	- 0.3i 0.3 + 0.3i 2.4i
	.
	....


void test_nag_complex_qr()
{	
	int i, j, m, n;

	matrix<complex> theta;
	theta.SetSize(1,10);					 //{0.0-0.3i,0.3+0.3i,0.0+2.4i}};
	matrix<complex> a;
	a.SetSize(20,10);
	a[0][0] = 0.0+0.5i;
	a[0][1] = -0.5+1.5i;
	a[0][2] = -1.0+1.0i;
	
	a[1][0] = 0.4+0.3i;
	a[1][1] = 0.9+1.3i;
	a[1][2] = 0.2+1.4i;
	
	a[2][0] = 0.4+0.0i;
	a[2][1] = -0.4+0.4i;
	a[2][2] = 1.8+0.0i;
	
	a[3][0] = 0.3-0.4i;
	a[3][1] = 0.1+0.7i;
	a[3][2] = 0.0+0.0i;
	
	a[4][0] = 0.0-0.3i;
	a[4][1] = 0.3+0.3i;
	a[4][2] = 0.0+2.4i;
	m = 5;
	n = 3;
	printf("A=\n");
	for (i=0; i<m; ++i)
	{
		for (j=0; j<n; ++j)
		{	
			complex temp = a[i][j];
			
			printf(" (%7.4f,%8.4f)", temp.m_re, temp.m_im);
			if(j%3==2 || j==n-1)
				printf("\n");
		}
	}
	nag_complex_qr(m, n, a, 10, theta);
	printf("QR factorization of A\n");
	printf("Vector THETA\n");
	for (i=0; i<n; ++i)
	{
		complex temp = theta[0][i];
		printf(" (%7.4f,%8.4f)", temp.m_re, temp.m_im);
		if(i%3==2 || i==n-1)
			printf("\n");
	}
	printf("\nMatrix A after factorization (upper triangular part is R)\n");
	for (i=0; i<m; ++i)
	{
		for (j=0; j<n; ++j)
		{	
			complex temp = a[i][j];
			
			printf(" (%7.4f,%8.4f)", temp.m_re, temp.m_im);
			if(j%3==2 || j==n-1)
				printf("\n");
		}
	}
}

	The output is following:
	A=
	 ( 0.0000,  0.5000) (-0.5000,  1.5000) (-1.0000,  1.0000)
     ( 0.4000,  0.3000) ( 0.9000,  1.3000) ( 0.2000,  1.4000)
     ( 0.4000,  0.0000) (-0.4000,  0.4000) ( 1.8000,  0.0000)
     ( 0.3000, -0.4000) ( 0.1000,  0.7000) ( 0.0000,  0.0000)
     ( 0.0000, -0.3000) ( 0.3000,  0.3000) ( 0.0000,  2.4000)

	
	QR factorization of A
	Vector THETA
	
	( 1.0000, 0.5000) ( 1.0954, -0.3333) ( 1.2649, 0.0000)
	
	Matrix A after factorization (upper triangular part is R)
	
	( 1.0000, 0.0000) 	( 1.0000, 1.0000)  ( 1.0000, 1.0000)
	(-0.2000, -0.4000)  (-2.0000, 0.0000)  (-1.0000, -1.0000)
	(-0.3200, -0.1600)  (-0.3505, 0.2629)  (-3.0000, 0.0000)
	(-0.4000, 0.2000)   ( 0.0000, 0.5477)  ( 0.0000, 0.0000)
	(-0.1200, 0.2400)   ( 0.1972, 0.2629)  ( 0.0000, 0.6325)


Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	17: On entry, m = _value_ while n = _value_. These parameters must satisfy m = n.  On entry, tda = _value_ while n = _value_. These parameters must satisfy tda = n.
	11: On entry, n must not be less than 0: n = _value_.
	
	successfully call of the nag_complex_qr function.

*/
	int  nag_complex_qr(
	int m, // the number of rows of A.
	int n, // the number of columns of A.
	complex a[], // Input: the leading m by n part of the array a must contain the matrix to be factorized.  Output: the n by n upper triangular part of a will contain the upper triangular matrix R, with the imaginary parts of the diagonal elements set to zero, and the m by n strictly lower triangular part of a will contain details of the factorization as described above.
	int tda, // the second dimension of the array a as declared in the function from which nag complex qr is called.
	complex theta[]  // the scalar theta_k for the kth transformation. 
	);

 
/**	f01rdc
		) performs one of the transformations B := QB or B := QHB, where B 
		is an m by ncolb complex matrix and Q is an m by m unitary matrix, 
		given as the product of Householder transformation matrices.


Example:
	To obtain the matrix QHB for the matrix B given by
	B =
	.
	....
	-0.55 + 1.05i 0.45 + 1.05i
	0.49 + 0.93i 1.09 + 0.13i
	0.56 - 0.16i 0.64 + 0.16i
	0.39 + 0.23i -0.39 - 0.23i
	1.1 3 +0.83i -1.1 3 +0.77i
	.
	....
	following the QR factorization of the 5 by 3 matrix A given by
	A =
	.
	....
	0.5i -0.5 + 1.5i -1.0+ 1.0i
	0.4 + 0.3i 0.9 + 1.3i 0.2+ 1.4i
	0.4 -0.4 + 0.4i 1.8
	0.3 - 0.4i 0.1 + 0.7i 0.0
	-0.3i 0.3 + 0.3i 2.4i
	.
	....

void test_nag_complex_apply_q()
{


	int i,j,m, n,ncolb;
	matrix<complex> a, b, theta;
	a.SetSize(20,10);
	b.SetSize(20,5);
	theta.SetSize(1,10);


	m = 5;
	n = 3;
	
	a[0][0] = 0.0+0.5i;
	a[0][1] = -0.5+1.5i;
	a[0][2] = -1.0+1.0i;
	
	a[1][0] = 0.4+0.3i;
	a[1][1] = 0.9+1.3i;
	a[1][2] = 0.2+1.4i;
	
	a[2][0] = 0.4+0.0i;
	a[2][1] = -0.4+0.4i;
	a[2][2] = 1.8+0.0i;
	
	a[3][0] = 0.3-0.4i;
	a[3][1] = 0.1+0.7i;
	a[3][2] = 0.0+0.0i;
	
	a[4][0] = 0.0-0.3i;
	a[4][1] = 0.3+0.3i;
	a[4][2] = 0.0+2.4i;
	ncolb = 2;
	b[0][0] =-0.55 +1.05i;
	b[0][1] =0.45+1.05i; 
	 
	b[1][0] = 0.49 + 0.93i;
	b[1][1] = 1.09+0.13i;
	
	b[2][0] =0.56 - 0.16i;
	b[2][1] = 0.64+0.16i;
	
	b[3][0] =0.39+0.23i;
	b[3][1] = -0.39-0.23i;
	
	b[4][0] =1.13+0.83i;
	b[4][1] = -1.13+0.77i;

	
	nag_complex_qr(m,n, a, 10, theta);
	
	nag_complex_apply_q(ConjugateTranspose, Nag_ElementsSeparate, m,n,a, 10,
						theta, ncolb, b, 5);

	printf("\nMatrix conjg( Q' )*B\n");
	for (i=0; i<m; ++i)
	{
		for (j=0; j<ncolb; ++j)
		{
			complex temp = b[i][j];
			printf(" (%7.4f,%8.4f)", temp.m_re,temp.m_im);
			if(j%2==1 || j==n-1)
				printf("\n");
			
		}
	}

}

	The output is following:
	
	Matrix conjg( Q' )*B
	
	( 1.0000,1.0000) 	( 1.0000, -1.0000)
	(-1.0000,0.0000) 	(-1.0000,0.0000)
	(-1.0000,1.0000) 	(-1.0000, -1.0000)
	(-0.0600,-0.0200) 	(-0.0400, 0.1200)
	( 0.0400,0.1200) 	(-0.0600,0.0200)


Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	70: On entry, parameter trans had an illegal value.  On entry, parameter wheret had an illegal value.
	17: On entry, m = _value_ while n = _value_. These parameters must satisfy m = n.  On entry, tda = _value_ while n = _value_. These parameters must satisfy tda = n.  On entry, tdb = _value_while ncolb = _value_. These parameters must satisfy tdb = ncolb.
	11: On entry, n must not be less than 0: n = _value_.  On entry, ncolb must not be less than 0: ncolb = _value_.
	73: Memory allocation failed.

	successfully call of the nag_complex_apply_q function.
*/
	int  nag_complex_apply_q(
	MatrixTranspose trans,
	Nag_WhereElements wheret,
	int m, // the number of rows of A.
	int n, // the number of columns of A.
	complex a[], // the leading m by n strictly lower triangular part of the array a must contain details of the matrix Q.  
	int tda, // the second dimension of the array a as declared in the function from which nag complex apply q is called.
	complex theta[], // contain the elements of theta. 
	int ncolb,  // Input: ncolb, the number of columns of B. When ncolb = 0 then an immediate return is effected.
	complex b[], // Input: the leading m by ncolb part of the array b must contain the matrix to be transformed.  Output: b is overwritten by the transformed matrix.
	int tdb // the second dimension of the array b as declared in the function from which nag complex apply q is called.
	);


/**	f01rec
		returns the .rst ncolq columns of the m by m unitary matrix Q, where
		Q is given as the product of Householder transformation matrices.


Example:
	To obtain the 5 by 5 unitary matrix Q following the QR factorization 
	of the 5 by 3 matrix A given by
	A =
	.
	....
	0.5i -0.5 + 1.5i -1.0+ 1.4i
	0.4 + 0.3i 0.9 + 1.3i 0.2+ 1.4i
	0.4 -0.4 + 0.4i 1.8
	0.3 - 0.4i 0.1 + 0.7i 0.0
	-0.3i 0.3 + 0.3i 2.4i
	.
	....
	

void test_nag_complex_form_q()
{
	int i, j, m, n, ncolq;
	matrix<complex>a, q, theta;
	a.SetSize(20,10);
	q.SetSize(20,20);
	theta.SetSize(1,10);
	m = 5;
	n = 3;
	a[0][0] = 0.0+0.5i;
	a[0][1] = -0.5+1.5i;
	a[0][2] = -1.0+1.0i;
	
	a[1][0] = 0.4+0.3i;
	a[1][1] = 0.9+1.3i;
	a[1][2] = 0.2+1.4i;
	
	a[2][0] = 0.4+0.0i;
	a[2][1] = -0.4+0.4i;
	a[2][2] = 1.8+0.0i;
	
	a[3][0] = 0.3-0.4i;
	a[3][1] = 0.1+0.7i;
	a[3][2] = 0.0+0.0i;
	
	a[4][0] = 0.0-0.3i;
	a[4][1] = 0.3+0.3i;
	a[4][2] = 0.0+2.4i;
	
	
	nag_complex_qr(m, n, a, 10, theta);
	
	for (j=0; j<n; ++j)
		for (i=0; i<m; ++i)
			q[i][j] = a[i][j]; 
	

	ncolq = m;
	nag_complex_form_q(Nag_ElementsSeparate, m, n, ncolq, q, 20,theta);

	printf("\nMatrix Q\n");
	for (i=0; i<m; ++i)
	{
		for (j=0; j<ncolq; ++j)
		{
			complex temp = q[i][j];
			printf(" (%5.2f,%5.2f)", temp.m_re, temp.m_im);
			if(j%5==4 || j==ncolq-1)
				printf("\n");
		}	
	}
}


	The output is following:
	
	Matrix Q
	
	( 0.00, 0.50) ( 0.00,-0.50) ( 0.00, 0.00) ( 0.50, 0.00) ( 0.40, 0.30)
	( 0.40, 0.30) (-0.40,-0.30) ( 0.00, 0.00) (-0.30, 0.40) (-0.48, 0.14)
	( 0.40, 0.00) ( 0.40, 0.00) (-0.60, 0.00) (-0.24,-0.32) ( 0.00, 0.40)
	( 0.30,-0.40) ( 0.30,-0.40) ( 0.00, 0.00) ( 0.50, 0.00) (-0.40,-0.30)
	( 0.00,-0.30) ( 0.00,-0.30) ( 0.00,-0.80) (-0.24, 0.18) ( 0.30, 0.00)


Parameters:

Return:
	This function returns NAG error code, 0 if no error.
	
	70: On entry, parameter wheret had an illegal value.
	17: On entry, m = _value_ while n = _value_. These parameters must satisfy m = n.  On entry, tda = _value_ while max(n,ncolq) = _value_. These parameters must satisfy tda = n.
	11: On entry, n must not be less than 0: n = _value_.  On entry, ncolq must not be less than 0: ncolq = _value_.
	19: On entry, ncolq = _value_ while m = _value_. These parameters must satisfy ncolq = m.
	73: Memory allocation failed.

	successfully call of the nag_complex_form_q function.
*/
	int  nag_complex_form_q(
	Nag_WhereElements wheret,
	int m, // the number of rows of A.
	int n, // the number of columns of A.
	int ncolq, // the required number of columns of Q.  
	complex a[],  // Input: the leading m by n strictly lower triangular part of the array a must contain details of the matrix Q.  Output: the .rst ncolq columns of the array a are overwritten by the first ncolq columns of the m by m unitary matrix Q. When n = 0 then the .rst ncolq columns of a are overwritten by the .rst ncolq columns of the unit matrix.
	int tda,  // the second dimension of the array a as declared in the function from which nag complex form q is called.
	complex theta[] // contain the elements of theta.
	);

/**	f02aac
		calculates all the eigenvalues of a real symmetric matrix.
		
Example:
	To calculate all the eigenvalues of the real symmetric matrix
	.
	..
	0.5 0.0 2.3 -2.6
	0.0 0.5 -1.4 -0.7
	2.3 -1.4 0.5 0.0
	-2.6 -0.7 0.0 0.5
	.
	..
	.

void test_nag_real_symm_eigenvalues()
{	
	int i, j, n;
	matrix a;
	a.SetSize(8,8);
	double r[8];

	n = 4;
	a[0][0] = 0.5;
	a[0][1] = 0.0;
	a[0][2] = 2.3;
	a[0][3] = -2.6;
	
	a[1][0] = 0.0;
	a[1][1] = 0.5;
	a[1][2] = -1.4;
	a[1][3] = -0.7;
	
	a[2][0] = 2.3;
	a[2][1] = -1.4;
	a[2][2] = 0.5;
	a[2][3] = 0.0;
	
	a[3][0] = -2.6;
	a[3][1] = -0.7;
	a[3][2] = 0.0;
	a[3][3] = 0.5;
		
	nag_real_symm_eigenvalues(n, a, 8, r);
	printf("Eigenvalues\n");
	for (i=0; i<n; i++)
	{
		printf("%9.4f",r[i]);
		if(i%8==7 || i==n-1)
			printf("\n");
	}	
}

	The output is following:
	
	Eigenvalues
	
	-3.0000 -1.0000 2.0000 4.0000

Parameters:

Return:
	This function returns NAG error code, 0 if no error.