您的位置:首页 > 其它

vtk学习笔记——利用vtk显示地震数据

2012-03-17 14:27 731 查看
多年来一直从事技术开发,越来越觉得自己知识太少,基本属于无知的行列,所以不敢写blog,心虚。。。

最近在研究vtk显示地震数据的问题,通过阅读之前其他网友所提供的vtk学习笔记,收获不少,特此表示感谢。

本着抛砖引玉的原则,把自己的一些心得写出来,期待大家拍砖。谢谢!

——————————————————————————————————————————————

地震的常见数据格式之一,segy,目前存在两个标准,1975年和2002年,后一版本尽量做到了向前兼容。

通过网上下载到了英文版的数据格式说明和中文版,感谢中文版作者的分享。

本人通过阅读数据格式规范,结合java,进行了数据的vtk显示。

由于不知道怎么发附件,只能把源码贴在下面,见谅各位。

核心类为SegyReader,它包括了卷头、道头、道数据的组织,道头和道数据都封装在了SegyTrace类里,这个类里实现了从IBM浮点数到java的转换(ps:这个费了偶不少劲,遍阅很多资料,反复若干次才弄对,所谓对是通过发布文章上的16进制表示的ibm浮点数与十进制的对应关系验证而得)。其中SimpleVTK2这个类里,封装了显示部分。很简单,呵呵,没做特殊处理,只是把一道数据做成一条线,抛砖引玉,请各位评判。

闲话不说了,上代码:

____________________________ SegyReader类 _________________________________________________________________________________

package org.sinoprobe.segy;

import java.awt.BorderLayout;

import java.io.DataInputStream;

import java.io.FileInputStream;

import java.io.IOException;

import java.util.ArrayList;

import javax.swing.JFrame;

import javax.swing.SwingUtilities;

import org.sinoprobe.util.Util4JavaBean;

/**

* segy数据文件读取类

* 它包括了对segy文件所包含全部数据对象的封装

*

* @author dingyi

*

*/

public class SegyReader {

private String filePath;

private SegyTextFileHeader txtHdr;

private SegyReelHdr reelHdr;

private ArrayList<SegyTrace> traceList=new ArrayList<SegyTrace>();

/**

* 构造Segy读取对象

* @param filePath segy文件名,包括全路径

*/

public SegyReader(String filePath) {

super();

this.filePath = filePath;

}

public void buildFromFile() throws IOException, SegyFileFormatException{

FileInputStream fis = new FileInputStream(filePath);

DataInputStream dis = new DataInputStream(fis);

// 首先读取前3200字节的文本头

byte[] buf = new byte[3200];

dis.readFully(buf);

byte[] b = CodeUtil.EBCDICToASCII(buf);

txtHdr = new SegyTextFileHeader(b);

//判断是否合法,即是否每一行均以C开头,如果不合法,那么,可能会是存在带头标签

if(!txtHdr.isValid()){

throw new SegyFileFormatException("卷头数据不是以C开头,可能存在标签对象,请尝试用另外方法读取!");

}

//继续读取卷头数据

buf = new byte[400];

dis.readFully(buf);

reelHdr = new SegyReelHdr(buf);

//判断是否有扩展卷头

if(reelHdr.has_extend_text_hdr_flag!=0){

//此时有3200字节的扩展原文文件头

//需要执行读取动作,此处暂空,以抛出异常为处理手段

throw new SegyFileFormatException("存在扩展原文文件头,需要读取!");

}

//继续读取道数据

//强制读取卷头所标记的

//for (int i = 0; i <reelHdr.traces_per_record; i++)

while(dis.available()>0) {

// 首先读取道头

buf = new byte[240];

dis.readFully(buf);

SegyTrace oneTrace = new SegyTrace(buf);

// 根据卷头和道头得出要读取的数据块的大小

int dataCount = oneTrace.samples_in_this_trace

* reelHdr.getMultiple();

//读取道数据

buf = new byte[dataCount];

dis.readFully(buf);

oneTrace.readDataInIBMFormat(buf);

//添加到道列表中

traceList.add(oneTrace);

}

System.out.println("读取道数据结束");

dis.close();

fis.close();

}

/**

* 得到原文文件头,文本格式

* @return

*/

public SegyTextFileHeader getTxtHdr() {

return txtHdr;

}

/**

* 得到卷(道集)头数据

* @return

*/

public SegyReelHdr getReelHdr() {

return reelHdr;

}

/**

* 得到道数据列表

* @return

*/

public ArrayList<SegyTrace> getTraceList() {

return traceList;

}

/**

* 得到测量到的最大最小值

* @return 顺序为min,Max

*/

public double[] getMinMaxValue(){

double[] minMax=new double[2];

minMax[0]=0;

minMax[1]=0;

for(SegyTrace t:traceList){

for(double d: t.data){

if(d>minMax[1])

minMax[1]=d;

if(d<minMax[0])

minMax[0]=d;

}

}

return minMax;

}

public static void main(String[] args) {

SegyReader reader=new SegyReader("f:\\3002.sgy");

try {

reader.buildFromFile();

System.out.println("数据读取完成!开始构造vtk显示");

//构造vtk界面

final SimpleVTK2 vtkPanel=new SimpleVTK2(reader);

SwingUtilities.invokeLater(new Runnable() {

@Override

public void run() {

JFrame frame = new JFrame("SimpleVTK");

frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);

frame.getContentPane().setLayout(new BorderLayout());

frame.getContentPane().add(vtkPanel);

frame.setSize(640, 480);

frame.setLocationRelativeTo(null);

frame.setVisible(true);

}

});

} catch (IOException e) {

// TODO Auto-generated catch block

e.printStackTrace();

} catch (SegyFileFormatException e) {

// TODO Auto-generated catch block

e.printStackTrace();

}

}

}

———————————————— SegyReelHdr 类————————————————————————————————————

package org.sinoprobe.segy;

/**

* 表2. 二进制文件头

400字节二进制文件头

字节 描述

3201-3204 作业标识号

3205-3208 测线号。对3-D叠后数据而言,它将典型地包含纵向测线(In-line)号

3209-3212 卷号

3213-32145 每个道集的数据道数。叠前数据强制要求

3215-32165 每个道集的辅助道数。叠前数据强制要求

3217-32186 微秒(us)形式的采样间隔。叠前数据强制要求

3219-3220 微秒(us)形式的原始野外记录采样间隔

3221-32226 数据道采样点数。叠前数据强制要求

注释:二进制文件头中的采样间隔和采样点数应当是文件中地震数据的首要一组参数

3223-3224 原始野外记录每道采样点数

3225-32266 数据采样格式编码。叠前数据强制要求

1=4字节IBM浮点数

2=4字节,两互补整数

3=2字节,两互补整数

4=4字节带增益定点数(过时,不再使用)

5=4字节IEEE浮点数

6=现在没有使用

7=现在没有使用

8=1字节,两互补整数

3227-32287 道集覆盖次数――每个数据集的期望数据道数(例如CMP覆盖次数)。强烈推荐所有类型的数据使用

3229-32307 道分选码(即集合类型):

-1=其他(应在用户扩展文件头文本段中解释)

0=未知

1=同记录(未分选)

2=CDP道集

3=单次覆盖连续剖面

4=水平叠加

5=共炮点

6=共接收点

7=共偏移距

8=共中心点

9=共转换点

强烈推荐所有类型的数据使用

3231-3232 垂直求和码:

1=不求和

2=两次求和



M=M-1求和(M=2到32767)

3233-3234 起始扫描频率(Hz)

3235-3236 终止扫描频率(Hz)

3237-3238 扫描长度(ms)

3239-3240 扫描类型码:

1=线性

2=抛物线

3=指数

4=其他

3241-3242 扫描信道的道数

3243-3244 有斜坡时,以毫秒表示的扫描道起始斜坡长度(斜坡从零时刻开始,对这个长度有效)

3245-3246 以毫秒表示的扫描道终止斜坡长度(斜坡终止始于扫描长度减去斜坡结尾处的长度)

3247-3248 斜坡类型:

1=线性

2=cos2

3=其他

3249-3250 相关数据道:

1=无相关

2=相关

3251-3252 二进制增益恢复:

1=恢复

2=未恢复

3253-3254 振幅恢复方法:

1=无

2=球面扩散

3=自动增益控制

4=其他

3255-32567 测量系统:强烈推荐所有类型的数据使用。如文件中包含位置数据文本段,这条必须与位置数据文本段一致。如不同,最后位置数据文本段有控制权。

1=米

2=英尺

3257-3258 脉冲极化码:

1=压力增大或检波器向上运动在磁带上记作负数

2=压力减小或检波器向下运动在磁带上记作正数

3259-3260 可控源极化码:

地震信号滞后引导信号:

1=337.5°-22.5°

2=22.5°-67.5°

3=67.5°-112.5°

4=112.5°-157.5°

5=157.5°-202.5°

6=202.5°-247.5°

7=247.5°-292.5°

8=292.5°-337.5°

3261-3500 未赋值

3501-3502 6 SEG Y格式修订版号。这是一个16比特无符号数值,在第一和第二字节间有Q点。例如SEG Y修订版1.0,如文档定义,它将记录为0100 。此字段对所有SEG Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

3503-3504 6 固定长度道标志。1表示SEG Y文件中所有道确保具有相同的采样间隔和采样点数,即在原文文件头中3217-3218和3221-3222字节。0表示文件中的道长可能变化,此时道头中115-116字节的采样点数必须用来确认各道的实际长度。此字段对所有SEG Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

3505-3506 6 3200字节扩展原文文件头记录在二进制头后。0表示没有扩展原文文件头记录(即此文件无扩展原文文件头)。-1表示扩展原文文件头记录数可变,并且扩展原文文件头结尾用最终记录的一个文本段((SEG: ENDText))表示。正值表示有很多原文文件头。注意虽然具体的扩展原文文件头数目是个有用的信息,但是在写二进制头时它并不是总知道也不是强制要求在此记录正值。此字段对所有SEG Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

3507-3600 未赋值

5 此信息对叠前数据强制要求。

6 此信息对所有数据类型强制要求。

7 强烈建议此信息一直被记录。

* @author dingyi

*

*/

public class SegyReelHdr {

// 3201-3204 作业标识号

int job_id_number;

// 3205-3208 测线号。对3-D叠后数据而言,它将典型地包含纵向测线(In-line)号

int line_number;

// 3209-3212 卷号

int reel_number;

// 3213-3214 每个道集的数据道数。叠前数据强制要求

short traces_per_record;

// 3215-3216 每个道集的辅助道数。叠前数据强制要求

short aux_traces_per_record;

// 3217-3218 微秒(us)形式的采样间隔。叠前数据强制要求

short sample_data_interval_ms;

// 3219-3220 微秒(us)形式的原始野外记录采样间隔

short original_data_interval_ms;

// 3221-3222 数据道采样点数。叠前数据强制要求

// 注释:二进制文件头中的采样间隔和采样点数应当是文件中地震数据的首要一组参数

short samples_per_trace;

// 3223-3224 原始野外记录每道采样点数

short original_samples_per_trace;

// 3225-3226 数据采样格式编码。叠前数据强制要求

/**

* 1=4字节IBM浮点数 2=4字节,两互补整数 3=2字节,两互补整数 4=4字节带增益定点数(过时,不再使用) 5=4字节IEEE浮点数

* 6=现在没有使用 7=现在没有使用 8=1字节,两互补整数

*/

short data_sample_format_code;// 数据采样格式编码

// 3227-32287 道集覆盖次数――每个数据集的期望数据道数(例如CMP覆盖次数)。强烈推荐所有类型的数据使用

short CDP_fold;

// 3229-32307 道分选码(即集合类型):

/**

* -1=其他(应在用户扩展文件头文本段中解释) 0=未知 1=同记录(未分选) 2=CDP道集 3=单次覆盖连续剖面 4=水平叠加 5=共炮点

* 6=共接收点 7=共偏移距 8=共中心点 9=共转换点 强烈推荐所有类型的数据使用

*/

short trace_sorting_code;// 3229-32307 道分选码(即集合类型):

// 3231-3232 垂直求和码:

/*

* 1=不求和 2=两次求和 … M=M-1求和(M=2到32767)

*/

short vertical_sum_code;

// 3233-3234 起始扫描频率(Hz)

short sweep_frequency_start_hz;

// 3235-3236 终止扫描频率(Hz)

short sweep_frequency_end_hz;

// 3237-3238 扫描长度(ms)

short sweep_length_ms;

// 3239-3240 扫描类型码:

/**

* 1=线性 2=抛物线 3=指数 4=其他

*/

short sweep_type_code;

// 3241-3242 扫描信道的道数

short trace_number_of_sweep_channel;

// 3243-3244 有斜坡时,以毫秒表示的扫描道起始斜坡长度(斜坡从零时刻开始,对这个长度有效)

short sweep_trace_taper_length_start_ms;

// 3245-3246 以毫秒表示的扫描道终止斜坡长度(斜坡终止始于扫描长度减去斜坡结尾处的长度)

short sweep_trace_taper_length_end_ms;

/**

* 3247-3248 斜坡类型: 1=线性 2=cos2 3=其他

*/

short taper_type_code;

// 3249-3250 相关数据道:

// 1=无相关

// 2=相关

short correlated_data_traces_flag;

/**

* 3251-3252 二进制增益恢复: 1=恢复 2=未恢复

**/

short binary_gain_recovered_flag;

/**

* 3253-3254 振幅恢复方法: 1=无 2=球面扩散 3=自动增益控制 4=其他*

*/

short amplitude_recovery_method_code;

/**

* 3255-32567

* 测量系统:强烈推荐所有类型的数据使用。如文件中包含位置数据文本段,这条必须与位置数据文本段一致。如不同,最后位置数据文本段有控制权。 1=米

* 2=英尺

**/

short measurement_system;

/**

* 3257-3258 脉冲极化码: 1=压力增大或检波器向上运动在磁带上记作负数 2=压力减小或检波器向下运动在磁带上记作正数

**/

short impulse_signal_polarity;

/**

* 3259-3260 可控源极化码: 地震信号滞后引导信号: 1=337.5°-22.5° 2=22.5°-67.5° 3=67.5°-112.5°

* 4=112.5°-157.5° 5=157.5°-202.5° 6=202.5°-247.5° 7=247.5°-292.5°

* 8=292.5°-337.5°

**/

short vibratory_polarity_code;

// 3261-3500 未赋值

// 3501-3502 6 SEG Y格式修订版号。这是一个16比特无符号数值,在第一和第二字节间有Q点。例如SEG

// Y修订版1.0,如文档定义,它将记录为0100 。此字段对所有SEG Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

short segy_version;

// 3503-3504 6 固定长度道标志。1表示SEG

// Y文件中所有道确保具有相同的采样间隔和采样点数,即在原文文件头中3217-3218和3221-3222字节。0表示文件中的道长可能变化,此时道头中115-116字节的采样点数必须用来确认各道的实际长度。此字段对所有SEG

// Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

short fix_length_trace_flag;

// 3505-3506 6

// 3200字节扩展原文文件头记录在二进制头后。0表示没有扩展原文文件头记录(即此文件无扩展原文文件头)。-1表示扩展原文文件头记录数可变,并且扩展原文文件头结尾用最终记录的一个文本段((SEG:

// ENDText))表示。正值表示有很多原文文件头。注意虽然具体的扩展原文文件头数目是个有用的信息,但是在写二进制头时它并不是总知道也不是强制要求在此记录正值。此字段对所有SEG

// Y版本强制要求,尽管零值表示遵从1975年标准的“传统”SEG Y

short has_extend_text_hdr_flag;// 如果不为0,则有扩展的3200个标头

// 3507-3600 未赋值

public SegyReelHdr(byte[] b) {

int i = 0;

job_id_number = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 4;

line_number = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 8;

reel_number = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 12;

traces_per_record = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 14;

aux_traces_per_record = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 16;

sample_data_interval_ms = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 18;

original_data_interval_ms = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 20;

samples_per_trace = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 22;

original_samples_per_trace = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 24;

data_sample_format_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 26;

CDP_fold = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 28;

trace_sorting_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 30;

vertical_sum_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 32;

sweep_frequency_start_hz = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 34;

sweep_frequency_end_hz = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 36;

sweep_length_ms = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 38;

sweep_type_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 40;

trace_number_of_sweep_channel = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 42;

sweep_trace_taper_length_start_ms = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 44;

sweep_trace_taper_length_end_ms = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 46;

taper_type_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 48;

correlated_data_traces_flag = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 50;

binary_gain_recovered_flag = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 52;

amplitude_recovery_method_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 54;

measurement_system = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 56;

impulse_signal_polarity = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 58;

vibratory_polarity_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

// 3261-3500 未赋值

i = 300;

segy_version = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 302;

fix_length_trace_flag = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 304;

has_extend_text_hdr_flag = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

// 3507-3600 未赋值

}

/**

* 得到编码格式对应的倍数

* 依据2002年segy数据规范

*

* 1=4字节IBM浮点数 2=4字节,两互补整数 3=2字节,两互补整数 4=4字节带增益定点数(过时,不再使用) 5=4字节IEEE浮点数

* 6=现在没有使用 7=现在没有使用 8=1字节,两互补整数

*/

public int getMultiple() throws SegyFileFormatException {

switch (this.data_sample_format_code) {

case 1:

case 2:

case 4:

case 5:

return 4;

case 3:

return 2;

case 8:

return 1;

}

throw new SegyFileFormatException("编码类型参数错误");

}

}

___________________________ SegyTrace类 __________________________________________________________

package org.sinoprobe.segy;

/**

* Segy道

* 包括道头和道数据

* 参照Segy 2002年的标准编写

* 下面是2002年标准中道头数据的说明

*

*

* 表3 道头

240字节道头

字节 描述

1-4 测线中道顺序号――若一条测线有若干SEG Y文件号数递增。强烈推荐所有类型的数据使用

5-8 SEG Y文件中道顺序号――每个文件以道顺序1开始。

9-12 野外原始记录号。强烈推荐所有类型的数据使用

13-16 野外原始记录的道号。强烈推荐所有类型的数据使用

17-20 震源点号――当在相同有效地表位置多于一个记录时使用。建议在道头197-202字节定义新的条目用于炮点号。

21-24 道集号(即CDP,CMP,CRP等)

25-28 道集的道数――每个道集从道号1开始

29-30 道识别码:

-1=其他

0=未知

1=地震数据

2=死道

3=哑道

4=时断

5=井口

6=扫描

7=定时

8=水断

9=近场枪信号

10=远场枪信号

11=地震压力传感器

12=多分量地震传感器――垂直分量

13=多分量地震传感器――横向分量

14=多分量地震传感器――纵向分量

15=旋转多分量地震传感器――垂直分量

16=旋转多分量地震传感器――切向分量

17=旋转多分量地震传感器――径向分量

18=可控源反应质量

19=可控源底盘

20=可控源估计地面力

21=可控源参考

22=时间速度对

23…N=选用,(最大N=32767)

强烈推荐所有类型的数据使用

31-32 产生该道的垂直叠加道数。(1是一道,2是两道求和,…)

33-34 产生该道的水平叠加道数。(1是一道,2是两道求和,…)

35-36 数据用途:

1=生产

2=试验

37-40 从震源中心点到检波器组中心的距离(若与炮激发线方向相反取负)

41-44 检波器组高程(所有基准以上高程为正,以下为负) 道头中69-70字节的因子应用于这些数值。在二进制文件头3255-3256字节指定单位为英尺或米。垂直基准应通过位置数据文本段定义(见D-1节)。

45-48 震源地表高程

49-52 震源距地表深度(正数)

53-56 检波器组基准高程

57-60 震源基准高程

61-64 震源水深

65-68 检波器组水深

69-70 应用于所有在道头41-68字节给定的真实高程和深度的因子。因子=1,+10,+100,或+1000。若为正,因子为乘数;若为负,因子为除数。

71-72 应用于所有在道头73-88字节和181-188字节给定的真实坐标值的因子。因子=1,+10,+100,或+1000。若为正,因子为乘数;若为负,因子为除数。

73-76 震源坐标――X 参考坐标系应通过扩展头位置数据文本段识别(见D-1节)。

若坐标单位为弧度秒、小数度或度/分/秒(DMS),X值代表经度,Y值代表纬度。正值代表格林威治子午线以东或赤道以北,负值代表南或西。

77-80 震源坐标――Y

81-84 检波器组坐标――X

85-88 检波器组坐标――Y

89-90 坐标单位:

1=长度(米或英尺)

2=弧度秒

3=小数度

4=度,分,秒(DMS)

注意:为编码±DDDMMSS 89-90字节等于±DDD*104+MM*102+SS,71-72字节设置为1;为编码±DDDMMSS 89-90字节等于±DDD*106+MM*104+SS*102,71-72字节设置为-100。

91-92 风化层速度(如二进制文件头3255-3256字节指明的ft/s或m/s)

93-94 风化层下速度(如二进制文件头3255-3256字节指明的ft/s或m/s)

接表3 道头

95-96 震源处井口时间(毫秒) 毫秒表示的时间需应用道头215-216字节指定的因子。

97-98 检波器组处井口时间(毫秒)

99-100 震源的静校正量(毫秒)

101-102 检波器组的校正量(毫秒)

103-104 应用的总静校正量(毫秒)(如没有应用静校正量为零)

105-106 延迟时间A――以毫秒表示的240字节道识别头的结束和时间断点之间的时间。当时间断点出现在头之后,该值为正;当时间断点出现在头之前,该值为负。时间断点是最初脉冲,它由辅助道记录或由其他记录系统指定。

107-108 延迟时间B――以毫秒表示的时间断点到能量源起爆时间之间的时间。可正可负。

109-110 记录延迟时间――以毫秒表示的能量源起爆时间到数据采样开始记录之间的时间。在SEG Y修订版0中本条用来表示深水作业,如果数据记录不从0时间开始。该条可为负值以适应负的起始时间(即数据记录在零时间之前,假设静校正量应用于数据道的结果)。若某非零值(正或负)记录在该条,它造成的影响的注释应出现在原文文件头。

111-112 起始切除时间(毫秒)

113-114 终止切除时间(毫秒)

115-1167 该道采样点数。强烈推荐所有类型的数据使用

117-1187 该道采样间隔(微秒)。

一道记录的字节数必须和写在道头中的采样点数一致。对所有记录介质都重要,但对正确处理磁盘文件中的SEG Y数据尤为关键(见附录C)。

若二进制文件头中的3503-3504字节设置了固定长度道标志,SEG Y文件每道的采样间隔和采样点数必须与二进制文件头所记录的值一致。若没有设置固定长道标志,采样间隔和采样点数可能每道变化。

强烈推荐所有类型的数据使用

119-120 野外仪器增益类型:

1=固定

2=二进制

3=浮点

4…N=选用

121-122 仪器增益常数(分贝)

123-124 仪器初始增益(分贝)

125-126 相关:

1=无

2=有

127-128 起始扫描频率(赫兹)

129-130 终止扫描频率(赫兹)

131-132 扫描长度(毫秒)

133-134 扫描类型:

1=线性

2=抛物线

3=指数

4=其他

135-136 扫描道斜坡起始长度(毫秒)

137-138 扫描道斜坡终止长度(毫秒)

139-140 斜坡类型:

1=线性

2=cos2

3=其他

141-142 假频滤波频率(赫兹),若使用

143-144 假频滤波坡度(分贝/倍频程)

145-146 陷波频率(赫兹),若使用

147-148 陷波坡度(分贝/倍频程)

149-150 低截频率(赫兹),若使用

151-152 高截频率(赫兹),若使用

153-154 低截坡度(分贝/倍频程)

155-156 高截坡度(分贝/倍频程)

157-158 数据记录的年――1975年标准没有说清应该记成2位还是4位还是两者都用。除SGE Y修订版0外,年份需记录成完整的4位罗马历法年(即2001年应记录为200110(7D116))。

159-160 日(以格林尼治标准时间和通用协调时间为基准的公元日)

161-162 时(24小时制)

163-164 分

165-166 秒

167-168 时间基准码:

1=当地

2=格林尼治标准时间

3=其他,应在扩展原文文件头的用户定义文本段解释

4=通用协调时间

169-170 道加权因子――最小有效位数定义为2-N伏。(N=0,1,…,32767)

171-172 滚动开关位置1的检波器组号

173-174 野外原始记录中道号1的检波器组号

175-176 野外原始记录中最后一道的检波器组号

177-178 间隔大小(滚动时甩掉的总检波器组数)

179-180 相对测线斜坡起始或终止点的移动

1=下(或后)

2=上(或前)

181-184 该道的道集(CDP)位置X坐标(应用道头71-72字节的因子)。参考坐标系应通过扩展头位置数据文本段识别(见D-1节)

185-188 该道的道集(CDP)位置Y坐标(应用道头71-72字节的因子)。参考坐标系应通过扩展头位置数据文本段识别(见D-1节)

189-192 对于3-D叠后数据,本字段用来填纵向线号(In-line)。若每个SEG Y文件记录一条纵向线,文件中所有道的该值应相同,并且同样的值将记录在二进制文件头的3205-3206字节中。

193-196 对于3-D叠后数据,本字段用来填横向线号(Cross-line)。它应与道头21-24字节中的道集(CDP)号的值一致,但这并不是实例。

197-200 炮点号――这可能只应用于2-D叠后数据。

注意在此假设炮点号相对于特定道最靠近叠加(CDP)位置震源位置。如果不是这种情况,应在原文文件头中有注释解释炮点实际参考点。

201-202 应用于道头中197-200字节中炮点号的因子,以得到实际数值。若为正,因子用作乘数;若为负,因子用作除数;若为零;炮点号不用于因子作用(即它是一个整数。典型的值是-10,允许炮点号小数点后有一位小数)。

203-204 道值测量单位:

-1=其他(应在数据采样测量单位文本段描述)

0=未知

1=帕斯卡(Pa)

2=伏特(V)

3=毫伏(mV)

4=安培(A)

5=米(m)

6=米每秒(m/s)

7=米每秒二次方(m/s2)

8=牛顿(N)

9=瓦特(W)

205-210 转换常数――该倍数用于将数据道采样转换成转换单位(道头211-212字节指定)。本常数以4字节编码,尾数是两互补整数(205-208字节)和2字节,十的指数幂是两互补整数(209-210字节)(即(205-208字节)×10**(209-210字节))。

211-212 转换单位――经乘以道头205-210字节中的转换常数后的数据道采样测量单位。

-1=其他(应在数据采样测量单位文本段36页描述)

0=未知

1=帕斯卡(Pa)

2=伏特(V)

3=毫伏(mV)

4=安培(A)

5=米(m)

6=米每秒(m/s)

7=米每秒二次方(m/s2)

8=牛顿(N)

9=瓦特(W)

213-214 设备/道标识――与数据道关联的单位号或设备号(即4368对应可控源号4368或20316对应2船3线16枪)。本字段允许道关联横跨独立于道号的道集(道头25-28字节)

215-216 在道头95-114字节给出的作用于时间的因子,以得到真实的毫秒表示的时间值。因子=1,+10,+100,+1000或+10000。若为正,因子用作乘数;若为负,因子用作除数。为零设定因子为一。

217-218 震源类型/方位――定义类型或能量源的方位。垂直项、横向项、纵向项作为正交坐标系的三个轴。坐标系轴的绝对角度方位可在面元网格定义文本段中定义(27页)

-1到-n=其他(应在震源类型/方位文本段38页描述)

0=未知

1=可控震源――垂直方位

2=可控震源――横向方位

3=可控震源――纵向方位

4=冲击源――垂直方位

5=冲击源――横向方位

6=冲击源――纵向方位

7=分布式冲击源――垂直方位

8=分布式冲击源――横向方位

9=分布式冲击源――纵向方位

219-224 相对震源方位的震源能量方向――正方位方向在道头217-218字节定义。能量方向以度数长度编码(即347.8°编码成3478)

225-230 震源测量――描述产生道的震源效应。测量可以简单,定量的测量如使用炸药总重量或气枪压力峰值或可控源振动次数和扫描周期时间。尽管这些简单的测量可接受,但最好使用真实的能量或工作测量单位。

本常数编码成4字节,尾数为两互补整数(225-228字节)和2字节,十的指数幂是两互补整数(209-230字节)(即(225-228字节)×10**(229-230字节))。

231-232 震源测量单位――用于震源测量、道头225-230字节的单位。

-1=其他(应在震源测量单位文本段39页描述)

0=未知

1=焦耳(J)

2=千瓦(kW)

3=帕斯卡(Pa)

4=巴(Bar)

4=巴-米(Bar-m)

5=牛顿(N)

6=千克(kg)

233-240 未赋值――为任选信息预留

* @author dingyi

*

*/

public class SegyTrace {

// 1-4 测线中道顺序号――若一条测线有若干SEG Y文件号数递增。强烈推荐所有类型的数据使用

int trace_sequence_number_within_line;

// 5-8 SEG Y文件中道顺序号――每个文件以道顺序1开始。

int trace_sequence_number_within_reel;

// 9-12 野外原始记录号。强烈推荐所有类型的数据使用

int original_field_record_number;

// 13-16 野外原始记录的道号。强烈推荐所有类型的数据使用

int trace_sequence_number_within_original_field_record;

// 17-20 震源点号――当在相同有效地表位置多于一个记录时使用。

// 建议在道头197-202字节定义新的条目用于炮点号。

int energy_source_point_number;

// 21-24 道集号(即CDP,CMP,CRP等)

int cdp_ensemble_number;

// 25-28 道集的道数――每个道集从道号1开始

int trace_sequence_number_within_cdp_ensemble;

/**

* 29-30 道识别码: -1=其他 0=未知 1=地震数据 2=死道 3=哑道 4=时断 5=井口 6=扫描 7=定时 8=水断 9=近场枪信号

* 10=远场枪信号 11=地震压力传感器 12=多分量地震传感器――垂直分量 13=多分量地震传感器――横向分量 14=多分量地震传感器――纵向分量

* 15=旋转多分量地震传感器――垂直分量 16=旋转多分量地震传感器――切向分量 17=旋转多分量地震传感器――径向分量 18=可控源反应质量

* 19=可控源底盘 20=可控源估计地面力 21=可控源参考 22=时间速度对 23…N=选用,(最大N=32767) 强烈推荐所有类型的数据使用

* *

*/

short trace_identification_code;// 道识别码

// 31-32 产生该道的垂直叠加道数。(1是一道,2是两道求和,…)

short number_of_vertically_summed_traces_yielding_this_trace;

// 33-34 产生该道的水平叠加道数。(1是一道,2是两道求和,…)

short number_of_horizontally_stacked_traced_yielding_this_trace;

/**

* 35-36 数据用途: 1=生产 2=试验*

*/

short data_use;

// 37-40 从震源中心点到检波器组中心的距离(若与炮激发线方向相反取负)

int distance_from_source_point_to_receiver_group;

// 41-44 检波器组高程(所有基准以上高程为正,以下为负)

// 道头中69-70字节的因子应用于这些数值。

//在二进制文件头3255-3256字节指定单位为英尺或米。

//垂直基准应通过位置数据文本段定义(见D-1节)。

int receiver_group_elevation;

// 45-48 震源地表高程

int surface_elevation_at_source;

// 49-52 震源距地表深度(正数)

int source_depth_below_surface;

// 53-56 检波器组基准高程

int datum_elevation_at_receiver_group;

// 57-60 震源基准高程

int datum_elevation_at_source;

// 61-64 震源水深

int water_depth_at_source;

// 65-68 检波器组水深

int water_depth_at_receiver_group;

// 69-70

// 应用于所有在道头41-68字节给定的真实高程和深度的因子。因子=1,+10,+100,或+1000。若为正,因子为乘数;若为负,因子为除数。

short scalar_for_elevations_and_depths;

// 71-72

// 应用于所有在道头73-88字节和181-188字节给定的真实坐标值的因子。因子=1,+10,+100,或+1000。若为正,因子为乘数;若为负,因子为除数。

short scalar_for_coordinates;

/**

* 73-76 震源坐标――X 参考坐标系应通过扩展头位置数据文本段识别(见D-1节)。

* 若坐标单位为弧度秒、小数度或度/分/秒(DMS),X值代表经度,

* Y值代表纬度。

* 正值代表格林威治子午线以东或赤道以北,

* 负值代表南或西。

*/

int x_source_coordinate;

// 77-80 震源坐标――Y

int y_source_coordinate;

// 81-84 检波器组坐标――X

int x_receiver_group_coordinate;

// 85-88 检波器组坐标――Y

int y_receiver_group_coordinate;

/**

* 89-90 坐标单位: 1=长度(米或英尺) 2=弧度秒 3=小数度 4=度,分,秒(DMS) 注意:为编码±DDDMMSS

* 89-90字节等于±DDD*104+MM*102+SS,71-72字节设置为1;为编码±DDDMMSS

* 89-90字节等于±DDD*106+MM*104+SS*102,71-72字节设置为-100。

*/

short coordinate_units;

// 91-92 风化层速度(如二进制文件头3255-3256字节指明的ft/s或m/s)

short weathering_velocity;

// 93-94 风化层下速度(如二进制文件头3255-3256字节指明的ft/s或m/s)

short subweathering_velocity;

// 95-96 震源处井口时间(毫秒) 毫秒表示的时间需应用道头215-216字节指定的因子。

short uphole_time_at_source;

// 97-98 检波器组处井口时间(毫秒)

short uphole_time_at_group;

// 99-100 震源的静校正量(毫秒)

short source_static_correction;

// 101-102 检波器组的校正量(毫秒)

short group_static_correction;

// 103-104 应用的总静校正量(毫秒)(如没有应用静校正量为零)

short total_static_applied;

// 105-106

// 延迟时间A――以毫秒表示的240字节道识别头的结束和时间断点之间的时间。当时间断点出现在头之后,该值为正;当时间断点出现在头之前,该值为负。时间断点是最初脉冲,它由辅助道记录或由其他记录系统指定。

short lag_time_a;

// 107-108 延迟时间B――以毫秒表示的时间断点到能量源起爆时间之间的时间。可正可负。

short lag_time_b;

// 109-110 记录延迟时间――以毫秒表示的能量源起爆时间到数据采样开始记录之间的时间。在SEG

// Y修订版0中本条用来表示深水作业,如果数据记录不从0时间开始。该条可为负值以适应负的起始时间(即数据记录在零时间之前,假设静校正量应用于数据道的结果)。若某非零值(正或负)记录在该条,它造成的影响的注释应出现在原文文件头。

short delay_according_time;

// 111-112 起始切除时间(毫秒)

short brute_time_start;

// 113-114 终止切除时间(毫秒)

short mute_time_end;

// 115-1167 该道采样点数。强烈推荐所有类型的数据使用

short samples_in_this_trace;

// 117-1187 该道采样间隔(微秒)。

// 一道记录的字节数必须和写在道头中的采样点数一致。对所有记录介质都重要,但对正确处理磁盘文件中的SEG Y数据尤为关键(见附录C)。

// 若二进制文件头中的3503-3504字节设置了固定长度道标志,SEG

// Y文件每道的采样间隔和采样点数必须与二进制文件头所记录的值一致。若没有设置固定长道标志,采样间隔和采样点数可能每道变化。

short sample_intervall;

/**

* 119-120 野外仪器增益类型: 1=固定 2=二进制 3=浮点 4…N=选用

*/

short gain_type_instruments;

// 121-122 仪器增益常数(分贝)

short instrument_gain_para;

// 123-124 仪器初始增益(分贝)

short instrument_initial_gain_para;

// 125-126 相关: 1=无 2=有

short has_relation;

// 127-128 起始扫描频率(赫兹)

short freq_start;

// 129-130 终止扫描频率(赫兹)

short freq_end;

// 131-132 扫描长度(毫秒)

short sweep_length;

// 133-134 扫描类型: 1=线性 2=抛物线 3=指数 4=其他

short sweep_type;

// 135-136 扫描道斜坡起始长度(毫秒)

short sweep_start_length;

// 137-138 扫描道斜坡终止长度(毫秒)

short sweep_end_length;

// 139-140 斜坡类型: 1=线性 2=cos2 3=其他

short sweep_gradient_type;

// 141-142 假频滤波频率(赫兹),若使用

short fake_filter_freq;

// _143-144 假频滤波坡度(分贝/倍频程)

short fake_filter_gradient;

// 145-146 陷波频率(赫兹),若使用

short trap_freq;

// 147-148 陷波坡度(分贝/倍频程)

short trap_gradient;

// 149-150 低截频率(赫兹),若使用

short low_freq;

// 151-152 高截频率(赫兹),若使用

short high_freq;

// 153-154 低截坡度(分贝/倍频程)

short low_gradient;

// 155-156 高截坡度(分贝/倍频程)

short hight_gradient;

// 157-158 数据记录的年――1975年标准没有说清应该记成2位还是4位还是两者都用。

// 除SGE Y修订版0外,年份需记录成完整的4位罗马历法年(即2001年应记录为2001(7D116))。

short year;

// 159-160 日(以格林尼治标准时间和通用协调时间为基准的公元日)

short day;

// 161-162 时(24小时制)

short hour;

// 163-164 分

short minute;

// 165-166 秒

short second;

// 167-168 时间基准码:

/**

* 1=当地 2=格林尼治标准时间 3=其他,应在扩展原文文件头的用户定义文本段解释 4=通用协调时间

*/

short datum_time;

// 169-170 道加权因子――最小有效位数定义为2-N伏。(N=0,1,…,32767)

short trace_factor;

// 171-172 滚动开关位置1的检波器组号

short switch_1_receiver_group;

// 173-174 野外原始记录中道号1的检波器组号

short switch_1_filed_receiver_group;

// 175-176 野外原始记录中最后一道的检波器组号

short last_trace_receiver_group;

// 177-178 间隔大小(滚动时甩掉的总检波器组数)

short interval;

// 179-180 相对测线斜坡起始或终止点的移动 1=下(或后) 2=上(或前)

short translate_start_or_end;

/**

* 81-184 该道的道集(CDP)位置X坐标(应用道头71-72字节的因子)。 参考坐标系应通过扩展头位置数据文本段识别(见D-1节)

*/

int cdp_x;

/**

* 185-188 该道的道集(CDP)位置Y坐标(应用道头71-72字节的因子)。 参考坐标系应通过扩展头位置数据文本段识别(见D-1节)

*/

int cdp_y;

/**

* 189-192 对于3-D叠后数据,本字段用来填纵向线号(In-line)。 若每个SEG Y文件记录一条纵向线,文件中所有道的该值应相同,

* 并且同样的值将记录在二进制文件头的3205-3206字节中。

*/

int in_line_num;

/**

* 193-196 对于3-D叠后数据,本字段用来填横向线号(Cross-line)。

* 它应与道头21-24字节中的道集(CDP)号的值一致,但这并不是实例。

*/

int cross_line_num;

/**

* 197-200 炮点号――这可能只应用于2-D叠后数据。

* 注意在此假设炮点号相对于特定道最靠近叠加(CDP)位置震源位置。如果不是这种情况,应在原文文件头中有注释解释炮点实际参考点。

*/

int source_power_num;

/**

* 201-202 应用于道头中197-200字节中炮点号的因子,以得到实际数值。 若为正,因子用作乘数;若为负,因子用作除数;

* 若为零;炮点号不用于因子作用(即它是一个整数。典型的值是-10, 允许炮点号小数点后有一位小数)。

*/

short source_power_factor;

/**

* 203-204 道值测量单位: -1=其他(应在数据采样测量单位文本段描述) 0=未知 1=帕斯卡(Pa) 2=伏特(V) 3=毫伏(mV)

* 4=安培(A) 5=米(m) 6=米每秒(m/s) 7=米每秒二次方(m/s2) 8=牛顿(N) 9=瓦特(W)

*/

short trace_measure_unit;

/**

* 205-210 转换常数――该倍数用于将数据道采样转换成转换单位 (道头211-212字节指定)。

* 本常数以4字节编码,尾数是两互补整数(205-208字节)和2字节, 十的指数幂是两互补整数(209-210字节)

* (即(205-208字节)×10**(209-210字节))。

*/

long trans_const_factor;

/**

* 211-212 转换单位――经乘以道头205-210字节中的转换常数后的 数据道采样测量单位。 -1=其他(应在数据采样测量单位文本段36页描述)

* 0=未知 1=帕斯卡(Pa) 2=伏特(V) 3=毫伏(mV) 4=安培(A) 5=米(m) 6=米每秒(m/s) 7=米每秒二次方(m/s2)

* 8=牛顿(N) 9=瓦特(W)

*/

short trans_unit;

/**

* 213-214 设备/道标识――与数据道关联的单位号或设备号 (即4368对应可控源号4368或20316对应2船3线16枪)。

* 本字段允许道关联横跨独立于道号的道集(道头25-28字节)

*/

short instrument_trace_id;

/**

* 215-216 在道头95-114字节给出的作用于时间的因子, 以得到真实的毫秒表示的时间值。

* 因子=1,+10,+100,+1000或+10000。 若为正,因子用作乘数;若为负,因子用作除数。为零设定因子为一。

*/

short time_factor;

/*

* 217-218

* 震源类型/方位――定义类型或能量源的方位。垂直项、横向项、纵向项作为正交坐标系的三个轴。坐标系轴的绝对角度方位可在面元网格定义文本段中定义

* (27页) -1到-n=其他(应在震源类型/方位文本段38页描述) 0=未知 1=可控震源――垂直方位 2=可控震源――横向方位

* 3=可控震源――纵向方位 4=冲击源――垂直方位 5=冲击源――横向方位 6=冲击源――纵向方位 7=分布式冲击源――垂直方位

* 8=分布式冲击源――横向方位 9=分布式冲击源――纵向方位

*/

short source_type_or_direction;

/**

* 219-224 相对震源方位的震源能量方向―― 正方位方向在道头217-218字节定义。 能量方向以度数长度编码(即347.8°编码成3478)

*/

long source_power_dirction;

/*

* 225-230 震源测量――描述产生道的震源效应。 测量可以简单,定量的测量如使用炸药总重量或气枪压力峰值 或可控源振动次数和扫描周期时间。

* 尽管这些简单的测量可接受,但最好使用真实的能量或工作测量单位。 本常数编码成4字节,尾数为两互补整数(225-228字节)和2字节,

* 十的指数幂是两互补整数(209-230字节)(即(225-228字节)×10**(229-230字节))。

*/

long source_power_measure;

/*

* 231-232 震源测量单位――用于震源测量、道头225-230字节的单位。 -1=其他(应在震源测量单位文本段39页描述) 0=未知

* 1=焦耳(J) 2=千瓦(kW) 3=帕斯卡(Pa) 4=巴(Bar) 4=巴-米(Bar-m) 5=牛顿(N) 6=千克(kg)

*/

short source_power_measure_unit;//

// 233-240 未赋值――为任选信息预留

public SegyTrace(byte[] b) {

int i = 0;

trace_sequence_number_within_line = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 4;

trace_sequence_number_within_reel = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 8;

original_field_record_number = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 12;

trace_sequence_number_within_original_field_record = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 16;

energy_source_point_number = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 20;

cdp_ensemble_number = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 24;

trace_sequence_number_within_cdp_ensemble = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 28;

trace_identification_code = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 30;

number_of_vertically_summed_traces_yielding_this_trace = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 32;

number_of_horizontally_stacked_traced_yielding_this_trace = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 34;

data_use = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 36;

distance_from_source_point_to_receiver_group = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 40;

receiver_group_elevation = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 44;

surface_elevation_at_source = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 48;

source_depth_below_surface = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 52;

datum_elevation_at_receiver_group = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 56;

datum_elevation_at_source = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 60;

water_depth_at_source = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 64;

water_depth_at_receiver_group = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 68;

scalar_for_elevations_and_depths = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 70;

scalar_for_coordinates = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 72;

x_source_coordinate = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 76;

y_source_coordinate = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 80;

x_receiver_group_coordinate = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 84;

y_receiver_group_coordinate = b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000)

| (b[i] << 24 & 0xff000000);

i = 88;

coordinate_units = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 90;

weathering_velocity = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 92;

subweathering_velocity = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 94;

uphole_time_at_source = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 96;

uphole_time_at_group = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 98;

source_static_correction = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 100;

group_static_correction = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 102;

total_static_applied = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 104;

lag_time_a = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 106;

lag_time_b = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 108;

delay_according_time = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 110;

brute_time_start = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 112;

mute_time_end = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 114;

samples_in_this_trace = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 116;

sample_intervall = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 118;

gain_type_instruments = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 120;

instrument_gain_para = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 122;

instrument_initial_gain_para = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 124;

has_relation = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 126;

freq_start = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 128;

freq_end = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 130;

sweep_length = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 132;

sweep_type = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 134;

sweep_start_length = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 136;

sweep_end_length = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 138;

sweep_gradient_type = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 140;

fake_filter_freq = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 142;

fake_filter_gradient = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 144;

trap_freq = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 146;

trap_gradient = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 148;

low_freq = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 150;

high_freq = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 152;

low_gradient = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 154;

hight_gradient = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 156;

year = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 158;

day = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 160;

hour = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 162;

minute = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 164;

second = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 166;

datum_time = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 168;

trace_factor = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 170;

switch_1_receiver_group = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 172;

switch_1_filed_receiver_group = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 174;

last_trace_receiver_group = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 176;

interval = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 178;

translate_start_or_end = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 180;

cdp_x = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 184;

cdp_y = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 188;

in_line_num = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 192;

cross_line_num = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 196;

source_power_num = b[i + 3] & 0xff | (b[i + 2] << 8 & 0xff00)

| (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000);

i = 200;

source_power_factor = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 202;

trace_measure_unit = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 204;

trans_const_factor = (long) ((b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000)) * Math

.pow(10, b[i + 5] | (b[i + 4] << 8 & 0xff00)));

i = 210;

trans_unit = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 212;

instrument_trace_id = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 214;

time_factor = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 216;

source_type_or_direction = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));

i = 218;

source_power_dirction = b[i + 5] | (b[i + 4] << 8 & 0xff00)

| (b[i + 3] << 16 & 0xff0000) | (b[i + 2] << 24 & 0xff000000)

| (b[i + 1] << 24 & 0xff000000) << 8

| (b[i] << 24 & 0xff000000) << 16;

i = 224;

source_power_measure = (long) ((b[i + 3] & 0xff

| (b[i + 2] << 8 & 0xff00) | (b[i + 1] << 16 & 0xff0000) | (b[i] << 24 & 0xff000000)) * Math

.pow(10, b[i + 5] | (b[i + 4] << 8 & 0xff00)));

i = 230;

source_power_measure_unit = (short) ((b[i + 1] & 0xff | (b[i] << 8 & 0xff00)));//

}

double[] data;// 存储道值

/**

* 按照IBM浮点数的方式读取数据,对应带头中的数据类型1

* @param buf

* @return

* @throws SegyFileFormatException

*/

public double[] readDataInIBMFormat(byte[] buf) throws SegyFileFormatException {

data = new double[samples_in_this_trace];

for (int i = 0; i < samples_in_this_trace; i++) {

// 按照字节位序一致的方式解读

int sign = (buf[4 * i] >> 7) & 0x01;

int exp = buf[4 * i] & 0x7f;

exp = exp - 64;

int frac = ((buf[4 * i + 1] << 16) & 0x00ff0000)

| ((buf[4 * i + 2] << 8) & 0x0000ff00)

| (buf[4 * i + 3] & 0x000000ff) & 0x00ffffff;

if (frac == 0) {

data[i] = 0;

continue;

}

double result = (1 - 2 * sign) * frac * Math.pow(2, 4 * exp - 24);

// 按照字节位序相反的方式解读,事实证明这样解读是不对的。

//

// byte b1=inverseBit(buf[4*i]);

// byte b2=inverseBit(buf[4*i+1]);

// byte b3=inverseBit(buf[4*i+2]);

// byte b4=inverseBit(buf[4*i+3]);

// int sign=(b1>>7 ) & 0x01;

//

// int exp=b1 & 0x7f;

// exp=exp-64;

// int frac=((b2<<16) & 0x00ff0000) |(( b3 <<8) & 0x0000ff00) | (b4

// & 0x000000ff) & 0x00ffffff;

// if(frac==0 ){

// data[i]=0;

// continue;

// }

// double result=(1-2*sign)*frac/Math.pow(2,24)*Math.pow(16, exp);

data[i] = result;

}

return data;

}

// 翻转字节序

private static byte inverseBit(byte b) {

String temp = Integer.toBinaryString(b & 0xff);

if (temp.length() >= 8) {

temp = temp.substring(temp.length() - 8);

} else {

// 前面补零

int len = temp.length();

for (int i = 0; i < 8 - len; i++) {

temp = "0" + temp;

}

}

System.out.println("转换前:" + temp);

// 把它翻转

String s = "";

for (int i = temp.length() - 1; i >= 0; i--) {

s = s + temp.charAt(i);

}

System.out.println("转换后:" + s);

return (byte) Integer.parseInt(s, 2);

}

public static void main(String[] args) {

System.out.println(Integer.toBinaryString(inverseBit((byte) -128)));

}

/**

* 得到震源地表高程

* 计算方法: (值-基准值)再乘或者除 高程因子,参看规范

* @return

*/

public int getSurface_elevation_at_source() {

if(scalar_for_elevations_and_depths>=0){

return (surface_elevation_at_source-datum_elevation_at_source)*scalar_for_elevations_and_depths;

}else

return -(surface_elevation_at_source-datum_elevation_at_source)/scalar_for_elevations_and_depths;

}

/**

* 得到检波器地表高程

* 计算方法: (值-基准值)再乘或者除 高程因子,参看规范

* @return

*/

public int getReceiver_group_elevation() {

if(scalar_for_elevations_and_depths>=0){

return (receiver_group_elevation-datum_elevation_at_receiver_group)*scalar_for_elevations_and_depths;

}else

return -(receiver_group_elevation-datum_elevation_at_receiver_group)/scalar_for_elevations_and_depths;

}

}

/**

* D-1

*

* 表5 位置数据文本段

文本段头和关键字 格式 注释

((SEG: 位置数据版本1.0)) 文本 文本段名称

以下关键字应用于所有参考坐标系(CRS)

CRS类型= 源于列举表:投影、地理、复合 投影=地图网格

地理=纬度、经度和3-DCRS情况下附加的椭圆高度

复合=准三维坐标系,由2D地理或带与重力相关的高度系统的投影CRS

CRS名称= 文本 参考坐标系名称

大地基准名称= 文本 大地基准名称

本初子午线名称= 文本 若不是“格林威治”,强制

注意:多数,但不是所有,参考坐标系用格林威治作本初子午线(PM)

PM格林威治经度= 实数 CRS的本初子午线相对于格林威治子午线的经度,当在格林威治的东边为正。若本初子午线名称=“格林威治”是不需要

PM格林威治经度单位名称= 文本 若本初子午线名称=“格林威治”是不需要

椭圆名称= 文本

椭圆半主轴= 实数

半主轴单位名称= 文本

椭圆扁率倒数= 实数

坐标轴1名称= 文本 道头73-76,81-84和181-184字节的坐标轴(CS)的名称或缩写。例如东,X,E,或经度

CS轴1方位= 文本 轴1的正方向。例如:“东”或“北”

坐标轴2名称= 文本 道头77-80,85-88和185-188字节的坐标轴(CS)的名称或缩写。例如北,Y,N,或纬度

CS轴2方位= 文本 轴2的正方向。例如:“北”或“东”

垂直基准名称= 文本 垂直基准名称。若使用椭圆高度不需要(多数高度和深度是与重力相关的,而不是椭圆)

坐标轴3名称= 文本 道头41-68字节的高程和深度的坐标轴的名称或缩写。例如:重力相关高度,椭圆高度

CS轴3方位= 文本 轴3的正方向。例如:“上”

以下关键字是投影坐标参考系额外需要的,它是当在字节89-90指明坐标单位是长度或当面元网格定义文本段或数据地理范围扩展文本段或覆盖区域文本段包含在扩展文件头中时

投影区名称= 文本

投影方法名称= 文本 例如:“横向墨卡托”,“朗伯正交圆锥投影(1SP)”,“朗伯正交圆锥投影(2SP)”。

数字和参数类型决定于地图投影方法。对于横向墨卡托和朗伯正交圆锥投影(1SP),需要的5个参数是:

• 原始自然纬度

• 原始自然经度

• 原始自然比例因子

• 东偏

• 北偏

地图投影方法为朗伯正交圆锥投影(2SP)需要的6个参数见下面的示例

投影参数1名称= 文本

投影参数1数值= 实数

投影单位1单位名称= 文本

投影参数2名称= 文本

投影参数2数值= 实数

投影单位2单位名称= 文本



投影参数7名称= 文本

投影参数7数值= 实数

投影单位7单位名称= 文本

D-1.2. 位置数据文本段示例

((SEG: 位置数据版本1.0))

CRS类型= 投影

CRS名称= NAD27/德克萨斯南中心

大地基准名称= 北美1927基准

椭圆名称= 克拉克1866

椭圆半主轴= 6378206.4

半主轴单位名称= 米

椭圆扁率倒数= 294.9786982

坐标轴1名称= Y

CS轴1方位= 北

坐标轴2名称= X

CS轴2方位= 东

投影区名称= 得克萨斯CS27南中心区

投影方法名称= 朗伯正交圆锥投影(2SP)

投影参数1名称= 原始纬度偏差

投影参数1数值= 27.5

投影单位1单位名称= DDD.MMSSsss

投影参数2名称= 原始经度偏差

投影参数2数值= -99

投影单位2单位名称= 度

投影参数3名称= 第一标准平行纬度

投影参数3数值= 28.23

投影单位3单位名称= DDD.MMSSsss

投影参数4名称= 第二标准平行纬度

投影参数4数值= 30.17

投影单位4单位名称= DDD.MMSSsss

投影参数5名称= 原始东偏

投影参数5数值= 2000000.0

投影单位5单位名称= US测量步

投影参数6名称= 原始北偏

投影参数6数值= 0.0

投影单位6单位名称= US测量步

*

*

*/

———————————————vtk显示类,用于测试和现实segy读取的结果—————————————————————————————————

package org.sinoprobe.segy;

import java.awt.BorderLayout;

import java.awt.event.ActionEvent;

import java.awt.event.ActionListener;

import java.io.UnsupportedEncodingException;

import javax.swing.JButton;

import javax.swing.JFrame;

import javax.swing.JPanel;

import javax.swing.SwingUtilities;

import vtk.vtkActor;

import vtk.vtkAxesActor;

import vtk.vtkBandedPolyDataContourFilter;

import vtk.vtkCamera;

import vtk.vtkCanvas;

import vtk.vtkCellArray;

import vtk.vtkConeSource;

import vtk.vtkDEMReader;

import vtk.vtkDoubleArray;

import vtk.vtkImageDataGeometryFilter;

import vtk.vtkImageShrink3D;

import vtk.vtkLODActor;

import vtk.vtkLookupTable;

import vtk.vtkNativeLibrary;

import vtk.vtkOrientationMarkerWidget;

import vtk.vtkPanel;

import vtk.vtkPoints;

import vtk.vtkPolyData;

import vtk.vtkPolyDataMapper;

import vtk.vtkPolyDataNormals;

import vtk.vtkRenderWindow;

import vtk.vtkRenderWindowInteractor;

import vtk.vtkRenderWindowPanel;

import vtk.vtkRenderer;

import vtk.vtkSTLReader;

import vtk.vtkShrinkPolyData;

import vtk.vtkWarpScalar;

/**

* An application that displays a 3D cone. The button allow to close the

* application

*/

public class SimpleVTK2 extends JPanel implements ActionListener {

private static final long serialVersionUID = 1L;

private vtkCanvas renPanel;

private JButton exitButton;

// -----------------------------------------------------------------

// Load VTK library and print which library was not properly loaded

static {

if (!vtkNativeLibrary.LoadAllNativeLibraries()) {

for (vtkNativeLibrary lib : vtkNativeLibrary.values()) {

if (!lib.IsLoaded()) {

System.out.println(lib.GetLibraryName() + " not loaded");

}

}

}

vtkNativeLibrary.DisableOutputWindow(null);

}

private vtkPolyData getGridPolyData(SegyReader segy){

int traceNums=segy.getTraceList().size();

double interval=segy.getReelHdr().sample_data_interval_ms/1000.0;//原来是微秒,转化为毫秒

int sampleNums=segy.getReelHdr().samples_per_trace;

double left=-1000;

double top=segy.getTraceList().get(0).getReceiver_group_elevation();//得到第一组的高程

int space=40;//缺省定义的道间距

int vSpace=2;

double[] minMax=segy.getMinMaxValue();

double dist=minMax[1]-minMax[0];

System.out.println("测量震幅的最大值:"+minMax[1]+" 最小值:"+minMax[0]+" 两者差值:"+dist);

//缺省认为最大最小值跨越 ? 个道间距

double amplitudeFactor=8*space/dist;

System.out.println("震幅计算因子:"+amplitudeFactor);

//此处以5:1的比例抽样,以加快速度,另外避免机器内存溢出

//如果为1,则依照原始数据绘图;

//为了避免偏振效应,需要按照奇数进行抽样,所以,不能为偶数

int sampleFactor=5;

vtkPolyData rtn=new vtkPolyData();

//生成所有的点

vtkPoints points=new vtkPoints();

//多边形

vtkCellArray cellArry=new vtkCellArray();

//标量值

vtkDoubleArray scalars=new vtkDoubleArray();

//开始构建数据

for(int i=0;i<traceNums;i++){

double x=left+i*space;

SegyTrace trace=segy.getTraceList().get(i);

cellArry.InsertNextCell(trace.samples_in_this_trace/sampleFactor + ((trace.samples_in_this_trace %sampleFactor >0) ? 1:0));

for(int j=0;j<trace.samples_in_this_trace;j=j+sampleFactor){

double z=top+j*vSpace;

int id=points.InsertNextPoint(x+trace.data[j]*amplitudeFactor,0,z);

scalars.InsertNextTuple1(trace.data[j]);

cellArry.InsertCellPoint(id);

}

}

//设置点

rtn.SetPoints(points);

//设置为线

rtn.SetLines(cellArry);

rtn.GetPointData().SetScalars(scalars);

System.out.println("vtk数据构造完成,开始显示");

return rtn;

}

// -----------------------------------------------------------------

public SimpleVTK2(SegyReader segy) {

super(new BorderLayout());

// build VTK Pipeline

//

renPanel = new vtkRenderWindowPanel();

vtkPolyData seg=getGridPolyData(segy);

vtkPolyDataMapper demMapper = new vtkPolyDataMapper();

demMapper.SetInput(seg);

demMapper.ScalarVisibilityOff();

vtkActor demActor = new vtkActor();

demActor.SetMapper(demMapper);

demActor.GetProperty().SetColor(0,0, 0);

// # Create the RenderWindow, Renderer and both Actors

vtkRenderer ren = renPanel.GetRenderer();

vtkRenderWindowInteractor iren = renPanel.getRenderWindowInteractor();

// # Add the actors to the renderer, set the background and size

ren.AddActor(demActor);

// 添加坐标系

vtkAxesActor axes = new vtkAxesActor();

vtkOrientationMarkerWidget widget = new vtkOrientationMarkerWidget();

widget.SetOutlineColor(0.9300, 0.5700, 0.1300);

widget.SetOrientationMarker(axes);

widget.SetInteractor(iren);

widget.SetViewport(0.7, 0.7, 1, 1);

widget.SetEnabled(1);

widget.InteractiveOff();

ren.SetBackground(1, 1, 1);

vtkCamera cam = new vtkCamera();

cam.SetPosition(0, 1000, 0);

cam.SetFocalPoint(0, 0, 0);

cam.SetViewUp(0, 0, -1);

ren.SetActiveCamera(cam);

ren.ResetCamera();

cam.Zoom(2);

add(renPanel, BorderLayout.CENTER);

}

/** An ActionListener that listens to the button. */

public void actionPerformed(ActionEvent e) {

if (e.getSource().equals(exitButton)) {

System.exit(0);

}

}

}

______________________________代码上传到此结束__________________________________________________________________________

几句题外话:

1,Segy数据读取类里只实现了IBM浮点数一种转换办法,其他的还没弄完,等弄完了再分享给大家;

2,在进行数据转换的时候,尽量多用移位和位操作,这样比计算加减乘除快很多;

3,对于增益、修正等参数,还没有做,只是对原始测量数据做了个显示;

4,vtk显示效果还不错,用java来使用vtk,是个好选择,建议大家多尝试。谁说java不能做这种大规模数据处理和运算的呢?谁说Java的显示效果不如c++呢?偶私下做了对比,这部分的执行效果,与网上下载的c++所写的segy显示的程序,执行速度要快很多,偶尝试了一个30M的单炮文件和100M左右的叠加后的文件,执行速度比c++的那个vsegyview快不少。而且,vtk可以进行任意的放大缩小漫游,以及三维转换,不是一个静态图。

5,配色和参数配置上,效果不理想,专业人员看了可能觉得差距太大,这部分有待于完善。欢迎大家交流。谢谢!

附上两张成图的效果:



上面是叠加后的

下面是单炮记录

内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: