用Python写Fortran非格式化文件

Question

我有一些用Fortran77写的单精度小端格式的未格式化数据文件。我正在用Python读取这些文件，使用的命令如下：

import numpy as np
original_data = np.dtype('float32')
f = open(file_name,'rb')                                                                                                 
original_data = np.fromfile(f,dtype='float32',count=-1)                                                                            
f.close()

在用Python对数据进行了一些处理后，我（想要）把它们以原来的格式写回去，使用的命令如下：

out_file = open(output_file,"wb")                                                                                             
s = struct.pack('f'*len(manipulated_data), *manipulated_data)                                                                     
out_file.write(s)
out_file.close()

但是似乎没有成功。有没有人知道用Python把数据写回原来的Fortran未格式化格式的正确方法是什么？

问题的详细情况：

我能够从Fortran读取到处理后的最终文件。但是，我想用一个软件（Paraview）来可视化这些数据。为此，我把未格式化的数据文件转换成*h5格式。我能够使用h5工具将原始数据和处理后的数据都转换成h5格式。但是，Paraview可以读取从原始数据创建的*h5文件，却无法读取从处理后的数据创建的*h5文件。我猜测在转换过程中可能丢失了一些信息。

这是我在Fortran中打开Python写入的文件的方式（单精度数据）：

open (in_file_id,FILE=in_file,form='unformatted',access='direct',recl=4*n*n*n)

这是我用Fortran写入原始未格式化数据的方式：

open(out_file_id,FILE=out_file,form="unformatted")

这些信息够吗？

Answer 1

这段内容是在讲如何创建一个没有格式的顺序访问文件：

open(out_file_id,FILE=out_file,form="unformatted")

假设你在写一个单一的数组 real a(n,n,n)，如果你直接用 write(out_file_id)a 来写，你会看到文件大小是 4*n^3+8 字节。多出来的 8 字节是一个 4 字节的整数（=4n^3），它在记录的开头和结尾各重复一次。

第二种形式：

open (in_file_id,FILE=in_file,form='unformatted',access='direct',recl=4*n*n*n)

打开的是直接访问，这种方式没有那些头信息。现在写入时，你需要用 write(unit,rec=1)a。如果你用直接访问去读取你的 顺序 访问文件，虽然不会出错，但你会把那个整数头当成浮点数（垃圾数据）读取到（1,1,1）的数组值里，然后其他数据都会错位。你说你可以用 Fortran 来读取，但你是想确认你读取到的真的是你期待的内容吗？

解决这个问题的最好办法是修改你原来的 Fortran 代码，让它在读写时都使用没有格式的直接访问。这样你就能得到一个“普通”的原始二进制文件，没有头信息。

另外，在你的 Python 代码中，你需要先读取那个 4 字节的整数，然后再读取数据。在输出时，你可以选择是否把整数头信息放回去，这取决于你的 Paraview 过滤器的要求。

---------- 下面是用 Python 来读取、修改和写入一个包含单个记录的没有格式的顺序 Fortran 文件的代码：

import struct
import numpy as np
f=open('infile','rb')
recl=struct.unpack('i',f.read(4))[0]
numval=recl/np.dtype('float32').itemsize
data=np.fromfile(f,dtype='float32',count=numval)
endrec=struct.unpack('i',f.read(4))[0]
if endrec is not recl: print "error unexpected end rec"
f.close()
f=open('outfile') 
f.write(struct.pack('i',recl))
for i in range(0,len(data)):data[i] = data[i]**2  #example data modification
data.tofile(f)
f.write(struct.pack('i',recl)

你只需循环处理多个记录……注意这里的数据是作为向量读取的，并假设都是浮点数。当然，你需要知道实际的数据类型才能正确使用它……另外，根据平台的不同，你可能还需要处理字节顺序的问题。

Answer 2

你有没有试过使用被处理过的数据数组的 .tofile 方法？这个方法可以按照C语言的顺序把数组写入文件，而且能够写成纯二进制格式。

关于 .tofile 的说明文档也提到，这个方法和下面的写法是一样的：

with open(outfile, 'wb') as fout:
    fout.write(manipulated_data.tostring())