回答此问题可获得 20 贡献值,回答如果被采纳可获得 50 分。
<p>我对Python还不太熟悉,我从头开始用Python编写可视化代码(以避免使用昂贵的专有程序,如IDL)。到目前为止,我一直使用IDL和gnuplot。我想做的是:</p>
<p>我使用fortran将二维数组写入未格式化的直接访问文件,我希望能够用python读取这些文件。下面给出了准确的测试代码。实际的代码是一个巨大的并行代码,但数据输出几乎是完全相同的格式。</p>
<pre><code>program binary_out
implicit none
integer :: i,j,t,rec_array
double precision, dimension(100,100) :: fn
double precision, parameter :: p=2*3.1415929
INQUIRE(IOLENGTH=rec_array) fn
open(unit=10,file='test',status='new',form='unformatted',access='direct',recl=rec_array)
fn=0
write(10,rec=1) fn
do t=1,3
do i=1,100
do j=1,100
fn(i,j)=sin(i*p*t/100)*cos(j*p*t/100)
enddo
enddo
write(10,rec=t+1) fn
enddo
close(10)
end program binary_out
</code></pre>
<p>对于t=1,程序应该给我零,对于t值的增加,应该给我增加“孤岛”的数量。但是当我使用下面给出的python代码阅读它时,我得到的只是零。如果删除第一个零的write语句,那么无论在python代码中使用的“time slice”值是什么,都只得到第一个时间片。我目前掌握的代码是:</p>
<pre><code>#!/usr/bin/env python
import scipy
import glob
import numpy as np
import matplotlib.pyplot as plt
import os, sys
from pylab import *
def readslice(inputfilename,field,nx,ny,timeslice):
f=open(inputfilename,'r')
f.seek(timeslice*nx*ny)
field=np.fromfile(inputfilename,dtype='d',count=nx*ny)
field=np.reshape(field,(nx,ny))
return field
a=np.dtype('d')
a=readslice('test',a,100,100,2)
im=plt.imshow(a)
plt.show()
</code></pre>
<p>我希望def readslice能够在时间间隔等于I的情况下在I处读取记录。为此,我尝试使用f.seek,但它似乎不起作用。numpy.fromfile似乎开始读取第一条记录本身。如何使numpy.from file从文件中的特定点读取?</p>
<p>我仍在努力适应Python风格并深入研究文档。任何帮助和指点都将不胜感激。</p>