SABER
model_wrf.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_wrf_coord
3 !> Get WRF coordinates
4 !----------------------------------------------------------------------
5 subroutine model_wrf_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) :: dx,dy
18 real(kind_real),allocatable :: lon(:,:),lat(:,:),pres(:)
19 character(len=1024),parameter :: subr = 'model_wrf_coord'
20 
21 ! Open file and get dimensions
22 model%ntile = 1
23 call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/grid.nc',nf90_share,ncid))
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))
31 
32 ! Allocation
33 call model%alloc
34 allocate(lon(model%nlon,model%nlat))
35 allocate(lat(model%nlon,model%nlat))
36 allocate(pres(model%nlev))
37 
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))
48 
49 ! Convert to radian
50 lon = lon*deg2rad
51 lat = lat*deg2rad
52 
53 ! Model grid
54 img = 0
55 do ilon=1,model%nlon
56  do ilat=1,model%nlat
57  img = img+1
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)
63  end do
64 end do
65 model%area = dx*dy/req**2
66 model%mask = .true.
67 
68 ! Vertical unit
69 do img=1,model%nmg
70  if (nam%logpres) then
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)
73  else
74  model%vunit(img,:) = real(nam%levs(1:model%nl0),kind_real)
75  end if
76 end do
77 
78 ! Release memory
79 deallocate(lon)
80 deallocate(lat)
81 deallocate(pres)
82 
83 end subroutine model_wrf_coord
84 
85 !----------------------------------------------------------------------
86 ! Subroutine: model_wrf_read
87 !> Read WRF field
88 !----------------------------------------------------------------------
89 subroutine model_wrf_read(model,mpl,nam,filename,fld)
90 
91 implicit none
92 
93 ! Passed variables
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
99 
100 ! Local variables
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'
106 
107 if (mpl%main) then
108  ! Allocation
109  allocate(fld_tmp(model%nlon,model%nlat,model%nl0))
110 
111  ! Open file
112  call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/'//trim(filename)//'.nc',nf90_nowrite,ncid))
113 end if
114 
115 do iv=1,nam%nv
116  if (mpl%main) then
117  ! Get variable id
118  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
119 
120  ! Check field size
121  call mpl%ncerr(subr,nf90_inquire_variable(ncid,fld_id,ndims=ndims))
122 
123  ! Read data
124  select case (ndims)
125  case (3)
126  ! 2D data
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/)))
132  end if
133  case (4)
134  ! 3D data
135  do il0=1,nam%nl
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)))
138  case ('U')
139  do ilat=1,model%nlat
140  do ilon=1,model%nlon
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)
143  end do
144  end do
145  case ('V')
146  do ilat=1,model%nlat
147  do ilon=1,model%nlon
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)
150  end do
151  end do
152  end select
153  end do
154  case default
155  call mpl%abort(subr,'wrong number of dimensions for variable '//trim(nam%variables(iv)))
156  end select
157 
158 
159  ! Read data
160  do il0=1,nam%nl
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)))
163  case ('U')
164  do ilat=1,model%nlat
165  do ilon=1,model%nlon
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)
168  end do
169  end do
170  case ('V')
171  do ilat=1,model%nlat
172  do ilon=1,model%nlon
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)
175  end do
176  end do
177  end select
178  end do
179 
180  ! Pack
181  do il0=1,model%nl0
182  do img=1,model%nmg
183  ilon = model%mg_to_lon(img)
184  ilat = model%mg_to_lat(img)
185  fld_mg(img,il0) = fld_tmp(ilon,ilat,il0)
186  end do
187  end do
188  end if
189  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
190 end do
191 
192 if (mpl%main) then
193  ! Close file
194  call mpl%ncerr(subr,nf90_close(ncid))
195 
196  ! Release memory
197  deallocate(fld_tmp)
198 end if
199 
200 end subroutine model_wrf_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::req
real(kind_real), parameter, public req
Definition: tools_const.F90:17
tools_const::deg2rad
real(kind_real), parameter, public deg2rad
Definition: tools_const.F90:15
tools_kinds::kind_real
integer, parameter, public kind_real
Definition: tools_kinds.F90:18