关于seed文件的读取方法

本文是我捣鼓了一个多月的总结,目的是自己的备忘录,也可供系统内有需求的同仁以参考。

这里的seed文件的读取方法只限于地震测震和强震数据,其他应用领域的数据没有测试过所以也不知道是否适用。

参考资料《地震波形数据交换格式》和 jopens软件程序

程序实例语言 :golang

这里的seed文件是指以seed为后缀名的文件,该格式广泛应用与数字地震波形数据存储和数据交换。其编码结构在《地震波形数据交换格式》中有详细说明,在这里不做赘述,有兴趣者可以下载《地震波形数据交换格式》自行研究。本文只介绍在seed文件读取解码过程中的关键点和方法步骤(只少我是认为这些是比较关键的)

seed文件一般是由两部分组成,一部分是头文件内容(使用ASCII格式存储),一部分是数据内容(使用二进制格式存储);头文件信息的正确读取直接决定着后面的数据内容是否可以正常读取。

头文件内容是由多个数据块组成,数据块长度不固定,其数据格式大体是“数据块类型”-》“数据块长度”-》“数据块存储的各类信息”     如下所示:

卷台站头段索引子块 [11]
0110021001JXH  000005
解析:
011   子块类型
0021  子块长度
001   台站数
JXH   台站标识码(这里是5位有两个空格)
000005  台站头端顺序号   (S标志)      
卷时间片索引子块 [12]
012006300012023,101,16:00:00.0000~2023,102,16:00:00.0000~000012
解析:
012   子块类型
0063  子块长度
0001  卷中时间片数量
2023,101,16:00:00.0000   时间片的开始时间
~     分隔符
2023,102,16:00:00.0000   时间片的结束时间
~     分隔符
000012                   时间片头段顺序号(T标志)

数据块和数据块之间一般没有间隔,数据块的识别我是通过文件顺序读取检索关键字的方式处理的。代码如下:

// 按照关键字查找返回信息字段起始位置和信息长度
// file 文件指针;p_start 游标开始位置; keyword 关键字(关键字后三位必须是子块类型码)
func (s *HeadInfo) GetInfoSL(file *os.File, p_start int, keyword string) (int, int, error) {
	scanner := bufio.NewScanner(file)

	file.Seek(int64(p_start), 0)
	str_length := 0
	p := 0
	buf := make([]byte, 0, 64*1024)
	scanner.Buffer(buf, 1024*1024)
	for scanner.Scan() {

		str := scanner.Text()
		ide_s := strings.Index(str, keyword)
		if ide_s != -1 {

			p += ide_s
			buffer := make([]byte, 4)
			n, err := file.ReadAt(buffer, int64(p_start+p+len(keyword)))
			if err != nil {

				//panic(err)
				return 0, 0, err

			}
			str_length, err = strconv.Atoi(string(buffer[:n]))
			if err != nil {
               log.Println("keyword:",keyword)
				//panic(err)
				return 0, 0, err

			}


			break
		} else {
			p += len(str)
		}
	}

	if err := scanner.Err(); err != nil {
        return -1, -1, err
	}

	return p + p_start + len(keyword) - 3, str_length, nil
}

这里可以同过顺序读取关键字来定位数据块的起始位置和数据块长度(其中数据块长度是包含数据块名称)

seed文件是以"000001V  "字符串为文件内容开头的,这个也可以作为辨别是否为seed文件格式的标志。

直接影响后面数据内容读取的seed头关键数据块信息由:

子块010:

获取逻辑卷存储长度(这个是后面数据内容中数据块的存储大小(大块));

获取逻辑卷存储时间(这个信息对后面的数据内容提取意义不大);

子块011:

获取台站标识码;

台站头顺序号(用于台站数据块信息识别);

子块012:

获取数据记录的开始时间和结束时间;

时间片头段顺序号(用于识别时间片数据块信息);

子块050:

获取台站信息,这里包括台站的一些经纬度、高程、台站代码、通道个数等信息(这里的一部分数据属于敏感数据,提取后处理时最好脱敏存储)

子块052:

获取通道信息,这里有多个通道需要循环读取。通道信息中比较重要的信息由:通道名称、采样率、存储长度(这里的数据n是指2^n)

子块070:

识别时间片标识数据块

子块074:

获取各通道的通道标识符、第一个数据的序号、最后一个记录的序号;(这里也是多个通道信息)

子块1000:

这个子块是个用二进制存储的数据块,它在固定头端区之后,一般在固定头端区最后一个信息:字段18 (第一个子块)(偏移量,一般是在固定头段区后48字节的位置)

该子块需要获取的内容由:编码格式、字序、数据记录长度等信息(这里有些信息在上面的数据块中已获取到)

另外固定头端区里还有一个比较重要的信息可以获取后用于后面数据解码后数据个数的验证,这个信息就是字段9样本数目。(固定头端区是由ASCII(20个字节)和二进制(28个字节)组成,在后面的每一个编码数据块前都有一个固定头端区。这里的编码数据块长度就是052中的存储长度。测震数据一般为512字节,强震数据为4096字节)

到这里影响后面二进制数据内容读取的信息基本结束。

读取信息的主要代码如下:

// 获取头信息
func (s *HeadInfo) GetHeadInfo(file *os.File) (headInfo *HeadInfo, err error) {
	buffer := make([]byte, 8)
	n, err := file.ReadAt(buffer, 0)
	if err != nil {
		panic(err)
	}
	if string(buffer[:n]) != "000001V " {
		return nil,errors.New("该文件不是seed格式,不能正常读取!!!")
	}
	//获取逻辑卷存储长度
	buffer = make([]byte, 2)
	n, err = file.ReadAt(buffer, 19)
	if err != nil {
		panic(err)
	}
	s.VolumeLength, _ = s.countSaveLen(string(buffer[:n]))
	//获取逻辑卷存储时间
	buffer = make([]byte, 22)
	n, err = file.ReadAt(buffer, 67)
	if err != nil {
		panic(err)
	}
	s.SaveTime = string(buffer[:n])
	//获取010长度
	pos_start, length, err := s.GetInfoSL(file, 0, "000001V 010")
	if err != nil {
		log.Println("GetInfoSL error: ", err)
		return nil,err
	}

	if length != 0 {
		//获取数据字典
		_,err = s.GetDataFormatSubInfo(file," 030")
		if err != nil {
			panic(err)
		}
		//获取台站头信息
		//台站标识码
		buffer = make([]byte, 5)
		n, err = file.ReadAt(buffer, int64(pos_start+length+10))
		if err != nil {
			panic(err)
		}
		s.StationMark = string(buffer[:n])
		//台站头端顺序号
		buffer = make([]byte, 6)
		n, err = file.ReadAt(buffer, int64(pos_start+length+15))
		if err != nil {
			panic(err)
		}
		s.StationDataMark = string(buffer[:n]) + "S"
		//获取时间信息
		//数据记录开始时间
		buffer = make([]byte, 22)
		n, err = file.ReadAt(buffer, int64(pos_start+length+32))
		if err != nil {
			panic(err)
		}
		s.TimeStart = string(buffer[:n])
		//数据记录结束时间
		n, err = file.ReadAt(buffer, int64(pos_start+length+55))
		if err != nil {
			panic(err)
		}
		s.TimeEnd = string(buffer[:n])
		//时间片头段顺序号
		buffer = make([]byte, 6)
		n, err = file.ReadAt(buffer, int64(pos_start+length+78))
		if err != nil {
			panic(err)
		}
		s.TimeDataMark = string(buffer[:n]) + "T"
	}

	//台站信息
	if len(s.StationDataMark) != 0 {
		pos_start, err = s.GetStationInfo(file, s.StationDataMark+" 050")
		if err != nil {
			log.Println("GetStationInfo error: ", err)
			return nil,err

		}

	}
	//时间头信息
	if len(s.TimeDataMark) != 0 && len(s.StationDataMark) != 0 && len(s.SInfo.ChannelNum) != 0 {

		channel_n, err := strconv.Atoi(s.SInfo.ChannelNum)
		log.Println("channel_n: ",channel_n)
		if err != nil {
			panic(err)
		}


		_, err = s.GetTimeInfo(file, pos_start, channel_n, s.TimeDataMark+" 070")
		if err != nil {
			log.Println("GetTimeInfo error: ", err)

		}
		
	}

	return &HeadInfo{
		VolumeLength:    s.VolumeLength,
		SaveTime:        s.SaveTime,
		StationDataMark: s.StationDataMark,
		StationMark:     s.StationMark,
		TimeStart:       s.TimeStart,
		TimeEnd:         s.TimeEnd,
		TimeDataMark:    s.TimeDataMark,
		SInfo:           s.SInfo,
		TInfo:           s.TInfo,
		DFInfo:          s.DFInfo,
	}, nil
}
// 获取台站信息信息,返回当前位置
func (s *HeadInfo) GetStationInfo(file *os.File, keyword string) (int, error) {
	pos_start, length, err := s.GetInfoSL(file, 0, keyword)
	if err != nil {
		log.Fatal("GetInfoSL error: ", err)
	}

	//获取字节交换顺序
	buffer_str := make([]byte,length)
	n, err := file.ReadAt(buffer_str, int64(pos_start))
	if err != nil {
		panic(err)
	}
	str := string(buffer_str[:n])

	tmp := strings.Split(str,"~")
	if len(tmp)>1 {
		s.SInfo.ByteOrder32 = tmp[1][3:7]
		s.SInfo.ByteOrder16 = tmp[1][7:9]
	}

	buffer := make([]byte, 10)
	pos_c := 0
	if length != 0 {
		//纬度
		n, err := file.ReadAt(buffer, int64(pos_start+12))
		if err != nil {
			panic(err)
		}
		s.SInfo.Lat = string(buffer[:n])
		//经度
		buffer = make([]byte, 11)
		n, err = file.ReadAt(buffer, int64(pos_start+22))
		if err != nil {
			panic(err)
		}
		s.SInfo.Lon = string(buffer[:n])
		//高程
		buffer = make([]byte, 7)
		n, err = file.ReadAt(buffer, int64(pos_start+33))
		if err != nil {
			panic(err)
		}
		s.SInfo.Altitude = string(buffer[:n])
		//通道数
		buffer = make([]byte, 4)
		n, err = file.ReadAt(buffer, int64(pos_start+40))
		if err != nil {
			panic(err)
		}
		s.SInfo.ChannelNum = string(buffer[:n])
		//台网代码
		buffer = make([]byte, 2)
		n, err = file.ReadAt(buffer, int64(pos_start+length-2))
		if err != nil {
			panic(err)
		}
		s.SInfo.StationCode = string(buffer[:n])
		//通道信息

		channel_n, err := strconv.Atoi(s.SInfo.ChannelNum)
		//log.Println("channel_n: ",channel_n)
		if err != nil {
			panic(err)
		}



		pos_c = pos_start + length

		str, err = s.getChannelStr(file, pos_c)
		//log.Println("getChannelStr====>>>>str: ",str)
		if err != nil {
			log.Println("getChannelStr error: ", err)
			return 0, err
		}
		if len(str)==0 {
			log.Println("getChannelStr error: null ")
			return 0, errors.New("通道信息为空!!!")
		}


		for i := 0; i < channel_n; i++ {

			//获取052通道信息
			str_a, err := s.getChannelInfo(str, "052")
			if err != nil {
				log.Println("getChannelInfo error: ", err)
			}
            str = str_a


		}

	}

	return pos_c, nil
}
// 获取时间头信息,返回结束位置
func (s *HeadInfo) GetTimeInfo(file *os.File, p_start int, channel_num int, keyword string) (int, error) {

	pos_start, length, err := s.GetInfoSL(file, p_start, keyword)
	if err != nil {
		log.Fatal("GetInfoSL error: ", err)
	}

	pos_c := pos_start + length

	if length != 0 {

		for i := 0; i < channel_num; i++ {
			pos_start, length, err = s.GetInfoSL(file, pos_c, "074")
			if err != nil {
				log.Fatal("GetInfoSL error: ", err)
			}


			//通道标识符
			buffer := make([]byte, 3)
			n, err := file.ReadAt(buffer, int64(pos_start+14))
			if err != nil {
				panic(err)
			}
			cm := string(buffer[:n])

			//第一个数据的序号
			buffer = make([]byte, 6)
			n, err = file.ReadAt(buffer, int64(pos_start+40))
			if err != nil {
				panic(err)
			}
			sm, err := strconv.Atoi(string(buffer[:n]))
			if err != nil {
				panic(err)
			}

			//最后一个记录的序号
			n, err = file.ReadAt(buffer, int64(pos_start+71))
			if err != nil {
				panic(err)
			}
			em, err := strconv.Atoi(string(buffer[:n]))
			if err != nil {
				panic(err)
			}

			s.TInfo = append(s.TInfo, TimeHeadInfo{
				ChannelMark: cm,
				StartMark:   sm,
				EndMark:     em,
			})
			pos_c = pos_start + length
		}

	}

	return pos_c, nil
}

因内容太长,本篇只介绍seed文件头信息读取,seed文件内容读取及steim格式解码,下篇后续。

  • 1
    点赞
  • 3
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
% Known encoding formats are the following FDSN codes: % 0: ASCII % 1: 16-bit integer % 2: 24-bit integer (untested) % 3: 32-bit integer % 4: IEEE float32 % 5: IEEE float64 % 10: Steim-1 % 11: Steim-2 % 12: GEOSCOPE 24-bit (untested) % 13: GEOSCOPE 16/3-bit gain ranged % 14: GEOSCOPE 16/4-bit gain ranged (untested) % 19: Steim-3 (alpha and untested) % % See also MKMSEED to export data in miniSEED format. % % % Author: Franois Beauducel % Institut de Physique du Globe de Paris % Created: 2010-09-17 % Updated: 2012-04-21 % % Acknowledgments: % Ljupco Jordanovski, Jean-Marie Saurel, Mohamed Boubacar, Jonathan Berger, % Shahid Ullah. % % References: % IRIS (2010), SEED Reference Manual: SEED Format Version 2.4, May 2010, % IFDSN/IRIS/USGS, http://www.iris.edu % Trabant C. (2010), libmseed: the Mini-SEED library, IRIS DMC. % Steim J.M. (1994), 'Steim' Compression, Quanterra Inc. % History: % [2012-04-21] % - Correct bug with Steim + little-endian coding % (thanks to Shahid Ullah) % [2012-03-21] % - Adds IDs for warning messages % [2011-11-10] % - Correct bug with multiple channel name length (thanks to % Jonathan Berger) % [2011-10-27] % - Add LocationIdentifier to X.ChannelFullName % [2011-10-24] % - Validation of IEEE double encoding (with PQL) % - Import/plot data even with file integrity problem (like PQL) % [2011-07-21] % - Validation of ASCII encoding format (logs) % - Blockettes are now stored in substructures below a single % field X.BLOCKETTES % - Add import of blockettes 500 and 2000 % - Accept multi-channel files with various data coding % [2010-10-16] % - Alpha-version of Steim-3 decoding... % - Extend output parameters with channel detection % - Add gaps and overlaps on plots % - Add possibility to force the plot % [2010-10-02] % - Add the input formats for GEOSCOPE multiplexed old data files % - Additional output argument with gap and overlap analysis % - Create a plot when no output argument are specified % - Optimize script coding (30 times faster STEIM decoding!) % % [2010-09-28] % - Correction of a problem with STEIM-1 nibble 3 decoding (one % 32-bit difference) % - Add reading of files without blockette 1000 with additional % input arguments (like Seismic Handler output files). % - Uses warning() function instead of fprintf().

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值