SABER
model_fv3.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_fv3_coord
3 !> Get FV3 coordinates
4 !----------------------------------------------------------------------
5 subroutine model_fv3_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,itile,il0
16 integer :: ncid,nlon_id,nlat_id,ntile_id,nlev_id,lon_id,lat_id,a_id,b_id
17 real(kind_real) :: sigmaup,sigmadn,sigma
18 real(kind_real),allocatable :: lon(:,:,:),lat(:,:,:),a(:),b(:)
19 character(len=1024),parameter :: subr = 'model_fv3_coord'
20 
21 ! Open file and get dimensions
22 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/grid.nc',nf90_share,ncid))
23 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'fxdim',nlon_id))
24 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'fydim',nlat_id))
25 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'ntile',ntile_id))
26 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlon_id,len=model%nlon))
27 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlat_id,len=model%nlat))
28 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,ntile_id,len=model%ntile))
29 model%nmg = model%nlon*model%nlat*model%ntile
30 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'xaxis_1',nlev_id))
31 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlev_id,len=model%nlev))
32 model%nlev = model%nlev-1
33 
34 ! Allocation
35 call model%alloc
36 allocate(lon(model%nlon,model%nlat,model%ntile))
37 allocate(lat(model%nlon,model%nlat,model%ntile))
38 allocate(a(model%nlev+1))
39 allocate(b(model%nlev+1))
40 
41 ! Read data and close file
42 call mpl%ncerr(subr,nf90_inq_varid(ncid,'flons',lon_id))
43 call mpl%ncerr(subr,nf90_inq_varid(ncid,'flats',lat_id))
44 call mpl%ncerr(subr,nf90_inq_varid(ncid,'ak',a_id))
45 call mpl%ncerr(subr,nf90_inq_varid(ncid,'bk',b_id))
46 call mpl%ncerr(subr,nf90_get_var(ncid,lon_id,lon))
47 call mpl%ncerr(subr,nf90_get_var(ncid,lat_id,lat))
48 call mpl%ncerr(subr,nf90_get_var(ncid,a_id,a,(/1,1/),(/model%nlev+1,1/)))
49 call mpl%ncerr(subr,nf90_get_var(ncid,b_id,b,(/1,1/),(/model%nlev+1,1/)))
50 call mpl%ncerr(subr,nf90_close(ncid))
51 
52 ! Model grid
53 img = 0
54 do itile=1,model%ntile
55  do ilon=1,model%nlon
56  do ilat=1,model%nlat
57  img = img+1
58  model%mg_to_tile(img) = itile
59  model%mg_to_lon(img) = ilon
60  model%mg_to_lat(img) = ilat
61  model%lon(img) = lon(ilon,ilat,itile)
62  model%lat(img) = lat(ilon,ilat,itile)
63  call lonlatmod(model%lon(img),model%lat(img))
64  end do
65  end do
66 end do
67 model%area = 4.0*pi/real(model%nmg,kind_real)
68 model%mask = .true.
69 
70 ! Vertical unit
71 do il0=1,model%nl0
72  if (il0<nam%nl) then
73  sigmaup = a(nam%levs(il0)+1)/ps+b(nam%levs(il0)+1)
74  sigmadn = a(nam%levs(il0))/ps+b(nam%levs(il0))
75  sigma = 0.5*(sigmaup+sigmadn)
76  else
77  sigma = 1.0
78  end if
79  if (nam%logpres) then
80  model%vunit(:,il0) = log(sigma*ps)
81  if (model%nl0>nam%nl) model%vunit(:,model%nl0) = log(ps)
82  else
83  model%vunit(:,il0) = sigma
84  end if
85 end do
86 
87 ! Release memory
88 deallocate(lon)
89 deallocate(lat)
90 deallocate(a)
91 deallocate(b)
92 
93 end subroutine model_fv3_coord
94 
95 !----------------------------------------------------------------------
96 ! Subroutine: model_fv3_read
97 !> Read FV3 field
98 !----------------------------------------------------------------------
99 subroutine model_fv3_read(model,mpl,nam,filename,fld)
100 
101 implicit none
102 
103 ! Passed variables
104 class(model_type),intent(inout) :: model !< Model
105 type(mpl_type),intent(inout) :: mpl !< MPI data
106 type(nam_type),intent(in) :: nam !< Namelist
107 character(len=*),intent(in) :: filename !< File name
108 real(kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
109 
110 ! Local variables
111 integer :: iv,il0,imgt,img,ilon,ilat
112 integer :: ncid,fld_id,dimids(4),zlen
113 real(kind_real) :: fld_mgt(model%nmgt,model%nl0)
114 real(kind_real),allocatable :: fld_tmp(:,:,:)
115 character(len=1) :: ctile
116 character(len=1024),parameter :: subr = 'model_fv3_read'
117 
118 if (model%ioproc(model%mytile)==mpl%myproc) then
119  ! Allocation
120  allocate(fld_tmp(model%nlon,model%nlat,model%nl0))
121 
122  ! Open file
123  write(ctile,'(i1.1)') model%mytile
124  call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/'//trim(filename)//'_tile'//ctile//'.nc',nf90_nowrite,ncid))
125 end if
126 
127 do iv=1,nam%nv
128  if (model%ioproc(model%mytile)==mpl%myproc) then
129  ! Get variable id
130  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
131 
132  ! Check field size
133  call mpl%ncerr(subr,nf90_inquire_variable(ncid,fld_id,dimids=dimids))
134  call mpl%ncerr(subr,nf90_inquire_dimension(ncid,dimids(3),len=zlen))
135 
136  ! Read data
137  if (zlen==1) then
138  ! 2D data
139  fld_tmp = mpl%msv%valr
140  if (trim(nam%lev2d)=='first') then
141  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,1),(/1,1,1,1/),(/model%nlon,model%nlat,1,1/)))
142  elseif (trim(nam%lev2d)=='last') then
143  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,model%nl0),(/1,1,1,1/),(/model%nlon,model%nlat,1,1/)))
144  end if
145  else
146  ! 3D data
147  do il0=1,nam%nl
148  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/)))
149  end do
150  end if
151 
152  ! Pack
153  do imgt=1,model%nmgt
154  img = model%mgt_to_mg(imgt)
155  ilon = model%mg_to_lon(img)
156  ilat = model%mg_to_lat(img)
157  fld_mgt(imgt,:) = fld_tmp(ilon,ilat,:)
158  end do
159  end if
160  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmgt,model%mga_to_mgt,fld_mgt,fld(:,:,iv), &
161  & model%ioproc(model%mytile),model%tilepool(:,model%mytile))
162 end do
163 
164 if (model%ioproc(model%mytile)==mpl%myproc) then
165  ! Close file
166  call mpl%ncerr(subr,nf90_close(ncid))
167 
168  ! Release memory
169  deallocate(fld_tmp)
170 end if
171 
172 end subroutine model_fv3_read
tools_const::ps
real(kind_real), parameter, public ps
Definition: tools_const.F90:19
tools_func::lonlatmod
subroutine, public lonlatmod(lon, lat)
Set latitude between -pi/2 and pi/2 and longitude between -pi and pi.
Definition: tools_func.F90:66
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
type_rng::a
integer(kind=int64), parameter a
Definition: type_rng.F90:18