source: branches/2015/dev_r5187_UKMO13_simplification/NEMOGCM/TOOLS/SCOORD_GEN/src/utils.F90 @ 5269

Last change on this file since 5269 was 5269, checked in by timgraham, 6 years ago

Removed use of 3D variables.
Compiles and runs with no errors but still need to check output

File size: 9.7 KB
Line 
1MODULE utils
2
3   USE netcdf
4
5   IMPLICIT NONE
6   PUBLIC             ! allows the acces to par_oce when dom_oce is used
7   !                  ! exception to coding rules... to be suppressed ???
8
9!   PUBLIC dom_oce_alloc
10
11   INTEGER, PARAMETER   :: dp=8 , sp=4, wp=dp
12
13   !! All coordinates
14   !! ---------------
15   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   gdep3w_0           !: depth of t-points (sum of e3w) (m)
16   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   gdept_0, gdepw_0   !: analytical (time invariant) depth at t-w  points (m)
17   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e3v_0  , e3f_0     !: analytical (time invariant) vertical scale factors at  v-f
18   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e3t_0  , e3u_0     !:                                      t-u  points (m)
19   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e3vw_0             !: analytical (time invariant) vertical scale factors at  vw
20   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e3w_0  , e3uw_0    !:                                      w-uw points (m)
21
22   !! s-coordinate and hybrid z-s-coordinate
23   !! =----------------======---------------
24   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)   ::   gsigt, gsigw       !: model level depth coefficient at t-, w-levels (analytic)
25   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)   ::   gsi3w              !: model level depth coefficient at w-level (sum of gsigw)
26   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)   ::   esigt, esigw       !: vertical scale factor coef. at t-, w-levels
27
28   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hbatv , hbatf      !: ocean depth at the vertical of  v--f
29   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hbatt , hbatu      !:                                 t--u points (m)
30   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   scosrf, scobot     !: ocean surface and bottom topographies
31   !                                                                           !  (if deviating from coordinate surfaces in HYBRID)
32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hifv  , hiff       !: interface depth between stretching at v--f
33   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hift  , hifu       !: and quasi-uniform spacing             t--u points (m)
34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   rx1                !: Maximum grid stiffness ratio
35
36   !!----------------------------------------------------------------------
37   !! masks, bathymetry
38   !! ---------------------------------------------------------------------
39   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mbathy             !: number of ocean level (=0, 1, ... , jpk-1)
40   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   bathy              !: ocean depth (meters - read from file)
41
42   !! Other variables needed by scoord_gen
43   INTEGER  ::   jpi, jpj, jpk            ! Size of the domain - read from bathy or namelist?
44   INTEGER  ::   ji, jj, jk, jl           ! dummy loop argument
45   INTEGER  ::   iip1, ijp1, iim1, ijm1   ! temporary integers
46   INTEGER  ::   ios                      ! Local integer output status for namelist read and allocation
47   INTEGER,PARAMETER  ::   numnam=8       ! File handle for namelist
48   REAL(wp) ::   zrmax, ztaper   ! temporary scalars
49   REAL(wp) ::   zrfact
50   !
51   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ztmpi1, ztmpi2, ztmpj1, ztmpj2
52   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zenv, ztmp, zmsk, zri, zrj, zhbat
53
54   !Namelist variables
55   REAL(wp) :: rn_jpk, rn_sbot_min, rn_sbot_max, rn_hc, rn_rmax, rn_theta
56   REAL(wp) :: rn_thetb, rn_bb, rn_alpha, rn_efold, rn_zs, rn_zb_a, rn_zb_b
57   LOGICAL :: ln_s_sh94, ln_s_sf12, ln_sigcrit
58
59   NAMELIST/namzgr_sco/rn_jpk, ln_s_sh94, ln_s_sf12, ln_sigcrit, rn_sbot_min, rn_sbot_max, rn_hc, rn_rmax,rn_theta, &
60                        rn_thetb, rn_bb, rn_alpha, rn_efold, rn_zs, rn_zb_a, rn_zb_b
61
62  ! IDs for output netcdf file
63  INTEGER :: id_x, id_y, id_z
64  INTEGER :: ncout
65  INTEGER, DIMENSION(11) :: var_ids  !Array to contain all variable IDs
66
67   CONTAINS
68
69   INTEGER FUNCTION dom_oce_alloc()
70      !!----------------------------------------------------------------------
71      INTEGER, DIMENSION(4) :: ierr
72      !!----------------------------------------------------------------------
73      ierr(:) = 0
74      !
75      ALLOCATE( zenv(jpi,jpj), ztmp(jpi,jpj), zmsk(jpi,jpj), zri(jpi,jpj), zrj(jpi,jpj), &
76         &      zhbat(jpi,jpj) , ztmpi1(jpi,jpj), ztmpi2(jpi,jpj), ztmpj1(jpi,jpj), ztmpj2(jpi,jpj), STAT=ierr(1) )
77         !
78      ALLOCATE( gdep3w_0(jpi,jpj) , e3v_0(jpi,jpj) , e3f_0(jpi,jpj) ,                         &
79         &      gdept_0(jpi,jpj) , e3t_0(jpi,jpj) , e3u_0 (jpi,jpj) ,                         &
80         &      gdepw_0(jpi,jpj) , e3w_0(jpi,jpj) , e3vw_0(jpi,jpj) , e3uw_0(jpi,jpj) , STAT=ierr(2) )
81         !
82         !
83         !
84      ALLOCATE( hbatv (jpi,jpj) , hbatf (jpi,jpj) ,     &
85         &      hbatt (jpi,jpj) , hbatu (jpi,jpj) ,     &
86         &      scosrf(jpi,jpj) , scobot(jpi,jpj) ,     &
87         &      hifv  (jpi,jpj) , hiff  (jpi,jpj) ,     &
88         &      hift  (jpi,jpj) , hifu  (jpi,jpj) , rx1 (jpi,jpj) , STAT=ierr(3) )
89
90      ALLOCATE( mbathy(jpi,jpj) , STAT=ierr(4) )
91     !
92      dom_oce_alloc = MAXVAL(ierr)
93      !
94   END FUNCTION dom_oce_alloc
95   
96
97   SUBROUTINE read_bathy()
98     !! Read bathymetry from input netcdf file
99     INTEGER :: var_id, ncin
100
101     CALL check_nf90( nf90_open('bathy.nc', NF90_NOWRITE, ncin), 'Error opening bathy.nc file' )
102
103     ! Find the size of the input bathymetry
104     CALL dimlen(ncin, 'lon', jpi)   
105     CALL dimlen(ncin, 'lat', jpj)   
106     
107     ALLOCATE( bathy(jpi, jpj) )
108     
109     ! Read the bathymetry variable from file
110     CALL check_nf90( nf90_inq_varid( ncin, 'Bathymetry', var_id ), 'Cannot get variable ID for bathymetry')
111     CALL check_nf90( nf90_get_var( ncin, var_id, bathy, (/ 1,1 /), (/ jpi, jpj /) ) )
112
113     CALL check_nf90( nf90_close(ncin), 'Error closing bathy.nc file' )
114
115   END SUBROUTINE read_bathy
116
117   SUBROUTINE dimlen( ncid, dimname, len )
118     ! Determine the length of dimension dimname
119     INTEGER, INTENT(in)          :: ncid
120     CHARACTER(LEN=*), INTENT(in) :: dimname
121     INTEGER, INTENT(out)         :: len
122     ! Local variables
123     INTEGER :: id_var, istatus
124
125     id_var = 1
126     CALL check_nf90( nf90_inq_dimid(ncid, dimname, id_var), 'Dimension not found in file')
127     CALL check_nf90( nf90_inquire_dimension(ncid,id_var,len=len))
128
129   END SUBROUTINE dimlen
130 
131 
132   SUBROUTINE make_coord_file()
133     ! Create new coordinates file and define dimensions and variables ready for
134     ! writing
135     
136
137     !Create the file
138     CALL check_nf90( nf90_create('coord_zgr.nc', NF90_CLOBBER, ncout), 'Could not create output file')
139     !
140     !Define dimensions
141     CALL check_nf90( nf90_def_dim(ncout, 'x', jpi, id_x) )
142     CALL check_nf90( nf90_def_dim(ncout, 'y', jpj, id_y) )
143     CALL check_nf90( nf90_def_dim(ncout, 'z', jpk, id_z) )
144     !
145     !Define variables
146     CALL check_nf90( nf90_def_var(ncout, 'gdept_0', nf90_double, (/id_x, id_y,id_z/), var_ids(1)) )
147     CALL check_nf90( nf90_def_var(ncout, 'gdepw_0', nf90_double, (/id_x, id_y,id_z/), var_ids(2)) )
148     CALL check_nf90( nf90_def_var(ncout, 'gdep3w_0', nf90_double, (/id_x, id_y,id_z/), var_ids(3)) )
149     CALL check_nf90( nf90_def_var(ncout, 'e3f_0', nf90_double, (/id_x, id_y,id_z/), var_ids(4)) )
150     CALL check_nf90( nf90_def_var(ncout, 'e3t_0', nf90_double, (/id_x, id_y,id_z/), var_ids(5)) )
151     CALL check_nf90( nf90_def_var(ncout, 'e3u_0', nf90_double, (/id_x, id_y,id_z/), var_ids(6)) )
152     CALL check_nf90( nf90_def_var(ncout, 'e3v_0', nf90_double, (/id_x, id_y,id_z/), var_ids(7)) )
153     CALL check_nf90( nf90_def_var(ncout, 'e3w_0', nf90_double, (/id_x, id_y,id_z/), var_ids(8)) )
154     CALL check_nf90( nf90_def_var(ncout, 'e3uw_0', nf90_double, (/id_x, id_y,id_z/), var_ids(9)) )
155     CALL check_nf90( nf90_def_var(ncout, 'e3vw_0', nf90_double, (/id_x, id_y,id_z/), var_ids(10)) )
156     CALL check_nf90( nf90_def_var(ncout, 'mbathy', nf90_double, (/id_x, id_y/), var_ids(11)) )
157     
158     ! End define mode
159     CALL check_nf90( nf90_enddef(ncout) )
160     
161     WRITE(*,*) 'Opened coord_zgr.nc file and defined variables'
162
163   END SUBROUTINE make_coord_file
164
165   SUBROUTINE write_netcdf_vars(kk)
166   ! Write  variables to the netcdf file at level kk
167     INTEGER, INTENT(in) :: kk
168
169     CALL check_nf90( nf90_put_var(ncout, var_ids(1), gdept_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
170     CALL check_nf90( nf90_put_var(ncout, var_ids(2), gdepw_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
171     CALL check_nf90( nf90_put_var(ncout, var_ids(3), gdep3w_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
172     CALL check_nf90( nf90_put_var(ncout, var_ids(4), e3f_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
173     CALL check_nf90( nf90_put_var(ncout, var_ids(5), e3t_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
174     CALL check_nf90( nf90_put_var(ncout, var_ids(6), e3u_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
175     CALL check_nf90( nf90_put_var(ncout, var_ids(7), e3v_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
176     CALL check_nf90( nf90_put_var(ncout, var_ids(8), e3w_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
177     CALL check_nf90( nf90_put_var(ncout, var_ids(9), e3uw_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
178     CALL check_nf90( nf90_put_var(ncout, var_ids(10), e3vw_0, (/ 1,1,kk /), (/ jpi, jpj,1 /) ) )
179
180   END SUBROUTINE write_netcdf_vars
181
182   SUBROUTINE check_nf90( istat, message )
183      !Check for netcdf errors
184      INTEGER, INTENT(in) :: istat
185      CHARACTER(LEN=*), INTENT(in), OPTIONAL :: message
186
187      IF (istat /= nf90_noerr) THEN
188         WRITE(*,*) 'ERROR! : '//TRIM(nf90_strerror(istat))
189         IF ( PRESENT(message) ) THEN ; WRITE(*,*) message ; ENDIF
190         STOP
191      ENDIF
192
193   END SUBROUTINE check_nf90
194
195
196END MODULE utils
Note: See TracBrowser for help on using the repository browser.