sacio.cpp

提供交互的方式选取合适的震相

#include <string.h>
#include <assert.h>

#include "SacIO.h"

using namespace std;

// colorref 0xaabbggrr
#define AXIS_PEN_CLR	0x00808080
#define IASP91_PEN_CLR	0x0080ff80

struct sac_header CSacFile::sac_null = {
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345., -12345., -12345., -12345., -12345.,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
-12345, -12345, -12345, -12345, -12345,
{ '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ',' ',' ',' ',' ',' ',' ',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }, { '-','1','2','3','4','5',' ',' ' },
{ '-','1','2','3','4','5',' ',' ' }
};

bool CSacFile::SetFileName( const char * in ) {

	assert( NULL != in );
	if( NULL == in )
		return false;

	strcpy( m_strFileName, in );
	
	return true;
}

const char * CSacFile::GetFileName( void) {

	return m_strFileName;

}

bool CSacFile::LoadData( const char * inFileName ) {

	char filename[FILENAME_MAX];
	size_t re;

	if( NULL == inFileName )
		strcpy( filename, m_strFileName );
	else {
		strcpy( filename, inFileName );
		strcpy( m_strFileName, inFileName );
	}

	if( !ReadHead( filename ) )
		return false;

	FILE * fp = fopen( filename, "rb" );
	if( NULL == fp )
		return false;

	fseek( fp, sizeof(struct sac_header), SEEK_SET );

	try {
		if( NULL != m_pData )
			delete []m_pData;
		m_pData = new float[m_header.npts];
	}
	catch(...) {
		printf("Failed to memory operation.\n");
		throw;
	}
	re = fread( m_pData, sizeof(float), m_header.npts, fp );

	assert( (size_t)m_header.npts == re );

	fclose(fp);

	m_eDrawIndex = m_header.npts - 1;

	return true;
}

bool CSacFile::ReadHead( const char * inFileName ) {

	char filename[FILENAME_MAX];
	size_t re;

	if( NULL == inFileName )
		strcpy( filename, m_strFileName );
	else {
		strcpy( filename, inFileName );
		strcpy( m_strFileName, inFileName );
	}

	char * pStr;
	pStr = strrchr( filename, '.' );
	assert( strlen(pStr) > 1 );
	switch(pStr[1]) {
	case 'e':
	case 'E':
		m_cmpnm = 'e';
		break;
	case 'n':
	case 'N':
		m_cmpnm = 'n';
		break;
	case 'z':
	case 'Z':
		m_cmpnm = 'z';
		break;
	}

	FILE * fp = fopen( filename, "rb" );
	if( NULL == fp )
		return false;

	re = fread( &m_header, sizeof(struct sac_header), 1, fp);

	assert( 1 == re );

	fclose(fp);

	m_npts = m_header.npts;
	m_delta = m_header.delta;
	m_StLa = m_header.stla;
	m_StLo = m_header.stlo;
	m_EvLa = m_header.evla;
	m_EvLo = m_header.evlo;
	m_gcArc = m_header.gcarc;
	if( m_header.evdp == sac_null.evdp )
		m_EvDp = 10;
	else
		m_EvDp = m_header.evdp;

	m_begin = m_header.b;
	m_end = m_header.e;
	if( m_header.nzyear == sac_null.nzyear )
		m_nzYear = 1998;
	else
		m_nzYear = m_header.nzyear;
	if( m_header.nzjday == sac_null.nzjday )
		m_nzjDay = 1;
	else
		m_nzjDay = m_header.nzjday;
	if( m_header.nzhour == sac_null.nzhour )
		m_nzHour = 1;
	else
		m_nzHour = m_header.nzhour;
	if( m_header.nzmin == sac_null.nzmin )
		m_nzMinute = 1;
	else
		m_nzMinute = m_header.nzmin;
	if( m_header.nzsec == sac_null.nzsec )
		m_nzSecond = 1;
	else
		m_nzSecond = m_header.nzsec;
	if( m_header.nzmsec == sac_null.nzmsec )
		m_nzMSec = 1;
	else
		m_nzMSec = m_header.nzmsec;

	string::size_type pos;
	string strEvtNm, strStaNm;
	strEvtNm = filename;
	pos = strEvtNm.rfind('/');
	strcpy( m_strEvtNm , strEvtNm.substr( pos+1,  strEvtNm.rfind('.') - pos-1).c_str() );

	strStaNm = filename;
	pos = strStaNm.rfind('/');
	pos = strStaNm.rfind('/',pos-1);
	strcpy(m_strStaNm , strStaNm.substr(pos+1,  strStaNm.rfind('/') - pos-1).c_str() );
	
	m_FirstPhaseTm = m_header.a;
	m_dataMean = m_header.depmen;
	m_dataMin = m_header.depmin;
	m_dataMax = m_header.depmax;

	return true;
}

#if defined(_AFXDLL)

void CSacFile::Draw(CDC &dc, CRect drawRc, bool bShowIASP91) {

	if( NULL == m_pData )
		return;
	if( m_npts <= 0 )
		return;

	CDC &memDC = dc;
	
//	memDC.CreateCompatibleDC(&dc);
	CPen *pOldPen;
	CPen drawPen;
	CPen axisPen;
	CPen iasp91Pen;
	LOGPEN logpen;

	drawPen.CreateStockObject(BLACK_PEN);
	logpen.lopnColor = AXIS_PEN_CLR;
	logpen.lopnStyle = PS_DOT;
	logpen.lopnWidth = CPoint(1,1);
	axisPen.CreatePenIndirect( &logpen );
	logpen.lopnColor = IASP91_PEN_CLR;
	iasp91Pen.CreatePenIndirect( &logpen );

	CPoint pt = drawRc.CenterPoint();
	const int iCenter = pt.y;
	pt.x = drawRc.left;

	//draw out the file name
	string filename = m_strFileName;
	int pos = filename.rfind('/');
	pos = filename.rfind('/', pos-1);
	filename = filename.substr( pos+1, filename.length());
	int nBkMode = memDC.SetBkMode( TRANSPARENT );
	memDC.TextOut( pt.x, 
				pt.y,
				filename.c_str(),
				filename.length() );
	memDC.SetBkMode( nBkMode);

	// draw axis
	pOldPen = memDC.SelectObject(&axisPen);
	memDC.MoveTo(drawRc.right, iCenter);
	memDC.LineTo(pt);

	// draw data
	double vScale;
	int bDraw = m_bDrawIndex;
	int eDraw = m_eDrawIndex;

	vScale = drawRc.Height() * m_vScale/ ( m_dataMax - m_dataMin );
	memDC.SelectObject(&drawPen);
	while(bDraw <= eDraw) {
		pt.x += long(m_hScale);
		pt.y = iCenter + long( (m_pData[bDraw] - m_dataMean) * vScale);
		memDC.LineTo(pt);
		if( pt.x > drawRc.right ) {
			break;
		}
		++bDraw;
	}

	// draw iasp91 marker
	if( bShowIASP91 ) {
		memDC.SelectObject(&iasp91Pen);
		if( m_IASP91.size() > 0 ) {
			map<float,string>::const_iterator iter = m_IASP91.begin();

			// indepth数据没有震源发震时刻数据，暂时用与首波的到时差来定位震相
			// draw first arrival phase name
			//转换到数据时间
			const long iasp91Begin = iter->first/m_delta; // 第一个震相的到时
			const long DtEvBegin = (m_FirstPhaseTm - m_begin)/m_delta;

			nBkMode = memDC.SetBkMode( TRANSPARENT );

			while( iter != m_IASP91.end() ) {
				pt.x = long( (-m_bDrawIndex + (iter->first/m_delta - iasp91Begin) + DtEvBegin ) * m_hScale );
				if( pt.x > drawRc.right ) {
					break;
				}
				else if( pt.x > drawRc.left ) {
					pt.y = drawRc.top;
					memDC.MoveTo(pt);
					pt.y = drawRc.bottom;
					memDC.LineTo(pt);
					//if( 0 == iter->first.compare("SKSac   ") ) {
					{
						memDC.TextOut( pt.x, 
							pt.y - 20,
							iter->second.c_str(),
							iter->second.size() );
					}
				}

				++iter;
			}
			memDC.SetBkMode( nBkMode);
		}
	}

	// blit the final draw to screen.
//	dc.BitBlt( drawRc.left, drawRc.top, drawRc.Width(), drawRc.Height(), &memDC,
//				drawRc.left, drawRc.top, SRCCOPY );

}

long CSacFile::GetPhaseDrawPos( int indexPhase ) {

	if( m_IASP91.size() > 0 ) {
		map<float, string>::const_iterator iter = m_IASP91.begin();
	
		//转换到数据时间
		const long iasp91Begin = iter->first/m_delta; // 第一个震相的到时
		const long DtEvBegin = (m_FirstPhaseTm - m_begin)/m_delta;
		while( iter != m_IASP91.end() ) {
			if( 0 == indexPhase-- ) {				
				return ((iter->first/m_delta - iasp91Begin) + DtEvBegin);
			}
			++iter;
		}

	}

	return -1;
}

CTime & CSacFile::GetTime( CTime & time, double & delta, int scrnOffset ) {
	double dIndex =  (m_bDrawIndex + scrnOffset/m_hScale) * m_delta;// 转换到秒
	dIndex += m_begin; // 加上SAC文件的自变量开始值的位移(秒)
	long bIndex = dIndex;
	delta = dIndex - bIndex;
	int sec = bIndex%60;
	int min = (bIndex/60)%60;
	int hour = (bIndex/3600)%24;
	long day = (bIndex/3600/24);
	CTimeSpan tmSpan(day, hour, min, sec);

	GetFileBaseTime(time);

	time += tmSpan;
	return time;
}
void CSacFile::Align2Time( const CTime & time, double delta ) {
	//find out the month and the day
	int mon,day;
	day = 0;
	for( mon = 1; mon<= 12; mon++) {
		day += daysInMonth( mon, m_nzYear );
		if( day >= m_nzjDay ) {
			day = m_nzjDay - (day - daysInMonth( mon, m_nzYear ));
			break;
		}
	}

	CTime myTime(m_nzYear, mon, day, m_nzHour, m_nzMinute, m_nzSecond);
	
	CTimeSpan tmSpan( time - myTime );
	// 需要减去SAC文件的自变量开始值的位移(秒)
	const double seconds = tmSpan.GetTotalSeconds() + delta - m_begin;

	m_bDrawIndex = seconds/m_delta; // 换算到数据的时间单位
	if( m_bDrawIndex < 0 || m_bDrawIndex > m_npts )
		m_bDrawIndex = m_npts - 1;
}

#endif