1 | MODULE 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(wp), DIMENSION(jpi,jpj), INTENT( inout ) :: presto |
---|
22 | REAL(wp), DIMENSION(jpi,jpj) :: zdct |
---|
23 | REAL(wp) :: zinfl = 1000.e3_wp ! 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_wp * ( 1._wp - 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(wp), 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(wp) :: zdate0 ! local scalar |
---|
65 | REAL(wp), POINTER, DIMENSION(:,:) :: zxt, zyt, zzt, zmask |
---|
66 | REAL(wp), 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._wp |
---|
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._wp ) |
---|
104 | llcotv(ji,jj) = ( tmask(ji,jj ) + tmask(ji ,jj+1) == 1._wp ) |
---|
105 | llcotf(ji,jj) = ( zmask(ji,jj) > 0._wp ) .AND. ( zmask(ji,jj) < 4._wp ) |
---|
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._wp ) THEN |
---|
199 | pdct(ji,jj) = 0._wp |
---|
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 | |
---|
220 | END MODULE coastdist |
---|