1 | !---------------------------------------------------------------------- |
---|
2 | ! NEMO system team, System and Interface for oceanic RElocable Nesting |
---|
3 | !---------------------------------------------------------------------- |
---|
4 | ! |
---|
5 | ! MODULE: interp |
---|
6 | ! |
---|
7 | ! DESCRIPTION: |
---|
8 | !> @brief |
---|
9 | !> This module manage linear interpolation on regular grid. |
---|
10 | !> |
---|
11 | !> @details |
---|
12 | !> to compute linear interpolation:<br/> |
---|
13 | !> @code |
---|
14 | !> CALL interp_linear_fill(dd_value, dd_fill, id_detect, id_rho, ld_even [,ld_discont] ) |
---|
15 | !> @endcode |
---|
16 | !> - dd_value is 2D array of variable value |
---|
17 | !> - dd_fill is the FillValue of variable |
---|
18 | !> - id_detect is 2D array of point to be interpolated (see interp module) |
---|
19 | !> - id_rho is array of refinment factor |
---|
20 | !> - ld_even indicates even refinment or not |
---|
21 | !> - ld_discont indicates longitudinal discontinuity (-180°/180°, 0°/360°) or not |
---|
22 | !> |
---|
23 | !> @author |
---|
24 | !> J.Paul |
---|
25 | ! REVISION HISTORY: |
---|
26 | !> @date September, 2014 -Initial version |
---|
27 | !> |
---|
28 | !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) |
---|
29 | !---------------------------------------------------------------------- |
---|
30 | |
---|
31 | MODULE interp_linear |
---|
32 | |
---|
33 | USE netcdf ! nf90 library |
---|
34 | USE global ! global variable |
---|
35 | USE kind ! F90 kind parameter |
---|
36 | USE logger ! log file manager |
---|
37 | USE fct ! basic useful function |
---|
38 | USE extrap ! extrapolation manager |
---|
39 | |
---|
40 | IMPLICIT NONE |
---|
41 | ! NOTE_avoid_public_variables_if_possible |
---|
42 | |
---|
43 | ! type and variable |
---|
44 | |
---|
45 | ! function and subroutine |
---|
46 | PUBLIC :: interp_linear_fill !< compute interpolation using linear method |
---|
47 | |
---|
48 | PRIVATE :: interp_linear__2D !< compute bilinear interpolation on 2D gid |
---|
49 | PRIVATE :: interp_linear__1D !< compute linear interpolation on 1D gid |
---|
50 | PRIVATE :: interp_linear__2D_coef !< compute coefficient for bilinear interpolation |
---|
51 | PRIVATE :: interp_linear__2D_fill !< fill value using bilinear interpolation |
---|
52 | PRIVATE :: interp_linear__1D_coef !< compute coefficient for linear interpolation |
---|
53 | PRIVATE :: interp_linear__1D_fill !< fill value using linear interpolation |
---|
54 | PRIVATE :: interp_linear__get_weight2D !< compute interpoaltion weight for 2D array. |
---|
55 | PRIVATE :: interp_linear__get_weight1D !< compute interpoaltion weight for 1D array. |
---|
56 | |
---|
57 | CONTAINS |
---|
58 | !------------------------------------------------------------------- |
---|
59 | !> @brief |
---|
60 | !> This subroutine compute horizontal linear interpolation on 4D array of value. |
---|
61 | !> |
---|
62 | !> @details |
---|
63 | !> |
---|
64 | !> @author J.Paul |
---|
65 | !> - September, 2014- Initial Version |
---|
66 | !> |
---|
67 | !> @param[inout] dd_value 2D array of variable value |
---|
68 | !> @param[in] dd_fill FillValue of variable |
---|
69 | !> @param[inout] id_detect 2D array of point to be interpolated |
---|
70 | !> @param[in] id_rho array of refinment factor |
---|
71 | !> @param[in] ld_even even refinment or not |
---|
72 | !> @param[in] ld_discont longitudinal discontinuity (-180°/180°, 0°/360°) or not |
---|
73 | !------------------------------------------------------------------- |
---|
74 | SUBROUTINE interp_linear_fill(dd_value, dd_fill, id_detect, & |
---|
75 | & id_rho, ld_even, ld_discont ) |
---|
76 | IMPLICIT NONE |
---|
77 | ! Argument |
---|
78 | REAL(dp) , DIMENSION(:,:,:,:), INTENT(INOUT) :: dd_value |
---|
79 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
80 | INTEGER(I4) , DIMENSION(:,:,:) , INTENT(INOUT) :: id_detect |
---|
81 | INTEGER(I4) , DIMENSION(:) , INTENT(IN ) :: id_rho |
---|
82 | LOGICAL , DIMENSION(:) , INTENT(IN ) :: ld_even |
---|
83 | LOGICAL , INTENT(IN ), OPTIONAL :: ld_discont |
---|
84 | |
---|
85 | ! local variable |
---|
86 | INTEGER(i4), DIMENSION(4) :: il_shape |
---|
87 | |
---|
88 | LOGICAL :: ll_discont |
---|
89 | |
---|
90 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_weight_IJ |
---|
91 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_weight_I |
---|
92 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_weight_J |
---|
93 | |
---|
94 | ! loop indices |
---|
95 | INTEGER(i4) :: ji |
---|
96 | INTEGER(i4) :: jj |
---|
97 | INTEGER(i4) :: jk |
---|
98 | INTEGER(i4) :: jl |
---|
99 | !---------------------------------------------------------------- |
---|
100 | ll_discont=.FALSE. |
---|
101 | IF( PRESENT(ld_discont) ) ll_discont=ld_discont |
---|
102 | |
---|
103 | il_shape(:)=SHAPE(dd_value) |
---|
104 | |
---|
105 | ! compute vect2D |
---|
106 | ALLOCATE(dl_weight_IJ(16,((id_rho(jp_I)+1)*(id_rho(jp_J)+1))) ) |
---|
107 | CALL interp_linear__get_weight2D(dl_weight_IJ(:,:), & |
---|
108 | & id_rho(:), ld_even(:)) |
---|
109 | |
---|
110 | ALLOCATE( dl_weight_I( 4,((id_rho(jp_I)+1) )) ) |
---|
111 | ALLOCATE( dl_weight_J( 4,( (id_rho(jp_J)+1))) ) |
---|
112 | CALL interp_linear__get_weight1D(dl_weight_I(:,:), & |
---|
113 | & id_rho(jp_I), ld_even(jp_I)) |
---|
114 | CALL interp_linear__get_weight1D(dl_weight_J(:,:), & |
---|
115 | & id_rho(jp_J), ld_even(jp_J)) |
---|
116 | |
---|
117 | DO jl=1,il_shape(4) |
---|
118 | ! loop on vertical level |
---|
119 | DO jk=1,il_shape(3) |
---|
120 | |
---|
121 | ! I-J plan |
---|
122 | CALL interp_linear__2D(dd_value(:,:,jk,jl), dd_fill,& |
---|
123 | & id_detect(:,:,jk), & |
---|
124 | & dl_weight_IJ(:,:), & |
---|
125 | & id_rho(jp_I), id_rho(jp_J), & |
---|
126 | & ll_discont) |
---|
127 | IF( ANY(id_detect(:,:,jk)==1) )THEN |
---|
128 | ! I direction |
---|
129 | DO jj=1,il_shape(2) |
---|
130 | CALL interp_linear__1D( dd_value(:,jj,jk,jl), dd_fill,& |
---|
131 | & id_detect(:,jj,jk), & |
---|
132 | & dl_weight_I(:,:), & |
---|
133 | & id_rho(jp_I), ll_discont ) |
---|
134 | ENDDO |
---|
135 | IF( ALL(id_detect(:,:,jk)==0) )THEN |
---|
136 | CYCLE |
---|
137 | ELSE |
---|
138 | ! J direction |
---|
139 | DO ji=1,il_shape(1) |
---|
140 | CALL interp_linear__1D( dd_value(ji,:,jk,jl), dd_fill,& |
---|
141 | & id_detect(ji,:,jk), & |
---|
142 | & dl_weight_J(:,:), & |
---|
143 | & id_rho(jp_J), ll_discont ) |
---|
144 | ENDDO |
---|
145 | ENDIF |
---|
146 | ENDIF |
---|
147 | |
---|
148 | ENDDO |
---|
149 | ENDDO |
---|
150 | |
---|
151 | DEALLOCATE(dl_weight_IJ) |
---|
152 | DEALLOCATE(dl_weight_I) |
---|
153 | DEALLOCATE(dl_weight_J) |
---|
154 | |
---|
155 | END SUBROUTINE interp_linear_fill |
---|
156 | !------------------------------------------------------------------- |
---|
157 | !> @brief |
---|
158 | !> This subroutine compute linear interpolation on 2D array of value. |
---|
159 | !> |
---|
160 | !> @details |
---|
161 | !> |
---|
162 | !> @author J.Paul |
---|
163 | !> - September, 2014- Initial Version |
---|
164 | !> |
---|
165 | !> @param[inout] dd_value 2D array of variable value |
---|
166 | !> @param[in] dd_fill FillValue of variable |
---|
167 | !> @param[inout] id_detect 2D array of point to be interpolated |
---|
168 | !> @param[in] id_rhoi refinment factor in i-direction |
---|
169 | !> @param[in] id_rhoj refinment factor in j-direction |
---|
170 | !> @param[in] id_rhok refinment factor in k-direction |
---|
171 | !> @param[in] ld_even even refinment or not |
---|
172 | !> @param[in] ld_discont longitudinal discontinuity (-180°/180°, 0°/360°) or not |
---|
173 | !------------------------------------------------------------------- |
---|
174 | SUBROUTINE interp_linear__2D( dd_value, dd_fill, & |
---|
175 | & id_detect, & |
---|
176 | & dd_weight, & |
---|
177 | & id_rhoi, id_rhoj, & |
---|
178 | & ld_discont ) |
---|
179 | |
---|
180 | IMPLICIT NONE |
---|
181 | ! Argument |
---|
182 | REAL(dp) , DIMENSION(:,:), INTENT(INOUT) :: dd_value |
---|
183 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
184 | INTEGER(I4) , DIMENSION(:,:), INTENT(INOUT) :: id_detect |
---|
185 | REAL(dp) , DIMENSION(:,:), INTENT(IN ) :: dd_weight |
---|
186 | INTEGER(I4) , INTENT(IN ) :: id_rhoi |
---|
187 | INTEGER(I4) , INTENT(IN ) :: id_rhoj |
---|
188 | LOGICAL , INTENT(IN ) :: ld_discont |
---|
189 | |
---|
190 | ! local variable |
---|
191 | INTEGER(I4) :: il_xextra |
---|
192 | INTEGER(I4) :: il_yextra |
---|
193 | INTEGER(i4), DIMENSION(2) :: il_shape |
---|
194 | INTEGER(i4), DIMENSION(2) :: il_dim |
---|
195 | |
---|
196 | REAL(dp) :: dl_min |
---|
197 | REAL(dp) :: dl_max |
---|
198 | REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_coef |
---|
199 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_coarse |
---|
200 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_tmp |
---|
201 | |
---|
202 | ! loop indices |
---|
203 | INTEGER(i4) :: ji |
---|
204 | INTEGER(i4) :: jj |
---|
205 | INTEGER(i4) :: ii |
---|
206 | INTEGER(i4) :: ij |
---|
207 | |
---|
208 | !---------------------------------------------------------------- |
---|
209 | |
---|
210 | IF( ANY(id_detect(:,:)==1) )THEN |
---|
211 | il_shape(:)=SHAPE(dd_value) |
---|
212 | |
---|
213 | ! compute coarse grid dimension |
---|
214 | il_xextra=id_rhoi-1 |
---|
215 | il_dim(1)=(il_shape(1)+il_xextra)/id_rhoi |
---|
216 | |
---|
217 | il_yextra=id_rhoj-1 |
---|
218 | il_dim(2)=(il_shape(2)+il_yextra)/id_rhoj |
---|
219 | |
---|
220 | ALLOCATE( dl_coarse(il_dim(1),il_dim(2)) ) |
---|
221 | |
---|
222 | ! value on coarse grid |
---|
223 | dl_coarse(:,:)=dd_value( 1:il_shape(1):id_rhoi, & |
---|
224 | & 1:il_shape(2):id_rhoj ) |
---|
225 | |
---|
226 | ALLOCATE( dl_tmp(2,2) ) |
---|
227 | ALLOCATE( dl_coef(4) ) |
---|
228 | |
---|
229 | DO jj=1,il_shape(2)-1,id_rhoj |
---|
230 | ij=((jj-1)/id_rhoj)+1 |
---|
231 | DO ji=1,il_shape(1)-1,id_rhoi |
---|
232 | ii=((ji-1)/id_rhoi)+1 |
---|
233 | |
---|
234 | ! check if point to be interpolated |
---|
235 | IF( ALL(id_detect(ji:ji+id_rhoi, & |
---|
236 | & jj:jj+id_rhoj)==0) ) CYCLE |
---|
237 | ! check data to needed to interpolate |
---|
238 | IF( ANY(dl_coarse(ii:ii+1,ij:ij+1)==dd_fill) ) CYCLE |
---|
239 | ! check longitude discontinuity |
---|
240 | dl_tmp(:,:)=dl_coarse(ii:ii+1,ij:ij+1) |
---|
241 | IF( ld_discont )THEN |
---|
242 | |
---|
243 | dl_min=MINVAL( dl_tmp(:,:), dl_tmp(:,:)/=dd_fill ) |
---|
244 | dl_max=MAXVAL( dl_tmp(:,:), dl_tmp(:,:)/=dd_fill ) |
---|
245 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
246 | WHERE( dl_tmp(:,:) < 0_dp ) |
---|
247 | dl_tmp(:,:) = dl_tmp(:,:)+360._dp |
---|
248 | END WHERE |
---|
249 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
250 | WHERE( dl_tmp(:,:) > 180_dp ) |
---|
251 | dl_tmp(:,:) = dl_tmp(:,:)-180._dp |
---|
252 | END WHERE |
---|
253 | ENDIF |
---|
254 | |
---|
255 | ENDIF |
---|
256 | |
---|
257 | ! compute bilinear coefficient |
---|
258 | dl_coef(:)=interp_linear__2D_coef(dl_tmp(:,:),& |
---|
259 | & dd_fill ) |
---|
260 | |
---|
261 | ! compute value on detetected point |
---|
262 | CALL interp_linear__2D_fill(dd_value( ji:ji+id_rhoi, & |
---|
263 | & jj:jj+id_rhoj ), & |
---|
264 | & id_detect(ji:ji+id_rhoi, & |
---|
265 | & jj:jj+id_rhoj ), & |
---|
266 | & dd_weight(:,:), dl_coef(:),& |
---|
267 | & dd_fill, id_rhoi, id_rhoj ) |
---|
268 | |
---|
269 | IF( ld_discont )THEN |
---|
270 | WHERE( dd_value( ji:ji+id_rhoi, & |
---|
271 | & jj:jj+id_rhoj ) >= 180._dp .AND. & |
---|
272 | & dd_value( ji:ji+id_rhoi, & |
---|
273 | & jj:jj+id_rhoj ) /= dd_fill ) |
---|
274 | dd_value( ji:ji+id_rhoi, & |
---|
275 | & jj:jj+id_rhoj ) = & |
---|
276 | & dd_value( ji:ji+id_rhoi, & |
---|
277 | & jj:jj+id_rhoj ) - 360._dp |
---|
278 | END WHERE |
---|
279 | ENDIF |
---|
280 | |
---|
281 | ENDDO |
---|
282 | ENDDO |
---|
283 | |
---|
284 | DEALLOCATE(dl_coef) |
---|
285 | DEALLOCATE(dl_tmp ) |
---|
286 | |
---|
287 | DEALLOCATE( dl_coarse ) |
---|
288 | ENDIF |
---|
289 | |
---|
290 | END SUBROUTINE interp_linear__2D |
---|
291 | !------------------------------------------------------------------- |
---|
292 | !> @brief |
---|
293 | !> This subroutine compute linear interpolation on 1D array of value. |
---|
294 | !> |
---|
295 | !> @details |
---|
296 | !> |
---|
297 | !> @author J.Paul |
---|
298 | !> - September, 2014- Initial Version |
---|
299 | !> |
---|
300 | !> @param[inout] dd_value 1D array of variable value |
---|
301 | !> @param[in] dd_fill FillValue of variable |
---|
302 | !> @param[inout] id_detect 1D array of point to be interpolated |
---|
303 | !> @param[in] id_rhoi refinment factor |
---|
304 | !> @param[in] ld_even even refinment or not |
---|
305 | !> @param[in] ld_discont longitudinal discontinuity (-180°/180°, 0°/360°) or not |
---|
306 | !------------------------------------------------------------------- |
---|
307 | SUBROUTINE interp_linear__1D( dd_value, dd_fill, & |
---|
308 | & id_detect, & |
---|
309 | & dd_weight, & |
---|
310 | & id_rhoi, & |
---|
311 | & ld_discont ) |
---|
312 | |
---|
313 | IMPLICIT NONE |
---|
314 | ! Argument |
---|
315 | REAL(dp) , DIMENSION(:) , INTENT(INOUT) :: dd_value |
---|
316 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
317 | INTEGER(I4) , DIMENSION(:) , INTENT(INOUT) :: id_detect |
---|
318 | REAL(dp) , DIMENSION(:,:), INTENT(IN ) :: dd_weight |
---|
319 | INTEGER(I4) , INTENT(IN ) :: id_rhoi |
---|
320 | LOGICAL , INTENT(IN ) :: ld_discont |
---|
321 | |
---|
322 | ! local variable |
---|
323 | INTEGER(I4) :: il_xextra |
---|
324 | INTEGER(i4), DIMENSION(1) :: il_shape |
---|
325 | INTEGER(i4), DIMENSION(1) :: il_dim |
---|
326 | |
---|
327 | REAL(dp) :: dl_min |
---|
328 | REAL(dp) :: dl_max |
---|
329 | REAL(dp) , DIMENSION(:), ALLOCATABLE :: dl_coef |
---|
330 | REAL(dp) , DIMENSION(:), ALLOCATABLE :: dl_coarse |
---|
331 | REAL(dp) , DIMENSION(:), ALLOCATABLE :: dl_tmp |
---|
332 | |
---|
333 | ! loop indices |
---|
334 | INTEGER(i4) :: ji |
---|
335 | INTEGER(i4) :: ii |
---|
336 | |
---|
337 | !---------------------------------------------------------------- |
---|
338 | |
---|
339 | IF( ANY(id_detect(:)==1) )THEN |
---|
340 | il_shape(:)=SHAPE(dd_value) |
---|
341 | |
---|
342 | ! compute coarse grid dimension |
---|
343 | il_xextra=id_rhoi-1 |
---|
344 | il_dim(1)=(il_shape(1)+il_xextra)/id_rhoi |
---|
345 | |
---|
346 | ALLOCATE( dl_coarse(il_dim(1)) ) |
---|
347 | |
---|
348 | ! value on coarse grid |
---|
349 | dl_coarse(:)=dd_value( 1:il_shape(1):id_rhoi ) |
---|
350 | |
---|
351 | ALLOCATE( dl_tmp(2) ) |
---|
352 | ALLOCATE( dl_coef(4) ) |
---|
353 | |
---|
354 | DO ji=1,il_shape(1)-1,id_rhoi |
---|
355 | ii=((ji-1)/id_rhoi)+1 |
---|
356 | |
---|
357 | ! check if point to be interpolated |
---|
358 | IF( ALL(id_detect(ji:ji+id_rhoi)==0) ) CYCLE |
---|
359 | ! check data needed to interpolate |
---|
360 | IF( ANY(dl_coarse(ii:ii+1)==dd_fill) ) CYCLE |
---|
361 | ! check longitude discontinuity |
---|
362 | dl_tmp(:)=dl_coarse(ii:ii+1) |
---|
363 | IF( ld_discont )THEN |
---|
364 | |
---|
365 | dl_min=MINVAL( dl_tmp(:), dl_tmp(:)/=dd_fill ) |
---|
366 | dl_max=MAXVAL( dl_tmp(:), dl_tmp(:)/=dd_fill ) |
---|
367 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
368 | WHERE( dl_tmp(:) < 0_dp ) |
---|
369 | dl_tmp(:) = dl_tmp(:)+360._dp |
---|
370 | END WHERE |
---|
371 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
372 | WHERE( dl_tmp(:) > 180_dp ) |
---|
373 | dl_tmp(:) = dl_tmp(:)-180._dp |
---|
374 | END WHERE |
---|
375 | ENDIF |
---|
376 | |
---|
377 | ENDIF |
---|
378 | |
---|
379 | ! compute bilinear coefficient |
---|
380 | dl_coef(:)=interp_linear__1D_coef(dl_tmp(:), & |
---|
381 | & dd_fill ) |
---|
382 | |
---|
383 | ! compute value on detetected point |
---|
384 | CALL interp_linear__1D_fill( dd_value( ji:ji+id_rhoi ), & |
---|
385 | & id_detect(ji:ji+id_rhoi ), & |
---|
386 | & dd_weight(:,:), dl_coef(:),& |
---|
387 | & dd_fill, id_rhoi ) |
---|
388 | |
---|
389 | IF( ld_discont )THEN |
---|
390 | WHERE( dd_value( ji:ji+id_rhoi ) >= 180._dp .AND. & |
---|
391 | & dd_value( ji:ji+id_rhoi ) /= dd_fill ) |
---|
392 | dd_value(ji:ji+id_rhoi) = dd_value(ji:ji+id_rhoi) - 360._dp |
---|
393 | END WHERE |
---|
394 | ENDIF |
---|
395 | |
---|
396 | ENDDO |
---|
397 | |
---|
398 | DEALLOCATE(dl_coef) |
---|
399 | DEALLOCATE(dl_tmp ) |
---|
400 | |
---|
401 | DEALLOCATE( dl_coarse ) |
---|
402 | ENDIF |
---|
403 | |
---|
404 | END SUBROUTINE interp_linear__1D |
---|
405 | !------------------------------------------------------------------- |
---|
406 | !> @brief |
---|
407 | !> This subroutine compute 2D array of coefficient for linear interpolation. |
---|
408 | !> |
---|
409 | !> @author J.Paul |
---|
410 | !> - September, 2014- Initial Version |
---|
411 | !> |
---|
412 | !> @param[in] dd_value 2D array of value |
---|
413 | !> @param[in] dd_fill FillValue of variable |
---|
414 | !------------------------------------------------------------------- |
---|
415 | FUNCTION interp_linear__2D_coef( dd_value, dd_fill ) |
---|
416 | IMPLICIT NONE |
---|
417 | ! Argument |
---|
418 | REAL(dp), DIMENSION(:,:) , INTENT(IN) :: dd_value |
---|
419 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
420 | |
---|
421 | ! function |
---|
422 | REAL(dp), DIMENSION(4) :: interp_linear__2D_coef |
---|
423 | |
---|
424 | ! local variable |
---|
425 | REAL(dp), DIMENSION(4,4), PARAMETER :: dl_matrix = RESHAPE( & |
---|
426 | & (/ 1 ,-1 ,-1 , 1 ,& |
---|
427 | 0 , 1 , 0 ,-1 ,& |
---|
428 | 0 , 0 , 1 ,-1 ,& |
---|
429 | 0 , 0 , 0 , 1 /), & |
---|
430 | & (/ 4, 4 /) ) |
---|
431 | |
---|
432 | REAL(dp), DIMENSION(4) :: dl_vect |
---|
433 | |
---|
434 | !---------------------------------------------------------------- |
---|
435 | ! init |
---|
436 | interp_linear__2D_coef(:)=dd_fill |
---|
437 | |
---|
438 | dl_vect( 1: 4)=PACK(dd_value(:,:),.TRUE. ) |
---|
439 | interp_linear__2D_coef(:)=MATMUL(dl_matrix(:,:),dl_vect(:)) |
---|
440 | |
---|
441 | END FUNCTION interp_linear__2D_coef |
---|
442 | !------------------------------------------------------------------- |
---|
443 | !> @brief |
---|
444 | !> This subroutine compute linear interpolation of a 2D array of value. |
---|
445 | !> |
---|
446 | !> @author J.Paul |
---|
447 | !> - September, 2014- Initial Version |
---|
448 | !> |
---|
449 | !> @param[inout] dd_value 2D array of mixed grid value |
---|
450 | !> @param[inout] id_detect 2D array of point to be interpolated |
---|
451 | !> @param[in] dd_coef 2D array of coefficient |
---|
452 | !> @param[in] dd_fill FillValue of variable |
---|
453 | !> @param[in] ld_even even refinment or not |
---|
454 | !> @param[in] id_rhoi refinement factor in i-direction |
---|
455 | !> @param[in] id_rhoj refinement factor in j-direction |
---|
456 | !------------------------------------------------------------------- |
---|
457 | SUBROUTINE interp_linear__2D_fill( dd_value, id_detect, & |
---|
458 | & dd_weight, dd_coef, & |
---|
459 | & dd_fill, id_rhoi, id_rhoj ) |
---|
460 | IMPLICIT NONE |
---|
461 | ! Argument |
---|
462 | REAL(dp) , DIMENSION(:,:), INTENT(INOUT) :: dd_value |
---|
463 | INTEGER(i4) , DIMENSION(:,:), INTENT(INOUT) :: id_detect |
---|
464 | REAL(dp) , DIMENSION(:,:), INTENT(IN ) :: dd_weight |
---|
465 | REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_coef |
---|
466 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
467 | INTEGER(I4) , INTENT(IN ) :: id_rhoi |
---|
468 | INTEGER(I4) , INTENT(IN ) :: id_rhoj |
---|
469 | |
---|
470 | ! local variable |
---|
471 | |
---|
472 | ! loop indices |
---|
473 | INTEGER(i4) :: ii |
---|
474 | |
---|
475 | INTEGER(i4) :: ji |
---|
476 | INTEGER(i4) :: jj |
---|
477 | !---------------------------------------------------------------- |
---|
478 | |
---|
479 | IF( ANY( dd_coef(:)==dd_fill ) )THEN |
---|
480 | CALL logger_error("INTERP LINEAR FILL: fill value detected in coef. "//& |
---|
481 | & "can not compute interpolation.") |
---|
482 | ELSE |
---|
483 | |
---|
484 | ii=0 |
---|
485 | DO jj=1,id_rhoj+1 |
---|
486 | DO ji=1,id_rhoi+1 |
---|
487 | |
---|
488 | ii=ii+1 |
---|
489 | IF(id_detect(ji,jj)==1)THEN |
---|
490 | |
---|
491 | dd_value(ji,jj)=DOT_PRODUCT(dd_coef(:),dd_weight(:,ii)) |
---|
492 | id_detect(ji,jj)=0 |
---|
493 | |
---|
494 | ENDIF |
---|
495 | |
---|
496 | ENDDO |
---|
497 | ENDDO |
---|
498 | |
---|
499 | ENDIF |
---|
500 | |
---|
501 | END SUBROUTINE interp_linear__2D_fill |
---|
502 | !------------------------------------------------------------------- |
---|
503 | !> @brief |
---|
504 | !> This subroutine compute 1D array of coefficient for linear interpolation. |
---|
505 | !> |
---|
506 | !> @author J.Paul |
---|
507 | !> - September, 2014- Initial Version |
---|
508 | !> |
---|
509 | !> @param[in] dd_value 1D array of value |
---|
510 | !> @param[in] dd_fill FillValue of variable |
---|
511 | !------------------------------------------------------------------- |
---|
512 | FUNCTION interp_linear__1D_coef( dd_value, dd_fill ) |
---|
513 | IMPLICIT NONE |
---|
514 | ! Argument |
---|
515 | REAL(dp), DIMENSION(:) , INTENT(IN) :: dd_value |
---|
516 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
517 | |
---|
518 | ! function |
---|
519 | REAL(dp), DIMENSION(2) :: interp_linear__1D_coef |
---|
520 | |
---|
521 | ! local variable |
---|
522 | REAL(dp), DIMENSION(2,2), PARAMETER :: dl_matrix = RESHAPE( & |
---|
523 | & (/ 1 ,-1 ,& |
---|
524 | 0 , 1 /), & |
---|
525 | & (/ 2, 2 /) ) |
---|
526 | |
---|
527 | REAL(dp), DIMENSION(2) :: dl_vect |
---|
528 | |
---|
529 | !---------------------------------------------------------------- |
---|
530 | ! init |
---|
531 | interp_linear__1D_coef(:)=dd_fill |
---|
532 | |
---|
533 | dl_vect( 1: 2)=PACK(dd_value(:),.TRUE. ) |
---|
534 | interp_linear__1D_coef(:)=MATMUL(dl_matrix(:,:),dl_vect(:)) |
---|
535 | |
---|
536 | END FUNCTION interp_linear__1D_coef |
---|
537 | !------------------------------------------------------------------- |
---|
538 | !> @brief |
---|
539 | !> This subroutine compute linear interpolation of a 1D array of value. |
---|
540 | !> |
---|
541 | !> @author J.Paul |
---|
542 | !> - September, 2014- Initial Version |
---|
543 | !> |
---|
544 | !> @param[inout] dd_value 1D array of mixed grid value |
---|
545 | !> @param[inout] id_detect 1D array of point to be interpolated |
---|
546 | !> @param[in] dd_coef 1D array of coefficient |
---|
547 | !> @param[in] dd_fill FillValue of variable |
---|
548 | !> @param[in] ld_even even refinment or not |
---|
549 | !> @param[in] id_rho refinement factor |
---|
550 | !------------------------------------------------------------------- |
---|
551 | SUBROUTINE interp_linear__1D_fill( dd_value, id_detect, & |
---|
552 | & dd_weight, dd_coef, & |
---|
553 | & dd_fill, id_rho ) |
---|
554 | IMPLICIT NONE |
---|
555 | ! Argument |
---|
556 | REAL(dp) , DIMENSION(:) , INTENT(INOUT) :: dd_value |
---|
557 | INTEGER(i4) , DIMENSION(:) , INTENT(INOUT) :: id_detect |
---|
558 | REAL(dp) , DIMENSION(:,:), INTENT(IN ) :: dd_weight |
---|
559 | REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_coef |
---|
560 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
561 | INTEGER(I4) , INTENT(IN ) :: id_rho |
---|
562 | |
---|
563 | ! local variable |
---|
564 | |
---|
565 | ! loop indices |
---|
566 | INTEGER(i4) :: ji |
---|
567 | !---------------------------------------------------------------- |
---|
568 | |
---|
569 | IF( ANY( dd_coef(:)==dd_fill ) )THEN |
---|
570 | CALL logger_error("INTERP LINEAR FILL: fill value detected. "//& |
---|
571 | & "can not compute interpolation") |
---|
572 | ELSE |
---|
573 | |
---|
574 | DO ji=1,id_rho+1 |
---|
575 | |
---|
576 | IF(id_detect(ji)==1)THEN |
---|
577 | |
---|
578 | dd_value(ji)=DOT_PRODUCT(dd_coef(:),dd_weight(:,ji)) |
---|
579 | id_detect(ji)=0 |
---|
580 | |
---|
581 | ENDIF |
---|
582 | |
---|
583 | ENDDO |
---|
584 | |
---|
585 | ENDIF |
---|
586 | |
---|
587 | END SUBROUTINE interp_linear__1D_fill |
---|
588 | !------------------------------------------------------------------- |
---|
589 | !> @brief |
---|
590 | !> This subroutine compute interpoaltion weight for 2D array. |
---|
591 | !> |
---|
592 | !> @author J.Paul |
---|
593 | !> - September, 2014- Initial Version |
---|
594 | !> |
---|
595 | !> @param[in] dd_weight interpolation weight of 2D array |
---|
596 | !> @param[in] ld_even even refinment or not |
---|
597 | !> @param[in] id_rho refinement factor |
---|
598 | !------------------------------------------------------------------- |
---|
599 | SUBROUTINE interp_linear__get_weight2D(dd_weight, & |
---|
600 | & id_rho, ld_even) |
---|
601 | IMPLICIT NONE |
---|
602 | ! Argument |
---|
603 | REAL(dp) , DIMENSION(:,:), INTENT(INOUT) :: dd_weight |
---|
604 | INTEGER(I4), DIMENSION(:) , INTENT(IN ) :: id_rho |
---|
605 | LOGICAL , DIMENSION(:) , INTENT(IN ) :: ld_even |
---|
606 | ! local variable |
---|
607 | REAL(dp) :: dl_dx |
---|
608 | REAL(dp) :: dl_x |
---|
609 | REAL(dp) :: dl_dy |
---|
610 | REAL(dp) :: dl_y |
---|
611 | |
---|
612 | ! loop indices |
---|
613 | INTEGER(i4) :: ii |
---|
614 | |
---|
615 | INTEGER(i4) :: ji |
---|
616 | INTEGER(i4) :: jj |
---|
617 | !---------------------------------------------------------------- |
---|
618 | |
---|
619 | IF( ld_even(jp_I) )THEN |
---|
620 | dl_dx=1./REAL(id_rho(jp_I)-1) |
---|
621 | ELSE ! odd refinement |
---|
622 | dl_dx=1./REAL(id_rho(jp_I)) |
---|
623 | ENDIF |
---|
624 | |
---|
625 | IF( ld_even(jp_J) )THEN |
---|
626 | dl_dy=1./REAL(id_rho(jp_J)-1) |
---|
627 | ELSE ! odd refinement |
---|
628 | dl_dy=1./REAL(id_rho(jp_J)) |
---|
629 | ENDIF |
---|
630 | |
---|
631 | ii=0 |
---|
632 | DO jj=1,id_rho(jp_J)+1 |
---|
633 | |
---|
634 | IF( ld_even(jp_J) )THEN |
---|
635 | dl_y=(jj-1)*dl_dy - dl_dy*0.5 |
---|
636 | ELSE ! odd refinement |
---|
637 | dl_y=(jj-1)*dl_dy |
---|
638 | ENDIF |
---|
639 | |
---|
640 | DO ji=1,id_rho(jp_I)+1 |
---|
641 | |
---|
642 | ! iter |
---|
643 | ii=ii+1 |
---|
644 | |
---|
645 | IF( ld_even(jp_I) )THEN |
---|
646 | dl_x=(ji-1)*dl_dx - dl_dx*0.5 |
---|
647 | ELSE ! odd refinement |
---|
648 | dl_x=(ji-1)*dl_dx |
---|
649 | ENDIF |
---|
650 | |
---|
651 | dd_weight(:,ii)=(/1._dp, dl_x, dl_y, dl_x*dl_y /) |
---|
652 | |
---|
653 | ENDDO |
---|
654 | ENDDO |
---|
655 | |
---|
656 | END SUBROUTINE interp_linear__get_weight2D |
---|
657 | !------------------------------------------------------------------- |
---|
658 | !> @brief |
---|
659 | !> This subroutine compute interpoaltion weight for 1D array. |
---|
660 | !> |
---|
661 | !> @author J.Paul |
---|
662 | !> - September, 2014- Initial Version |
---|
663 | !> |
---|
664 | !> @param[in] dd_weight interpolation weight of 1D array |
---|
665 | !> @param[in] ld_even even refinment or not |
---|
666 | !> @param[in] id_rho refinement factor |
---|
667 | !------------------------------------------------------------------- |
---|
668 | SUBROUTINE interp_linear__get_weight1D(dd_weight, & |
---|
669 | & id_rho, ld_even) |
---|
670 | IMPLICIT NONE |
---|
671 | ! Argument |
---|
672 | REAL(dp) , DIMENSION(:,:), INTENT(INOUT) :: dd_weight |
---|
673 | INTEGER(I4) , INTENT(IN ) :: id_rho |
---|
674 | LOGICAL , INTENT(IN ) :: ld_even |
---|
675 | |
---|
676 | ! local variable |
---|
677 | REAL(dp) :: dl_dx |
---|
678 | REAL(dp) :: dl_x |
---|
679 | |
---|
680 | ! loop indices |
---|
681 | INTEGER(i4) :: ji |
---|
682 | !---------------------------------------------------------------- |
---|
683 | |
---|
684 | IF( ld_even )THEN |
---|
685 | dl_dx=1./REAL(id_rho-1) |
---|
686 | ELSE ! odd refinement |
---|
687 | dl_dx=1./REAL(id_rho) |
---|
688 | ENDIF |
---|
689 | |
---|
690 | DO ji=1,id_rho+1 |
---|
691 | IF( ld_even )THEN |
---|
692 | dl_x=(ji-1)*dl_dx - dl_dx*0.5 |
---|
693 | ELSE ! odd refinement |
---|
694 | dl_x=(ji-1)*dl_dx |
---|
695 | ENDIF |
---|
696 | |
---|
697 | dd_weight(:,ji)=(/1._dp, dl_x /) |
---|
698 | ENDDO |
---|
699 | |
---|
700 | END SUBROUTINE interp_linear__get_weight1D |
---|
701 | END MODULE interp_linear |
---|