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