1 !----------------------------------------------------------------------
2 ! Subroutine: model_wrf_coord
4 !----------------------------------------------------------------------
5 subroutine model_wrf_coord(model,mpl,nam)
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
15 integer :: img,ilon,ilat
16 integer :: ncid,nlon_id,nlat_id,nlev_id,lon_id,lat_id,pres_id
18 real(
kind_real),allocatable :: lon(:,:),lat(:,:),pres(:)
19 character(len=1024),parameter :: subr =
'model_wrf_coord'
21 ! Open file and get dimensions
23 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)
24 call mpl%ncerr(subr,nf90_inq_dimid(ncid,
'west_east',nlon_id))
25 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlon_id,len=model%nlon))
26 call mpl%ncerr(subr,nf90_inq_dimid(ncid,
'south_north',nlat_id))
27 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlat_id,len=model%nlat))
28 model%nmg = model%nlon*model%nlat
29 call mpl%ncerr(subr,nf90_inq_dimid(ncid,
'bottom_top',nlev_id))
30 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlev_id,len=model%nlev))
34 allocate(lon(model%nlon,model%nlat))
35 allocate(lat(model%nlon,model%nlat))
36 allocate(pres(model%nlev))
38 ! Read data and close file
39 call mpl%ncerr(subr,nf90_inq_varid(ncid,
'XLONG',lon_id))
40 call mpl%ncerr(subr,nf90_inq_varid(ncid,
'XLAT',lat_id))
41 call mpl%ncerr(subr,nf90_inq_varid(ncid,
'PB',pres_id))
42 call mpl%ncerr(subr,nf90_get_var(ncid,lon_id,lon,(/1,1,1/),(/model%nlon,model%nlat,1/)))
43 call mpl%ncerr(subr,nf90_get_var(ncid,lat_id,lat,(/1,1,1/),(/model%nlon,model%nlat,1/)))
44 call mpl%ncerr(subr,nf90_get_var(ncid,pres_id,pres))
45 call mpl%ncerr(subr,nf90_get_att(ncid,nf90_global,
'DX',dx))
46 call mpl%ncerr(subr,nf90_get_att(ncid,nf90_global,
'DY',dy))
47 call mpl%ncerr(subr,nf90_close(ncid))
58 model%mg_to_tile(img) = 1
59 model%mg_to_lon(img) = ilon
60 model%mg_to_lat(img) = ilat
61 model%lon(img) = lon(ilon,ilat)
62 model%lat(img) = lat(ilon,ilat)
65 model%area = dx*dy/
req**2
71 model%vunit(img,1:nam%nl) = log(pres(nam%levs(1:nam%nl)))
72 if (model%nl0>nam%nl) model%vunit(img,model%nl0) = log(
ps)
74 model%vunit(img,:) = real(nam%levs(1:model%nl0),
kind_real)
83 end subroutine model_wrf_coord
85 !----------------------------------------------------------------------
86 ! Subroutine: model_wrf_read
88 !----------------------------------------------------------------------
89 subroutine model_wrf_read(model,mpl,nam,filename,fld)
94 class(model_type),intent(inout) :: model !< Model
95 type(mpl_type),intent(inout) :: mpl !< MPI data
96 type(nam_type),intent(in) :: nam !< Namelist
97 character(len=*),intent(in) :: filename !< File name
98 real(
kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
101 integer :: iv,il0,img,ilon,ilat,ndims
102 integer :: ncid,fld_id
103 real(
kind_real) :: fld_tmp2,fld_mg(model%nmg,model%nl0)
104 real(
kind_real),allocatable :: fld_tmp(:,:,:)
105 character(len=1024),parameter :: subr =
'model_wrf_read'
109 allocate(fld_tmp(model%nlon,model%nlat,model%nl0))
112 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)
118 call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
121 call mpl%ncerr(subr,nf90_inquire_variable(ncid,fld_id,ndims=ndims))
127 fld_tmp = mpl%msv%valr
128 if (trim(nam%lev2d)==
'first') then
129 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,1),(/1,1,1/),(/model%nlon,model%nlat,1/)))
130 elseif (trim(nam%lev2d)==
'last') then
131 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,model%nl0),(/1,1,1/),(/model%nlon,model%nlat,1/)))
136 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/)))
137 select case (trim(nam%variables(iv)))
141 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp2,(/ilon+1,ilat,nam%levs(il0),1/)))
142 fld_tmp(ilon,ilat,il0) = 0.5*(fld_tmp(ilon,ilat,il0)+fld_tmp2)
148 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp2,(/ilon,ilat+1,nam%levs(il0),1/)))
149 fld_tmp(ilon,ilat,il0) = 0.5*(fld_tmp(ilon,ilat,il0)+fld_tmp2)
155 call mpl%abort(subr,
'wrong number of dimensions for variable '
161 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/)))
162 select case (trim(nam%variables(iv)))
166 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp2,(/ilon+1,ilat,nam%levs(il0),1/)))
167 fld_tmp(ilon,ilat,il0) = 0.5*(fld_tmp(ilon,ilat,il0)+fld_tmp2)
173 call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp2,(/ilon,ilat+1,nam%levs(il0),1/)))
174 fld_tmp(ilon,ilat,il0) = 0.5*(fld_tmp(ilon,ilat,il0)+fld_tmp2)
183 ilon = model%mg_to_lon(img)
184 ilat = model%mg_to_lat(img)
185 fld_mg(img,il0) = fld_tmp(ilon,ilat,il0)
189 call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
194 call mpl%ncerr(subr,nf90_close(ncid))
200 end subroutine model_wrf_read