SABER
model_arp.inc
Go to the documentation of this file.
1 !----------------------------------------------------------------------
2 ! Subroutine: model_arp_coord
3 !> Get ARPEGE coordinates
4 !----------------------------------------------------------------------
5 subroutine model_arp_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,a_id,b_id
17 real(kind_real),allocatable :: lon(:,:),lat(:,:),a(:),b(:)
18 character(len=1024),parameter :: subr = 'model_arp_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,'longitude',nlon_id))
24 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'latitude',nlat_id))
25 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlon_id,len=model%nlon))
26 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlat_id,len=model%nlat))
27 call mpl%ncerr(subr,nf90_inq_dimid(ncid,'Z',nlev_id))
28 call mpl%ncerr(subr,nf90_inquire_dimension(ncid,nlev_id,len=model%nlev))
29 
30 ! Allocation
31 call model%alloc
32 allocate(lon(model%nlon,model%nlat))
33 allocate(lat(model%nlon,model%nlat))
34 allocate(a(model%nlev+1))
35 allocate(b(model%nlev+1))
36 
37 ! Read data and close file
38 call mpl%ncerr(subr,nf90_inq_varid(ncid,'longitude',lon_id))
39 call mpl%ncerr(subr,nf90_inq_varid(ncid,'latitude',lat_id))
40 call mpl%ncerr(subr,nf90_inq_varid(ncid,'hybrid_coef_A',a_id))
41 call mpl%ncerr(subr,nf90_inq_varid(ncid,'hybrid_coef_B',b_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,a_id,a))
45 call mpl%ncerr(subr,nf90_get_var(ncid,b_id,b))
46 call mpl%ncerr(subr,nf90_close(ncid))
47 
48 ! Grid size
49 model%nmg = count(lon>-1000.0)
50 
51 ! Convert to radian
52 lon = lon*deg2rad
53 lat = lat*deg2rad
54 
55 ! Model grid
56 img = 0
57 do ilon=1,model%nlon
58  do ilat=1,model%nlat
59  if (lon(ilon,ilat)>-1000.0) then
60  img = img+1
61  model%mg_to_tile(img) = 1
62  model%mg_to_lon(img) = ilon
63  model%mg_to_lat(img) = ilat
64  model%lon(img) = lon(ilon,ilat)
65  model%lat(img) = lat(ilon,ilat)
66  end if
67  end do
68 end do
69 model%area = 4.0*pi/real(model%nmg,kind_real)
70 model%mask = .true.
71 
72 ! Vertical unit
73 do img=1,model%nmg
74  if (nam%logpres) then
75  model%vunit(img,1:nam%nl) = log(0.5*(a(nam%levs(1:nam%nl))+a(nam%levs(1:nam%nl)+1)) &
76  & +0.5*(b(nam%levs(1:nam%nl))+b(nam%levs(1:nam%nl)+1))*ps)
77  if (model%nl0>nam%nl) model%vunit(img,model%nl0) = log(ps)
78  else
79  model%vunit(img,:) = real(nam%levs(1:model%nl0))
80  end if
81 end do
82 
83 ! Release memory
84 deallocate(lon)
85 deallocate(lat)
86 deallocate(a)
87 deallocate(b)
88 
89 end subroutine model_arp_coord
90 
91 !----------------------------------------------------------------------
92 ! Subroutine: model_arp_read
93 !> Read ARPEGE field
94 !----------------------------------------------------------------------
95 subroutine model_arp_read(model,mpl,nam,filename,fld)
96 
97 implicit none
98 
99 ! Passed variables
100 class(model_type),intent(inout) :: model !< Model
101 type(mpl_type),intent(inout) :: mpl !< MPI data
102 type(nam_type),intent(in) :: nam !< Namelist
103 character(len=*),intent(in) :: filename !< File name
104 real(kind_real),intent(out) :: fld(model%nmga,model%nl0,nam%nv) !< Field
105 
106 ! Local variables
107 integer :: iv,il0,img,ilon,ilat,info
108 integer :: ncid,fld_id
109 real(kind_real) :: fld_mg(model%nmg,model%nl0)
110 real(kind_real),allocatable :: fld_tmp(:,:,:)
111 character(len=3) :: ilchar
112 character(len=1024),parameter :: subr = 'model_arp_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  do il0=1,nam%nl
125  ! Get id
126  write(ilchar,'(i3.3)') nam%levs(il0)
127  info = nf90_inq_varid(ncid,'S'//ilchar//trim(nam%variables(iv)),fld_id)
128 
129  ! Read data
130  if (info==nf90_noerr) then
131  ! 3D data
132  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,il0)))
133  else
134  ! Check if the variable exists as 2D data
135  call mpl%ncerr(subr,nf90_inq_varid(ncid,nam%variables(iv),fld_id))
136 
137  ! 2D data
138  if (((trim(nam%lev2d)=='first').and.(il0==1)).or.((trim(nam%lev2d)=='last').and.(il0==model%nl0))) then
139  fld_tmp = mpl%msv%valr
140  call mpl%ncerr(subr,nf90_get_var(ncid,fld_id,fld_tmp(:,:,il0)))
141  end if
142  end if
143  end do
144 
145  ! Pack
146  do il0=1,model%nl0
147  do img=1,model%nmg
148  ilon = model%mg_to_lon(img)
149  ilat = model%mg_to_lat(img)
150  fld_mg(img,il0) = fld_tmp(ilon,ilat,il0)
151  end do
152  end do
153  end if
154  call mpl%glb_to_loc(model%nl0,model%nmga,model%nmg,model%mga_to_mg,fld_mg,fld(:,:,iv))
155 end do
156 
157 if (mpl%main) then
158  ! Close file
159  call mpl%ncerr(subr,nf90_close(ncid))
160 
161  ! Release memory
162  deallocate(fld_tmp)
163 end if
164 
165 end subroutine model_arp_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::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
type_rng::a
integer(kind=int64), parameter a
Definition: type_rng.F90:18