fortran协同工作

例子:TELSEM的使用

Liting Mai (麦李婷)

Fortran的协同工作办法有:

方法一:将多个.f90文件一块编译

gfortran test1.f90 test2.f90 -o exe

方法二:建立静态链接库.lib或动态链接库.dll。需要被链接的库的rwx权限

可参考: 臭石头和雪球的教程:静态库和动态链接库, fortran-hdf库, fortran-netcdf库

方法一.例子.TELSEM的使用

TELSEM由CNRS的F.Aires和C. Priegent 在《A Tool to Estimate Land-Surface Emissivities at Microwave frequencies (TELSEM) for use in numerical weather prediction》 中提出。 TELSEM基于SSM/I地表发射率数据集(C. Prigent, 1997)的19V&H, 22V, 37V&H, 85V&H等通道的Emissivity进行频率插值,得到月平均的,具有一定空间分辨率的(或固定经纬度)的全球陆面微波发射率, 其适用频段为10 – 700 GHz. (cited from胡继恒)。

paper1: Microwave land surface emissivities estimated from SSM/I observations. J. Geophys. Res. Space Phys. 1997

paper2: Land Surface Microwave Emissivities over the Globe for a Decade. Bull. Am. Meteorol. Soc. 2006: 19-85 GHz

paper3: Evaluation of modeled microwave land surface emissivities with satellite-based estimates. J. Geophys. Res. Atmos. 2015: less than 40GHz (JGR-Atmosphere). https://agupubs.onlinelibrary.wiley.com/doi/full/10.1002/2014JD021817

TELSEM在RTTOV中被用于计算地表发射率,也可以单独使用。

TELSEM官网下载包, 解压:

$ cd /home/hjh/rttov13/emis_data
$ tar -xvf telsem2_mw_atlas.tar.bz2

解压后的文件中,mod_mwatlas_m2.F90提供了使用telsem的接口,test_telsem2.F90提供了示例脚本。

胡继恒师兄的目录 /home/hjh/rttov13/emis_data中已解压完毕。我在使用时,只能修改自己的目录文件/data04/1/mlt/run/TELSEM/TELSEM166.f90

不妨直接联合编译和运行:

# 在目录/data04/1/mlt/run/TELSEM下
$h5fc TELSEM166.f90 /home/hjh/rttov13/emis_data/telsem2/mod_mwatlas_m2.F90 -o test.exe
$./test.exe

附: TELSEM166.f90源代码(根据示例test_telsum2.F90修改,模拟了89VH和166VH的地表反射率)

PROGRAM globe_mlse_telsem
	! call TELSEM2 atlas and interpolator to generate global 0.25 degree land surface emissivity on monthly basis
	! By Qingyang Liu, Aug, 2023, USTC
	! adapted from test_telsem.F90
	
		USE mod_mwatlas_m2
		use HDF5
	
		IMPLICIT NONE
		INTEGER       :: m,n
		CHARACTER(len=130) fout
		CHARACTER(len=2) mm
	
		INTEGER(jpim) :: error_status
	
		LOGICAL(jplm) :: verbose   ! For atlas reading subroutine
		INTEGER(jpim) :: verb      !=1 for TRUE and 0 for FALSE - for emissivity routines
	
		INTEGER(jpim) :: month     ! (1->12)
		CHARACTER(LEN=256) :: dir  ! directory of emis database
	
		REAL(jprb)    :: resol     ! horizontal resolution for the user
		REAL(jprb)    :: lat       ! (-90->90)
		REAL(jprb)    :: lon       ! (0->360)
		REAL(jprb)    :: theta     ! (0->60�)
	
		! For individual freq interpolations
		REAL(jprb)    :: ev, eh, stdv, stdh, covvh
		REAL(jprb)    :: freq      ! (19->85GHz)
		INTEGER(jpim) :: i
	
		! For multiple freq interpolations
		INTEGER(jpim), PARAMETER :: nchan = 2
		REAL(jprb)    :: ev2(nchan), eh2(nchan), std(2*nchan,2*nchan)
		REAL(jprb)    :: freq2(nchan)
	
		! Structure containing atlas data
		TYPE(telsem2_atlas_data) :: atlas

	!   !MLSE output
	! integer(HID_T) :: fid,dset_id,space,dcpl
	! integer :: error
	! integer(hsize_t),parameter :: dimo5(3)=(/4,100,50/)    ! ouput: Emiss& Tb
	!   real(jprb) :: MLSEout(4,100,50)!1440,720)

		!--- End of header ----------------------------------
		verbose = .TRUE.  ! Verbose output for reading subroutine
		verb = 0          ! No verbose output for emissivity subroutines
		!====================================================
		! Read the atlas
		!====================================================
		dir = '/home/hjh/rttov13/emis_data/'      ! liting 20240112
		month = 6                  !! modify or do loop to simulate other months
	
		! WRITE(0,'(a,i3)') 'Reading atlas for month ',month
		CALL rttov_readmw_atlas(TRIM(dir), month, atlas, verbose, error_status)
		IF (error_status /= 0) THEN
		WRITE(0,'(a)') 'Error reading atlas'
		STOP 1
		ENDIF
	
		!====================================================
		! Calculate emissivities
		!====================================================
		! Multiple frequencies, spatially averaged, return emissivities only
	!   write(mm,"%0.2i") month !!  to the form '06'
	!   yue = (/"01","02","03","04","05","06","07","08","09","10","11","12"/)
		fout     ='TELSEM_89GHz_166GHz_Jun.txt'
		resol    = 0.25
		theta    = 52.8  !53.1
		freq     = 10.65
	!   freq2(:) = (/10.65, 18.7, 23.8, 36.5, 89.0/) ! For multiple frequency experiment
		freq2(:) = (/89.0,166.5/) ! For multiple frequency experiment
		i        = 0                          ! Index of freq in freq2 array

		open(15,file= trim(adjustl(fout)))       !!final product
		lat      = -90.+0.125
		do m     = 1,720
		lon  = 0.125
		do n = 1,1440
			CALL emis_interp_int_mult(lat, lon, resol, theta, freq2, nchan, atlas, ev2, eh2, verb = verb)
			WRITE(15,100) lat,lon,(ev2(i),i=1,nchan),(eh2(i),i=1,nchan)
			lon = lon+0.25
		end do
		lat = lat+0.25
		! write(15,100)lat,lon,(ev2(i),i=1,nchan),(eh2(i),i=1,nchan)
		print*,m
		! MLSEout(:2,n,m)=ev2(:)
		! MLSEout(3:,n,m)=eh2(:)
		end do
		CLOSE(15)

		!====================================================
		! Close the atlas
		!====================================================
		CALL rttov_closemw_atlas(atlas)
	100 FORMAT(12(1x,f14.6))

	!   call h5fcreate_f("/data04/1/mlt/run/TELSEM/mlse.hdf",H5F_ACC_TRUNC_F,fid,error)
	!   call h5screate_simple_f(3, dimo5, space, error)
	! PRINT*,'BEGIN'
	!   call h5pcreate_f(H5P_DATASET_CREATE_F, dcpl, error)
	! PRINT*,'BEGIN2'
	! call h5dcreate_f(fid, "MLSE", H5T_IEEE_F32LE, space, dset_id, error, dcpl)
	! PRINT*,'BEGIN3'
	! call h5dwrite_f(dset_id, H5T_IEEE_F32LE, MLSEout, dimo5, error)
	! PRINT*,'BEGIN4'
	! CALL h5pclose_f(dcpl , error)
	!   CALL h5dclose_f(dset_id , error)
	!   CALL h5sclose_f(space, error)
	
	END PROGRAM globe_mlse_telsem