SABER
model_ifs.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_ifs_coord
3 !> Get IFS coordinates
4 !----------------------------------------------------------------------
5 subroutine model_ifs_coord(model,mpl,nam)
6 
7 implicit none
8 
9 ! Passed variables
10 class(model_type),intent(inout) :: model !< Model
11 type(mpl_type),intent(inout) :: mpl !< MPI data
12 type(nam_type),intent(in) :: nam !< Namelist
13 
14 ! Local variables
15 integer :: img,ilon,ilat
16 integer :: ncid,nlon_id,nlat_id,nlev_id,lon_id,lat_id,pres_id
17 real(kind_real),allocatable :: lon(:),lat(:),pres(:)
18 character(len=1024),parameter :: subr = 'model_ifs_coord'
19 
20 ! Open file and get dimensions
21 model%ntile = 1
22 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/grid.nc',nf90_share,ncid))
23 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'longitude',nlon_id))
24 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'latitude',nlat_id))
25 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlon_id,len=model%nlon))
26 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlat_id,len=model%nlat))
27 model%nmg = model%nlon*model%nlat
28 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'level',nlev_id))
29 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlev_id,len=model%nlev))
30 
31 ! Allocation
32 call model%alloc
33 allocate(lon(model%nlon))
34 allocate(lat(model%nlat))
35 allocate(pres(model%nlev))
36 
37 ! Read data and close file
38 call mpl%ncerr(subr,nf90_inq_varid(ncid,'longitude',lon_id))
39 call mpl%ncerr(subr,nf90_inq_varid(ncid,'latitude',lat_id))
40 call mpl%ncerr(subr,nf90_inq_varid(ncid,'pf',pres_id))
41 call mpl%ncerr(subr,nf90_get_var(ncid,lon_id,lon))
42 call mpl%ncerr(subr,nf90_get_var(ncid,lat_id,lat))
43 call mpl%ncerr(subr,nf90_get_var(ncid,pres_id,pres))
44 call mpl%ncerr(subr,nf90_close(ncid))
45 
46 ! Convert to radian
47 lon = lon*deg2rad
48 lat = lat*deg2rad
49 
50 ! Model grid
51 img = 0
52 do ilon=1,model%nlon
53  do ilat=1,model%nlat
54  img = img+1
55  model%mg_to_tile(img) = 1
56  model%mg_to_lon(img) = ilon
57  model%mg_to_lat(img) = ilat
58  model%lon(img) = lon(ilon)
59  model%lat(img) = lat(ilat)
60  end do
61 end do
62 model%area = 4.0*pi/real(model%nmg,kind_real)
63 model%mask = .true.
64 
65 ! Vertical unit
66 do img=1,model%nmg
67  if (nam%logpres) then
68  model%vunit(img,1:nam%nl) = log(pres(nam%levs(1:nam%nl)))
69  if (model%nl0>nam%nl) model%vunit(img,model%nl0) = log(ps)
70  else
71  model%vunit(img,:) = real(nam%levs(1:model%nl0),kind_real)
72  end if
73 end do
74 ! Release memory
75 deallocate(lon)
76 deallocate(lat)
77 
78 end subroutine model_ifs_coord
79 
80 !----------------------------------------------------------------------
81 ! Subroutine: model_ifs_read
82 !> Read IFS field
83 !----------------------------------------------------------------------
84 subroutine model_ifs_read(model,mpl,nam,filename,fld)
85 
86 implicit none
87 
88 ! Passed variables
89 class(model_type),intent(inout) :: model !< Model
90 type(mpl_type),intent(inout) :: mpl !< MPI data
91 type(nam_type),intent(in) :: nam !< Namelist
92 character(len=*),intent(in) :: filename !< File name
93 real(kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
94 
95 ! Local variables
96 integer :: iv,il0,img,ilon,ilat
97 integer :: ncid,fld_id
98 real(kind_real) :: fld_mg(model%nmg,model%nl0)
99 real(kind_real),allocatable :: fld_tmp(:,:,:)
100 character(len=1024),parameter :: subr = 'model_ifs_read'
101 
102 if (mpl%main) then
103  ! Allocation
104  allocate(fld_tmp(model%nlon,model%nlat,model%nl0))
105 
106  ! Open file
107  call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/'//trim(filename)//'.nc',nf90_nowrite,ncid))
108 end if
109 
110 do iv=1,nam%nv
111  if (mpl%main) then
112  ! Get variable id
113  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
114 
115  ! Read data
116  do il0=1,nam%nl
117  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,il0),(/1,1,nam%levs(il0),1/),(/model%nlon,model%nlat,1,1/)))
118  end do
119 
120  ! Pack
121  do il0=1,model%nl0
122  do img=1,model%nmg
123  ilon = model%mg_to_lon(img)
124  ilat = model%mg_to_lat(img)
125  fld_mg(img,il0) = fld_tmp(ilon,ilat,il0)
126  end do
127  end do
128  end if
129  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
130 end do
131 
132 if (mpl%main) then
133  ! Close file
134  call mpl%ncerr(subr,nf90_close(ncid))
135 
136  ! Release memory
137  deallocate(fld_tmp)
138 end if
139 
140 end subroutine model_ifs_read
main
int main(int argc, char **argv)
Definition: bump_main.cc:17
tools_const::ps
real(kind_real), parameter, public ps
Definition: tools_const.F90:19
tools_const::deg2rad
real(kind_real), parameter, public deg2rad
Definition: tools_const.F90:15
tools_const::pi
real(kind_real), parameter, public pi
Definition: tools_const.F90:14
tools_kinds::kind_real
integer, parameter, public kind_real
Definition: tools_kinds.F90:18