SABER
model_qg.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_qg_coord
3 !> Get QG coordinates
4 !----------------------------------------------------------------------
5 subroutine model_qg_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,il0
16 integer :: ncid,nx_id,ny_id,nz_id,lon_id,lat_id,z_id,area_id,lmask_id
17 real(kind_real),allocatable :: lon(:,:),lat(:,:),z(:),area(:,:),lmask(:,:,:)
18 character(len=1024),parameter :: subr = 'model_qg_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,'nx',nx_id))
24 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'ny',ny_id))
25 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nx_id,len=model%nlon))
26 model%nlon = model%nlon+1
27 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,ny_id,len=model%nlat))
28 model%nmg = model%nlon*model%nlat
29 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'nz',nz_id))
30 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nz_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(z(model%nlev))
37 allocate(area(model%nlon,model%nlat))
38 allocate(lmask(model%nlon,model%nlat,model%nlev))
39 
40 ! Read data and close file
41 call mpl%ncerr(subr,nf90_inq_varid(ncid,'lon',lon_id))
42 call mpl%ncerr(subr,nf90_inq_varid(ncid,'lat',lat_id))
43 call mpl%ncerr(subr,nf90_inq_varid(ncid,'z',z_id))
44 call mpl%ncerr(subr,nf90_inq_varid(ncid,'area',area_id))
45 call mpl%ncerr(subr,nf90_inq_varid(ncid,'lmask',lmask_id))
46 call mpl%ncerr(subr,nf90_get_var(ncid,lon_id,lon(1:model%nlon-1,:)))
47 call mpl%ncerr(subr,nf90_get_var(ncid,lat_id,lat(1:model%nlon-1,:)))
48 call mpl%ncerr(subr,nf90_get_var(ncid,z_id,z))
49 call mpl%ncerr(subr,nf90_get_var(ncid,area_id,area(1:model%nlon-1,:)))
50 call mpl%ncerr(subr,nf90_get_var(ncid,lmask_id,lmask(1:model%nlon-1,:,:)))
51 call mpl%ncerr(subr,nf90_close(ncid))
52 
53 ! Add redundant longitude for tests
54 lon(model%nlon,:) = lon(1,:)
55 lat(model%nlon,:) = lat(1,:)
56 area(model%nlon,:) = area(1,:)
57 lmask(model%nlon,:,:) = lmask(1,:,:)
58 
59 ! Convert to radian
60 lon = lon*deg2rad
61 lat = lat*deg2rad
62 
63 ! Model grid
64 img = 0
65 do ilon=1,model%nlon
66  do ilat=1,model%nlat
67  img = img+1
68  model%mg_to_tile(img) = 1
69  model%mg_to_lon(img) = ilon
70  model%mg_to_lat(img) = ilat
71  model%lon(img) = lon(ilon,ilat)
72  model%lat(img) = lat(ilon,ilat)
73  model%area(img) = area(ilon,ilat)/req**2
74  do il0=1,model%nl0
75  model%mask(img,il0) = (lmask(ilon,ilat,nam%levs(il0))>0.5)
76  end do
77  end do
78 end do
79 
80 ! Vertical unit
81 do il0=1,model%nl0
82  model%vunit(1:model%nmg,il0) = z(nam%levs(il0))
83 end do
84 
85 ! Release memory
86 deallocate(lon)
87 deallocate(lat)
88 deallocate(z)
89 deallocate(area)
90 deallocate(lmask)
91 
92 end subroutine model_qg_coord
93 
94 !----------------------------------------------------------------------
95 ! Subroutine: model_qg_read
96 !> Read QG field
97 !----------------------------------------------------------------------
98 subroutine model_qg_read(model,mpl,nam,filename,fld)
99 
100 implicit none
101 
102 ! Passed variables
103 class(model_type),intent(inout) :: model !< Model
104 type(mpl_type),intent(inout) :: mpl !< MPI data
105 type(nam_type),intent(in) :: nam !< Namelist
106 character(len=*),intent(in) :: filename !< File name
107 real(kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
108 
109 ! Local variables
110 integer :: iv,il0,img,ilon,ilat,ndims
111 integer :: ncid,fld_id
112 real(kind_real) :: fld_mg(model%nmg,model%nl0)
113 real(kind_real),allocatable :: fld_tmp(:,:,:)
114 character(len=1024),parameter :: subr = 'model_qg_read'
115 
116 if (mpl%main) then
117  ! Allocation
118  allocate(fld_tmp(model%nlon,model%nlat,model%nl0))
119 
120  ! Open file
121  call mpl%ncerr(subr,nf90_open(trim(nam%datadir)//'/'//trim(filename)//'.nc',nf90_nowrite,ncid))
122 end if
123 
124 do iv=1,nam%nv
125  if (mpl%main) then
126  ! Get variable id
127  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
128 
129  ! Check field size
130  call mpl%ncerr(subr,nf90_inquire_variable(ncid,fld_id,ndims=ndims))
131 
132  ! Read data
133  select case (ndims)
134  case (2)
135  ! 2D data
136  fld_tmp = mpl%msv%valr
137  if (trim(nam%lev2d)=='first') then
138  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(1:model%nlon-1,:,1),(/1,1/), &
139  & (/model%nlon-1,model%nlat/)))
140  fld_tmp(model%nlon,:,1) = fld_tmp(1,:,1)
141  elseif (trim(nam%lev2d)=='last') then
142  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(1:model%nlon-1,:,model%nl0),(/1,1/), &
143  & (/model%nlon-1,model%nlat/)))
144  fld_tmp(model%nlon,:,model%nl0) = fld_tmp(1,:,model%nl0)
145  end if
146  case (3)
147  ! 3D data
148  do il0=1,nam%nl
149  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(1:model%nlon-1,:,il0),(/1,1,nam%levs(il0)/), &
150  & (/model%nlon-1,model%nlat,1/)))
151  fld_tmp(model%nlon,:,il0) = fld_tmp(1,:,il0)
152  end do
153  case default
154  call mpl%abort(subr,'wrong number of dimensions for variable '//trim(nam%variables(iv)))
155  end select
156 
157  ! Pack
158  do il0=1,model%nl0
159  do img=1,model%nmg
160  ilon = model%mg_to_lon(img)
161  ilat = model%mg_to_lat(img)
162  fld_mg(img,il0) = fld_tmp(ilon,ilat,il0)
163  end do
164  end do
165  end if
166  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
167 end do
168 
169 if (mpl%main) then
170  ! Close file
171  call mpl%ncerr(subr,nf90_close(ncid))
172 
173  ! Release memory
174  deallocate(fld_tmp)
175 end if
176 
177 end subroutine model_qg_read
main
int main(int argc, char **argv)
Definition: bump_main.cc:17
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