welch.cpp

笔者在前人研究的基础上

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "msp.h"

float psdr[4096];
complex x[128];
float xr[128];
float xi[128];

//#define DEBUG

//void mperpsd(complex x[],float psdr[], int n,int nshift,int nsamp,int iwidow,float ts);
int readData()
{ 
	FILE *fpr,*fpi; 
	int i=0; 
	
	fpr=fopen("F:\\myprojects\\Welch\\datar.txt","r");
	fpi=fopen("F:\\myprojects\\Welch\\datai.txt","r");
    
	if((fpr==NULL)||(fpi==NULL))
	{ 
		printf("错误,文件打不开,请检查文件名\n"); 
	} 
	else
	{ 
		while(!feof(fpr))
		{ 
			fscanf(fpr,"   %f",&x[i].real); 			
			i++; 
		} 	
		
		i=0;		
		while(!feof(fpi))
		{ 		
			fscanf(fpi,"   %f",&x[i].imag); 
			i++; 
		}
		fclose(fpr);
		fclose(fpi);

		return(1); 
	} 
	return(0); 
} 

int writeData(float *src,int len)
{ 
	FILE *fpr; 
	int i=0; 
	
	fpr=fopen("F:\\myprojects\\Welch\\data.txt","w");
	    
	if(fpr==NULL)
	{ 
		printf("错误,文件打不开,请检查文件名\n"); 
	} 
	else
	{ 
		while(i<len)
		{ 
#ifdef DEBUG
			printf("   %f\n",src[i]);
#endif
			fprintf(fpr,"   %f\n",src[i]);
			i++; 
		} 
		fclose(fpr);
		return(1); 
	} 
	return(0); 
} 

int main()
{
   readData();
   for(int i=0;i<128;i++)
   {
     xr[i]=x[i].real;
     xi[i]=x[i].imag;   
#ifdef DEBUG
	printf(" %f",x[i].real);
#endif
   }
   

#ifdef DEBUG
   mrelfft(xr,xi,128,-1);
   writeData(xi,128);    
#endif      
   mperpsd(x,128,16,33,-1,0.2);
   
   return 0;
}