关于模式中fortran读.a和.b文件不匹配的问题
本帖最后由 taotao0718 于 2014-12-2 09:18 编辑首先，先说一下，我是在学习一个模式，叫HYCOM，如果论坛里有HYCOM达人，那真是最好不过 了。
这里 .a和.b文件是网格构造生成的两个地形文件，.a文件里是ASCII码，.b戳开就能看，
.b文件的内容是：
1440 'idm '
960 'jdm '
-1 'mapflg'
plon:min,max = -179.99986 179.99980
plat:min,max = -79.00232 89.91841
qlon:min,max = -179.99979 179.99986
qlat:min,max = -79.02532 89.95042
ulon:min,max = -179.99977 179.99992
ulat:min,max = -79.00232 89.92188
vlon:min,max = -179.99983 179.99965
vlat:min,max = -79.02532 89.94513
pang:min,max = -179.95934 179.86636
scpx:min,max = 5120.85742 39728.21484
scpy:min,max = 5120.86133 39728.18359
scqx:min,max = 5110.66846 39728.27734
scqy:min,max = 5110.67236 39728.24609
scux:min,max = 5120.85791 39728.21875
scuy:min,max = 5120.86133 39728.18750
scvx:min,max = 5110.66846 39728.27734
scvy:min,max = 5110.67236 39728.24609
cori:min,max = -0.0001427841 0.0001454440
pasp:min,max = 0.99998 1.00004
.a文件里就是详细的经纬度的分布，可以用matlab的程序读出。
我在构造网格后，检查了.a文件的最大最小值，和.b文件匹配，
但是在最后一步提交作业时，出现了一个错误，如下
reading grid file from regional.grid.
-1 'mapflg'
plon:min,max = -179.99986 179.99980 %%（标注） q1
error - .a and .b files not consistent:
.a,.b min =-3.400760E+38 -1.799999E+02 -3.400760E+38 %%（标注） q2
.a,.b max = 3.401695E+381.799998E+023.401695E+38
就是说.a .b的最大最小值不匹配。
因此我查了这一部分的相关F文件
if (k.eq.1) then
call zaiord(plon, ip,.false., hmina,hmaxa, 9) %%这里 对应 q1 即 hmina=-179.999986hmaxa=179.99980(hmina hmaxa即 .a文件中plon的最大最小值)
......
if (abs(hmina-hminb).gt.abs(hminb)*1.e-4 .or.
& abs(hmaxa-hmaxb).gt.abs(hmaxb)*1.e-4 ) then
if (mnproc.eq.1) then
write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
& 'error - .a and .b files not consistent:',
& '.a,.b min = ',hmina,hminb,hmina-hminb, %% 这里对应q2 即hmina= -3.400760E+38， 可是之前读出来是和q1行对应的-179.999986啊，这里真的搞不明白了
& '.a,.b max = ',hmaxa,hmaxb,hmaxa-hmaxb
endif
call xcstop('(geopar)')
stop '(geopar)'
endif
红字部分就是我的疑问了，我对模式还不是很熟悉，想不太明白为什么会出现这种问题，真心请教各位，谢谢！
第一，我不知道什么 HYCOM，也不知道什么模式。
第二，世界上就2种文件：ASCII文件（也叫文本文件），二进制文件。
没有“十进制”文件的说法。从你的内容来看，b 文件也是ASCII文件。
第三，你的问题比较诡异，需要检查
if (k.eq.1) then
call zaiord(plon, ip,.false., hmina,hmaxa, 9)
这里至
'.a,.b min = ',hmina,hminb,hmina-hminb,
这里之间的代码过程，例如函数传递是否正确、中间是否改变了 hmina 和 hminb。
如果你检查起来有困难，请给出全部代码（如有特殊输入文件，请一并给出，如果HYCOM需要特殊的环境，恐怕就很难办了。这个环境允许单步Debug吗？） :-L对很多术语都不是很清楚，.b文件反正就是戳开就能看了的。
大概了解您说的单步debug说的是啥意思了，一定要这么做一定是有办法做到，就是比较麻烦，因为里面还调用了很多别的F文件里的函数
这一段的程序如下：
implicit none
c
include 'common_blocks.h'
c
real dp0kf,dpm,dpms,ds0kf,dsm,dsms
real hmina,hminb,hmaxa,hmaxb
integer i,ios,j,k,ktr,l
character preambl(5)*79,cline*80
#if defined(USE_CCSM3)
real plinei(itdm),plinej(jtdm)
save plinei,plinej
#endif
c
real aspmax
parameter (aspmax=2.0)! maximum grid aspect ratio for diffusion
* parameter (aspmax=1.0)! ignoregrid aspect ratio indiffusion
c
c --- read grid location,spacing,coriolis arrays
c
if (mnproc.eq.1) then! .b file from 1st tile only
write (lp,'(3a)') ' reading grid file from ',
& trim(flnmgrd),'.'
open (unit=uoff+9,file=trim(flnmgrd)//'.b',
& status='old')
endif
call xcsync(flush_lp)
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
read(cline,*) i
c
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
read (cline,*) j
c
if (i.ne.itdm .or. j.ne.jtdm) then
if (mnproc.eq.1) then
write(lp,'(/ a /)')
& 'error - wrong array size in grid file'
endif
call xcstop('(geopar)')
stop '(geopar)'
endif
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
if (mnproc.eq.1) then
write (lp,'(a)') trim(cline)
endif
read (cline,*) mapflg
c
call zaiopf(trim(flnmgrd)//'.a','old', 9)
c
do k= 1,15
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
i = index(cline,'=')
read (cline(i+1:),*) hminb,hmaxb
if (mnproc.eq.1) then
write (lp,'(a)') trim(cline)
endif
call xcsync(flush_lp)
c
if (k.eq.1) then
call zaiord(plon, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.2) then
call zaiord(plat, ip,.false., hmina,hmaxa, 9)
do i= 1,2!skip qlon,qlat
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',
& uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
call zaiosk(9)
enddo
elseif (k.eq.3) then
call zaiord(ulon, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.4) then
call zaiord(ulat, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.5) then
call zaiord(vlon, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.6) then
call zaiord(vlat, ip,.false., hmina,hmaxa, 9)
call zagetc(cline,ios, uoff+9)
if (ios.ne.0) then
if (mnproc.eq.1) then
write(lp,'(/ a,i4,i9 /)')
& 'geopar: I/O error from zagetc, iunit,ios = ',uoff+9,ios
endif !1st tile
call xcstop('(geopar)')
stop '(geopar)'
endif
#if defined(USE_CCSM3)
c pang in ANGLET
i = index(cline,'=')
read (cline(i+1:),*) hminb,hmaxb
if (mnproc.eq.1) then
write (lp,'(a)') trim(cline)
endif
call xcsync(flush_lp)
call zaiord(ANGLET, ip,.false., hmina,hmaxa, 9)
#else
c skip pang
call zaiosk(9)
#endif
elseif (k.eq.7) then
call zaiord(scpx, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.8) then
call zaiord(scpy, ip,.false., hmina,hmaxa, 9)
elseif (k.eq.9) then
call zaiord(scqx, iq,.false., hmina,hmaxa, 9)
elseif (k.eq.10) then
call zaiord(scqy, iq,.false., hmina,hmaxa, 9)
elseif (k.eq.11) then
call zaiord(scux, iu,.false., hmina,hmaxa, 9)
elseif (k.eq.12) then
call zaiord(scuy, iu,.false., hmina,hmaxa, 9)
elseif (k.eq.13) then
call zaiord(scvx, iv,.false., hmina,hmaxa, 9)
elseif (k.eq.14) then
call zaiord(scvy, iv,.false., hmina,hmaxa, 9)
else
call zaiord(corio,iq,.false., hmina,hmaxa, 9)
endif
c
if (abs(hmina-hminb).gt.abs(hminb)*1.e-4 .or.
& abs(hmaxa-hmaxb).gt.abs(hmaxb)*1.e-4 ) then
if (mnproc.eq.1) then
write(lp,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
& 'error - .a and .b files not consistent:',
& '.a,.b min = ',hmina,hminb,hmina-hminb,
& '.a,.b max = ',hmaxa,hmaxb,hmaxa-hmaxb
endif
call xcstop('(geopar)')
stop '(geopar)'
endif
enddo
c
call zaiocl(9)
if (mnproc.eq.1) then! .b file from 1st tile only
close(unit=uoff+9)
endif
c
if (itest.gt.0 .and. jtest.gt.0) then
i=itest
j=jtest
write (lp,'(/ a,2i5,a,f8.3,a,f12.9,2f10.2/)')
& ' i,j=',i+i0,j+j0,
& ' plat=',plat(i,j),
& ' corio,scux,vy=',corio(i,j),scux(i,j),scvy(i,j)
endif
call xcsync(flush_lp)
#if defined(USE_CCSM3)
c --- printout similar to ccsm ice model
call xclget(plinei,itdm, plon, 1,1, +1, 0, 1)
call xclget(plinej,jtdm, plat, 1,1,0,+1, 1)
if (mnproc.eq.1) then
write (lp,*)
write (lp,'(a,4f9.3,a,4f9.3)')
& '(domain) plon(:,1): ',
& plinei(1:4),' ...', plinei(itdm-2:itdm)
write (lp,'(a,4f9.3,a,4f9.3)')
& '(domain) plat(1,:): ',
& plinej(1:4),' ...', plinej(jtdm-2:jtdm)
write (lp,*)
endif
call xcsync(flush_lp)
#endif
vvt 发表于 2014-12-2 08:57
第一，我不知道什么 HYCOM，也不知道什么模式。
第二，世界上就2种文件：ASCII文件（也叫文本文件），二进 ...
全部代码贴上了，但是没选在回复里，请您看看，谢谢:-loveliness: 凡事打开能看懂的，就是 ASCII 文件，也叫文本文件。二进制文件一般是看不懂的。（不要说十进制文件，没有这个说法）
对于你的问题，单步调试是最简单、有效的办法。
求助往往意义不大。因为我没有你的包含文件 include 'common_blocks.h' 等等，没有你的输入数据，甚至我不确定 HYCOM 是否需要特定的环境。
从代码来看，两端之间还有很多次的 call zaiord，这个函数也可能会改变 hmina 和 hminb 的值。
单步debug并不麻烦，点几下鼠标的事儿。如果你确实有困难，可以尝试一步一步把 hmina 和 hminb 都 write 出来，看是哪一步改变了，然后再研究这一步做了什么，该不该做 vvt 发表于 2014-12-2 10:51
凡事打开能看懂的，就是 ASCII 文件，也叫文本文件。二进制文件一般是看不懂的。（不要说十进制文件，没有 ...
谢谢，认识到这种问题求助也没啥用了，也谢谢您对于文件说法的纠正。
页:
[1]