SURFEX v8.1
General documentation of Surfex
bilin_gridin.F90
Go to the documentation of this file.
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SFX_LIC for details. version 1.
5 ! #########
6  SUBROUTINE bilin_gridin (PFIELD1,PFIELD_X,PFIELD_Y,PFIELD_XY)
7 ! #########################################################################
8 !
9 !!**** *BILIN_GRIDINEAR * - subroutine to interpolate surface FIELD
10 !!
11 !! PURPOSE
12 !! -------
13 !!
14 !!
15 !!** METHOD
16 !! ------
17 !!
18 !! Interpolation is bilinear, and uses 9 grid points, located in the
19 !! center of model 1 grid mesh, and at the boundaries of this grid
20 !! mesh (2 X limits, 2 Y limits and 4 corners).
21 !! This implies that the grid mesh values located around the model 1
22 !! grid mesh are not used directly. The values at the boundaries of the
23 !! grid mesh are defined by the average between the middle point
24 !! (this grid mesh value), and the one in the considered direction.
25 !! So the eight grid meshes around the considered grid mesh are used
26 !! equally.
27 !! This is important to note that these average values are erased
28 !! and replaced by zero if they are at the limit of any grid
29 !! mesh with the zero value. This allows to insure zero value in model 2
30 !! grid meshes where there was not the considered class in corresponding
31 !! model 1 grid mesh, and to insure continuity of the FIELD type
32 !! at such boundaries.
33 !!
34 !!
35 !! The arrays and array index are defined on the following (model1) grid:
36 !!
37 !!
38 !! XFIELD XFIELD XFIELD
39 !! * * *
40 !! i-1,j+1 i,j+1 i+1,j+1
41 !!
42 !!
43 !!
44 !! PFIELD_XY PFIELD_Y PFIELD_XY
45 !! * * *
46 !! i,j+1 i,j+1 i+1,j+1
47 !!
48 !!
49 !!
50 !! XFIELD PFIELD_X XFIELD PFIELD_X XFIELD
51 !! * * * * *
52 !! i-1,j i,j i,j i+1,j i+1,j
53 !!
54 !!
55 !!
56 !! PFIELD_XY PFIELD_Y PFIELD_XY
57 !! * * *
58 !! i,j i,j i+1,j
59 !!
60 !!
61 !!
62 !! XFIELD XFIELD XFIELD
63 !! * * *
64 !! i-1,j-1 i,j-1 i+1,j-1
65 !!
66 !!
67 !!
68 !!
69 !!
70 !! AUTHOR
71 !! ------
72 !!
73 !! V. Masson * METEO-FRANCE *
74 !!
75 !! MODIFICATIONS
76 !! -------------
77 !!
78 !! Original 01/2004
79 ! TD&DD: added OpenMP directives
80 
81 !-------------------------------------------------------------------------------
82 !
83 !* 0. DECLARATIONS
84 ! ------------
85 !
86 USE modd_surfex_mpi, ONLY : nrank
87 USE modd_surf_par, ONLY : xundef
88 !
89 USE yomhook ,ONLY : lhook, dr_hook
90 USE parkind1 ,ONLY : jprb
91 !
92 IMPLICIT NONE
93 !
94 !* 0.1 Declarations of dummy arguments :
95 !
96 REAL, DIMENSION(:,:), INTENT(INOUT) :: PFIELD1 ! FIELD on regular input grid
97 !
98 !* 0.2 Declarations of local variables for print on FM file
99 !
100 REAL, DIMENSION (:,:), INTENT(OUT) :: PFIELD_X ! FIELD at mesh interface
101 REAL, DIMENSION (:,:), INTENT(OUT) :: PFIELD_Y ! FIELD at mesh interface
102 REAL, DIMENSION (:,:), INTENT(OUT) :: PFIELD_XY! FIELD at mesh corner
103 !
104 REAL, DIMENSION(:,:),ALLOCATABLE :: ZW
105 INTEGER :: IIU ! model 1 X size
106 INTEGER :: IJU ! model 1 Y size
107 !
108 INTEGER :: JI ! grid 1 x index
109 !
110 INTEGER :: JL ! grid 2 index
111 !
112 REAL :: ZEPS=1.e-3
113 REAL(KIND=JPRB) :: ZHOOK_HANDLE, ZHOOK_HANDLE_OMP
114 !-------------------------------------------------------------------------------
115 !
116 IF (lhook) CALL dr_hook('BILIN_GRIDIN_1',0,zhook_handle)
117 
118 iiu=SIZE(pfield1,1)
119 iju=SIZE(pfield1,2)
120 !
121 ALLOCATE(zw(iiu,iju))
122 WHERE (pfield1/=xundef)
123  zw=1.
124 ELSEWHERE
125  zw=0.
126 END WHERE
127 !
128 IF (lhook) CALL dr_hook('BILIN_GRIDIN_1',1,zhook_handle)
129 !-------------------------------------------------------------------------------
130 !
131 !* 1. FIELD type at grid mesh interfaces (in X directions)
132 ! ----------------------------------
133 !
134 pfield_x(:,:) = 0.
135 !
136 pfield_x(1,:) = zw(1,:) * pfield1(1,:)
137 pfield_x(iiu+1,:) = zw(iiu,:) * pfield1(iiu,:)
138 !
139 IF (iiu>1) THEN
140 !$OMP PARALLEL PRIVATE(ZHOOK_HANDLE_OMP)
141  IF (lhook) CALL dr_hook('BILIN_GRIDIN_2',0,zhook_handle_omp)
142 !$OMP DO SCHEDULE(DYNAMIC,1) PRIVATE(JI,JL)
143  DO jl = 1,iju
144  DO ji = 2,iiu
145  IF (pfield1(ji-1,jl)/=0. .AND. pfield1(ji,jl)/=0.) THEN
146  pfield_x(ji,jl) = (zw(ji-1,jl)*pfield1(ji-1,jl) + zw(ji,jl)*pfield1(ji,jl)) / &
147  max(1.,(zw(ji-1,jl) + zw(ji,jl)) )
148  ENDIF
149  ENDDO
150  ENDDO
151 !$OMP END DO
152  IF (lhook) CALL dr_hook('BILIN_GRIDIN_2',1,zhook_handle_omp)
153 !$OMP END PARALLEL
154 ENDIF
155 !
156 !-------------------------------------------------------------------------------
157 !
158 !* 2. FIELD type at grid mesh interfaces (in X directions)
159 ! ----------------------------------
160 !
161 pfield_y(:,:) = 0.
162 !
163 pfield_y(:,1) = zw(:,1) * pfield1(:,1)
164 pfield_y(:,iju+1) = zw(:,iju) * pfield1(:,iju)
165 !
166 IF (iju>1) THEN
167 !$OMP PARALLEL PRIVATE(ZHOOK_HANDLE_OMP)
168  IF (lhook) CALL dr_hook('BILIN_GRIDIN_3',0,zhook_handle_omp)
169 !$OMP DO SCHEDULE(DYNAMIC,1) PRIVATE(JI,JL)
170  DO jl = 2,iju
171  DO ji = 1,iiu
172  IF (pfield1(ji,jl-1)/=0. .AND. pfield1(ji,jl)/=0.) THEN
173  pfield_y(ji,jl) = (zw(ji,jl-1)*pfield1(ji,jl-1) + zw(ji,jl)*pfield1(ji,jl)) / &
174  max(1.,(zw(ji,jl-1) + zw(ji,jl)) )
175  ENDIF
176  ENDDO
177  ENDDO
178 !$OMP END DO
179  IF (lhook) CALL dr_hook('BILIN_GRIDIN_3',1,zhook_handle_omp)
180 !$OMP END PARALLEL
181 ENDIF
182 !
183 !-------------------------------------------------------------------------------
184 !
185 !* 3. FIELD type at grid mesh corners
186 ! -------------------------------
187 !
188 pfield_xy(:,:) = 0.
189 !
190 IF (iiu>1 .AND. iju>1) THEN
191  !
192 !$OMP PARALLEL PRIVATE(ZHOOK_HANDLE_OMP)
193  IF (lhook) CALL dr_hook('BILIN_GRIDIN_4',0,zhook_handle_omp)
194 !$OMP DO SCHEDULE(DYNAMIC,1) PRIVATE(JI,JL)
195  DO jl = 2,iju
196  DO ji = 2,iiu
197  IF (pfield1(ji-1,jl-1)/=0. .AND. pfield1(ji-1,jl)/=0. .AND. &
198  pfield1(ji ,jl-1)/=0. .AND. pfield1(ji ,jl)/=0. ) THEN
199  pfield_xy(ji,jl) = (zw(ji-1,jl-1)*pfield1(ji-1,jl-1) + zw(ji-1,jl)*pfield1(ji-1,jl) + &
200  zw(ji, jl-1)*pfield1(ji, jl-1) + zw(ji, jl)*pfield1(ji, jl))/&
201  max(1.,zw(ji-1,jl-1) + zw(ji-1,jl) + zw(ji,jl-1) + zw(ji,jl))
202  ENDIF
203  ENDDO
204  ENDDO
205 !$OMP END DO
206  IF (lhook) CALL dr_hook('BILIN_GRIDIN_4',1,zhook_handle_omp)
207 !$OMP END PARALLEL
208  !
209 ENDIF
210 !
211 IF (lhook) CALL dr_hook('BILIN_GRIDIN_5',0,zhook_handle)
212 !
213 IF (iju>1) THEN
214  DO jl = 2,iju
215  IF (pfield1(1, jl-1)/=0. .AND. pfield1(1, jl)/=0.) THEN
216  pfield_xy(1 ,jl) = (zw(1,jl-1)*pfield1(1,jl-1) + zw(1,jl)*pfield1(1,jl)) / &
217  max(1.,zw(1,jl-1) + zw(1,jl))
218  ENDIF
219  IF (pfield1(iiu,jl-1)/=0. .AND. pfield1(iiu,jl)/=0.) THEN
220  pfield_xy(iiu+1,jl) = (zw(iiu,jl-1)*pfield1(iiu,jl-1) + zw(iiu,jl)*pfield1(iiu,jl)) / &
221  max(1.,zw(iiu,jl-1) + zw(iiu,jl))
222  ENDIF
223  ENDDO
224 ENDIF
225 !
226 IF (iiu>1) THEN
227  DO ji = 2,iiu
228  IF (pfield1(ji-1,1 )/=0. .AND. pfield1(ji,1 )/=0.) THEN
229  pfield_xy(ji,1 ) = (zw(ji-1,1)*pfield1(ji-1,1) + zw(ji,1)*pfield1(ji,1)) / &
230  max(1.,zw(ji-1,1) + zw(ji,1))
231  ENDIF
232  IF (pfield1(ji-1,iju)/=0. .AND. pfield1(ji,iju)/=0.) THEN
233  pfield_xy(ji,iju+1) = (zw(ji-1,iju)*pfield1(ji-1,iju) + zw(ji,iju)*pfield1(ji,iju)) / &
234  max(1.,zw(ji-1,iju) + zw(ji,iju))
235  ENDIF
236  ENDDO
237 ENDIF
238 !
239 pfield_xy(1 ,1 ) = pfield1(1 ,1 )
240 pfield_xy(iiu+1,1 ) = pfield1(iiu,1 )
241 pfield_xy(1 ,iju+1) = pfield1(1 ,iju)
242 pfield_xy(iiu+1,iju+1) = pfield1(iiu,iju)
243 !
244 IF (lhook) CALL dr_hook('BILIN_GRIDIN_5',1,zhook_handle)
245 !
246 DEALLOCATE(zw)
247 !
248 !-------------------------------------------------------------------------------
249 !
250 END SUBROUTINE bilin_gridin
real, parameter xundef
integer, parameter jprb
Definition: parkind1.F90:32
subroutine bilin_gridin(PFIELD1, PFIELD_X, PFIELD_Y, PFIELD_XY)
Definition: bilin_gridin.F90:7
logical lhook
Definition: yomhook.F90:15