6

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)
Vladimir F Героям слава
  • 57,977
  • 4
  • 76
  • 119
marcelorodrigues
  • 919
  • 5
  • 12
  • 22
  • Possible duplicate of [Reading a direct access fortran unformatted file in Python](http://stackoverflow.com/questions/10475839/reading-a-direct-access-fortran-unformatted-file-in-python) – Alexander Vogt May 30 '16 at 22:16
  • Possible duplicate of [Python: Reading Fortran Binary file using numpy or scipy](http://stackoverflow.com/questions/30307305/python-reading-fortran-binary-file-using-numpy-or-scipy) – glls May 30 '16 at 22:19
  • Thank you guys but I try both solutions. In my fortran code there a 3d variable and real*4 precision. – marcelorodrigues May 30 '16 at 22:33

3 Answers3

5

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).

roygvib
  • 7,218
  • 2
  • 19
  • 36
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.]]
Serenity
  • 35,289
  • 20
  • 120
  • 115
  • 3
    Your answer appears to start "you should change the Fortran program to use direct, rather than sequential, access in the output". If that's the intention (which is a significant change), perhaps it would be better to be more explicit about that? – francescalus May 31 '16 at 07:40
  • Yes, @francescalus was right, most of the time changing the Fortran code to create direct-access files is simply not an option. And yes, this numpy solution does NOT work for classical unformatted Fortran binaries. – CoolKoon Dec 08 '22 at 00:45
-1

Easiest way will be to use data_py package. To install enter pip install data-py

Example use

from data_py import datafile

NoOfLines=0   
lineNumber=2  # Line number to be read (Excluding lines starting with '#')
df1=datafile("C:/Folder/SubFolder/data-file-name.txt")
df1.separator=","  # No need to specify if separator is space(" "). For 'tab' separated values use '\t'
NoOfLines=df1.lines  # Total number of lines in the data file (Excluding lines starting with '#')
Col=["Null"]*5  # This will create 5 column variables with an intial string 'Null'. 
                # Number of column variables (here 5) should not be greater than number of columns in data file.
df1.read(Col,lineNumber) # Will read first five columns from the data file at the line number given, and stores in Col.
print(Col)

For details visit: https://www.respt.in/p/python-package-datapy.html

Dwaipayan Deb
  • 11
  • 1