SABER
model_mpas.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_mpas_coord
3 !> Get MPAS coordinates
4 !----------------------------------------------------------------------
5 subroutine model_mpas_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
16 integer :: ncid,nmg_id,nlev_id,lon_id,lat_id,pres_id
17 real(kind=4),allocatable :: lon(:),lat(:),pres(:)
18 character(len=1024),parameter :: subr = 'model_mpas_coord'
19 
20 ! Open file and get dimensions
21 model%ntile = 1
22 model%nlon = mpl%msv%vali
23 model%nlat = mpl%msv%vali
24 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/grid.nc',nf90_share,ncid))
25 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'nCells',nmg_id))
26 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nmg_id,len=model%nmg))
27 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'nVertLevels',nlev_id))
28 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlev_id,len=model%nlev))
29 
30 ! Allocation
31 call model%alloc
32 allocate(lon(model%nmg))
33 allocate(lat(model%nmg))
34 allocate(pres(model%nlev))
35 
36 ! Read data and close file
37 call mpl%ncerr(subr,nf90_inq_varid(ncid,'lonCell',lon_id))
38 call mpl%ncerr(subr,nf90_inq_varid(ncid,'latCell',lat_id))
39 call mpl%ncerr(subr,nf90_inq_varid(ncid,'pressure_base',pres_id))
40 call mpl%ncerr(subr,nf90_get_var(ncid,lon_id,lon))
41 call mpl%ncerr(subr,nf90_get_var(ncid,lat_id,lat))
42 call mpl%ncerr(subr,nf90_get_var(ncid,pres_id,pres))
43 call mpl%ncerr(subr,nf90_close(ncid))
44 
45 ! Model grid
46 model%mg_to_tile = 1
47 model%lon = lon
48 model%lat = lat
49 model%area = 4.0*pi/real(model%nmg,kind_real)
50 model%mask = .true.
51 
52 ! Vertical unit
53 do img=1,model%nmg
54  if (nam%logpres) then
55  model%vunit(img,1:nam%nl) = log(pres(nam%levs(1:nam%nl)))
56  if (model%nl0>nam%nl) model%vunit(img,model%nl0) = log(ps)
57  else
58  model%vunit(img,:) = real(nam%levs(1:model%nl0),kind_real)
59  end if
60 end do
61 
62 ! Release memory
63 deallocate(lon)
64 deallocate(lat)
65 deallocate(pres)
66 
67 end subroutine model_mpas_coord
68 
69 !----------------------------------------------------------------------
70 ! Subroutine: model_mpas_read
71 !> Read MPAS field
72 !----------------------------------------------------------------------
73 subroutine model_mpas_read(model,mpl,nam,filename,fld)
74 
75 implicit none
76 
77 ! Passed variables
78 class(model_type),intent(inout) :: model !< Model
79 type(mpl_type),intent(inout) :: mpl !< MPI data
80 type(nam_type),intent(in) :: nam !< Namelist
81 character(len=*),intent(in) :: filename !< File name
82 real(kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
83 
84 ! Local variables
85 integer :: iv,il0,ndims
86 integer :: ncid,fld_id
87 real(kind_real) :: fld_mg(model%nmg,model%nl0)
88 character(len=1024),parameter :: subr = 'model_mpas_read'
89 
90 ! Open file
91 if (mpl%main) call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/'//trim(filename)//'.nc',nf90_nowrite,ncid))
92 
93 do iv=1,nam%nv
94  if (mpl%main) then
95  ! Get variable id
96  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
97 
98  ! Check field size
99  call mpl%ncerr(subr,nf90_inquire_variable(ncid,fld_id,ndims=ndims))
100 
101  ! Read data
102  select case (ndims)
103  case (2)
104  ! 2D data
105  fld_mg = mpl%msv%valr
106  if (trim(nam%lev2d)=='first') then
107  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_mg(:,1),(/1,1/),(/model%nmg,1/)))
108  elseif (trim(nam%lev2d)=='last') then
109  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_mg(:,model%nl0),(/1,1/),(/model%nmg,1/)))
110  end if
111  case (3)
112  ! 3D data
113  do il0=1,nam%nl
114  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_mg(:,il0),(/nam%levs(il0),1,1/),(/1,model%nmg,1/)))
115  end do
116  case default
117  call mpl%abort(subr,'wrong number of dimensions for variable '//trim(nam%variables(iv)))
118  end select
119  end if
120  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
121 end do
122 
123 ! Close file
124 if (mpl%main) call mpl%ncerr(subr,nf90_close(ncid))
125 
126 end subroutine model_mpas_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::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