使用Python读SEGY道头

 

自己遇到的一个问题, WestGeco的3D地震数据, 提取导航数据:

 

 

#!/bin/env python

import sys

import struct

 

try:

        f=open(sys.argv[1],'rb')

except (IOError,Exception):

    print '''usage:

        scriptname segyfilename

'''

    sys.exit(1)

 

#skip EBCDIC header

try:

    f.seek(3200)

except Exception:

    print 'Oops! your file is broken..'

 

#read binary header

binhead = f.read(400)

ns = struct.unpack('>h',binhead[20:22])[0]

if ns < 0:

    print 'file read error'

    sys.exit(1)

 

#read trace header

while True:

    trchead = f.read(240)

    if trchead == '':

        break

    nav = trchead[224:228]+/

          trchead[228: