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