Как прочитать неформатированный двоичный файл fortran 77 в python

У меня есть скрипт в IDL, который читает неформатированный двоичный файл (F77) и выводит его как файл .sav, но я хочу преобразовать этот скрипт в python и сохранить как файл .npz, и у меня возникают проблемы с чтением строки.

IDL-код:

;create save file for QBO model output

;---------------------------------------------------------------------
;input data here (Only adjust this info)
;---------------------------------------------------------------------

  time=7300  ;duration from fortran code
  tstep=5   ;time step from fortran code
  inputfile='/home/cwilliams/Metr_205/qbo/uvwtom.'; binary file from fortran
  outputsavefile='~cwilliams/Metr_205/qbo/qbo_FULLX.sav'; name of save file with     variables
;--------------------------------------------------------------------------------
;--------------------------------------------------------------------------------





;==============================================================================
;   Do Not Adjust the Following Code
;==============================================================================
fname='dat'
nstep=time/tstep ;time/timestep
nm=nstep-1

um=fltarr(17,34)
vm=fltarr(17,34)
wm=fltarr(17,34)
tm=fltarr(17,34)
pm=fltarr(17,34)
om=fltarr(17,34)
vpm=fltarr(17,34)
vdm=fltarr(17,34)


u=fltarr(17,34,nstep)
v=fltarr(17,34,nstep)
w=fltarr(17,34,nstep)
temp=fltarr(17,34,nstep)
pressure=fltarr(17,34,nstep)
ozone=fltarr(17,34,nstep)

date=findgen(nstep)*(2.*tstep/365.)
ulevels2=[-30, -25, -20, -15, -10, -5, 5, 10, 15,  20, 25, 30]
ulevels3=findgen(15)*3.-12
lab=findgen(20)*0+1
lat=findgen(17)*2.7-1.35

ht=findgen(34)*0.5+17.


openr,37,inputfile+fname,/f77_unformatted & title='base'
for n=0,nm do begin
    readu,37,um,vm,wm,tm,om,pm
    u(*,*,n)=um/100.
    v(*,*,n)=vm
    w(*,*,n)=wm
    temp(*,*,n)=tm*700000./2.87e6
    pressure(*,*,n)=pm/10000.
    ozone(*,*,n)=om
endfor
close,37

save,filename=outputsavefile,u,v,w,temp,pressure,ozone,date,ulevels2,ulevels3,lab,$
lat,ht
stop
end

Код Python: #создать файл сохранения для вывода модели QBO

#---------------------------------------------------------------------
#input data here (Only adjust this info)
#---------------------------------------------------------------------

time=7300  #duration from fortran code
tstep=5   #time step from fortran code
inputfile='/home/cwilliams/metr51/lab12/uvwtom.dat'# binary file from fortran
outputsavefile='~cwilliams/metr51/lab12/qbo_FULLX.sav'# name of save file with variables
#--------------------------------------------------------------------------------
#--------------------------------------------------------------------------------


import numpy as n


#==============================================================================
#   Do Not Adjust the Following Code
#==============================================================================
nstep=time/tstep #time/timestep
nm=nstep-1

um=n.empty((17,34))
vm=n.empty((17,34))
wm=n.empty((17,34))
tm=n.empty((17,34))
pm=n.empty((17,34))
om=n.empty((17,34))
vpm=n.empty((17,34))
vdm=n.empty((17,34))


u=n.empty((17,34,nstep))
v=n.empty((17,34,nstep))
w=n.empty((17,34,nstep))
temp=n.empty((17,34,nstep))
pressure=n.empty((17,34,nstep))
ozone=n.empty((17,34,nstep))

date=n.arange(nstep)*(2.*tstep/365.)
ulevels2=n.array([-30, -25, -20, -15, -10, -5, 5, 10, 15,  20, 25, 30])
ulevels3=n.arange(15)*3.-12
lab=n.arange(20)*0+1
lat=n.arange(17)*2.7-1.35

ht=n.arange(34)*0.5+17.

Вот где мне нужна помощь:

f=open(inputfile,'rb')     
data=f.read()
for i in range(nm+1):
#    readu,37,um,vm,wm,tm,om,pm
#    u(:,:,i)=um/100.
#    v(:,:,i)=vm
#    w(:,:,i)=wm
#    temp(:,:,i)=tm*700000./2.87e6
#    pressure(:,:,i)=pm/10000.
#    ozone(*,*,n)=om
#
#n.savez(filename=outputsavefile,u=u,v=v,w=w,temp=temp,pressure=pressure,ozone=ozone,date=date,ulevels2=ulevels2,ulevels3=ulevels3,lab=lab,lat=lat,ht=ht)

Я знаю, что есть проблема с порядком строк/столбцов между IDL и Python, но я думаю, что смогу исправить это, как только прочитаю код.


python-2.7 fortran idl-programming-language
person CRogers    schedule 29.04.2014    source источник
comment
stackoverflow.com/questions/14985311/   -  person agentp    schedule 30.04.2014
comment
Какая у тебя беда?   -  person Vladimir F    schedule 30.04.2014
comment
@george Я видел этот пост, и он не помог, так как я читаю несколько переменных, которые имеют несколько измерений.   -  person CRogers    schedule 30.04.2014
comment
@VladimirF Мне нужна помощь с последним разделом (кода), преобразующим его из IDL в python.   -  person CRogers    schedule 30.04.2014
comment
код в моем ответе очень близок к тому, что вам нужно. Вам нужно знать типы данных (в примере float32) и выполнять некоторые циклы. Я не знаю IDL, поэтому я не могу точно сказать, что он делает.   -  person agentp    schedule 30.04.2014

Ответы (1)


arrow_upward
2
arrow_downward

Во-первых, я должен предположить, что вы написали исходный файл Fortran, используя последовательный доступ, а не прямой доступ, это очень важное отличие. Я знаю IDL, и команда чтения, которую вы настроили в своем примере, согласуется с последовательным доступом.

Это важно, потому что последовательный доступ Fortran добавляет «маркеры записи» в файл для каждого вызова WRITE(), поэтому вам придется учитывать это в своей процедуре чтения Python, как и в ответе, который Джордж связал в комментариях.

Что касается порядка строк/столбцов, то на самом деле имеет значение только то, какое измерение быстрее итерируется в памяти. Двумерный массив Fortran, например, будет записываться в память с первым измерением массива как с самой быстрой итерацией. Когда вы читаете то же самое в 2D-переменной Python, это будет последнее измерение — имейте это в виду при изменении формы плоского массива, и вам придется транспонировать его, если вы хотите сохранить свою индексацию.

Я думаю, что следующий код будет делать примерно то, что вы хотите, для краткости я просто включил ваши переменные um и vm:

import numpy as np

um=np.empty((34,17), dtype='float32') # Make these dimensions "backwards" for easier reshaping
vm=np.empty((34,17), dtype='float32') # Also watch out, I believe the default type is float64

f = open(inputfile,'rb')
recl = np.zeros(1,dtype=np.uint32)
for i in range(nm+1):
    recl = np.fromfile(f, dtype='uint32', count=1)
    tmpu = np.fromfile(f, dtype='float32', count=um.size) # These arrays will be flat
    tmpv = np.fromfile(f, dtype='float32', count=vm.size)
    recl = np.fromfile(f, dtype='uint32', count=1)

    um = np.transpose(np.reshape(tmpu, um.shape))
    vm = np.transpose(np.reshape(tmpv, vm.shape))

Затем вы также можете делать все, что пожелаете, в цикле с другими вашими переменными, и вы должны иметь возможность индексировать um и vm, как вы привыкли в IDL или Fortran. Это может быть не самый простой или лучший способ сделать это, но я думаю, что это, по крайней мере, поможет вам начать.

person Ajean    schedule 01.05.2014