读取Fortran无格式文件时记录标记不一致

2 投票
3 回答
2282 浏览
提问于 2025-04-17 20:12

我正在尝试用Python读取一个非常大的Fortran未格式化的二进制文件。这个文件里有2的30次方个整数。

我发现记录标记有点让人困惑(第一个标记是-2147483639),不过我已经成功找到了数据结构(这些想要的整数都是相似的,因此和记录标记不同),并写了下面的代码(在这里得到了帮助)。

但是,我们可以看到每个记录开始和结束的标记并不相同。这是为什么呢?

是因为数据的大小太长了(536870910 = (2的30次方 - 4) / 2)吗?但是(2的31次方 - 1)/ 4 = 536870911大于536870910。

还是说只是数据文件的作者犯了一些错误?

还有一个问题,记录开始的标记是什么类型,是整数还是无符号整数?

fp = open(file_path, "rb")

rec_len1, = struct.unpack( '>i', fp.read(4) )
data1 = np.fromfile( fp, '>i', 536870910)
rec_end1, = struct.unpack( '>i', fp.read(4) )

rec_len2, = struct.unpack( '>i', fp.read(4) )
data2 = np.fromfile( fp, '>i', 536870910)
rec_end2, = struct.unpack( '>i', fp.read(4) )

rec_len3, = struct.unpack( '>i', fp.read(4) )
data3 = np.fromfile( fp, '>i', 4)
rec_end3, = struct.unpack( '>i', fp.read(4) )
data = np.concatenate([data1, data2, data3])

(rec_len1,rec_end1,rec_len2,rec_end2,rec_len3,rec_end3)

这里是读取到的记录长度的值,如上所示:

(-2147483639, -2176, 2406, 589824, 1227787, -18)

3 个回答

0

因为这个问题似乎经常出现,所以这里有一段Python工具代码,用来扫描一个二进制文件,判断它是否可能是一个Fortran的非格式化顺序访问文件。这个工具通过尝试几种不同的文件头格式来工作。当然,由于“非格式化”这种格式并没有统一标准,可能会有其他变种,但这段代码应该能覆盖到最常见的几种情况。

注意,左括号被转义了,所以如果你复制这段内容,可能需要把& #060;改回为小于号('<')。


def scanfbinary(hformat,file,fsize):
 """ scan a file to see if it has the simple structure typical of
     an unformatted sequential access fortran binary:
     recl1,<data of length recl1 bytes>,recl1,recl2,<data of length recl2 bytes>,recl2 ...
 """
 import struct
 print 'scan type',hformat,
 if 'qQ'.find(hformat[1])>=0:  hsize=8
 elif 'iIlL'.find(hformat[1])>=0:  hsize=4
 if hformat[0] == '<':  endian='little'
 elif hformat[0] == '>':  endian='big'
 print '(',endian,'endian',hsize,'byte header)',
 f.seek(0)
 nrec = 0
 while fsize > 0:
  h0=struct.unpack(hformat,f.read(hsize))[0]
  if h0 < 0 :   print 'invalid integer ',h0; return 1
  if h0 > fsize - 2*hsize:
   print 'invalid header size ',h0,' exceeds file size ',fsize
   if nrec > 0:print 'odd perhaps a corrupe file?'
   return 2
# to read the data replace the next line with code to read h0 bytes..
# eg 
#  import numpy
#  dtype = numpy.dtype('<i')
#  record=numpy.fromfile(f,dtype,h0/dtype.itemsize) 
  f.seek(h0,1)   
  h=struct.unpack(hformat,f.read(hsize))[0]
  if h0!=h :  print 'unmatched header';   return 3
  nrec+=1
  if nrec == 1:print
  if nrec < 10:print 'read record',nrec,'size',h
  fsize-=(h+2*hsize)
 print 'successfully read ',nrec,' records with unformatted fortran header type',hformat
 return 0
f=open('binaryfilename','r')
f.seek(0,2)
fsize=f.tell()
res=[scanfbinary(hformat,f,fsize) for hformat in ('<q','>q','<i','>i')]
if res.count(0)==0:
 print 'no match found, file size ',fsize, 'starts..'
 f.seek(0)
 for i in range(0,12): print f.read(2).encode('hex_codec'),
 print 

0

我发现使用 f2py 是让 Python 访问 Fortran 数据的一个更方便的方法。不过,记录标记的奇怪行为还是个问题。至少我们可以避免深入了解(有时让人困惑的)Fortran 非格式化文件结构。而且它和 numpy 结合得很好。

F2PY 用户指南和参考手册可以在 这里 找到。这里有一个 Fortran 源文件的例子,展示了如何打开和关闭文件,读取一维整数数组和二维浮点数组。注意以 !f2py 开头的注释,它们能帮助 f2py 更“聪明”。

要使用它,你需要把它封装成一个模块,然后导入到 Python 会话中。这样你就可以像调用 Python 函数一样调用这些函数了。

!ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
!cc                                                         cc
!cc      FORTRAN MODULE for PYTHON PROGRAM CALLING          cc
!cc                                                         cc
!ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

!Usage: 
!   Compile:   f2py -c fortio.f90 -m fortio
!   Import:    from fortio import *
!       or     import fortio
!Note:
!   Big endian: 1; Little endian: 0


!cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
SUBROUTINE open_fortran_file(fileUnit, fileName, endian, error)
  implicit none

  character(len=256) :: fileName
  integer*4 :: fileUnit, error, endian
  !f2py integer*4 optional, intent(in) :: endian=1
  !f2py integer*4 intent(out) :: error

  if(endian .NE. 0) then
     open(unit=fileUnit, FILE=fileName, form='unformatted', status='old', &
          iostat=error, convert='big_endian')
  else
     open(unit=fileUnit, FILE=fileName, form='unformatted', status='old', &
          iostat=error)
  endif
END SUBROUTINE 

!cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
SUBROUTINE read_fortran_integer4(fileUnit, arr, leng)
  implicit none

  integer*4 :: fileUnit, leng
  integer*4 :: arr(leng)
  !f2py integer*4 intent(in) :: fileUnit, leng 
  !f2py integer*4 intent(out), dimension(leng), depend(leng) :: arr(leng)

  read(fileUnit) arr
END SUBROUTINE

!cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
SUBROUTINE read_fortran_real4(fileUnit, arr, row, col)
  implicit none

  integer*4 :: fileUnit, row, col
  real*4 :: arr(row,col)
  !f2py integer*4 intent(in):: fileUnit, row, col
  !f2py real*4 intent(out), dimension(row, col), depend(row, col) :: arr(row,col)

  read(fileUnit) arr
END SUBROUTINE

!cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
SUBROUTINE close_fortran_file(fileUnit, error)
  implicit none

  integer*4 :: fileUnit, error
  !f2py integer*4 intent(in) :: fileUnit
  !f2py integer*4 intent(out) :: error

  close(fileUnit, iostat=error)
END SUBROUTINE 
2

最后,事情似乎变得更清楚了。

这里有一个英特尔Fortran编译器的用户和参考指南, 可以查看记录类型:可变长度记录这一部分。

如果记录的长度超过2,147,483,639字节,那么这个记录会被分成几个子记录。每个子记录的长度可以是1到2,147,483,639字节之间的任意值,包括这两个值。

开头的长度字段的符号位表示这个记录是否继续。结尾的长度字段的符号位表示前面是否有子记录。符号位的位置取决于文件的字节序。

一个继续的子记录的开头长度字段的符号位值为1。构成一个记录的最后一个子记录的开头长度字段的符号位值为0。一个有前面子记录的子记录的结尾长度字段的符号位值为1。构成一个记录的第一个子记录的结尾长度字段的符号位值为0。如果符号位的值是1,那么记录的长度是用二进制补码表示法存储的。

经过多次尝试,我意识到我被二进制补码表示法误导了,记录标记只是根据上述规则改变符号,而不是在符号位为1时转换为二进制补码表示法。不过,也有可能我的数据是用不同的编译器创建的。

下面是解决方案。

数据超过2GB,所以被分成了几个子记录。 我们看到第一个记录的开始标记是-2147483639, 所以第一个记录的长度是2147483639,这正好是子记录的最大长度,而不是我之前认为的2147483640,也不是-2147483639的二进制补码表示法2147483638。

如果我们跳过2147483639字节去读取记录的结束标记,你会得到2147483639, 因为这是第一个子记录,其结束标记是正数。

下面是检查记录标记的代码:

fp = open(file_path, "rb")
while 1:
    prefix, = struct.unpack( '>i', fp.read(4) )
    fp.seek(abs(prefix), 1)    #or read |prefix| bytes data as you want
    suffix, = struct.unpack( '>i', fp.read(4) )
    print prefix, suffix
    if abs(suffix) - abs(prefix): 
        print "suffix != prefix!"
        break
    if prefix > 0: break

还有屏幕打印的结果

-2147483639 2147483639
-2147483639 -2147483639
18 -18

我们可以看到记录的开始标记和结束标记总是相同,除了符号位不同。 三个记录的长度分别是2147483639、2147483639和18字节,不一定是4的倍数。所以第一个记录以某个整数的前3个字节结束,第二个记录以剩下的1个字节开始。

撰写回答