fits_sci.cpp

来自「本源码是一个用于卫星数据处理的程序」· C++ 代码 · 共 729 行 · 第 1/3 页

			}
			
			
		}
		
		*fptr = fp;
		return 0;
}

/*****************************************************
** 函数名： write_fits_hx_sci
** 输入：*fptr，data，time，n，ptime 
**        fptr    --- fits文件结构（structure used to store basic HDU information）
**		  data    --- 数组，存放能量和通道号
**		  time    --- 时间结构
**        n
**        ptime
** 
** 输出： 
**        返回0函数执行成功
** 调用说明：	
**				FindZT() .......得到对应的姿态数据

** 功能描述： 按照需求把需要的字段写入fits文件中
****************************************************************************************/
int write_fits_hx_sci( fitsfile *fptr, const char *data, hx_t &time, int n,  const timespec_t *ptime )
{
	if( fptr == NULL || ptime == NULL || data == 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;
		const ZT_t * Minzt;
		ZT_t zt;
		Minzt = FindZT(get_sci_data(), ptime);
		ASSERT( Minzt != NULL );
		if( Minzt->utc.time < ptime->time-100 || Minzt->utc.time > ptime->time+100 )//事例时间不在所给姿态文件中记录的时间范围 
		{
			Minzt = &zt;//初始值已经设置成0，见ZT_T的定义
		}
		
		double time_b = ptime->time;//事例的秒以上部分的对应时间
		double time_b_usec = ptime->usec / 1000000.0;//事例对应的微妙时间
		time_b = time_b + time_b_usec;
		
		unsigned char var[10];
		memset( var, 0, sizeof(var) );
		u_long t = htonl(time.time);
		memcpy( var, &t, 4 );
		t = htonl(time.usec);
		memcpy( var+4, &t, 4 );
		
		/***********************逐列写入每一行***************************/
		fits_write_col(fptr, TBYTE, 1, firstrow, firstelem, 7, (void *)var, &status);
		fits_write_col(fptr, TDOUBLE, 2, firstrow, firstelem, 1, (void *)&time_b, &status);
		fits_write_col(fptr, TDOUBLE, 3, firstrow, firstelem, 1, (void *)&Minzt->trans_lon, &status);
		fits_write_col(fptr, TDOUBLE, 4, firstrow, firstelem, 1, (void *)&Minzt->trans_lon_angle, &status);
		fits_write_col(fptr, TDOUBLE, 5, firstrow, firstelem, 1, (void *)&Minzt->trans_lat, &status);
		fits_write_col(fptr, TDOUBLE, 6, firstrow, firstelem, 1, (void *)&Minzt->trans_lat_angle, &status);
		fits_write_col(fptr, TBYTE, 7, firstrow, firstelem, 1, (void *)data, &status);
		fits_write_col(fptr, TBYTE, 8, firstrow, firstelem, 1, (void *)(data+1), &status);

	}
	
	return 0;
}

/*****************************************************
** 函数名： open_fits_sx_sci
** 输入： hxmtfile
**        hxmtfile    --- 文件名结构 
** 输出： **fptr
** **fptr --- 通过**fptr返回一个结构(structure used to store basic HDU information)
** 调用说明：	
**				get_fits_filename() .......得到生成fits文件的保存路径
**				DirExist1()................判断路径是否存在
**				CreateFilePath()...........建立文件路径
**				time2strUTC_fitshead().....把秒时间转换成对应的年月日时分秒格式
** 功能描述： 输入一个文件名结构变量，生成fits头并写入相关信息
****************************************************************************************/
/*******************函数注释参照上面的高能X射线**********************************/
int open_fits_sx_sci( 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;
		
		timespec_t t;
		time( &t.time );
		t.usec = 0;
		char buf[40];
		time2strUTC_fitshead( buf, sizeof(buf), &t );
		
		/* 写入需要的关键字 */
		if( modetype ){
		
		if ( fits_update_key(fp, TSTRING , "HDUNAME", "PRIMARY", "Name of the data", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "ORIGIN", "HDPS", "Source of FITS file", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "CREATOR", "Fits_mfc.exe", "tool that create this output", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "CHECKSUM", "2", "HUD checksum", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "FILEDATE", buf, "Date and time of file creaton", &status) ) 
			return status;		
		if ( fits_update_key(fp, TSTRING , "DATE-OBS", buf, "Date and time of observation start", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "DATE-END", buf, "Date and time of observation stop", &status) )
			 return status;		
		if ( fits_update_key(fp, TSTRING , "OBTSTART", "asda", "Start orbit number", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "OBT-END", "adad", "End orbit number", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "MJDREF", "55197", "MJD zero point for times", &status) ) 
			return status;
		if ( fits_update_key(fp, TSTRING , "TSTART", "adas", "Observation start time", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "TSTOP", "asdas", "Observation end time", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "TELESCOP", "HXMT", "telescope", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "INSTRUME", "HXT", "instrument", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "OBS_MODE", "SURVEY", "Observation mode", &status) ) 		
			 return status;
		if ( fits_update_key(fp, TSTRING , "OBS_ID", "", "Observation id", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "SEQ_NUM", "", "Sequence number", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "FILENAME", (void *)fits_name.c_str(), "Name of the current data file", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "FILESIZE", "2", "Size of the current data file", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "REVISION", "1", "Edition of the data package", &status) ) 		
			 return status;

		status = 0;
		int tfields   = 8;       /* table will have 2 columns */
		long nrows    = 0;
		char extname[] = "SCID_Binary";           /* extension name */
		char *ttype[] = { "time", "time_d", "trans_lon", "trans_lon_angle", "trans_lat", "trans_lat_angle", "energy", "detector" };
		char *tform[] = { "7B", "1D", "1D", "1D", "1D", "1D", "1B", "1B" };
		char *tunit[] = { "\0", "s", "\0", "\0", "\0", "\0", "\0", "\0" };
		if ( fits_create_tbl( fp, BINARY_TBL, nrows, tfields, ttype, tform, tunit, extname, &status) ){
			 return status;
		}
		if ( fits_update_key(fp, TSTRING , "CHECKSUM", "0", "HUD checksum updated", &status) )
			return status;
		if ( fits_update_key(fp, TSTRING , "DATASUM", "0", "data unit chechsum", &status) )
			return status;
		if ( fits_update_key(fp, TSTRING , "FILEDATE", buf, "Date and time of file creaton", &status) ) 
			return status;		
		if ( fits_update_key(fp, TSTRING , "DATE-OBS", buf, "Date and time of observation start", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "DATE-END", buf, "Date and time of observation stop", &status) )
			 return status;		
		if ( fits_update_key(fp, TSTRING , "TSTART", "asda", "observation start time", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "TSTOP", "adad", "observation stop time", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "MISSION", "HXMT", "Mission", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "TELESCOP", "HXMT", "telescope", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "INSTRUME", "HXT", "instrument", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "SEQ_NUM", "", "Sequence number", &status) )
			 return status;

		}
		else{
		
		if ( fits_update_key(fp, TSTRING , "HDUNAME", "PRIMARY", "Name of the data", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "ORIGIN", "HDPS", "Source of FITS file", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "CREATOR", "Fits_mfc.exe", "tool that create this output", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "CHECKSUM", "0", "HUD checksum", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "FILEDATE", buf, "Date and time of file creaton", &status) ) 
			return status;		
		if ( fits_update_key(fp, TSTRING , "DATE-OBS", buf, "Date and time of observation start", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "DATE-END", buf, "Date and time of observation stop", &status) )
			 return status;		
		if ( fits_update_key(fp, TSTRING , "MJDREF", "55197", "MJD zero point for times", &status) ) 
			return status;
		if ( fits_update_key(fp, TSTRING , "TSTART", "adas", "Observation start time", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "TSTOP", "asdas", "Observation end time", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "DURATION", "0", "observation length", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "PI_ID", "", "", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "PI_NAME", "", "", &status) ) 		
			 return status;
		if ( fits_update_key(fp, TSTRING , "RA_CEN", "", "", &status) )
			 return status; 
		if ( fits_update_key(fp, TSTRING , "DEC_CED", "", "", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "OBJECT", "", "", &status) )
			 return status;
		if ( fits_update_key(fp, TSTRING , "AREA", "", "", &status) )