登录后查看更多精彩内容~
您需要 登录 才可以下载或查看，没有帐号？立即注册
x
为了对linda耐心学习表示鼓励，今晚写了一个fortran处理160站资料并进行格式转换的比较详细的参考，希望对需要学习fortran读写站点文件的朋友们有用，由于我的能力有限，源码中如果出现了让你看不顺的语句请多多提出批评指正！
由于时间限制，今晚先贴出转为Micaps 3格式的程序，并对程序做较详细的说明
本程序使用的fortran平台：fortran power station 4.0
fortran语法版本：f90
所用资料：气候中心160站数据（免费注册后即可下载），station.txt
注意：资料下载后请按此命名：
R1601,R1602...R1612,T1601,T1602...T1612
无需后缀，如果加了后缀，请在程序中打开文件处增加.txt后缀
program main
integer,parameter::sta=160,year=61
real lon(sta),lat(sta)
integer temp(year,sta)
integer stid(sta),h(sta),lab
integer imon,iyea
character yea*4,mon*2,kind*6,bt*30,filename*10,next*1
!************************变量说明**********************************
! sta：站点数（固定值），year：年份，固定参数，每年修改一次即可
! lon/lat：经纬度数组
! temp：读取出来的某月160数据，本例为提供批处理方式，如有需要，请自行改写
! stid：站号，h：参见micaps手册
! lab：标记要转换的是降水还是温度，kind：温度降水字符串
! imon,iyea：整型的月份和年份
! mon,yea：字符型的月份和年份（用于生成文件名）
! filename：打开的文件名称（由输入的年份和月份而定）
! bt：生成的文件名称
! next：记录是否继续转换
!************************************************************************************
!选择转换类型,今天只写了Micaps 3类格式的转换，因此略去选择转换类型的提示
!******************
!******************
!读取站点经纬数据，只需执行一次
open(1,file='station.txt',status='old') !station.txt 存放和程序放在同一目录下
do i=1,sta
read(1,*)stid(i),lon(i),lat(i),h(i) !循环读取
enddo
close(1)
!数组赋初值（年份不够的赋值为缺测，可以测试一下读取今年还没有资料的月份）
do i=1,year
do j=1,sta
temp(i,j)=9999
enddo
enddo
!if 转为micaps 3类格式 then
!*****************micaps数据格式*************
10 continue !用于循环操作
print*,"请输入年、月，格式为：2010 08 ："
print*,""
read*,iyea,imon
if(imon<1.or.imon>12.or.iyea<1951.or.iyea>1950+year)then !简单的判断输入是否合法
print*,"输入有误，请重新输入："
goto 10
endif
!调用子程序对整型的年份和月份数据进行转换
call i2c(imon,mon,iyea,yea)
!生成文件名的后半部分，如输入了8，则生成"1608"
filename="16"//mon
11 continue
print*,"降水请选1，温度请选2：" !选择转换的要素，并生成完整的要打开的文件名（如果你的数据后面有.txt，
print*,"" !请自行修改）
read*,lab
if(lab==1)then !根据输入的标识判断要转换的要素
kind="降水量"
filename="R"//trim(filename)
elseif(lab==2)then
kind="温度"
filename="T"//trim(filename)
else
print*,"输入有误，请重新输入："
goto 11
endif
!读取数据文件
call read160(filename,temp)
!此处若采用子函数的形式，会使读取站点经纬度不方便（要么重复劳动，要么传入参数太多，定义全局变量也不方便）
!生成输出文件的标题，可自行修改
bt=yea//"年"//mon//"月"//trim(kind)//"值" !//表示字符串连接，trim()函数去除首尾空格
open(1,file=''//trim(bt)//'.txt',status='replace')
write(1,20)"diamond","3",bt !写入Micaps第三类数据格式标题
20 format(1x,A,1x,A,6x,A)
write(1,*)"0 0 0 0 0 0 1 0 0 1 1000" !这个我随便写的，实际使用可参考Micaps手册自行修改
do i=1,sta
write(1,30)stid(i),lon(i),lat(i),h(i),temp(iyea-1950,i) !输出数据
enddo
30 format(1x,i5,2(3x,f8.2),3x,i3,3x,i5)
close(1)
print*,"转换成功！继续转换请输入 y"
read*,next !判断是否要继续执行
if(next.eq."y")then
goto 10
endif
end !程序结束
!读取160站的某个月份的资料到数组 temp，只要传入文件名即可
subroutine read160(filename,temp)
integer,parameter::sta=160,year=61 !此处定义需和主程序一致
integer temp(year,sta)
character filename*10
open(1,file=''//trim(filename)//'',status='old',err=101) !如果读取出错则停止执行
!也许有人会问为什么不需要知道要读取的行数，请参考end标签
read(1,*,end=100)((temp(i,j),j=1,sta),i=1,year) !当读取到文件末尾时跳到指定行执行
goto 100
101 print*,"文件",trim(filename),"读取出错！程序终止！"
close(1)
stop !意外终止，停止执行
close(1)
100 continue !fortran按地址传递，因此，只要传入后在子程序中改变temp的数据即可
end
!整型转字符型函数
subroutine i2c(imon,mon,iyea,yea)
integer imon,iyea
character mon*2,yea*4 !此处定义最好和主程序中完全一致，并且按实际需要的长短定义
write(yea,'(i4)')iyea !四位数整型年份转字符串
if(imon<=9)then !如果月份小于9，则写入一位整型到字符串
write(mon,'(i1)')imon
mon="0"//trim(mon) !在字符串前接0，如整型 8，处理后为"08"
else
write(mon,'(i2)')imon
endif
end
读取160站资料是非常基本的东西，如果只是纯读取，也就几句话而已，我之所以写这么多，是想新手了解一下一个完整的处理过程，并且判断是否该使用子程序来完成功能，如果想要程序执行的时候方便，写的时候就要费点心，如果怕写稍微复杂的代码，那么执行的时候就要手动改很多东西。
程序及资料下载：
station.txt
(3.84 KB, 下载次数: 541)
|
[复制链接]
窥视卡
雷达卡
发表于 2011-6-27 21:59:08
提升卡
置顶卡
沉默卡
喧嚣卡
变色卡
抢沙发
千斤顶
显身卡
发表于 2011-6-27 23:50:20
楼主