如何将未格式化的fortran文件(modflow输出)转换为numpy数组

2024-05-16 04:13:39 发布

您现在位置:Python中文网/ 问答频道 /正文

我有一个扩展名为hds的modflow输出文件。a file的谷歌驱动器链接。它是一个未格式化的fortran文件。我需要把它转换成numpy数组,我试过了:

floattype = 'f4'
a = np.fromfile("lake_example.hds", np.dtype([('kstp','i4'),('kper','i4'),('pertim',floattype),('totim',floattype),('text','a16'),('ncol','i4'),('nrow','i4'),('ilay','i4')]))
print a
print a.shape

{git链接^代码

我是一个尝试从this link开始的教程。因为我在linux上,所以不能使用flopy的方法从文件中获取输出数组。所以我试着用np.fromfile文件,但我在获取输出时遇到问题。在

我的输出是这样的:

^{pr2}$

我只包含了几行输出。在

关于标题信息,您可以参考它们的源代码:https://github.com/modflowpy/flopy/blob/master/flopy/utils/binaryfile.py#L30g


Tags: 文件numpy链接np数组file驱动器print
2条回答

代码与数据文件的结构不匹配:

00000000  2c 00 00 00 01 00 00 00  01 00 00 00 00 00 80 3f  |,..............?|
00000010  00 00 80 3f 20 20 20 20  20 20 20 20 20 20 20 20  |...?            |
00000020  48 45 41 44 0b 00 00 00  0b 00 00 00 01 00 00 00  |HEAD............|
00000030  2c 00 00 00 e4 01 00 00  00 00 c8 42 00 00 c8 42  |,..........B...B|
00000040  00 00 c8 42 00 00 c8 42  00 00 c8 42 00 00 c8 42  |...B...B...B...B|
00000050  00 00 c8 42 00 00 c8 42  00 00 c8 42 00 00 c8 42  |...B...B...B...B|

每个数据块都有自己的56字节头,包括: 3个整数(i4),2个浮点值(f4),16个字符,还有5个整数(i4):

^{pr2}$

然后数据块如下(11x11浮点值):

100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0
100.0 99.87313842773438 99.74627685546875 99.6254653930664 ...

我不确定这是否可以直接导入numpy数组。在

以下示例代码将在整个文件中循环,并为每个块提取头和数据:

#!/usr/bin/python

import struct
import numpy as np

infile = open("lake_example.hds","rb")

blockdata = []

while infile.read(1):
    infile.seek(-1,1)
    data = infile.read(56)
    n = struct.unpack('<3i4', data[0:12])
#    print n[0], n[1], n[2]
    n = struct.unpack('<2f4', data[12:20])
#    print n[0], n[1]
#    print data[20:36]
    n = struct.unpack('<5i4', data[36:56])
#    print n[0], n[1], n[2], n[3], n[4]
    ncol = n[0]
    nrow = n[1]
    a = np.fromfile(infile,dtype='f4',count=ncol*nrow).reshape((ncol,nrow))
    blockdata.append(a)
    data = infile.read(4)
    n = struct.unpack('<i4', data)
#    print n[0]

for block in blockdata:
    print block

很可能您还需要块头中的一些信息(请参阅print语句)。在

另请参见'flopy.utils.binaryfile文件模块': http://modflowpy.github.io/flopydoc/binaryfile.html

查看Flopy-3教程2(无侧限瞬态流模型),在绘图部分:

发件人:http://modflowpy.github.io/flopydoc/tutorial2.html

“headobj”首先定义为:

headobj = bf.HeadFile(modelname+'.hds')

头部提取如下:

^{pr2}$

在Debian上运行

相关问题 更多 >