Python Read Fortran Binary File

Question

I'm trying to read a binary file output from Fortran code below, but the results are not the same from output file.

Fortran 77 code:

    program test
    implicit none
    integer i,j,k,l
    real*4       pcp(2,3,4)
    open(10, file='pcp.bin', form='unformatted')
    l = 0
    do i=1,2
      do j=1,2
        do k=1,2
          print*,k+l*2
          pcp(i,j,k)=k+l*2
          l = l + 1
        enddo
      enddo
    enddo
    do k=1,4
       write(10)pcp(:,:,k)
    enddo
    close(10)
    stop
    end

I'm trying to use the Python code below:

from scipy.io import FortranFile
f = FortranFile('pcp.bin', 'r')
a = f.read_reals(dtype=float)
print(a)

Answer 1

Because you are writing real*4 data on a sequential file, simply try replacing dtype=float to dtype='float32' (or dtype=np.float32) in read_reals():

>>> from scipy.io import FortranFile
>>> f = FortranFile( 'pcp.bin', 'r' )
>>> print( f.read_reals( dtype='float32' ) )
[  1.   9.   5.  13.   0.   0.]
>>> print( f.read_reals( dtype='float32' ) )
[  4.  12.   8.  16.   0.   0.]
>>> print( f.read_reals( dtype='float32' ) )
[ 0.  0.  0.  0.  0.  0.]
>>> print( f.read_reals( dtype='float32' ) )
[ 0.  0.  0.  0.  0.  0.]

The obtained data correspond to each pcp(:,:,k) in Fortran, as verified by

do k=1,4
   print "(6f8.3)", pcp(:,:,k)
enddo

which gives (with pcp initialized to zero)

   1.0   9.0   5.0  13.0   0.0   0.0
   4.0  12.0   8.0  16.0   0.0   0.0
   0.0   0.0   0.0   0.0   0.0   0.0
   0.0   0.0   0.0   0.0   0.0   0.0

But because >>> help( FortranFile ) says

An example of an unformatted sequential file in Fortran would be written as::

OPEN(1, FILE=myfilename, FORM='unformatted')

WRITE(1) myvariable

Since this is a non-standard file format, whose contents depend on the compiler and the endianness of the machine, caution is advised. Files from gfortran 4.8.0 and gfortran 4.1.2 on x86_64 are known to work.

Consider using Fortran direct-access files or files from the newer Stream I/O, which can be easily read by numpy.fromfile.

it may be simpler to use numpy.fromfile() depending on cases (as shown in StanleyR's answer).

Answer 2

Use nupy.fromfile (http://docs.scipy.org/doc/numpy/reference/generated/numpy.fromfile.html)

I guess you have missed something in fortran code, to write to binary file apply this code:

program test
implicit none
integer i,j,k,l, reclen
real*4       pcp(2,3,4)

inquire(iolength=reclen)pcp(:,:,1)
open(10, file='pcp.bin', form='unformatted', access = 'direct', recl = reclen)
pcp = 0
l = 0
do i=1,2
do j=1,2
do k=1,2
   print*,i,j,k,k+l*2
   pcp(i,j,k)=k+l*2
   l = l + 1
enddo
enddo
enddo
do k=1,4
   write(10, rec=k)pcp(:,:,k)
enddo
close(10)
end

To read file by python:

import numpy as np
with open('pcp.bin','rb') as f:
    for k in xrange(4):
        data = np.fromfile(f, dtype=np.float32, count = 2*3)
        print np.reshape(data,(2,3))

Output:

[[  1.   9.   5.]
 [ 13.   0.   0.]]
[[  4.  12.   8.]
 [ 16.   0.   0.]]
[[ 0.  0.  0.]
 [ 0.  0.  0.]]
[[ 0.  0.  0.]
 [ 0.  0.  0.]]