process_data.c

解吸SEED格式的源代码

        
        	}
		else
		dump_LOG_rec(input_data_ptr + (input_data_hdr->bod - 8), 
				LRECL - input_data_hdr->bod,
				&this_time,
				this_station,
				this_channel,
				network_code,
				this_location);

		return;

	}


	this_nsamples = input_data_hdr->nsamples;

	/* if number of samples is zero, do nothing unless we are
	 * outputing miniseed or seed
	 */

	if ((this_nsamples == 0) && ((!mini_flag) && (!Seed_flag)))
		return;

	if ((this_nsamples == 0) && (mini_flag || Seed_flag))
		if (strip_flag)
			return;

	/* figure out current block ending time */
    	if (this_sample_rate != 0)
        	duration = ((double)this_nsamples*TIME_PRECISION)/this_sample_rate;
    	else duration = 0;
 
    	newtime = timeadd_double(this_time, duration);

    	/* if block time span is not overlapping a selected span,
     	 *  forget this block
     	*/ 

	/* if times needed are from summary file, match up the 
	 * stn, chn with times
	 */

	if (!read_summary_flag)
	{
		/* else use the normal entry (from keyboard) */
    		if (!chk_time (this_time, newtime))
		{
        		return;
		}
	}
	else
	{
		if (!chk_summary(this_time, newtime))
			return;
	}


/*                 +=======================================+                 */
/*=================| accumulate various descriptive info   |=================*/
/*                 +=======================================+                 */

	/* figure out current ending time of seismogram */
	if (last_sample_rate != 0)
		duration = (((double)last_nsamples * TIME_PRECISION) / last_sample_rate);

	else duration = 0;

	newtime = timeadd_double(last_time, duration);

/*                 +=======================================+                 */
/*=================|          determine continuity         |=================*/
/*                 +=======================================+                 */

	/* determine whether this is station-, channel-, and time-continuous */
	/* with previous data block */

	/* if no samples and we are this far, user is asking for mini or
	 * seed file output. Don't clip on zero sample record
	 */

	if (this_nsamples == 0)
		continuous = TRUE;
	else
	if ((strcmp (data_hdr->station, this_station) == 0) &&
		(strcmp (data_hdr->channel, this_channel) == 0) &&
		(strcmp (data_hdr->location, this_location) == 0) &&
		(timetol (newtime, this_time, last_nsamples, this_sample_rate) == 0))
	{
		continuous = TRUE;
	}
	else 
	{
		continuous = FALSE;
		output_overflow = FALSE;

	}

	/*
	printf ("%d\n", continuous);
	*/

	/* save information for next time around */
	last_time.year    = (int) this_time.year;
	last_time.day     = (int) this_time.day;
	last_time.hour    = (int) this_time.hour;
	last_time.minute  = (int) this_time.minute;
	last_time.second  = (int) this_time.second;
	last_time.fracsec = (int) this_time.fracsec;
	last_nsamples     = this_nsamples;
	last_sample_rate  = this_sample_rate;

	/*
	printf ("%f ", last_sample_rate);
	*/


/*                 +=======================================+                 */
/*=================|if continuous, append to current seisgm|=================*/
/*                 +=======================================+                 */

	/* DEBUG */
	/*
	if (continuous)
		fprintf (stderr, "%d: %d + %d = %d < %d\n",
		input.recordnumber, data_hdr->nsamples, this_nsamples,
		data_hdr->nsamples + this_nsamples, MAX_DATA_LENGTH);
	*/

	if (continuous)
	{
		/* check size of data; if length of seismogram is too large,
		 * fake a station/component/time discontinuity
		 * and add the current recordnumber to the desired seismogram list,
		 * which will cause (partial) seismogram files to be written */
		if ((data_hdr->nsamples + this_nsamples) > (seis_buffer_length)/data_hdr->num_mux_chan)
		{
			continuous = FALSE;
			output_overflow = TRUE;
			fprintf (stderr, "WARNING (process_data):  ");
			fprintf (stderr, "too much data (burp!).\n");
			fprintf (stderr, 
				"\tMaximum data length exceeded for station %s, channel %s.\n",
				this_station, this_channel);
			fprintf (stderr, "\tDue to indigestion, ");
			fprintf (stderr, "this seismogram will be split into sections.\n");
			fprintf (stderr, "\tTry increasing the Sample Buffer Length\n");
	
		}

		else 
		{
			output_overflow = FALSE;
			next_index = data_hdr->nsamples;

			/* update certain values in data_hdr */
			data_hdr->nsamples           += this_nsamples;
			data_hdr->activity_flags     |= input_data_hdr->activity_flags;
			data_hdr->io_flags           |= input_data_hdr->io_flags;
			data_hdr->data_quality_flags |= input_data_hdr->data_quality_flags;

			data_hdr->number_time_corrections += input_data_hdr->number_time_corrections;

			data_hdr->bod = input_data_hdr->bod;
			data_hdr->bofb = input_data_hdr->bofb;

		}
	}

/*                 +=======================================+                 */
/*=================| if not continuous, write out seisgm   |=================*/
/*                 +=======================================+                 */

	if ((!continuous) && (output) && (data_hdr->time.year != 0))
	{
		dump_seismic_buffer();    

	}

	/* data_hdr->time.year !=0 avoids spurious output at beginning of tape */

	current_channel = new_current_channel;
	current_station = new_current_station;


/*                 +=======================================+                 */
/*=================|  if not continuous, init new seisgm   |=================*/
/*                 +=======================================+                 */

	/* initialize for new seismogram */
	if (!continuous)
	{
		/* initialize for new seismogram */
		data_hdr->nsamples = 0;
		next_index = 0;
		seismic_data_ptr = seismic_data;
		start_record = input.recordnumber;

		/* find out if this seismogram is wanted */
		output = TRUE;

/*                 +=======================================+                 */
/*=================| build minimal header from known info  |=================*/
/*                 +=======================================+                 */

		/* build minimal data header from SEED fixed section of data header */
		convert_seedhdr (input_data_hdr);
		data_hdr->sample_rate = this_sample_rate;

 /* Added 4-30-92 (WAN) so that data_hdr has corrected start time */
		data_hdr->time.year    = this_time.year;
		data_hdr->time.day     = this_time.day;
		data_hdr->time.hour    = this_time.hour;
		data_hdr->time.minute  = this_time.minute;
		data_hdr->time.second  = this_time.second;
		data_hdr->time.fracsec = this_time.fracsec;

		memset((char *)&mini_data_hdr, 0, sizeof(struct mini_data_hdr));

		memcpy((char *)&mini_data_hdr.hdr, 
				(char *)(lrecord_ptr + 8), 
					sizeof(mini_data_hdr.hdr));

		mini_data_hdr.hdr.time.year		= data_hdr->time.year;
		mini_data_hdr.hdr.time.day		= data_hdr->time.day;
		mini_data_hdr.hdr.time.hour 	= data_hdr->time.hour;
		mini_data_hdr.hdr.time.minute 	= data_hdr->time.minute;
		mini_data_hdr.hdr.time.second 	= data_hdr->time.second; 
		mini_data_hdr.hdr.time.fracsec  = data_hdr->time.fracsec;

		if (mini_data_hdr.hdr.bofb)	
		{
			if (mini_data_hdr.hdr.bod == 0)
		/* move optional blockettes into buffer for later output */
				memcpy(mini_data_hdr.blockettes, 
				lrecord_ptr + mini_data_hdr.hdr.bofb,
					256 - mini_data_hdr.hdr.bofb);
			else
				if (mini_data_hdr.hdr.bod < 256)
		/* move optional blockettes into buffer for later output */
				memcpy(mini_data_hdr.blockettes, 
				lrecord_ptr + mini_data_hdr.hdr.bofb,
			mini_data_hdr.hdr.bod - mini_data_hdr.hdr.bofb);

		}


		data_hdr->num_mux_chan = 1;
		for (i=0;i<16;i++) data_hdr->mux_chan_name[i] = '1'+ i;


		/* find the current channel within the current station */
		if (current_station != NULL)
		{
			int n = 3;

			for (new_current_channel = current_station->type52_head;
				new_current_channel != NULL;
				new_current_channel = new_current_channel->next)
			{
				i = new_current_channel->subchannel - 1; 
				if ((i<0) || (i>9)) 
					i = 0;
				else /* multiplexed data */
					n = 2;

				if (strncmp(this_channel, new_current_channel->channel, n) == 0)
				{
					data_hdr->mux_chan_name[i] =
						determine_orient_code(new_current_channel);
					i++;
				}
			}
		} 
	}

/*                 +=======================================+                 */
/*=================|       read and decode the data        |=================*/
/*                 +=======================================+                 */

	if (output)
	{
/*                 +=======================================+                 */
/*=================| find compression/decompression format |=================*/
/*                 +=======================================+                 */

		/* find compression format */

		/* check for the presense fo a blockette 1000, 
		 * use the data type if found 
		 */
	
		strcpy(decode, "");
	
		if (data_hdr->bofb)
		{
			struct data_blk_1000 *p;

			/* since the input_data_hdr struct doesn't include the	
			 * 1st 8 bytes, whereas the bofb does, subtract out 
			 */
            		p = scan_for_blk_1000((char *)input_data_hdr + 
					      	 input_data_hdr->bofb - 8, 
							input_data_ptr - 8);


			if (p)
			{
				strcpy(decode, code_to_english(p->encoding_fmt));
			}

		}

		if (strlen(decode) == 0)
		{
			strncpy (decode, "     ", 5);
			for (this_type30 = type30_head; 
				this_type30 != NULL; this_type30 = this_type30->next)
				if (this_type30->code == current_channel->format_code) break;
			if (this_type30 != NULL)
			{
				strcpy (decode, this_type30->name);
			}
		}
		else
		{
			/* must still cycle to the proper type30 for
			 * geoscopian data
			 */
			for (this_type30 = type30_head;
                                this_type30 != NULL; this_type30 = this_type30->next)
                                	if (this_type30->code == current_channel->format_code) 
						break;
		}

		if (this_type30 == NULL)
		{
			fprintf (stderr, "WARNING (process_data):  ");
			fprintf (stderr, "data format description %d ",
			current_channel->format_code);
			fprintf (stderr, "not found in data format table.\n");
			fprintf (stderr, "\tData block Skipped.\n");
			return;
		}
	
		/* point to first sample of this data block */

		input_data_ptr = lrecord_ptr + input_data_hdr->bod;

		/* decode the data */
		strupr(decode);	    /* make all upper case comparisons */

		if (Seed_flag)
		{
			output_seed_data_file(lrecord_ptr + 8);

			return;

		}
	
		if ((this_nsamples == 0) && (mini_flag))
		{

			struct data_hdr d_hdr;

                        memset((char *)&d_hdr, 0, sizeof(d_hdr));

                        d_hdr.station = this_station;
                        d_hdr.channel = this_channel;
			d_hdr.location = this_location;


                        strncpy(d_hdr.network,
                                 network_code, strlen(network_code));

                        d_hdr.time = this_time;

			d_hdr.num_mux_chan = data_hdr->num_mux_chan;

			for (i = 0; i < sizeof(data_hdr->mux_chan_name); i++)
				d_hdr.mux_chan_name[i] = data_hdr->mux_chan_name[i];

			output_mini(&d_hdr, 0);

		}

		if (strcmp(decode, "") == 0)
			return;

		if (strstr(decode, "16-BIT"))
			decode_16bit(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "ASRO"))
			decode_asro(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "CDSN"))
			decode_cdsn(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "GRAEF"))
			decode_graef(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "DWWSS"))
			decode_dwwssn(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "STEIM2")) 
			decode_steim2(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "STEIM")) /* MUST come after STEIM2 */
			decode_steim(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "SRO G"))
			decode_sro(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "ECHY"))
			decode_echery(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "RSTN"))
			decode_rstn(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "32-BI"))
			decode_32bit(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "SUN I"))
			decode_ieeefloat(input_data_ptr, this_nsamples, next_index);
		else if (strstr(decode, "GEOSC"))
			decode_geoscope (input_data_ptr, 
					 this_nsamples, 
					 next_index, 
					 this_type30->name);
		else 
		if (strstr(decode, "USNSN"))
			decode_usnsn(input_data_ptr, input_data_hdr);
		else
			decode_ddl(input_data_ptr, this_nsamples, next_index, this_type30);

	}

/*                 +=======================================+                 */
/*=================|              clean up                 |=================*/
/*                 +=======================================+                 */

}


/*=======================================================================*/
/* SEED reader     |             chk_time              |    process_data */
/*=======================================================================*/
int chk_time(blk_start, blk_end)
struct time blk_start;
struct time blk_end;
{
	int i;

	/* if no start times, accept all blocks */
	if ((start_time_count == 0) && (end_time_count == 0)) 
		return(1);

	/* for each user entered time span */
	for (i=0;i<start_time_count;i++)
	{
		/* is start of block within a time span */
		if (timecmp (blk_start, start_time_point[i]) == 0) 
			return(1);

		if (timecmp (blk_start, start_time_point[i]) > 0)
		{
			if (i+1 > end_time_count) 
				return(1);
			else 
			if (timecmp(blk_start, end_time_point[i]) <= 0) 
				return(1);
		}
		/* is end of this block within a time span */
		else
		{
			if (timecmp (blk_end, start_time_point[i]) > 0) 
				return(1);
		}
	}

	/* no overlapping time spans found */

	return(0);
}

/* ------------------------------------------------------------------------ */
#define WITHIN(t, s, e) ((timecmp(t, s) >= 0) && (timecmp(t, e) <= 0))
	
/* ------------------------------------------------------------------------ */
struct type50 *get_station_rec(s, n, rec_time)
char *s;
char *n;
struct time *rec_time;

{
	struct type50 *p, *p2 = NULL;
	
	struct time start, end; 	/* response effective times */

	char stn[10];  			/* if we ever get more than 10 
					 * for station and network sizes,
					 * we need to upgrade 
					 */
	char net[10];

	strncpy(stn, s, 10);
	strncpy(net, n, 10);

	trim(stn);
	trim(net);

	for (p = type50_head; p != NULL;p = p->next)
	{
		if (strcmp(stn, p->station) != 0) 
			continue;

		if ((type10.version >= 2.3) && (!ignore_net_codes))  
			if (strcmp(net, p->network_code) != 0)
				continue;

		timecvt(&start, p->start);
		timecvt(&end, p->end);

		if (p->end == NULL) 
			end.year = 9999;

		/* if no time is given, just use the one with the channel recs attached */
		if (rec_time == NULL)
		{
			if (p->type52_head != NULL)
				p2 = p;
		}
		else	
			/* else, do time check */
			if (WITHIN((*rec_time), start, end))	
			{
				/* make check for channel records there */
				if (p->type52_head != NULL) 
					p2 = p;

			}
	}

	return p2;

}
/* ------------------------------------------------------------------------ */
struct type50 *find_station_rec(s, n, beg_time, end_time)
char *s, *n;
struct time *beg_time;
struct time *end_time;