爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 3517|回复: 4

[求助] 用Fortran读取nc数据出现问题

[复制链接]

新浪微博达人勋

发表于 2015-8-24 20:16:50 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
本人用fortran读取并处理nc数据,遇到了问题,求助各位大神,程序是在Linux下用ifort编译连接,并生成了可执行文件main.exe,但是执行main.exe后,并没有如愿生成主程序中28行生成的输入结果的文件,而且在读取地理信息数据的子程序的第17行后,也没有输出status,但是前面打开nc文件,查询变量id号都可以输出status为零,说明打开nc数据和查询变量id号都没有问题,但是读取纬度这一变量就不行,导致主程序中没有生成任何输出结果的文件,不知道为什么,我实在不知道自己哪里写的不对,之前这么读nc数据就没有问题,这回只是把读取过程写成了一个子程序,就一直出现这种情况,哪位前辈可以指点一下,谢谢!!!
部分程序如下:
Module constant
    Implicit None
    Save
    Integer,Parameter::nx=152,ny=110,nz=1,nfile=19,nt=72*nfile,nvar=5
    Real(Kind=4),Parameter::pi=3.1415926,r=6371.004
End Module constant
1Program Grapes_Cmaq
2   Use constant
3   Implicit None
4   !Include "/public/home/dengtao/atmosphere/netcdf/include/netcdf.inc"
5  Integer::i,j,k,l,m,n
6    Integer::ii_gz,jj_gz,ii_py,jj_py
7    !Integer,Parameter::nx=152,ny=110,nz=1,nt=72*19,nvar=5
8    !Real(Kind=4),Parameter::pi=3.1415926
9    Real(Kind=4)::lat_obs,lon_obs
10    Real(Kind=4)::lat(nx,ny,nz,1),lon(nx,ny,nz,1) !存储nc数据的经纬度,实际只有两维
11    Real(Kind=4)::p(nx,ny,nz,nt),t(nx,ny,nz,nt),q(nx,ny,nz,nt),ws(nx,ny,nz,nt),wd(nx,ny,nz,nt),rh(nx,ny,nz,nt) !存储nc数据中的气压、温度、混合比、风向、风速值
12    Real(Kind=4)::p_py(1,1,1,nt),t_py(1,1,1,nt),ws_py(1,1,1,nt),wd_py(1,1,1,nt),rh_py(1,1,1,nt)存储观测值
13    Real(Kind=4)::p_gz(1,1,1,nt),t_gz(1,1,1,nt),ws_gz(1,1,1,nt),wd_gz(1,1,1,nt),rh_gz(1,1,1,nt)
14   Real(Kind=4),Dimension(nvar)::meanerror_py,rel_py,rmse_py,meanerror_gz,rel_gz,rmse_gz !存储计算得到的观测地的模拟值与观测值的平均误差、相对湿度、均方根误差(2013年12月至2014年1月共57天)
15    Character(Len=100)::filename,str !存储nc文件名
16    Character,Dimension(nfile)::date
17    Real relative_humidity,average,sd,cc,meanerror,rmse !声明函数类型 此处出现错误:如果使用了implicit none 语句,则在函数和调用程序中都要声明函数的类型
18    date=(/'131204','131207','131210','131213','131216','131219','131222','131225','131228','131231',&
19    &'140103','140106','140109','140112','140115','140118','140121','140124','140127'/)
20    !读取地理信息数据并计算离观测站最近的点
21    filename='GRIDCRO2D_2013120400'
22    !print*,filename
23    Call readdata_from_nc(filename,lat,lon,p,t,q,ws,wd)
24    !open(30,file='lat.txt')
25    !write(30,'(f7.3)') ((lat(i,j,nz,1),i=1,nx),j=1,ny)
26    lat_obs=23.217;lon_obs=113.483
27    Call distance_of_twopoints(lat,lon,lat_obs,lon_obs,ii_gz,jj_gz)
28    Open(10,file='lat_lon.txt',status='new',access='append')
29    Write(*,'(1x,i3,1x,i3,1x,f7.3,1x,f7.3)') '离观测点1最近点的经纬度为:',ii_gz,jj_gz,lat(ii_gz,jj_gz,nz,1),lon(ii_gz,jj_gz,nz,1)
        Close(10)
30    lat_obs=22.933;lon_obs=113.317
31    Call distance_of_twopoints(lat,lon,lat_obs,lon_obs,ii_py,jj_py)
32    Open(10,file='lat_lon.txt',status='old',access='append')
33    Write(*,'(1x,i3,1x,i3,1x,f7.3,1x,f7.3)') '离观测站点2最近点的经纬度为:',ii_py,jj_py,lat(ii_py,jj_py,nz,1),lon(ii_py,jj_py,nz,1)
34   Close(10)
35End
读取nc地理数据的子程序
1Subroutine readgeodata_from_nc(filename,lat,lon)
2   Use constant
3   Implicit None
4   Include "/public/home/dengtao/atmosphere/netcdf/include/netcdf.inc"
5    Character(Len=100),Intent(In)::filename
6    Real(Kind=4),Intent(Out)::lat(:,:,:,:),lon(:,:,:,:)
7    Integer,Dimension(4)::start,count
8    Integer::ncid,latid,lonid,status
9    character(len=4)::str
10    data  start /1,1,1,1/
11    data count1 /nx,ny,nz,1/
12        status=nf_open(filename,nf_nowrite,ncid)
13        !print*,status,ncid
14        status=nf_inq_varid(ncid,'LAT',latid)
15       status=nf_inq_varid(ncid,'LON',lonid)
16        !print*,status
17        status=nf_get_var_real(ncid,latid,start,count,lat)!纬度
18       print*,status
19       status=nf_get_var_real(ncid,lonid,start,count,lon)!经度
20       status=nf_close(ncid)
  21  End Subroutine readgeodata_from_nc
读取nc气象数据的子程序(为了查错,只看读取地理数据的部分,因而这一部分注释掉)
!Subroutine readdata_from_nc(filename,p,t,q,ws,wd)
    !Use constant
    !Implicit None
    !Include "/public/home/dengtao/atmosphere/netcdf/include/netcdf.inc"
    !Character(Len=100),Intent(In)::filename
    !Real(Kind=4),Intent(Out)::p(:,:,:,:),t(:,:,:,:),q(:,:,:,:),ws(:,:,:,:),wd(:,:,:,:)
    !Integer,Dimension(4)::start,count
    !Integer::ncid,pid,tid,qid,wsid,wdid,status
    !data start /1,1,1,1/
    !data count /nx,ny,nz,nt/
    !    status=nf_open(filename_d,nf_nowrite,ncid)
    !    status=nf_inq_varid(ncid,'PRSFC',pid)    !气压
    !    status=nf_inq_varid(ncid,'TEMP2',tid)    !温度
    !    status=nf_inq_varid(ncid,'Q2',qid)       !混合比
    !    status=nf_inq_varid(ncid,'WSPD10',wsid)  !风速
    !    status=nf_inq_varid(ncid,'WDIR10',wdid)  !风向
    !    status=nf_get_vara_real(ncid,pid,start,count2,p_d)
    !    status=nf_get_vara_real(ncid,tid,start,count2,t_d)
    !    status=nf_get_vara_real(ncid,qid,start,count2,q_d)
    !    status=nf_get_vara_real(ncid,wsid,start,count2,ws_d)
    !    status=nf_get_vara_real(ncid,wdid,start,count2,wd_d)
    !    status=nf_close(ncid)
    !End Subroutine readdata_from_nc
计算离观测站点最近的点的子程序

        
Subroutine distance_of_twopoints(lat_d1,lon_d1,lat_obs_d,lon_obs_d,ii_d,jj_d)
    Use constant
    Implicit None
    Integer::i,j
    Integer,Intent(Out)::ii_d,jj_d
    Real(Kind=4),Intent(In)::lat_obs_d,lon_obs_d
    Real(Kind=4),Intent(In)::lat_d1(:,:,:,:),lon_d1(:,:,:,:)
    Real(Kind=4)::distance0,c0,distance,c
    c0=sin((90-lat_d1(1,1,1,1))*pi/180)*sin((90-lat_obs_d)*pi/180)*cos((lon_d1(1,1,1,1)-lon_obs_d)*pi/180)+cos((90-lat_d1(1,1,1,1))*pi/180)*cos((90-lat_obs_d)*pi/180)
    distance0=r*acos(c0)*pi/180
    Do j=1,ny
        Do i=1,nx
        c=sin((90-lat_d1(i,j,nz,1))*pi/180)*sin((90-lat_obs_d)*pi/180)*cos((lon_d1(i,j,nz,1)-lon_obs_d)*pi/180)+cos((90-lat_d1(i,j,nz,1))*pi/180)*cos((90-lat_obs_d)*pi/180)
        distance=r*acos(c)*pi/180
        if(distance<distance0)then
            distance0=distance
            ii_d=i
            jj_d=j
        endif
        enddo
    enddo
End Subroutine distance_of_twopoints


密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2015-8-24 21:27:41 | 显示全部楼层
能不能丢个文件上来,我用别的软件帮你看看?!
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2015-8-24 22:57:14 | 显示全部楼层
又是那隻貓 发表于 2015-8-24 21:27
能不能丢个文件上来,我用别的软件帮你看看?!

谢谢前辈,这个问题已经解决
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2017-12-2 15:57:50 | 显示全部楼层
xiuyuanyang 发表于 2015-8-24 22:57
谢谢前辈,这个问题已经解决

请问读取纬度这一变量不行是怎么解决的,我只是一个最简单的读取程序也出现这个问题了,读出来的经纬度不在原文件范围内
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2017-12-6 16:48:37 | 显示全部楼层
cxyykl 发表于 2017-12-2 15:57
请问读取纬度这一变量不行是怎么解决的,我只是一个最简单的读取程序也出现这个问题了,读出来的经纬度不 ...

能说的具体点吗,比如程序什么样,读取什么数据,什么叫读出来的“经纬度不在源文件范围内”
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表