New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
interp_nearest.f90 in branches/UKMO/r5518_INGV1_WAVE-coupling/NEMOGCM/TOOLS/SIREN/src – NEMO

source: branches/UKMO/r5518_INGV1_WAVE-coupling/NEMOGCM/TOOLS/SIREN/src/interp_nearest.f90 @ 7152

Last change on this file since 7152 was 7152, checked in by jcastill, 7 years ago

Initial implementation of wave coupling branch - INGV wave branch + UKMO wave coupling branch

File size: 11.5 KB
Line 
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 nearest interpolation on regular grid.
10!>
11!> @details
12!> to compute nearest interpolation:<br/>
13!> @code
14!> CALL interp_nearest_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!----------------------------------------------------------------------
30MODULE interp_nearest
31
32   USE netcdf                          ! nf90 library
33   USE global                          ! global variable
34   USE kind                            ! F90 kind parameter
35   USE logger                          ! log file manager
36   USE fct                             ! basic useful function
37
38   IMPLICIT NONE
39   ! NOTE_avoid_public_variables_if_possible
40
41   ! type and variable
42
43   ! function and subroutine
44   PUBLIC :: interp_nearest_fill  !< compute interpolation using nearest method
45
46   PRIVATE :: interp_nearest__2D       !< compute binearest interpolation on 2D gid
47   PRIVATE :: interp_nearest__1D       !< compute   nearest interpolation on 1D gid
48   PRIVATE :: interp_nearest__2D_fill  !< fill value using binearest interpolation
49   PRIVATE :: interp_nearest__1D_fill  !< fill value using   nearest interpolation
50
51CONTAINS   
52   !-------------------------------------------------------------------
53   !> @brief
54   !> This subroutine compute horizontal nearest interpolation on 4D array of value.
55   !>
56   !> @author J.Paul
57   !> @date September, 2014 - Initial Version
58   !>
59   !> @param[inout] dd_value  2D array of variable value
60   !> @param[inout] id_detect 2D array of point to be interpolated
61   !> @param[in]    id_rho    array of refinment factor
62   !-------------------------------------------------------------------
63   SUBROUTINE interp_nearest_fill(dd_value, id_detect, id_rho )
64      IMPLICIT NONE
65      ! Argument
66      REAL(dp)        , DIMENSION(:,:,:,:), INTENT(INOUT) :: dd_value 
67      INTEGER(I4)     , DIMENSION(:,:,:)  , INTENT(INOUT) :: id_detect
68      INTEGER(I4)     , DIMENSION(:)      , INTENT(IN   ) :: id_rho
69
70      ! local variable
71      INTEGER(i4), DIMENSION(4)                  :: il_shape
72
73      INTEGER(I4), DIMENSION(:,:,:), ALLOCATABLE :: il_detect
74
75      ! loop indices
76      INTEGER(i4) :: ji
77      INTEGER(i4) :: jj
78      INTEGER(i4) :: jk
79      INTEGER(i4) :: jl
80      !----------------------------------------------------------------
81
82      il_shape(:)=SHAPE(dd_value)
83
84      ALLOCATE(il_detect(il_shape(1),il_shape(2),il_shape(3)))
85      DO jl=1,il_shape(4)
86         il_detect(:,:,:)=id_detect(:,:,:)
87         ! loop on vertical level
88         DO jk=1,il_shape(3)
89
90            ! I-J plan
91            CALL interp_nearest__2D(dd_value(:,:,jk,jl),&
92            &                       il_detect(:,:,jk),  &
93            &                       id_rho(jp_I), id_rho(jp_J) )           
94            IF( ANY(il_detect(:,:,jk)==1) )THEN
95               ! I direction
96               DO jj=1,il_shape(2)
97                  CALL interp_nearest__1D( dd_value(:,jj,jk,jl),&
98                  &                        il_detect(:,jj,jk),  &
99                  &                        id_rho(jp_I) )
100               ENDDO
101               IF( ALL(il_detect(:,:,jk)==0) )THEN
102                  CYCLE
103               ELSE
104                  ! J direction
105                  DO ji=1,il_shape(1)
106                     CALL interp_nearest__1D( dd_value(ji,:,jk,jl),&
107                     &                        il_detect(ji,:,jk),  &
108                     &                        id_rho(jp_J) )
109                  ENDDO
110               ENDIF
111            ENDIF
112
113         ENDDO
114      ENDDO
115
116      id_detect(:,:,:)=il_detect(:,:,:)
117      DEALLOCATE(il_detect)
118
119   END SUBROUTINE interp_nearest_fill
120   !-------------------------------------------------------------------
121   !> @brief
122   !> This subroutine compute nearest interpolation on 2D array of value.
123   !>
124   !> @author J.Paul
125   !> @date September, 2014 - Initial Version
126   !>
127   !> @param[inout] dd_value  2D array of variable value
128   !> @param[inout] id_detect 2D array of point to be interpolated
129   !> @param[in] id_rhoi      refinment factor in i-direction
130   !> @param[in] id_rhoj      refinment factor in j-direction
131   !> @param[in] id_rhok      refinment factor in k-direction
132   !-------------------------------------------------------------------
133   SUBROUTINE interp_nearest__2D( dd_value, id_detect, &
134      &                           id_rhoi, id_rhoj )
135
136      IMPLICIT NONE
137      ! Argument
138      REAL(dp)        , DIMENSION(:,:), INTENT(INOUT) :: dd_value 
139      INTEGER(I4)     , DIMENSION(:,:), INTENT(INOUT) :: id_detect
140      INTEGER(I4)                     , INTENT(IN   ) :: id_rhoi
141      INTEGER(I4)                     , INTENT(IN   ) :: id_rhoj
142
143      ! local variable
144      INTEGER(i4), DIMENSION(2)                :: il_shape
145
146      ! loop indices
147      INTEGER(i4) :: ji
148      INTEGER(i4) :: jj
149
150      !----------------------------------------------------------------
151
152      IF( ANY(id_detect(:,:)==1) )THEN
153
154         il_shape(:)=SHAPE(dd_value)
155
156         DO jj=1,il_shape(2)-1,id_rhoj
157            DO ji=1,il_shape(1)-1,id_rhoi
158         
159               ! check if point to be interpolated
160               IF( ALL(id_detect(ji:ji+id_rhoi,   &
161               &                 jj:jj+id_rhoj)==0) ) CYCLE
162
163               ! compute value on detetected point
164               CALL interp_nearest__2D_fill(dd_value( ji:ji+id_rhoi,   &
165               &                                      jj:jj+id_rhoj ), &
166               &                            id_detect(ji:ji+id_rhoi,   &
167               &                                      jj:jj+id_rhoj ) )
168
169            ENDDO
170         ENDDO
171
172      ENDIF
173
174   END SUBROUTINE interp_nearest__2D
175   !-------------------------------------------------------------------
176   !> @brief
177   !> This subroutine compute nearest interpolation on 1D array of value.
178   !>
179   !> @author J.Paul
180   !> @date September, 2014 - Initial Version
181   !>
182   !> @param[inout] dd_value  1D array of variable value
183   !> @param[inout] id_detect 1D array of point to be interpolated
184   !> @param[in]    id_rhoi   refinment factor
185   !-------------------------------------------------------------------
186   SUBROUTINE interp_nearest__1D( dd_value,  id_detect, &
187      &                           id_rhoi )
188
189      IMPLICIT NONE
190      ! Argument
191      REAL(dp)        , DIMENSION(:), INTENT(INOUT) :: dd_value 
192      INTEGER(I4)     , DIMENSION(:), INTENT(INOUT) :: id_detect
193      INTEGER(I4)                   , INTENT(IN   ) :: id_rhoi
194
195      ! local variable
196      INTEGER(i4), DIMENSION(1)              :: il_shape
197
198      ! loop indices
199      INTEGER(i4) :: ji
200
201      !----------------------------------------------------------------
202
203      IF( ANY(id_detect(:)==1) )THEN
204         il_shape(:)=SHAPE(dd_value)
205
206         DO ji=1,il_shape(1)-1,id_rhoi
207         
208            ! check if point to be interpolated
209            IF( ALL(id_detect(ji:ji+id_rhoi)==0) ) CYCLE
210
211            ! compute value on detetected point
212            CALL interp_nearest__1D_fill( dd_value( ji:ji+id_rhoi ), &
213            &                             id_detect(ji:ji+id_rhoi ) )
214
215         ENDDO
216
217      ENDIF
218
219   END SUBROUTINE interp_nearest__1D
220   !-------------------------------------------------------------------
221   !> @brief
222   !> This subroutine compute nearest interpolation of a 2D array of value.
223   !>
224   !> @author J.Paul
225   !> @date September, 2014 - Initial Version
226   !>
227   !> @param[inout] dd_value  2D array of mixed grid value
228   !> @param[inout] id_detect 2D array of point to be interpolated
229   !-------------------------------------------------------------------
230   SUBROUTINE interp_nearest__2D_fill( dd_value, id_detect )
231      IMPLICIT NONE
232      ! Argument
233      REAL(dp)   , DIMENSION(:,:)  , INTENT(INOUT) :: dd_value 
234      INTEGER(i4), DIMENSION(:,:)  , INTENT(INOUT) :: id_detect
235
236      ! local variable
237      INTEGER(i4), DIMENSION(2) :: il_shape
238
239      INTEGER(i4) :: il_i1
240      INTEGER(i4) :: il_i2
241      INTEGER(i4) :: il_j1
242      INTEGER(i4) :: il_j2
243
244      INTEGER(i4) :: il_half1
245      INTEGER(i4) :: il_half2
246
247      ! loop indices
248      INTEGER(i4) :: ji
249      INTEGER(i4) :: jj
250      !----------------------------------------------------------------
251
252      il_shape(:)=SHAPE(dd_value(:,:))
253
254      il_i1=1
255      il_i2=il_shape(1)
256
257      il_j1=1
258      il_j2=il_shape(2)
259
260      il_half1=CEILING(il_shape(1)*0.5)
261      il_half2=CEILING(il_shape(2)*0.5)
262
263      DO jj=1,il_half2
264
265         DO ji=1,il_half1
266           
267            ! lower left point
268            IF(id_detect(ji,jj)==1)THEN
269
270               dd_value( ji,jj)=dd_value(il_i1,il_j1)
271               id_detect(ji,jj)=0
272
273            ENDIF
274
275            ! lower right point
276            IF(id_detect(il_shape(1)-ji+1,jj)==1)THEN
277
278               dd_value( il_shape(1)-ji+1,jj)=dd_value(il_i2,il_j1)
279               id_detect(il_shape(1)-ji+1,jj)=0
280
281            ENDIF
282
283            ! upper left point
284            IF(id_detect(ji,il_shape(2)-jj+1)==1)THEN
285
286               dd_value( ji,il_shape(2)-jj+1)=dd_value(il_i1,il_j2)
287               id_detect(ji,il_shape(2)-jj+1)=0
288
289            ENDIF           
290
291            ! upper right point
292            IF(id_detect(il_shape(1)-ji+1,il_shape(2)-jj+1)==1)THEN
293
294               dd_value( il_shape(1)-ji+1,il_shape(2)-jj+1)=dd_value(il_i2,il_j2)
295               id_detect(il_shape(1)-ji+1,il_shape(2)-jj+1)=0
296
297            ENDIF           
298
299         ENDDO
300
301      ENDDO
302
303   END SUBROUTINE interp_nearest__2D_fill
304   !-------------------------------------------------------------------
305   !> @brief
306   !> This subroutine compute nearest interpolation of a 1D array of value.
307   !>
308   !> @author J.Paul
309   !> @date September, 2014 - Initial Version
310   !>
311   !> @param[inout] dd_value  1D array of mixed grid value
312   !> @param[inout] id_detect 1D array of point to be interpolated
313   !-------------------------------------------------------------------
314   SUBROUTINE interp_nearest__1D_fill( dd_value, id_detect )
315      IMPLICIT NONE
316      ! Argument
317      REAL(dp)   , DIMENSION(:), INTENT(INOUT) :: dd_value 
318      INTEGER(i4), DIMENSION(:), INTENT(INOUT) :: id_detect
319
320      ! local variable
321      INTEGER(i4), DIMENSION(1) :: il_shape
322
323      INTEGER(i4) :: il_i1
324      INTEGER(i4) :: il_i2
325
326      INTEGER(i4) :: il_half1
327     
328      ! loop indices
329      INTEGER(i4) :: ji
330      !----------------------------------------------------------------
331
332      il_shape(:)=SHAPE(dd_value)
333
334      il_i1=1
335      il_i2=il_shape(1)
336
337      il_half1=CEILING(il_shape(1)*0.5)
338     
339      DO ji=1,il_half1
340         
341         ! lower left point
342         IF(id_detect(ji)==1)THEN
343
344            dd_value( ji)=dd_value(il_i1)
345            id_detect(ji)=0
346
347         ENDIF
348
349         ! lower right point
350         IF(id_detect(il_shape(1)-ji+1)==1)THEN
351
352            dd_value( il_shape(1)-ji+1)=dd_value(il_i2)
353            id_detect(il_shape(1)-ji+1)=0
354
355         ENDIF
356
357      ENDDO
358
359   END SUBROUTINE interp_nearest__1D_fill
360END MODULE interp_nearest
Note: See TracBrowser for help on using the repository browser.