nrutil.cpp

针对基于DSP的使用C++编译的数字信号的高低频滤波的编程源代码

#include   <stdio.h>
#include   <stddef.h>
#include   <stdlib.h>

//报告错误
void nrerror(char *error_text)
{
	printf("Numerical Recipes run-time error...\n");
	printf("%s\n",error_text);
	printf("...now exiting to system...\n");
	exit(1);
}

/****************************************************************************
*****函数名称:  vector
*****功    能:  申请一个实数数组
*****参    数:  
                nl -----数组的最低下标
				nh -----数组的最高下标
*****返    回:  数组地址
*/

float  *fvector(long nl,long nh)
{
	float *v;

	v=(float *)calloc((unsigned)(nh-nl+1),sizeof(float));
	if(!v)
		nrerror("数组分配失败!");
	return v-nl;
}

/*
double*	vector (long nl, long nh)             //申请一个实数数组
{
	double* 	v;
	v = (double*)malloc((size_t)((nh-nl+1+NR_END)*sizeof(double)));
	if (!v)
		nrerror("allocation failure in vector()");
	return	v-nl+NR_END;
}
*/
/****************************************************************************
*****函数名称:  ivector
*****功    能:  申请一个整型数组
*****参    数:  
                nl -----数组的最低下标
				nh -----数组的最高下标
*****返    回:  数组地址
*/

int  *ivector(long nl,long nh)
{
	int *v;

	v=(int *)calloc((unsigned)(nh-nl+1),sizeof(int));
	if(!v)
		nrerror("数组分配失败!");
	return v-nl;
}

/****************************************************************************
*****函数名称:  cvector
*****功    能:  申请一个无符号字符数组
*****参    数:  
                nl -----数组的最低下标
				nh -----数组的最高下标
*****返    回:  数组地址
*/

unsigned char  *cvector(long nl,long nh)
{
	unsigned char *v;

	v=(unsigned char *)calloc((unsigned)(nh-nl+1),sizeof(unsigned char));
	if(!v)
		nrerror("数组分配失败!");
	return v-nl;
}

/****************************************************************************
*****函数名称:  lvector
*****功    能:  申请一个长整型数组
*****参    数:  
                nl -----数组的最低下标
				nh -----数组的最高下标
*****返    回:  数组地址
*/

unsigned long  *lvector(long nl,long nh)
{
	unsigned long *v;

	v=(unsigned long *)calloc((unsigned)(nh-nl+1),sizeof(unsigned long));
	if(!v)
		nrerror("数组分配失败!");
	return v-nl;
}

/****************************************************************************
*****函数名称:  dvector
*****功    能:  申请一个双精度实数数组
*****参    数:  
                nl -----数组的最低下标
				nh -----数组的最高下标
*****返    回:  数组地址
*/

double  *dvector(long nl,long nh)
{
	double *v;

	v=(double *)calloc((unsigned)(nh-nl+1),sizeof(double));
	if(!v)
		nrerror("数组分配失败!");
	return v-nl;
}

/****************************************************************************
*****函数名称:  matrix
*****功    能:  申请一个实数二维数组(矩阵)
*****参    数:  
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****返    回:  数组地址
*/

float  **matrix(long nrl,long nrh,long ncl,long nch)
{
	int      i;
	float    **m;

	m=(float **)calloc((unsigned)(nrh-nrl+1),sizeof(float *));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=nrl;

	for(i=nrl;i<=nrh;i++)
	{
		m[i]=(float *)calloc((unsigned)(nch-ncl+1),sizeof(float));
		if(!m[i])
			nrerror("数组的列分配失败!");
		m[i]-=ncl;
	}
	return m;
}

/****************************************************************************
*****函数名称:  dmatrix
*****功    能:  申请一个双精度实数二维数组(矩阵)
*****参    数:  
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****返    回:  数组地址
*/

double  **dmatrix(long nrl,long nrh,long ncl,long nch)
{
	int      i;
	double   **m;

	m=(double **)calloc((unsigned)(nrh-nrl+1),sizeof(double *));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=nrl;

	for(i=nrl;i<=nrh;i++)
	{
		m[i]=(double *)calloc((unsigned)(nch-ncl+1),sizeof(double));
		if(!m[i])
			nrerror("数组的列分配失败!");
		m[i]-=ncl;
	}
	return m;
}

/****************************************************************************
*****函数名称:  imatrix
*****功    能:  申请一个整数二维数组(矩阵)
*****参    数:  
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****返    回:  数组地址
*/

int  **imatrix(long nrl,long nrh,long ncl,long nch)
{
	int      i;
	int    **m;

	m=(int **)calloc((unsigned)(nrh-nrl+1),sizeof(int*));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=nrl;

	for(i=nrl;i<=nrh;i++)
	{
		m[i]=(int *)calloc((unsigned)(nch-ncl+1),sizeof(int));
		if(!m[i])
			nrerror("数组的列分配失败!");
		m[i]-=ncl;
	}
	return m;
}

/****************************************************************************
*****函数名称:  cmatrix
*****功    能:  申请一个无符号字符型二维数组(矩阵)
*****参    数:  
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****返    回:  数组地址
*/

unsigned char  **cmatrix(long nrl,long nrh,long ncl,long nch)
{
	int      i;
	unsigned char    **m;

	m=(unsigned char **)calloc((unsigned)(nrh-nrl+1),sizeof(unsigned char *));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=nrl;

	for(i=nrl;i<=nrh;i++)
	{
		m[i]=(unsigned char *)calloc((unsigned)(nch-ncl+1),sizeof(unsigned char));
		if(!m[i])
			nrerror("数组的列分配失败!");
		m[i]-=ncl;
	}
	return m;
}

/****************************************************************************
*****函数名称:  submatrix
*****功    能:  求一个二维矩阵的子矩阵
*****参    数:  
                a     -----二维矩阵的首地址
                oldrl -----数组的行最低下标
				oldrh -----数组的行最高下标
                oldcl -----数组的列最低下标
				oldch -----数组的列最高下标
				newrl -----新数组的行最低下标
				newcl -----新数组的列最低下标
*****返    回:  数组地址
*/

float  **submatrix(float **a,long oldrl,long oldrh,long oldcl,long oldch,long newrl,long newcl)
{
	int      i,j;
	float    **m;

	m=(float **)calloc((unsigned)(oldrh-oldrl+1),sizeof(float *));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=newrl;

	for(i=oldrl,j=newrl;i<=oldrh;i++,j++)
	{
		m[j]=a[i]+oldcl-newcl;
	}
	return m;
}

/****************************************************************************
*****函数名称:  convert_matrix
*****功    能:  求一个实数一维数组转化为二维数组(矩阵)
*****参    数:  
                a     -----原一维数组
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标				
*****返    回:  数组地址
*/

float  **convert_matrix(float *a,long nrl,long nrh,long ncl,long nch)
{
	int      i,j,nrow,ncol;
	float    **m;

	nrow=nrh-nrl+1;
	ncol=nch-ncl+1;
	m=(float **)calloc((unsigned)(nrow),sizeof(float *));
	if(!m)
		nrerror("转换数组分配失败!");
	m-=nrl;

	for(i=0,j=nrl;i<=nrow-1;i++,j++)
	{
		m[j]=a+ncol*i-ncl;
	}
	return m;
}


/****************************************************************************
*****函数名称:  f3tensor
*****功    能:  申请一个实数二维数组(矩阵)
*****参    数:  
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
                ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****返    回:  数组地址
*/

float ***f3tensor(long nrl,long nrh,long ncl,long nch,long ndl,long ndh)
{
	int      i,j;
	float    ***m;

	m=(float ***)calloc((unsigned)(nrh-nrl+1),sizeof(float **));
	if(!m)
		nrerror("数组的行分配失败!");
	m-=nrl;

	for(i=nrl;i<=nrh;i++)
	{
		m[i]=(float **)calloc((unsigned)(nch-ncl+1),sizeof(float*));
		if(!m[i])
			nrerror("数组的列分配失败!");
		m[i]-=ncl;
		for(j=ncl;j<=nch;j++)
		{
			m[i][j]=(float *)calloc((unsigned)(nch-ncl+1),sizeof(float));
		    if(!m[i][j])
			    nrerror("数组的列分配失败!");
		    m[i][j]-=ncl;
		}

	}
	return m;
}
/////////////////////////////////////////////////////////////////////////////
/////////////////////////////////////////////////////////////////////////////

/****************************************************************************
*****函数名称:  free_vector
*****功    能:  释放一个实数一维数组
*****参    数:  
                nl -----数组的行最低下标
				nh -----数组的列最低下标
*****无 返 回 值 
*/

void  free_fvector(float *v,long nl,long nh)
{
	free((char *)(v+nl));
}

/****************************************************************************
*****函数名称:  free_ivector
*****功    能:  释放一个整型一维数组
*****参    数:  
                nl -----数组的行最低下标
				nh -----数组的列最低下标
*****无 返 回 值 
*/

void  free_ivector(int *v,long nl,long nh)
{
	free((char *)(v+nl));
}

/****************************************************************************
*****函数名称:  free_cvector
*****功    能:  释放一个无符号字符型一维数组
*****参    数:  
                nl -----数组的行最低下标
				nh -----数组的列最低下标
*****无 返 回 值 
*/

void  free_vector(unsigned char *v,long nl,long nh)
{
	free((char *)(v+nl));
}

/****************************************************************************
*****函数名称:  free_lvector
*****功    能:  释放一个长整型一维数组
*****参    数:  
                nl -----数组的行最低下标
				nh -----数组的列最低下标
*****无 返 回 值 
*/

void  free_lvector(unsigned long *v,long nl,long nh)
{
	free((char *)(v+nl));
}

/****************************************************************************
*****函数名称:  free_dvector
*****功    能:  释放一个双精度一维数组
*****参    数:  
                nl -----数组的行最低下标
				nh -----数组的列最低下标
*****无 返 回 值 
*/

void  free_dvector(double *v,long nl,long nh)
{
	free((char *)(v+nl));
}

/****************************************************************************
*****函数名称:  free_matrix
*****功    能:  释放一个实数二维数组(矩阵)
*****参    数:  
                m  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_matrix(float **m,long nrl,long nrh,long ncl,long nch)
{
	int i;

	for(i=nrh;i>=nrl;i--)
		free((char *)(m[i]+ncl));
	free((char *)(m+nrl));
}

/****************************************************************************
*****函数名称:  free_dmatrix
*****功    能:  释放一个双精度二维数组(矩阵)
*****参    数:  
                m  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_dmatrix(double **m,long nrl,long nrh,long ncl,long nch)
{
	int i;

	for(i=nrh;i>=nrl;i--)
		free((char *)(m[i]+ncl));
	free((char *)(m+nrl));
}

/****************************************************************************
*****函数名称:  free_imatrix
*****功    能:  释放一个整型二维数组(矩阵)
*****参    数:  
                m  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_imatrix(int **m,long nrl,long nrh,long ncl,long nch)
{
	int i;

	for(i=nrh;i>=nrl;i--)
		free((char *)(m[i]+ncl));
	free((char *)(m+nrl));
}

/****************************************************************************
*****函数名称:  free_cmatrix
*****功    能:  释放一个无符号字符型二维数组(矩阵)
*****参    数:  
                m  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_cmatrix(unsigned char **m,long nrl,long nrh,long ncl,long nch)
{
	int i;

	for(i=nrh;i>=nrl;i--)
		free((char *)(m[i]+ncl));
	free((char *)(m+nrl));
}

/****************************************************************************
*****函数名称:  free_submatrix
*****功    能:  释放一个子二维数组(矩阵)
*****参    数:  
                b  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_submatrix(float **b,long nrl,long nrh,long ncl,long nch)
{
	free((char *)(b+nrl));
}

/****************************************************************************
*****函数名称:  free_convert_submatrix
*****功    能:  释放一个转化后的二维数组(矩阵)
*****参    数:  
                b  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
*****无 返 回 值 
*/

void  free_convert_matrix(float **b,long nrl,long nrh,long ncl,long nch)
{
	free((char *)(b+nrl));
}

/****************************************************************************
*****函数名称:  free_convert_submatrix
*****功    能:  释放一个转化后的二维数组(矩阵)
*****参    数:  
                b  -----原数组地址
                nrl -----数组的行最低下标
				nrh -----数组的行最高下标
				ncl -----数组的列最低下标
				nch -----数组的列最高下标
				ndl -----数组的第三维最低下标
				ndh -----数组的第三维最高下标
*****无 返 回 值 
*/
void  free_f3tensor(float ***b,long nrl,long nrh,long ncl,long nch,long ndl,long ndh)
{
	free((char *)(b+nrl));
}