fits_eph.cpp

本源码是一个用于卫星数据处理的程序

/***************************************************************************
    Fits_eph.cpp:星历文件一级产品生成所调用的函数和方法

	LoadEPH      －－ 从文件中读取相应值并赋给所定义的结构体中的各个变量
    open_fits_eph   －－ 写 FITS binary table extension
	write_fits_eph －－ 在 bianry table 后插入一行并写入相应的数据值
	close_fits_eph       －－关闭Fits 文件
	
***************************************************************************/

/********此文件中函数前都有详细说明，函数中具体实现参照fits_asf.cpp中的注释******/

#include <stdio.h>
#include <vector>
#include "stdafx.h"
#include "ccsds_pro.h"
#include "fits_eph.h"
/*****************************************************
** 函数名： LoadEPH
** 输入： hxmtfile, s
**        hxmtfile    --- 文件名结构 
**        s           --- 存放多个容器的结构体
** 输出： 
**           函数执行成功返回0
** 调用说明：	
**			strUTC2time() .......把年月日是时分秒格式的时间字符串转换成对应的秒			
** 功能描述：输入一个星历文件，把文件中的每行数据的各个字段读入容器中
****************************************************************************************/
int LoadEPH( const hxmt_file_t *hxmtfile , struct sci_intern *s )
{
	FILE *fp = NULL;
	char *buffer = NULL;
	if( !hxmtfile || (fp = fopen(hxmtfile->filename.c_str(), "r")) == NULL  )
		return -1;

	if( (buffer = (char *)malloc(FREAD_BUFFER_SIZE)) == NULL )	{
		fclose(fp); 
		return -ERROR_NOT_ENOUGH_MEMORY;
	}

	while( fgets(buffer, FREAD_BUFFER_SIZE, fp) )
	{
		EPH_t eph;
		float f;
		char *p = NULL;
		if( (p = strchr(buffer, '\t')) == NULL )
			continue;
		*p = '\0';
		eph.sid = buffer; 

		sscanf( ++p, "%d", &eph.cir );
		if( (p = strchr( p, '\t')) == NULL )	continue;
		strUTC2time( ++p, &eph.utc );

		if( (p = strchr( p, '\t')) == NULL )	continue;
		strUTC2time( ++p, &eph.localtime );

		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &eph.longitude );

		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.latitude = f;

		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.x = f;
		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.y = f;
		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.z = f;

		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.vx = f;
		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.vy = f;
		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.vz = f;

		if( (p = strchr( p, '\t')) == NULL )	continue;
		sscanf( ++p, "%f", &f );
		eph.height = f;

		s->eph_table.push_back( eph );			
	}
	
	free( buffer );
	fclose( fp );

	return 0;	
}
/*****************************************************
** 函数名： process_eph
** 输入： hxmtfile, 
**        hxmtfile    --- 文件名结构 
** 输出： 
**           函数执行成功返回0
** 调用说明：	
**			get_sci_data() .......
**			get_fits_filename()...	
**			close_fits_eph() ......
**			open_fits_eph() .......
**			write_fits_eph() ......	
** 功能描述: 输入一个星历文件，处理星历文件中的数据，并生成相应的fits文件
****************************************************************************************/
int process_eph( const hxmt_file_t *hxmtfile )
{
	int modetype;
	std::string fits_name;
	fitsfile *fptr = NULL;
	std::string pathfilename;
	hxmt_file_t hxmt = *hxmtfile;
	std::string prev_fitsname;
	
	std::vector<EPH_t>::const_iterator end = get_sci_data()->eph_table.end();
	for( std::vector<EPH_t>::const_iterator ite = get_sci_data()->eph_table.begin();
			ite != end; ++ite )
	{
		hxmt.start.time = ite->utc.time;
		hxmt.start.usec = ite->utc.usec;
		{
			std::string new_fitsname;
		
			get_fits_filename( &hxmt, new_fitsname, pathfilename, modetype, fits_name );
			if( new_fitsname != prev_fitsname ) {
				prev_fitsname = new_fitsname;
				close_fits_eph( fptr );
				fptr = NULL;
				if( 0 != open_fits_eph( &hxmt, &fptr) ) {
					return -3;
				}
			}
		}

		write_fits_eph( fptr, ite, 1 );
	}
	
	close_fits_eph( fptr );
	return 0;
}
/*****************************************************
** 函数名： open_fits_eph
** 输入： hxmtfile, 
**        hxmtfile    --- 文件名结构 
** 输出： *fptr
**        *fptr       ---
** 调用说明：	
**			get_fits_filename() .......得到文件名和保存路径以及模式
**			DirExist1()         .......判断路径是否存在
**			CreateFilePath()    .......创建一个新路径
**			time2strUTC_fitshead() ....秒时间转换成对应的年月日时分秒格式的字符串
** 功能描述: 输入一个星历文件，生成fits文件，并写入fits头中的关键字，通过fptr返回生成fits文件的结构
****************************************************************************************/
int open_fits_eph( const hxmt_file_t *hxmtfile, fitsfile **fptr )
{
	fitsfile *fp;
	int status = 0;
	std::string filename;
	std::string pathname;
	int modetype;
	std::string fits_name;
	if( get_fits_filename(hxmtfile, filename, pathname, modetype, fits_name) )
		return -1;
	{
		extern std::string product_dest_dir;
		std::string path = product_dest_dir + pathname + "\\";
		if( !DirExist1(path.c_str()) )		{
			extern FILE *arch_filep;
			pathname += "\n";
			fwrite(pathname.c_str(), sizeof(char), pathname.length(), arch_filep);
		}
		
	}
	if( !CreateFilePath(filename.c_str()) )
		return -2;
	if( fits_open_file(&fp, filename.c_str(), READWRITE, &status) )
	{
		status = 0;
		if( fits_create_file(&fp, filename.c_str(), &status) )
			return status;

		if( fits_create_img(fp,  BYTE_IMG, 0, NULL, &status) )
			return status;
//		ffpcom(*fptr, "ALL TIME IN THIS FILE IS UTC", &status); /*增加必要的注释*/
		
		timespec_t t;
		time( &t.time );
		t.usec = 0;
		char buf[40];
		time2strUTC_fitshead( buf, sizeof(buf), &t );
		/* 写入需要的关键字 */
		if ( fits_update_key(fp, TSTRING , "DATE_CRT", buf, "FITS file creation date", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "TELESCOP", "HXMT", "Telescope or mission name", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "ORIGIN", "HDPS", "Origin of FITS file", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "HDPSLEVL", "RAW", "HDPS level of data processing", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "CREATOR", "POAC", "program that create this FITS file", &status) )                                                                                                                                   
			 return status;

		status = 0;
		int tfields   = 13;       /* table will have 13 columns */
		long nrows    = 0;
		char extname[] = "ORBIT_Binary";           /* extension name */
		char *ttype[] = { "SID", "CIR", "UTC_TIME", "LOCAL_TIME", "LONG", "LAT", "X", "Y", "Z", "VX", "VY", "VZ", "HIGH" };
		char *tform[] = { "2A", "1J", "27A", "27A", "1E", "1D", "1D", "1D", "1D", "1D", "1D", "1D", "1D" };
		char *tunit[] = { "\0", "\0", "\0", "\0", "deg", "deg", "\0", "\0", "\0", "\0", "\0", "\0", "m" };
		if ( fits_create_tbl( fp, BINARY_TBL, nrows, tfields, ttype, tform,	tunit, extname, &status) ) {
			fits_close_file( fp, &status );
			return status;
		}
	}
	
	*fptr = fp;
	return 0;
}
/*****************************************************
** 函数名： write_fits_eph
** 输入： fptr, eph，n
**        hxmtfile    --- 文件名结构 
**        eph          --- 星历数据结构
**		  n           ---
** 输出： 
**        函数执行成功返回0
** 调用说明：	
**			time2strUTC_fitshead() ....秒时间转换成对应的年月日时分秒格式的字符串
** 功能描述: 输入先前生成的fits文件结构和星历数据，把星历数据写入fits文件中
****************************************************************************************/
int write_fits_eph( fitsfile *fptr, const EPH_t *eph, int n )
{
	if( fptr == NULL )
		return -1;

	int status = 0;
	unsigned int nrows = 0;
	int hdutype;
    if ( fits_movabs_hdu(fptr, 2, &hdutype, &status) ) 
		 return status;
	if( fits_read_key(fptr, TINT, "NAXIS2", &nrows, NULL, &status) )
		return status;

//	for( int i = 0; i < n; i++, nrows++ )
	{
		if (fits_insert_rows(fptr, nrows, 1, &status))
			  return status;
    
		int firstrow = nrows + 1;
		int firstelem = 1;

		void *var[1];
		var[0] = (void *)(eph->sid.c_str());
		status = 0;
		fits_write_col(fptr, TSTRING, 1, firstrow, firstelem, n, var, &status);
		fits_write_col(fptr, TINT, 2, firstrow, firstelem, n, (void *)&eph->cir, &status);
		char strUTC[40];
		time2strUTC( strUTC, sizeof(strUTC), &eph->utc );
		var[0] = (void *)strUTC;
		fits_write_col(fptr, TSTRING, 3, firstrow, firstelem, n, var, &status);
		time2strUTC( strUTC, sizeof(strUTC), &eph->localtime );
		fits_write_col(fptr, TSTRING, 4, firstrow, firstelem, n, var, &status);
		fits_write_col(fptr, TFLOAT, 5, firstrow, firstelem, n, (void *)&eph->longitude, &status);
		fits_write_col(fptr, TDOUBLE, 6, firstrow, firstelem, n, (void *)&eph->latitude, &status);
		fits_write_col(fptr, TDOUBLE, 7, firstrow, firstelem, n, (void *)&eph->x, &status);
		fits_write_col(fptr, TDOUBLE, 8, firstrow, firstelem, n, (void *)&eph->y, &status);
		fits_write_col(fptr, TDOUBLE, 9, firstrow, firstelem, n, (void *)&eph->z, &status);
		fits_write_col(fptr, TDOUBLE, 10, firstrow, firstelem, n, (void *)&eph->vx, &status);
		fits_write_col(fptr, TDOUBLE, 11, firstrow, firstelem, n, (void *)&eph->vy, &status);
		fits_write_col(fptr, TDOUBLE, 12, firstrow, firstelem, n, (void *)&eph->vz, &status);
		fits_write_col(fptr, TDOUBLE, 13, firstrow, firstelem, n, (void *)&eph->height, &status);
		
	}
	return 0;
}

int close_fits_eph( fitsfile *fptr )//此函数用于关闭fits文件
{
	int status = 0;
	if( fits_close_file( fptr, &status) )
		return status;
	return 0;
}


int get_eph_orbit_number( const timespec_t *time )
{
	time_t t = time->time;
	time_t off = 0xfffffff;

	std::vector<EPH_t>::const_iterator end = get_sci_data()->eph_table.end();
	std::vector<EPH_t>::const_iterator dest = end;
	for( std::vector<EPH_t>::const_iterator ite = get_sci_data()->eph_table.begin();
			ite != end; ++ite )
	{
		if( t - ite->utc.time < off ) {
			off = t - ite->utc.time;
			dest = ite;
		}
	}
	if( dest == end )
		return -1;
	return dest->cir;
}