source: branches/2014/dev_r4650_UKMO3_masked_damping/NEMOGCM/TOOLS/DMP_TOOLS/src/coast_dist.F90 @ 4738

Last change on this file since 4738 was 4738, checked in by timgraham, 7 years ago

Modified tra_dmp module to read in restoration coefficient from a netcdf file

Added a tool to create the netcdf file - this replaces all of the hard coded resolution dependencies in tra_dmp_init

File size: 9.1 KB
Line 
1MODULE coastdist
2
3   USE utils
4   USE netcdf
5
6   IMPLICIT NONE
7   PUBLIC
8
9   CONTAINS
10
11   SUBROUTINE coast_dist_weight( presto ) 
12      !!----------------------------------------------------------------------
13      !!                 *** ROUTINE coast_dist_weight ***
14      !!
15      !! ** Purpose: Weight restoration coefficient by distance to coast
16      !!
17      !! ** Method: 1) Calculate distance to coast
18      !!            2) Reduce resto with 1000km of coast
19      !!
20      IMPLICIT NONE
21      REAL(8), DIMENSION(jpi,jpj), INTENT( inout ) :: presto
22      REAL(8), DIMENSION(jpi,jpj) :: zdct
23      REAL(8) :: zinfl = 1000.e3  ! Distance of influence of coast line (could be
24                                  ! a namelist setting)
25      INTEGER :: jj, ji           ! dummy loop indices
26     
27
28      CALL cofdis( zdct )
29      DO jj = 1, jpj
30         DO ji = 1, jpi
31            zdct(ji,jj) = MIN( zinfl, zdct(ji,jj) )
32            presto(ji,jj) = presto(ji, jj) * 0.5 * (  1. - COS( rpi*zdct(ji,jj)/zinfl) )
33         END DO
34      END DO
35
36   END SUBROUTINE coast_dist_weight
37
38
39   SUBROUTINE cofdis( pdct )
40      !!----------------------------------------------------------------------
41      !!                 ***  ROUTINE cofdis  ***
42      !!
43      !! ** Purpose :   Compute the distance between ocean T-points and the
44      !!      ocean model coastlines.
45      !!
46      !! ** Method  :   For each model level, the distance-to-coast is
47      !!      computed as follows :
48      !!       - The coastline is defined as the serie of U-,V-,F-points
49      !!      that are at the ocean-land bound.
50      !!       - For each ocean T-point, the distance-to-coast is then
51      !!      computed as the smallest distance (on the sphere) between the
52      !!      T-point and all the coastline points.
53      !!       - For land T-points, the distance-to-coast is set to zero.
54      !!
55      !! ** Action  : - pdct, distance to the coastline (argument)
56      !!              - NetCDF file 'dist.coast.nc'
57      !!----------------------------------------------------------------------
58      !!
59      REAL(8), DIMENSION(jpi,jpj), INTENT( out ) ::   pdct   ! distance to the coastline
60      !!
61      INTEGER ::   ji, jj, jl   ! dummy loop indices
62      INTEGER ::   iju, ijt, icoast, itime, ierr, icot   ! local integers
63      CHARACTER (len=32) ::   clname                     ! local name
64      REAL(8) ::   zdate0                               ! local scalar
65      REAL(8), POINTER, DIMENSION(:,:) ::  zxt, zyt, zzt, zmask
66      REAL(8), POINTER, DIMENSION(:  ) ::  zxc, zyc, zzc, zdis    ! temporary workspace
67      LOGICAL , ALLOCATABLE, DIMENSION(:,:) ::  llcotu, llcotv, llcotf   ! 2D logical workspace
68
69      !!----------------------------------------------------------------------
70      !
71      ALLOCATE( zxt(jpi,jpj) , zyt(jpi,jpj) , zzt(jpi,jpj) , zmask(jpi,jpj)    )
72      ALLOCATE(zxc(3*jpi*jpj), zyc(3*jpi*jpj), zzc(3*jpi*jpj), zdis(3*jpi*jpj)     )
73      ALLOCATE( llcotu(jpi,jpj), llcotv(jpi,jpj), llcotf(jpi,jpj)  )
74      ALLOCATE( gphiu(jpi,jpj), gphiv(jpi,jpj), gphif(jpi,jpj)  )
75      ALLOCATE( glamu(jpi,jpj), glamv(jpi,jpj), glamf(jpi,jpj), glamt(jpi,jpj)  )
76      ALLOCATE( umask(jpi,jpj), vmask(jpi,jpj), fmask(jpi,jpj)  )
77      !
78
79      CALL check_nf90( nf90_get_var( ncin, gphit_id, gphit, (/ 1,1 /), (/ jpi, jpj /) ) )
80      CALL check_nf90( nf90_get_var( ncin, gphiu_id, gphiu, (/ 1,1 /), (/ jpi, jpj /) ) )
81      CALL check_nf90( nf90_get_var( ncin, gphiv_id, gphiv, (/ 1,1 /), (/ jpi, jpj /) ) )
82      CALL check_nf90( nf90_get_var( ncin, gphif_id, gphif, (/ 1,1 /), (/ jpi, jpj /) ) )
83      CALL check_nf90( nf90_get_var( ncin, glamt_id, glamt, (/ 1,1 /), (/ jpi, jpj /) ) )
84      CALL check_nf90( nf90_get_var( ncin, glamu_id, glamu, (/ 1,1 /), (/ jpi, jpj /) ) )
85      CALL check_nf90( nf90_get_var( ncin, glamv_id, glamv, (/ 1,1 /), (/ jpi, jpj /) ) )
86      CALL check_nf90( nf90_get_var( ncin, glamf_id, glamf, (/ 1,1 /), (/ jpi, jpj /) ) )
87      CALL check_nf90( nf90_get_var( ncin, tmask_id, tmask, (/ 1,1 /), (/ jpi, jpj /) ) )
88      CALL check_nf90( nf90_get_var( ncin, umask_id, umask, (/ 1,1 /), (/ jpi, jpj /) ) )
89      CALL check_nf90( nf90_get_var( ncin, vmask_id, vmask, (/ 1,1 /), (/ jpi, jpj /) ) )
90      CALL check_nf90( nf90_get_var( ncin, fmask_id, fmask, (/ 1,1 /), (/ jpi, jpj /) ) )
91
92      pdct(:,:) = 0.
93      zxt(:,:) = COS( rad * gphit(:,:) ) * COS( rad * glamt(:,:) )
94      zyt(:,:) = COS( rad * gphit(:,:) ) * SIN( rad * glamt(:,:) )
95      zzt(:,:) = SIN( rad * gphit(:,:) )
96
97
98         ! Define the coastline points (U, V and F)
99         DO jj = 2, jpj-1
100            DO ji = 2, jpi-1
101               zmask(ji,jj) =  ( tmask(ji,jj+1) + tmask(ji+1,jj+1) &
102                   &           + tmask(ji,jj  ) + tmask(ji+1,jj  ) )
103               llcotu(ji,jj) = ( tmask(ji,jj ) + tmask(ji+1,jj  ) == 1. ) 
104               llcotv(ji,jj) = ( tmask(ji,jj  ) + tmask(ji  ,jj+1) == 1. ) 
105               llcotf(ji,jj) = ( zmask(ji,jj) > 0. ) .AND. ( zmask(ji,jj) < 4. )
106            END DO
107         END DO
108
109         ! Lateral boundaries conditions
110         llcotu(:, 1 ) = umask(:,  2  ) == 1
111         llcotu(:,jpj) = umask(:,jpj-1) == 1
112         llcotv(:, 1 ) = vmask(:,  2  ) == 1
113         llcotv(:,jpj) = vmask(:,jpj-1) == 1
114         llcotf(:, 1 ) = fmask(:,  2  ) == 1
115         llcotf(:,jpj) = fmask(:,jpj-1) == 1
116
117         IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) THEN
118            llcotu( 1 ,:) = llcotu(jpi-1,:)
119            llcotu(jpi,:) = llcotu(  2  ,:)
120            llcotv( 1 ,:) = llcotv(jpi-1,:)
121            llcotv(jpi,:) = llcotv(  2  ,:)
122            llcotf( 1 ,:) = llcotf(jpi-1,:)
123            llcotf(jpi,:) = llcotf(  2  ,:)
124         ELSE
125            llcotu( 1 ,:) = umask(  2  ,:) == 1
126            llcotu(jpi,:) = umask(jpi-1,:) == 1
127            llcotv( 1 ,:) = vmask(  2  ,:) == 1
128            llcotv(jpi,:) = vmask(jpi-1,:) == 1
129            llcotf( 1 ,:) = fmask(  2  ,:) == 1
130            llcotf(jpi,:) = fmask(jpi-1,:) == 1
131         ENDIF
132         IF( jperio == 3 .OR. jperio == 4 ) THEN
133            DO ji = 1, jpi-1
134               iju = jpi - ji + 1
135               llcotu(ji,jpj  ) = llcotu(iju,jpj-2)
136               llcotf(ji,jpj-1) = llcotf(iju,jpj-2)
137               llcotf(ji,jpj  ) = llcotf(iju,jpj-3)
138            END DO
139            DO ji = jpi/2, jpi-1
140               iju = jpi - ji + 1
141               llcotu(ji,jpj-1) = llcotu(iju,jpj-1)
142            END DO
143            DO ji = 2, jpi
144               ijt = jpi - ji + 2
145               llcotv(ji,jpj-1) = llcotv(ijt,jpj-2)
146               llcotv(ji,jpj  ) = llcotv(ijt,jpj-3)
147            END DO
148         ENDIF
149         IF( jperio == 5 .OR. jperio == 6 ) THEN
150            DO ji = 1, jpi-1
151               iju = jpi - ji
152               llcotu(ji,jpj  ) = llcotu(iju,jpj-1)
153               llcotf(ji,jpj  ) = llcotf(iju,jpj-2)
154            END DO
155            DO ji = jpi/2, jpi-1
156               iju = jpi - ji
157               llcotf(ji,jpj-1) = llcotf(iju,jpj-1)
158            END DO
159            DO ji = 1, jpi
160               ijt = jpi - ji + 1
161               llcotv(ji,jpj  ) = llcotv(ijt,jpj-1)
162            END DO
163            DO ji = jpi/2+1, jpi
164               ijt = jpi - ji + 1
165               llcotv(ji,jpj-1) = llcotv(ijt,jpj-1)
166            END DO
167         ENDIF
168
169         ! Compute cartesian coordinates of coastline points
170         ! and the number of coastline points
171         icoast = 0
172         DO jj = 1, jpj
173            DO ji = 1, jpi
174               IF( llcotf(ji,jj) ) THEN
175                  icoast = icoast + 1
176                  zxc(icoast) = COS( rad*gphif(ji,jj) ) * COS( rad*glamf(ji,jj) )
177                  zyc(icoast) = COS( rad*gphif(ji,jj) ) * SIN( rad*glamf(ji,jj) )
178                  zzc(icoast) = SIN( rad*gphif(ji,jj) )
179               ENDIF
180               IF( llcotu(ji,jj) ) THEN
181                  icoast = icoast+1
182                  zxc(icoast) = COS( rad*gphiu(ji,jj) ) * COS( rad*glamu(ji,jj) )
183                  zyc(icoast) = COS( rad*gphiu(ji,jj) ) * SIN( rad*glamu(ji,jj) )
184                  zzc(icoast) = SIN( rad*gphiu(ji,jj) )
185               ENDIF
186               IF( llcotv(ji,jj) ) THEN
187                  icoast = icoast+1
188                  zxc(icoast) = COS( rad*gphiv(ji,jj) ) * COS( rad*glamv(ji,jj) )
189                  zyc(icoast) = COS( rad*gphiv(ji,jj) ) * SIN( rad*glamv(ji,jj) )
190                  zzc(icoast) = SIN( rad*gphiv(ji,jj) )
191               ENDIF
192            END DO
193         END DO
194
195         ! Distance for the T-points
196         DO jj = 1, jpj
197            DO ji = 1, jpi
198               IF( tmask(ji,jj) == 0. ) THEN
199                  pdct(ji,jj) = 0.
200               ELSE
201                  DO jl = 1, icoast
202                     zdis(jl) = ( zxt(ji,jj) - zxc(jl) )**2   &
203                        &     + ( zyt(ji,jj) - zyc(jl) )**2   &
204                        &     + ( zzt(ji,jj) - zzc(jl) )**2
205                  END DO
206                  pdct(ji,jj) = ra * SQRT( MINVAL( zdis(1:icoast) ) )
207               ENDIF
208            END DO
209         END DO
210
211      DEALLOCATE( zxt , zyt , zzt , zmask    )
212      DEALLOCATE(zxc, zyc, zzc, zdis    )
213      DEALLOCATE( llcotu, llcotv, llcotf  )
214      DEALLOCATE( gphiu, gphiv, gphif  )
215      DEALLOCATE( glamu, glamv, glamf, glamt  )
216      DEALLOCATE( umask, vmask, fmask  )
217
218   END SUBROUTINE cofdis
219
220END MODULE coastdist
Note: See TracBrowser for help on using the repository browser.