[11201] | 1 | MODULE domclo |
---|
| 2 | !!====================================================================== |
---|
| 3 | !! *** MODULE domclo *** |
---|
| 4 | !! domclo : definition of closed sea mask needed by NEMO |
---|
| 5 | !!====================================================================== |
---|
| 6 | !! History : NEMO 4.0 ! 07-2019 (P. Mathiot) original version |
---|
| 7 | !!---------------------------------------------------------------------- |
---|
| 8 | !!---------------------------------------------------------------------- |
---|
| 9 | !! dom_clo : definition of closed sea mask needed by NEMO |
---|
| 10 | !!---------------------------------------------------------------------- |
---|
| 11 | USE dom_oce ! ocean space and time domain |
---|
| 12 | USE domngb ! closest point algorithm |
---|
| 13 | USE phycst ! rpi, rad, ra |
---|
| 14 | USE domutil ! flood filling algorithm (fill_pool) |
---|
| 15 | USE in_out_manager ! I/O manager |
---|
| 16 | USE lbclnk ! lateral boundary condition - MPP exchanges |
---|
| 17 | USE lib_mpp |
---|
| 18 | USE lib_fortran |
---|
| 19 | |
---|
| 20 | IMPLICIT NONE |
---|
| 21 | PRIVATE |
---|
| 22 | |
---|
| 23 | PUBLIC dom_clo ! routine called by domain module |
---|
| 24 | |
---|
| 25 | !! * Substitutions |
---|
| 26 | # include "vectopt_loop_substitute.h90" |
---|
| 27 | !!---------------------------------------------------------------------- |
---|
| 28 | !! NEMO/OPA 3.3 , NEMO Consortium (2010) |
---|
| 29 | !! $Id: closea.F90 7213 2016-11-09 09:07:42Z timgraham $ |
---|
| 30 | !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) |
---|
| 31 | !!---------------------------------------------------------------------- |
---|
| 32 | CONTAINS |
---|
| 33 | |
---|
| 34 | SUBROUTINE dom_clo |
---|
| 35 | !!--------------------------------------------------------------------- |
---|
| 36 | !! *** ROUTINE dom_clo *** |
---|
| 37 | !! |
---|
| 38 | !! ** Purpose : Closed sea mask definition |
---|
| 39 | !! |
---|
| 40 | !! ** Method : A flood filling algorithm is used to detect lake and open ocean |
---|
| 41 | !! Namelist provide seed for each lake and condition on how to |
---|
| 42 | !! manage net evap/runoff. |
---|
| 43 | !! |
---|
| 44 | !! ** Action : - compute open ocean mask |
---|
| 45 | !! - loop over all the lake in namelist |
---|
| 46 | !! - define lake msk (msk_csglo/emp/rnf) |
---|
| 47 | !! - define river mouth mask for each lake or group of lake with the same river mouth |
---|
| 48 | !! (msk_csgrpglo/emp/rnf) |
---|
| 49 | !! |
---|
| 50 | !!---------------------------------------------------------------------- |
---|
| 51 | |
---|
| 52 | TYPE closea |
---|
| 53 | CHARACTER(256) :: cname ! name |
---|
| 54 | REAL(wp) :: rlonsrc ! seed src location (lon) |
---|
| 55 | REAL(wp) :: rlatsrc ! seed src location (lat) |
---|
| 56 | REAL(wp) :: rlontrg ! seed trg location (lon) |
---|
| 57 | REAL(wp) :: rlattrg ! seed trg location (lat) |
---|
| 58 | CHARACTER(256) :: cloctrg ! where water is spread |
---|
| 59 | CHARACTER(256) :: cschtrg ! how water is spread |
---|
| 60 | REAL(wp) :: radtrg ! radius of closed sea river mouth |
---|
| 61 | INTEGER :: idtrg ! target id in case multiple lakes for the same river mouth |
---|
| 62 | END TYPE |
---|
| 63 | |
---|
| 64 | INTEGER, PARAMETER :: jp_closea = 98 ! number maximal of closea |
---|
| 65 | INTEGER :: ji, jj, jcs, jsch ! loop indexes |
---|
| 66 | INTEGER :: jglo, jemp, jrnf ! closed sea indexes for global, emp, rnf spreading |
---|
| 67 | INTEGER :: jiseed, jjseed ! seed indexes |
---|
| 68 | INTEGER :: jglo, jloc, jcoast ! counter |
---|
| 69 | INTEGER :: ios |
---|
| 70 | INTEGER :: nn_closea ! number of closed seas |
---|
| 71 | |
---|
| 72 | REAL(wp) :: zdistseed ! distance to seed |
---|
| 73 | REAL(wp) :: zarea ! river mouth area |
---|
| 74 | REAL(wp) :: rn_lon_opnsea, rn_lat_opnsea ! open sea seed |
---|
| 75 | REAL(wp), DIMENSION(1) :: zchk, zradtrg |
---|
| 76 | REAL(wp), DIMENSION(jpi,jpj) :: zmsksrc, zmsktrg, zmsk_coastline, zdist |
---|
| 77 | |
---|
| 78 | CHARACTER(256) :: csch, cloc ! scheme name for water spreading (glo, rnf, emp) |
---|
| 79 | TYPE(closea) , DIMENSION(jp_closea) :: sn_lake ! lake properties |
---|
| 80 | |
---|
| 81 | LOGICAL :: lskip ! flag in case lake seed on land or already filled (...) |
---|
| 82 | |
---|
| 83 | !!---------------------------------------------------------------------- |
---|
| 84 | !! 0 : Read namelist for closed sea definition |
---|
| 85 | !!---------------------------------------------------------------------- |
---|
| 86 | IF(lwp) WRITE(numout,*) |
---|
| 87 | IF(lwp) WRITE(numout,*)'dom_clo : closed seas ' |
---|
| 88 | IF(lwp) WRITE(numout,*)'~~~~~~~' |
---|
| 89 | |
---|
| 90 | NAMELIST/namclo/ rn_lon_opnsea, rn_lat_opnsea, nn_closea, sn_lake |
---|
| 91 | !!--------------------------------------------------------------------- |
---|
| 92 | PRINT *, rn_lon_opnsea, rn_lat_opnsea, nn_closea |
---|
| 93 | !! |
---|
| 94 | REWIND( numnam_ref ) ! Namelist namlbc in reference namelist : Lateral momentum boundary condition |
---|
| 95 | READ ( numnam_ref, namclo, IOSTAT = ios, ERR = 901 ) |
---|
| 96 | 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namclo in reference namelist', lwp ) |
---|
| 97 | REWIND( numnam_cfg ) ! Namelist namlbc in configuration namelist : Lateral momentum boundary condition |
---|
| 98 | READ ( numnam_cfg, namclo, IOSTAT = ios, ERR = 902 ) |
---|
| 99 | 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namclo in configuration namelist', lwp ) |
---|
| 100 | IF(lwm) WRITE ( numond, namclo ) |
---|
| 101 | |
---|
| 102 | !!---------------------------------------------------------------------- |
---|
| 103 | !! 1 : Define open ocean and closed sea masks |
---|
| 104 | !!---------------------------------------------------------------------- |
---|
| 105 | |
---|
| 106 | !! 1.1 get closest point to the lat/lon seed |
---|
| 107 | msk_opnsea(:,:) = ssmask(:,:) |
---|
| 108 | CALL dom_ngb(rn_lon_opnsea, rn_lat_opnsea, jiseed, jjseed, zdistseed, 'T') |
---|
| 109 | |
---|
| 110 | !! 1.2 fill select cell to -1 |
---|
| 111 | CALL fill_pool( jiseed, jjseed, msk_opnsea, -1._wp ) |
---|
| 112 | |
---|
| 113 | !! sanity check on the msk value (closea mask should be 99 otherwise, lake already processed) |
---|
| 114 | !! check if some lake are not connected |
---|
| 115 | zchk = 0._wp |
---|
| 116 | IF (mi0(jiseed) == mi1(jiseed) .AND. mj0(jjseed) == mj1(jjseed)) zchk = ssmask(mi0(jiseed),mj0(jjseed)) |
---|
| 117 | CALL mpp_max('domclo',zchk) |
---|
| 118 | IF (zchk(1) == 0._wp) CALL ctl_stop( 'STOP', 'open sea seed is on land, please update namelist (rn_lon_opnsea,rn_lat_opnsea)' ) |
---|
| 119 | |
---|
| 120 | !! 1.3 set to 0 everything >0 and revert mask |
---|
| 121 | IF (lwp) THEN |
---|
| 122 | WRITE(numout,*) |
---|
| 123 | WRITE(numout,*)'open ocean: ' |
---|
| 124 | WRITE(numout,'(a,2f7.2)')' lon/lat seed to detect main ocean is: ', rn_lon_opnsea, rn_lat_opnsea |
---|
| 125 | WRITE(numout,'(a,2i7)' )' i/j seed to detect main ocean is: ', jiseed, jjseed |
---|
| 126 | WRITE(numout,'(a,f8.0)' )' distance (m) between namelist seed location and model seed location: ', zdistseed |
---|
| 127 | WRITE(numout,*) |
---|
| 128 | END IF |
---|
| 129 | |
---|
| 130 | WHERE (msk_opnsea(:,:) > 0) msk_opnsea(:,:) = 0 ! mask all closed seas |
---|
| 131 | WHERE (msk_opnsea(:,:) < 0) msk_opnsea(:,:) = 1 ! restore mask value |
---|
| 132 | |
---|
| 133 | !! 1.4 Define closed sea mask (all of them, ie defined in the namelist or not) |
---|
| 134 | !! needed to remove the undefined closed seas at the end |
---|
| 135 | msk_closea = (ssmask - msk_opnsea) * 99 |
---|
| 136 | |
---|
| 137 | !!---------------------------------------------------------------------- |
---|
| 138 | !! 2 : compute closed sea mask for global, rnf and emp cases |
---|
| 139 | !!---------------------------------------------------------------------- |
---|
| 140 | |
---|
| 141 | !! 2.1 : initialisation of masks and lake indexes |
---|
| 142 | jglo = 1 ; jrnf = 1 ; jemp = 1 |
---|
| 143 | !! mask used to group lake by net evap/precip distribution technics |
---|
| 144 | msk_glo = msk_closea |
---|
| 145 | msk_rnf = msk_closea |
---|
| 146 | msk_emp = msk_closea |
---|
| 147 | |
---|
| 148 | !! mask used to group multiple lake with the same river mouth (great lake for example) |
---|
| 149 | msk_gloid = 0.0_wp |
---|
| 150 | msk_rnfid = 0.0_wp |
---|
| 151 | msk_empid = 0.0_wp |
---|
| 152 | |
---|
| 153 | IF (lwp) WRITE(numout,*)'closed seas: ' |
---|
| 154 | |
---|
| 155 | DO jcs = 1,nn_closea |
---|
| 156 | |
---|
| 157 | !! 2.2 : print out and sanity check |
---|
| 158 | !! flag changed in sanity checks |
---|
| 159 | lskip = .FALSE. |
---|
| 160 | |
---|
| 161 | !! define seed to extract the closed sea jcs |
---|
| 162 | CALL dom_ngb(sn_lake(jcs)%rlonsrc, sn_lake(jcs)%rlatsrc, jiseed, jjseed, zdistseed, 'T') |
---|
| 163 | |
---|
| 164 | !! how the excess of the closed seas is spread out: |
---|
| 165 | IF (lwp) THEN |
---|
| 166 | WRITE(numout,*) |
---|
| 167 | WRITE(numout,'(a,a10,a,2f7.2,a)')' processing lake ',TRIM(sn_lake(jcs)%cname) & |
---|
| 168 | & ,' ( lat/lon : ',sn_lake(jcs)%rlatsrc, sn_lake(jcs)%rlonsrc,' )' |
---|
| 169 | END IF |
---|
| 170 | cloc = sn_lake(jcs)%cloctrg |
---|
| 171 | csch = sn_lake(jcs)%cschtrg |
---|
| 172 | |
---|
| 173 | !! set up indexes and mask ! SELECT CASE ... |
---|
| 174 | SELECT CASE (csch) |
---|
| 175 | CASE('glo') |
---|
| 176 | jsch = jglo |
---|
| 177 | zmsksrc = msk_glo |
---|
| 178 | zmsktrg = msk_gloid |
---|
| 179 | IF (lwp) WRITE(numout,*)' net evap/precip will be spread globally' |
---|
| 180 | CASE('rnf') |
---|
| 181 | jsch = jrnf |
---|
| 182 | zmsksrc = msk_rnf |
---|
| 183 | zmsktrg = msk_rnfid |
---|
| 184 | IF (lwp) WRITE(numout,*)' net precip will be spread locally and net evap globally' |
---|
| 185 | CASE('emp') |
---|
| 186 | jsch = jemp |
---|
| 187 | zmsksrc = msk_emp |
---|
| 188 | zmsktrg = msk_empid |
---|
| 189 | IF (lwp) WRITE(numout,*)' net precip will be spread locally' |
---|
| 190 | CASE DEFAULT |
---|
| 191 | CALL ctl_stop( 'STOP', 'domclo: unknown target type for lake (should be glo, emp, rnf)' ) |
---|
| 192 | END SELECT |
---|
| 193 | |
---|
| 194 | !! sanity check on the msk value (closea mask should be 99 otherwise, lake already processed) |
---|
| 195 | !! check if some lake are not connected |
---|
| 196 | zchk = 0._wp |
---|
| 197 | IF (mi0(jiseed) == mi1(jiseed) .AND. mj0(jjseed) == mj1(jjseed)) zchk = msk_closea(mi0(jiseed),mj0(jjseed)) |
---|
| 198 | CALL mpp_max('domclo',zchk) |
---|
| 199 | |
---|
| 200 | IF (zchk(1) /= 99._wp) THEN |
---|
| 201 | IF (lwp) THEN |
---|
| 202 | WRITE(numout,*) |
---|
| 203 | WRITE(numout,*)'W A R N I N G: lake ',TRIM(sn_lake(jcs)%cname),' already processed or seed on land, we skip it.' |
---|
| 204 | WRITE(numout,*)'It is possible this lake is connected to another one already processed.' |
---|
| 205 | WRITE(numout,*)'It is possible the seed is on land because the lake is not present in the bathymetry file' |
---|
| 206 | WRITE(numout,*)'If concerned about it, check input file and namelist definition of this lake.' |
---|
| 207 | END IF |
---|
| 208 | lskip = .TRUE. |
---|
| 209 | END IF |
---|
| 210 | |
---|
| 211 | !! 2.3 : compute mask for global, rnf, emp case |
---|
| 212 | IF (.NOT. lskip) THEN |
---|
| 213 | !! fill close sea mask with counter value |
---|
| 214 | CALL fill_pool( jiseed, jjseed, zmsksrc, REAL(jsch ,8)) |
---|
| 215 | WHERE (zmsksrc == REAL(jsch,8)) |
---|
| 216 | msk_closea = 0._wp |
---|
| 217 | zmsktrg = sn_lake(jcs)%idtrg |
---|
| 218 | END WHERE |
---|
| 219 | |
---|
| 220 | !! compute location of river mouth and distance to river mouth |
---|
| 221 | IF (cloc /= 'global') THEN |
---|
| 222 | !! set a minimum value for radius of the river influence |
---|
| 223 | zradtrg = 0._wp |
---|
| 224 | IF (mi0(jiseed) == mi1(jiseed) .AND. mj0(jjseed) == mj1(jjseed)) THEN |
---|
| 225 | zradtrg(1) = MAX( sn_lake(jcs)%radtrg, e1t( mi0(jiseed), mj0(jjseed) ), e2t( mi0(jiseed), mj0(jjseed) ) ) |
---|
| 226 | END IF |
---|
| 227 | CALL mpp_max('domclo',zradtrg) |
---|
| 228 | |
---|
| 229 | !! compute seed location for print |
---|
| 230 | CALL dom_ngb(sn_lake(jcs)%rlontrg, sn_lake(jcs)%rlattrg, jiseed, jjseed, zdistseed, 'T') |
---|
| 231 | !! compute distance to river mouth |
---|
| 232 | zdist(:,:) = dist(sn_lake(jcs)%rlontrg, sn_lake(jcs)%rlattrg, glamt, gphit) |
---|
| 233 | END IF |
---|
| 234 | |
---|
| 235 | !! define estuary |
---|
| 236 | !! deal with local/coast cases |
---|
| 237 | SELECT CASE (cloc) |
---|
| 238 | CASE ('global') |
---|
| 239 | WHERE (msk_opnsea(:,:) == 1) zmsktrg = sn_lake(jcs)%idtrg |
---|
| 240 | |
---|
| 241 | CASE ('local') |
---|
| 242 | !! compute mask |
---|
| 243 | WHERE (zdist(:,:) < zradtrg(1) .AND. msk_opnsea(:,:) == 1 ) zmsktrg = sn_lake(jcs)%idtrg |
---|
| 244 | !! print |
---|
| 245 | IF (lwp) WRITE(numout,'(a,f7.0,a,2f7.2,a,2i7,a)')' river mouth area is defined by points within ',zradtrg(1) & |
---|
| 246 | & ,' m of lat/lon ', sn_lake(jcs)%rlontrg, sn_lake(jcs)%rlattrg & |
---|
| 247 | & ,' (closest point is i/j ',jiseed, jjseed,')' |
---|
| 248 | |
---|
| 249 | CASE ('coast') |
---|
| 250 | !! define coastline mask |
---|
| 251 | zmsk_coastline = 0._wp |
---|
| 252 | DO jj=2,jpj-1 |
---|
| 253 | DO ji=2,jpi-1 |
---|
| 254 | IF ( ssmask(ji,jj) == 1._wp .AND. SUM(ssmask(ji-1:ji+1,jj-1:jj+1)) < 9 ) zmsk_coastline(ji,jj) = 1._wp |
---|
| 255 | END DO |
---|
| 256 | END DO |
---|
| 257 | CALL lbc_lnk('domclo', zmsk_coastline,'T',1._wp) |
---|
| 258 | |
---|
| 259 | !! compute mask |
---|
| 260 | WHERE ( zdist(:,:) < zradtrg(1) .AND. zmsk_coastline(:,:) == 1 .AND. msk_opnsea(:,:) == 1 ) zmsktrg = sn_lake(jcs)%idtrg |
---|
| 261 | !! print |
---|
| 262 | IF (lwp) WRITE(numout,'(a,f7.0,a,2f7.2,a,2i7,a)')' river mouth area is defined by coastal points within ',zradtrg(1) & |
---|
| 263 | & ,' m of lat/lon ', sn_lake(jcs)%rlontrg, sn_lake(jcs)%rlattrg & |
---|
| 264 | & ,' (closest point is i/j ',jiseed, jjseed,')' |
---|
| 265 | |
---|
| 266 | CASE DEFAULT |
---|
| 267 | CALL ctl_stop( 'STOP', 'domclo: unknown target type for lake (should be global, coast or local)' ) |
---|
| 268 | |
---|
| 269 | END SELECT |
---|
| 270 | |
---|
| 271 | !! sanity check |
---|
| 272 | zarea = glob_sum('domclo',zmsktrg * msk_opnsea) |
---|
| 273 | IF (zarea == 0._wp) CALL ctl_stop( 'STOP', 'river mouth area is 0, tune lon/lat trg or radtrg for this lake') |
---|
| 274 | |
---|
| 275 | !! set up indexes and mask |
---|
| 276 | SELECT CASE (csch) |
---|
| 277 | CASE ('glo') |
---|
| 278 | jglo = jglo + 1 |
---|
| 279 | msk_glo = zmsksrc |
---|
| 280 | msk_gloid = zmsktrg |
---|
| 281 | CASE ('rnf') |
---|
| 282 | jrnf = jrnf + 1 |
---|
| 283 | msk_rnf = zmsksrc |
---|
| 284 | msk_rnfid = zmsktrg |
---|
| 285 | CASE ('emp') |
---|
| 286 | jemp = jemp + 1 |
---|
| 287 | msk_emp = zmsksrc |
---|
| 288 | msk_empid = zmsktrg |
---|
| 289 | END SELECT |
---|
| 290 | |
---|
| 291 | END IF ! lskip |
---|
| 292 | |
---|
| 293 | END DO ! nn_closea |
---|
| 294 | |
---|
| 295 | !!---------------------------------------------------------------------- |
---|
| 296 | !! 3 : clean the masks of possible remaining undefined closed seas |
---|
| 297 | !!---------------------------------------------------------------------- |
---|
| 298 | |
---|
| 299 | !! mask all the cells not defined as closed sea |
---|
| 300 | WHERE ( msk_glo == 99 ) msk_glo = 0 |
---|
| 301 | WHERE ( msk_rnf == 99 ) msk_rnf = 0 |
---|
| 302 | WHERE ( msk_emp == 99 ) msk_emp = 0 |
---|
| 303 | |
---|
| 304 | !! non defined closed sea |
---|
| 305 | WHERE (msk_closea > 0) msk_closea = 1 |
---|
| 306 | |
---|
| 307 | END SUBROUTINE dom_clo |
---|
| 308 | |
---|
| 309 | !!====================================================================== |
---|
| 310 | |
---|
| 311 | END MODULE domclo |
---|
| 312 | |
---|