database.cpp

新一代天气雷达基础数据读取原程序,具体参数请参看头文件.

#include "DataBase.h"
#include <windows.h>
#include <iostream.h>
#include <stdio.h>
//************************************************************************//
//* 1)在主函数中设置基数据文件名,并调用ReadBaseData(),SavedataIntoFiles  
//*   两个函数,分别进行读数据和写数据到文件中的操作                      
//* 2)定义了一个变量和4个数组:                                           
//*       NumValidCuts I*4,  记录当前体扫的有效仰角数(指角度不同的仰角)  
//*       Elevation[], R*4,  存放各层仰角
//*       RData[][][], R*4,  存放反射率体扫数据
//*       VData[][][], R*4,  存放速度体扫数据
//*       WData[][][], R*4,  存放谱宽体扫数据
//* 3)将三种基数据分别保存到3个文件中,以二进制的格式存放.
//*       RefArray.dat  按顺序保存下列数据
//*                     NumValidCuts I4
//*                     Elevation[]  NumValidCuts*R4 BYTES
//*                     RData[][][]  NumValidCuts*MaxRads*RGates*R4 BYTES
//*       VelArray.dat  按顺序保存下列数据
//*                     NumValidCuts I4
//*                     Elevation[]  NumValidCuts*R4 BYTES
//*                     VData[][][]  NumValidCuts*MaxRads*VGates*R4 BYTES
//*       SpwArray.dat  按顺序保存下列数据
//*                     NumValidCuts I4
//*                     Elevation[]  NumValidCuts*R4 BYTES
//*                     WData[][][]  NumValidCuts*MaxRads*WGates*R4 BYTES     
//* 4)对本程序作很小的改动，就可以分别读取C波段雷达基数据和S波段雷达基数据，
//*   在DataBase.h头文件中，有如下几行代码，用来定义距离库长：
//*     ... ...
//*		SA，SB雷达的距离库数
//*		#define RGates  460				//S BAND 反射率距离库数
//*		#define VGates  920				//S BAND 速度距离库数
//*		#define WGates  920				//S BAND 谱宽距离库数

//*		CB 雷达的距离库数
//*		#define RGates  800				//C BAND 反射率距离库数
//*		#define VGates  1600			//C BAND 速度距离库数
//*		#define WGates  1600			//C BAND 谱宽距离库数   
//*     ... ...       
//*   默认情况下可读S波段(SA,SB)雷达基数据, 定义C波段距离库长的代码被注释,
//*   如果需要读C波段(CB)雷达基数据,只需使定义C波段距离库长的代码有效,
//*   注释定义S波段距离库长的代码
//*************************************************************************//
//* MAIN函数，
//* 1)在其中设置基数据文件名，调用ReadBaseData函数，将基数据读出到数组中。
//* 2)调用SavedataIntoFiles函数，将基数据数组中的数据分别保存到3个文件中。
void main()
{
	char filename[]=/*"archive2.001";*/"2006062510.43A";//"ARCHIVE2.189";//
	bool res;

	printf("Begin reading...\n");
	res=ReadBaseData(filename);

	printf("Begin saving...\n");
	if(res)
		SavedataIntoFiles();

	printf("Complete!\n");

	return;
}

//**************************************************************//
//* 名称：ReadBaseData(char filename[80])
//* 类型：布耳型，若读文件操作失败（体扫不完整等），返回FALSE。
//* 功能：从基数据文件中读出仰角，反射率，速度，谱宽数据，保存到
//*       数组中。
//* 参数：filename，基数据文件名
//**************************************************************//
bool ReadBaseData(char filename[80])
{
	FILE *fp=0;

	int FstBin,LstBin,BinNum;
	float CurAz,CurEl;
	int ElIndex,AzIndex,BnIndex;
	int ptrPos;

	size_t readSize;
	int fileEndFlag;

	bool VolBeg=false;
	bool VolEnd=false;

	bool RFlag,VFlag,WFlag;

	fp = fopen(filename,"rb");
	if(fp==0) return false;
	
	pOneRadial = 0;
	pOneRadial = (RADIALDATA*)malloc(sizeof(RADIALDATA));
	if(pOneRadial==0) return false;

	//Initialize array
	for(ElIndex=0; ElIndex<MaxCuts; ElIndex++)
	{
		for(AzIndex=0; AzIndex<MaxRads; AzIndex++)
		{
			for(BnIndex=0; BnIndex<RGates; BnIndex++)
				RData[ElIndex][AzIndex][BnIndex] = VALUE_INVALID;
			for(BnIndex=0; BnIndex<VGates; BnIndex++)
				VData[ElIndex][AzIndex][BnIndex] = VALUE_INVALID;
			for(BnIndex=0; BnIndex<WGates; BnIndex++)
				WData[ElIndex][AzIndex][BnIndex] = VALUE_INVALID;
		}//end az

		Elevation[ElIndex]=VALUE_INVALID;
	}//end el

	do
	{
		//Initialize flags
		RFlag=VFlag=WFlag=false;

		readSize = fread(pOneRadial,sizeof(RADIALDATA),1,fp);
		fileEndFlag = feof(fp);

		//Start a volume scan
		if(pOneRadial->RadialStatus == VOL_BEG) 
		{
			ElIndex=0;
			CurEl = float((pOneRadial->El/8.)*(180./4096.));
			Elevation[ElIndex] = CurEl;

			VolBeg = true;
			//output text information
//			TRACE("VCP number is %3d\n",pOneRadial->VcpNumber);
			printf("VCP number is %3d\n",pOneRadial->VcpNumber);

//			TRACE("Elevation %3d  (%5.2f Degree) start...\n",ElIndex+1,CurEl);
			printf("Elevation %3d  (%5.2f Degree) start...\n",ElIndex+1,CurEl);
		}
		//Find the beginning of the volume scan
		if(!VolBeg) continue;		

		//Start an elevation
		if(pOneRadial->RadialStatus == ELV_BEG)
		{
			CurEl = float((pOneRadial->El/8.)*(180./4096.));
			if(CurEl-Elevation[ElIndex] > 0.4)
			{//different elevation angle
				ElIndex++;
				Elevation[ElIndex] = CurEl;
			}
//			TRACE("Elevation %3d  (%5.2f Degree) start...\n",ElIndex+1,CurEl);
			printf("Elevation %3d  (%5.2f Degree) start...\n",ElIndex+1,CurEl);
		}

		//Start an elevation
		if(pOneRadial->RadialStatus == VOL_END)
			VolEnd=true;

//		if(pOneRadial->RadialStatus == ELV_BEG || pOneRadial->RadialStatus == VOL_BEG)


		//Calculate azimuth angle and Azimuth Index
		CurAz = float((pOneRadial->Az/8.)*(180./4096.));
		if(CurAz >= 360.) CurAz = CurAz-360.;
		AzIndex = int(CurAz+0.5);

		//what kind of data in this cut
		if(pOneRadial->PtrOfReflectivity !=0)	RFlag=true;
		if(pOneRadial->PtrOfVelocity !=0)		VFlag=true;
		if(pOneRadial->PtrOfSpectrumWidth !=0)	WFlag=true;

		//Save reflectivity data into the array
		if(RFlag)
		{
			//Get first bin, last bin, and number of bins
			FstBin = int(pOneRadial->RangeToFirstGateOfRef/pOneRadial->GateSizeOfReflectivity+0.5);
			BinNum = pOneRadial->GatesNumberOfReflectivity;
			if(FstBin<0)
			{
				BinNum = FstBin+BinNum;
				FstBin = -1*FstBin;
			}
			LstBin = FstBin + BinNum;
			ptrPos = pOneRadial->PtrOfReflectivity;  
			//Save data
			for(BnIndex=FstBin; BnIndex<LstBin; BnIndex++)
				RData[ElIndex][AzIndex][BnIndex] = DecodeRef(pOneRadial->Echodata[ptrPos+BnIndex]);
		}

		//Save velocity data into the array
		if(VFlag)
		{
			//Get first bin, last bin, and number of bins
			FstBin = int(pOneRadial->RangeToFirstGateOfDop/pOneRadial->GateSizeOfDoppler+0.5);
			BinNum = pOneRadial->GatesNumberOfDoppler;
			if(FstBin<0)
			{
				BinNum = FstBin+BinNum;
				FstBin = -1*FstBin;
			}
			LstBin = FstBin + BinNum;
			ptrPos = pOneRadial->PtrOfVelocity;   
			//Save data
			for(BnIndex=FstBin; BnIndex<LstBin; BnIndex++)
				VData[ElIndex][AzIndex][BnIndex] = 
				DecodeVel(pOneRadial->Echodata[ptrPos+BnIndex],pOneRadial->ResolutionOfVelocity);
		}

		//Save spectrum width data into the array
		if(WFlag)
		{
			//Get first bin, last bin, and number of bins
			FstBin = int(pOneRadial->RangeToFirstGateOfDop/pOneRadial->GateSizeOfDoppler+0.5);
			BinNum = pOneRadial->GatesNumberOfDoppler;
			if(FstBin<0)
			{
				BinNum = FstBin+BinNum;
				FstBin = -1*FstBin;
			}
			LstBin = FstBin + BinNum;
			ptrPos = pOneRadial->PtrOfSpectrumWidth;  
			//Save data
			for(BnIndex=FstBin; BnIndex<LstBin; BnIndex++)
				WData[ElIndex][AzIndex][BnIndex] = DecodeSpw(pOneRadial->Echodata[ptrPos+BnIndex]);
		}
	}while(fileEndFlag==0 && !VolEnd && !(readSize<1));
	
	if(!VolEnd && VolBeg)
	{
		fclose(fp);
		free(pOneRadial);
		printf("Error! Incomplete Volume Scan\n");
//		TRACE("Error! Incomplete Volume Scan\n");
		return false;
	}

	NumValidCuts = ElIndex+1;
	fclose(fp);
	free(pOneRadial);
	return true;
}

//**************************************************************//
//* 名称：DecodeRef(unsigned char code)
//* 类型：实型，返回解码后的反射率（DBZ）。
//* 功能：将读出的反射率数据解码。
//* 参数：code，无符号的BYTE型，读出的反射率编码值。
//**************************************************************//
float DecodeRef(unsigned char code)
{
	if(code==CODE_INVALID)		return VALUE_INVALID;
	else if(code==CODE_RANFOLD)	return VALUE_RANFOLD;
	else						
		return (float((code-2.)/2.-32.5));
}

//**************************************************************//
//* 名称：DecodeVel(unsigned char code)
//* 类型：实型，返回解码后的径向速度（M/S）。
//* 功能：将读出的径向速度数据解码。
//* 参数：code，无符号的BYTE型，读出的径向速度编码值；
//*       ResType，2字节整型，速度精度标记，精度不同，解码方式不同
//*            ResType = 2：精度为0.5 M/S
//*            ResType = 4：精度为1.0 M/S
//**************************************************************//
float DecodeVel(unsigned char code, short ResType)
{
	if(code==CODE_INVALID)		return VALUE_INVALID;
	else if(code==CODE_RANFOLD)	return VALUE_RANFOLD;
	else						
	{
		if(ResType==RES_POINT_FIVE) //0.5 m/s 
			return (float((code-2.)/2.-63.5));
		else
			return (float((code-2)-127.));
	}
}

//**************************************************************//
//* 名称：DecodeSpw(unsigned char code)
//* 类型：实型，返回解码后的谱宽（M/S）。
//* 功能：将读出的谱宽数据解码。
//* 参数：code，无符号的BYTE型，读出的谱宽编码值。
//**************************************************************//
float DecodeSpw(unsigned char code)
{
	if(code==CODE_INVALID)		return VALUE_INVALID;
	else if(code==CODE_RANFOLD)	return VALUE_RANFOLD;
	else						
		return (float((code-2.)/2.-63.5));
}

//**************************************************************//
//* 名称：SavedataIntoFiles()
//* 类型：布耳型，保存数组到文件中的操作如果出错，返回FALSE。
//* 功能：将反射率数组，速度数组，谱宽数组，仰角数组分别保存到3个文件中。
//*       RefArray.dat中保存仰角数，仰角，和体扫反射率数组的数据
//*       VelArray.dat中保存仰角数，仰角，和体扫速度数组的数据
//*       spwArray.dat中保存仰角数，仰角，和体扫谱宽数组的数据
//*       数据以二进制的方式保存，存放顺序：
//*       1）当前仰角层数，类型为4字节整型，长度为4字节
//*       2）有效层数的仰角，类型为4字节实型，长度为4*仰角层数 字节
//*       3）R，V或W体扫数据，类型为4字节实型，长度为4*MaxRads*Gates**仰角层数 字节
//* 参数：code，无符号的BYTE型，读出的谱宽编码值。
//**************************************************************//
bool SavedataIntoFiles()
{
	FILE *fpR=0,*fpV=0,*fpW=0;
	char fileR[] = "RefArray.dat";
	char fileV[] = "VelArray.dat";
	char fileW[] = "spwArray.dat";

	//Open Ref. File for Saving
	fpR = fopen(fileR,"wb");
	if(fpR==0) return false;
	//Save Ref. Array
	fwrite(&NumValidCuts,sizeof(NumValidCuts),1,fpR);
	fwrite(Elevation,sizeof(float),NumValidCuts,fpR);
	fwrite(RData,sizeof(float),NumValidCuts*MaxRads*RGates,fpR);
	fclose(fpR);

	//Open Ref. File for Saving
	fpV = fopen(fileV,"wb");
	if(fpV==0) return false;
	//Save Ref. Array
	fwrite(&NumValidCuts,sizeof(NumValidCuts),1,fpV);
	fwrite(Elevation,sizeof(float),NumValidCuts,fpV);
	fwrite(VData,sizeof(float),NumValidCuts*MaxRads*VGates,fpV);
	fclose(fpV);

	//Open Ref. File for Saving
	fpW = fopen(fileW,"wb");
	if(fpW==0) return false;
	//Save Ref. Array
	fwrite(&NumValidCuts,sizeof(NumValidCuts),1,fpW);
	fwrite(Elevation,sizeof(float),NumValidCuts,fpW);
	fwrite(WData,sizeof(float),NumValidCuts*MaxRads*WGates,fpW);
	fclose(fpW);
	return true;
}