[4213] | 1 | !---------------------------------------------------------------------- |
---|
| 2 | ! NEMO system team, System and Interface for oceanic RElocable Nesting |
---|
| 3 | !---------------------------------------------------------------------- |
---|
| 4 | ! |
---|
| 5 | ! MODULE: extrap |
---|
| 6 | ! |
---|
| 7 | ! DESCRIPTION: |
---|
| 8 | !> @brief |
---|
[5832] | 9 | !> This module manage extrapolation. |
---|
[4213] | 10 | !> |
---|
| 11 | !> @details |
---|
[5832] | 12 | !> Extrapolation method to be used is specify inside variable |
---|
| 13 | !> strcuture, as array of string character.<br/> |
---|
| 14 | !> - td_var\%c_extrap(1) string character is the interpolation name choose between: |
---|
| 15 | !> - 'dist_weight' |
---|
| 16 | !> - 'min_error' |
---|
[4213] | 17 | !> |
---|
[5832] | 18 | !> @note Extrapolation method could be specify for each variable in namelist _namvar_, |
---|
| 19 | !> defining string character _cn\_varinfo_. By default _dist_weight_.<br/> |
---|
| 20 | !> Example: |
---|
| 21 | !> - cn_varinfo='varname1:dist_weight', 'varname2:min_error' |
---|
| 22 | !> |
---|
| 23 | !> to detect point to be extrapolated:<br/> |
---|
| 24 | !> @code |
---|
| 25 | !> il_detect(:,:,:)=extrap_detect(td_var, [td_level], [id_offset,] [id_rho,] [id_ext]) |
---|
| 26 | !> @endcode |
---|
| 27 | !> - il_detect(:,:,:) is 3D array of point to be extrapolated |
---|
| 28 | !> - td_var is coarse grid variable to be extrapolated |
---|
| 29 | !> - td_level is fine grid array of level (see vgrid_get_level) [optional] |
---|
| 30 | !> - id_offset is array of offset between fine and coarse grid [optional] |
---|
| 31 | !> - id_rho is array of refinment factor [optional] |
---|
| 32 | !> - id_ext is array of number of points to be extrapolated [optional] |
---|
| 33 | !> |
---|
| 34 | !> to extrapolate variable:<br/> |
---|
| 35 | !> @code |
---|
| 36 | !> CALL extrap_fill_value( td_var, [td_level], [id_offset], [id_rho], [id_iext], [id_jext], [id_kext], [id_radius], [id_maxiter]) |
---|
| 37 | !> @endcode |
---|
| 38 | !> - td_var is coarse grid variable to be extrapolated |
---|
| 39 | !> - td_level is fine grid array of level (see vgrid_get_level) [optional] |
---|
| 40 | !> - id_offset is array of offset between fine and coarse grid [optional] |
---|
| 41 | !> - id_rho is array of refinment factor [optional] |
---|
| 42 | !> - id_iext is number of points to be extrapolated in i-direction [optional] |
---|
| 43 | !> - id_jext is number of points to be extrapolated in j-direction [optional] |
---|
| 44 | !> - id_kext is number of points to be extrapolated in k-direction [optional] |
---|
| 45 | !> - id_radius is radius of the halo used to compute extrapolation [optional] |
---|
| 46 | !> - id_maxiter is maximum number of iteration [optional] |
---|
| 47 | !> |
---|
| 48 | !> to add extraband to the variable (to be extrapolated):<br/> |
---|
| 49 | !> @code |
---|
| 50 | !> CALL extrap_add_extrabands(td_var, [id_isize,] [id_jsize] ) |
---|
| 51 | !> @endcode |
---|
| 52 | !> - td_var is variable structure |
---|
| 53 | !> - id_isize : i-direction size of extra bands [optional] |
---|
| 54 | !> - id_jsize : j-direction size of extra bands [optional] |
---|
| 55 | !> |
---|
| 56 | !> to delete extraband of a variable:<br/> |
---|
| 57 | !> @code |
---|
| 58 | !> CALL extrap_del_extrabands(td_var, [id_isize,] [id_jsize] ) |
---|
| 59 | !> @endcode |
---|
| 60 | !> - td_var is variable structure |
---|
| 61 | !> - id_isize : i-direction size of extra bands [optional] |
---|
| 62 | !> - id_jsize : j-direction size of extra bands [optional] |
---|
| 63 | !> |
---|
| 64 | !> to compute first derivative of 1D array:<br/> |
---|
| 65 | !> @code |
---|
| 66 | !> dl_value(:)=extrap_deriv_1D( dd_value(:), dd_fill, [ld_discont] ) |
---|
| 67 | !> @endcode |
---|
| 68 | !> - dd_value is 1D array of variable |
---|
| 69 | !> - dd_fill is FillValue of variable |
---|
| 70 | !> - ld_discont is logical to take into account longitudinal East-West discontinuity [optional] |
---|
| 71 | !> |
---|
| 72 | !> to compute first derivative of 2D array:<br/> |
---|
| 73 | !> @code |
---|
| 74 | !> dl_value(:,:)=extrap_deriv_2D( dd_value(:,:), dd_fill, cd_dim, [ld_discont] ) |
---|
| 75 | !> @endcode |
---|
| 76 | !> - dd_value is 2D array of variable |
---|
| 77 | !> - dd_fill is FillValue of variable |
---|
| 78 | !> - cd_dim is character to compute derivative on first (I) or second (J) dimension |
---|
| 79 | !> - ld_discont is logical to take into account longitudinal East-West discontinuity [optional] |
---|
| 80 | !> |
---|
| 81 | !> to compute first derivative of 3D array:<br/> |
---|
| 82 | !> @code |
---|
| 83 | !> dl_value(:,:,:)=extrap_deriv_3D( dd_value(:,:,:), dd_fill, cd_dim, [ld_discont] ) |
---|
| 84 | !> @endcode |
---|
| 85 | !> - dd_value is 3D array of variable |
---|
| 86 | !> - dd_fill is FillValue of variable |
---|
| 87 | !> - cd_dim is character to compute derivative on first (I), second (J), or third (K) dimension |
---|
| 88 | !> - ld_discont is logical to take into account longitudinal East-West discontinuity [optional] |
---|
| 89 | !> |
---|
| 90 | !> @warning _FillValue must not be zero (use var_chg_FillValue()) |
---|
| 91 | !> |
---|
[4213] | 92 | !> @author |
---|
| 93 | !> J.Paul |
---|
| 94 | ! REVISION HISTORY: |
---|
| 95 | !> @date Nov, 2013 - Initial Version |
---|
[5832] | 96 | !> @date September, 2014 |
---|
| 97 | !> - add header |
---|
[4213] | 98 | !> |
---|
[5832] | 99 | !> @todo |
---|
| 100 | !> - create module for each extrapolation method |
---|
[4213] | 101 | !> |
---|
| 102 | !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) |
---|
| 103 | !---------------------------------------------------------------------- |
---|
| 104 | MODULE extrap |
---|
| 105 | USE netcdf ! nf90 library |
---|
[5832] | 106 | USE kind ! F90 kind parameter |
---|
| 107 | USE phycst ! physical constant |
---|
| 108 | USE global ! global variable |
---|
| 109 | USE fct ! basic useful function |
---|
| 110 | USE date ! date manager |
---|
| 111 | USE logger ! log file manager |
---|
| 112 | USE att ! attribute manager |
---|
| 113 | USE dim ! dimension manager |
---|
| 114 | USE var ! variable manager |
---|
[4213] | 115 | |
---|
| 116 | IMPLICIT NONE |
---|
| 117 | ! NOTE_avoid_public_variables_if_possible |
---|
| 118 | |
---|
| 119 | ! type and variable |
---|
[5832] | 120 | PRIVATE :: im_maxiter !< default maximum number of iteration |
---|
| 121 | PRIVATE :: im_minext !< default minumum number of point to extrapolate |
---|
| 122 | PRIVATE :: im_mincubic !< default minumum number of point to extrapolate for cubic interpolation |
---|
[4213] | 123 | |
---|
| 124 | ! function and subroutine |
---|
| 125 | PUBLIC :: extrap_detect !< detected point to be extrapolated |
---|
| 126 | PUBLIC :: extrap_fill_value !< extrapolate value over detected point |
---|
[5832] | 127 | PUBLIC :: extrap_add_extrabands !< add extraband to the variable (to be extrapolated) |
---|
| 128 | PUBLIC :: extrap_del_extrabands !< delete extraband of the variable |
---|
| 129 | PUBLIC :: extrap_deriv_1D !< compute first derivative of 1D array |
---|
| 130 | PUBLIC :: extrap_deriv_2D !< compute first derivative of 2D array |
---|
| 131 | PUBLIC :: extrap_deriv_3D !< compute first derivative of 3D array |
---|
[4213] | 132 | |
---|
[5832] | 133 | PRIVATE :: extrap__detect_wrapper ! detected point to be extrapolated wrapper |
---|
| 134 | PRIVATE :: extrap__detect ! detected point to be extrapolated |
---|
| 135 | PRIVATE :: extrap__fill_value_wrapper ! extrapolate value over detected point wrapper |
---|
| 136 | PRIVATE :: extrap__fill_value ! extrapolate value over detected point |
---|
| 137 | PRIVATE :: extrap__3D ! |
---|
| 138 | PRIVATE :: extrap__3D_min_error_coef ! |
---|
| 139 | PRIVATE :: extrap__3D_min_error_fill ! |
---|
| 140 | PRIVATE :: extrap__3D_dist_weight_coef ! |
---|
| 141 | PRIVATE :: extrap__3D_dist_weight_fill ! |
---|
[4213] | 142 | |
---|
| 143 | INTEGER(i4), PARAMETER :: im_maxiter = 10 !< default maximum number of iteration |
---|
| 144 | INTEGER(i4), PARAMETER :: im_minext = 2 !< default minumum number of point to extrapolate |
---|
| 145 | INTEGER(i4), PARAMETER :: im_mincubic= 4 !< default minumum number of point to extrapolate for cubic interpolation |
---|
| 146 | |
---|
| 147 | INTERFACE extrap_detect |
---|
[5832] | 148 | MODULE PROCEDURE extrap__detect_wrapper !< detected point to be extrapolated |
---|
[4213] | 149 | END INTERFACE extrap_detect |
---|
| 150 | |
---|
| 151 | INTERFACE extrap_fill_value |
---|
| 152 | MODULE PROCEDURE extrap__fill_value_wrapper !< detected point to be interpolated |
---|
| 153 | END INTERFACE extrap_fill_value |
---|
| 154 | |
---|
| 155 | CONTAINS |
---|
| 156 | !------------------------------------------------------------------- |
---|
| 157 | !> @brief |
---|
[5832] | 158 | !> This function detected point to be extrapolated, given variable structure. |
---|
[4213] | 159 | !> |
---|
| 160 | !> @details |
---|
[5832] | 161 | !> optionaly, you could sepcify fine grid level, refinment factor (default 1), |
---|
| 162 | !> offset between fine and coarse grid (default compute from refinment factor |
---|
| 163 | !> as offset=(rho-1)/2), number of point to be extrapolated in each direction |
---|
| 164 | !> (default im_minext).<br/> |
---|
[4213] | 165 | !> |
---|
[5832] | 166 | !> First coarsening fine grid level, if need be, then select point near |
---|
| 167 | !> grid point already inform. |
---|
| 168 | !> |
---|
| 169 | !> @note point to be extrapolated are selected using FillValue, |
---|
| 170 | !> so to avoid mistake FillValue should not be zero (use var_chg_FillValue) |
---|
| 171 | !> |
---|
[4213] | 172 | !> @author J.Paul |
---|
[5832] | 173 | !> - November, 2013- Initial Version |
---|
[4213] | 174 | ! |
---|
[5832] | 175 | !> @param[in] td_var0 coarse grid variable to extrapolate |
---|
| 176 | !> @param[in] td_level1 fine grid array of level |
---|
| 177 | !> @param[in] id_offset array of offset between fine and coarse grid |
---|
| 178 | !> @param[in] id_rho array of refinment factor |
---|
| 179 | !> @param[in] id_ext array of number of points to be extrapolated |
---|
| 180 | !> @return array of point to be extrapolated |
---|
[4213] | 181 | !------------------------------------------------------------------- |
---|
| 182 | FUNCTION extrap__detect( td_var0, td_level1, & |
---|
| 183 | & id_offset, id_rho, id_ext ) |
---|
| 184 | IMPLICIT NONE |
---|
| 185 | ! Argument |
---|
[5832] | 186 | TYPE(TVAR) , INTENT(IN ) :: td_var0 |
---|
[4213] | 187 | TYPE(TVAR) , DIMENSION(:) , INTENT(IN ), OPTIONAL :: td_level1 |
---|
| 188 | INTEGER(i4), DIMENSION(:,:), INTENT(IN ), OPTIONAL :: id_offset |
---|
| 189 | INTEGER(i4), DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_rho |
---|
| 190 | INTEGER(i4), DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_ext |
---|
| 191 | |
---|
| 192 | ! function |
---|
| 193 | INTEGER(i4), DIMENSION(td_var0%t_dim(1)%i_len,& |
---|
| 194 | & td_var0%t_dim(2)%i_len,& |
---|
| 195 | & td_var0%t_dim(3)%i_len ) :: extrap__detect |
---|
| 196 | |
---|
| 197 | ! local variable |
---|
| 198 | CHARACTER(LEN=lc) :: cl_level |
---|
| 199 | |
---|
[5832] | 200 | INTEGER(i4) :: il_ind |
---|
[4213] | 201 | INTEGER(i4) , DIMENSION(:,:,:), ALLOCATABLE :: il_detect |
---|
[5832] | 202 | INTEGER(i4) , DIMENSION(:,:,:), ALLOCATABLE :: il_tmp |
---|
[4213] | 203 | INTEGER(i4) , DIMENSION(:,:) , ALLOCATABLE :: il_offset |
---|
| 204 | INTEGER(i4) , DIMENSION(:,:) , ALLOCATABLE :: il_level1 |
---|
| 205 | INTEGER(i4) , DIMENSION(:,:) , ALLOCATABLE :: il_level1_G0 |
---|
| 206 | INTEGER(i4) , DIMENSION(:,:) , ALLOCATABLE :: il_extra |
---|
| 207 | INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_ext |
---|
| 208 | INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_rho |
---|
| 209 | INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_dim0 |
---|
| 210 | |
---|
| 211 | TYPE(TVAR) :: tl_var1 |
---|
| 212 | |
---|
| 213 | ! loop indices |
---|
| 214 | INTEGER(i4) :: ji0 |
---|
| 215 | INTEGER(i4) :: jj0 |
---|
| 216 | INTEGER(i4) :: jk0 |
---|
| 217 | INTEGER(i4) :: ji1 |
---|
| 218 | INTEGER(i4) :: jj1 |
---|
| 219 | INTEGER(i4) :: ji1m |
---|
| 220 | INTEGER(i4) :: jj1m |
---|
| 221 | INTEGER(i4) :: ji1p |
---|
| 222 | INTEGER(i4) :: jj1p |
---|
| 223 | !---------------------------------------------------------------- |
---|
| 224 | |
---|
| 225 | ! init |
---|
| 226 | extrap__detect(:,:,:)=0 |
---|
| 227 | |
---|
| 228 | ALLOCATE( il_dim0(3) ) |
---|
| 229 | il_dim0(:)=td_var0%t_dim(1:3)%i_len |
---|
| 230 | |
---|
| 231 | ! optional argument |
---|
| 232 | ALLOCATE( il_rho(ip_maxdim) ) |
---|
| 233 | il_rho(:)=1 |
---|
| 234 | IF( PRESENT(id_rho) ) il_rho(1:SIZE(id_rho(:)))=id_rho(:) |
---|
| 235 | |
---|
| 236 | ALLOCATE( il_offset(ip_maxdim,2) ) |
---|
| 237 | il_offset(:,:)=0 |
---|
| 238 | IF( PRESENT(id_offset) )THEN |
---|
| 239 | il_offset(1:SIZE(id_offset(:,:),DIM=1),& |
---|
| 240 | & 1:SIZE(id_offset(:,:),DIM=2) )= id_offset(:,:) |
---|
[5832] | 241 | ELSE |
---|
| 242 | il_offset(jp_I,:)=FLOOR(REAL(il_rho(jp_I)-1,dp)*0.5) |
---|
| 243 | il_offset(jp_J,:)=FLOOR(REAL(il_rho(jp_J)-1,dp)*0.5) |
---|
[4213] | 244 | ENDIF |
---|
| 245 | |
---|
| 246 | ALLOCATE( il_ext(ip_maxdim) ) |
---|
| 247 | il_ext(:)=im_minext |
---|
| 248 | IF( PRESENT(id_ext) ) il_ext(1:SIZE(id_ext(:)))=id_ext(:) |
---|
| 249 | |
---|
| 250 | ALLOCATE( il_detect(il_dim0(1),& |
---|
| 251 | & il_dim0(2),& |
---|
| 252 | & il_dim0(3)) ) |
---|
| 253 | il_detect(:,:,:)=0 |
---|
| 254 | |
---|
| 255 | ! select point already inform |
---|
[5832] | 256 | DO jk0=1,td_var0%t_dim(3)%i_len |
---|
| 257 | DO jj0=1,td_var0%t_dim(2)%i_len |
---|
| 258 | DO ji0=1,td_var0%t_dim(1)%i_len |
---|
| 259 | IF( td_var0%d_value(ji0,jj0,jk0,1) /= td_var0%d_fill ) il_detect(ji0,jj0,jk0)=1 |
---|
| 260 | ENDDO |
---|
| 261 | ENDDO |
---|
| 262 | ENDDO |
---|
[4213] | 263 | |
---|
| 264 | IF( PRESENT(td_level1) )THEN |
---|
| 265 | SELECT CASE(TRIM(td_var0%c_point)) |
---|
[5832] | 266 | CASE DEFAULT !'T' |
---|
| 267 | cl_level='tlevel' |
---|
| 268 | CASE('U') |
---|
| 269 | cl_level='ulevel' |
---|
| 270 | CASE('V') |
---|
| 271 | cl_level='vlevel' |
---|
| 272 | CASE('F') |
---|
| 273 | cl_level='flevel' |
---|
[4213] | 274 | END SELECT |
---|
[5832] | 275 | |
---|
| 276 | il_ind=var_get_index(td_level1(:),TRIM(cl_level)) |
---|
| 277 | IF( il_ind == 0 )THEN |
---|
[4213] | 278 | CALL logger_error("EXTRAP DETECT: can not compute point to be "//& |
---|
| 279 | & "extrapolated for variable "//TRIM(td_var0%c_name)//& |
---|
| 280 | & ". can not find "//& |
---|
| 281 | & "level for variable point "//TRIM(TRIM(td_var0%c_point))) |
---|
| 282 | ELSE |
---|
[5832] | 283 | tl_var1=var_copy(td_level1(il_ind)) |
---|
[4213] | 284 | |
---|
[5832] | 285 | ALLOCATE( il_level1_G0( il_dim0(1), il_dim0(2)) ) |
---|
| 286 | IF( ALL(tl_var1%t_dim(1:2)%i_len == il_dim0(1:2)) )THEN |
---|
[4213] | 287 | |
---|
[5832] | 288 | ! variable to be extrapolated use same resolution than level |
---|
| 289 | il_level1_G0(:,:)=INT(tl_var1%d_value(:,:,1,1),i4) |
---|
| 290 | |
---|
| 291 | ELSE |
---|
| 292 | ! variable to be extrapolated do not use same resolution than level |
---|
| 293 | ALLOCATE( il_level1(tl_var1%t_dim(1)%i_len, & |
---|
| 294 | & tl_var1%t_dim(2)%i_len) ) |
---|
| 295 | ! match fine grid vertical level with coarse grid |
---|
| 296 | il_level1(:,:)=INT(tl_var1%d_value(:,:,1,1),i4)/il_rho(jp_K) |
---|
[4213] | 297 | |
---|
[5832] | 298 | ALLOCATE( il_extra(ip_maxdim,2) ) |
---|
| 299 | ! coarsening fine grid level |
---|
| 300 | il_extra(jp_I,1)=CEILING(REAL(il_rho(jp_I)-1,dp)*0.5_dp) |
---|
| 301 | il_extra(jp_I,2)=FLOOR(REAL(il_rho(jp_I)-1,dp)*0.5_dp) |
---|
[4213] | 302 | |
---|
[5832] | 303 | il_extra(jp_J,1)=CEILING(REAL(il_rho(jp_J)-1,dp)*0.5_dp) |
---|
| 304 | il_extra(jp_J,2)=FLOOR(REAL(il_rho(jp_J)-1,dp)*0.5_dp) |
---|
[4213] | 305 | |
---|
[5832] | 306 | DO jj0=1,td_var0%t_dim(2)%i_len |
---|
| 307 | |
---|
| 308 | jj1=(jj0-1)*il_rho(jp_J)+1-il_offset(jp_J,1) |
---|
[4213] | 309 | |
---|
[5832] | 310 | jj1m=MAX( jj1-il_extra(jp_J,1), 1 ) |
---|
| 311 | jj1p=MIN( jj1+il_extra(jp_J,2), & |
---|
| 312 | & tl_var1%t_dim(2)%i_len-il_offset(jp_J,2) ) |
---|
| 313 | |
---|
| 314 | DO ji0=1,td_var0%t_dim(1)%i_len |
---|
| 315 | |
---|
| 316 | ji1=(ji0-1)*il_rho(jp_I)+1-id_offset(jp_I,1) |
---|
| 317 | |
---|
| 318 | ji1m=MAX( ji1-il_extra(jp_I,1), 1 ) |
---|
| 319 | ji1p=MIN( ji1+il_extra(jp_I,2), & |
---|
| 320 | & tl_var1%t_dim(1)%i_len-id_offset(jp_I,2) ) |
---|
[4213] | 321 | |
---|
[5832] | 322 | il_level1_G0(ji0,jj0)=MAXVAL(il_level1(ji1m:ji1p,jj1m:jj1p)) |
---|
[4213] | 323 | |
---|
[5832] | 324 | ENDDO |
---|
[4213] | 325 | ENDDO |
---|
| 326 | |
---|
[5832] | 327 | ! clean |
---|
| 328 | DEALLOCATE( il_extra ) |
---|
| 329 | DEALLOCATE( il_level1 ) |
---|
[4213] | 330 | |
---|
[5832] | 331 | ENDIF |
---|
[4213] | 332 | |
---|
[5832] | 333 | ! look for sea point |
---|
| 334 | DO jk0=1,td_var0%t_dim(3)%i_len |
---|
| 335 | WHERE( il_level1_G0(:,:) >= jk0) |
---|
| 336 | il_detect(:,:,jk0)=1 |
---|
| 337 | END WHERE |
---|
| 338 | ENDDO |
---|
[4213] | 339 | |
---|
| 340 | ! clean |
---|
[5832] | 341 | DEALLOCATE( il_level1_G0 ) |
---|
| 342 | CALL var_clean(tl_var1) |
---|
[4213] | 343 | |
---|
| 344 | ENDIF |
---|
| 345 | ENDIF |
---|
| 346 | |
---|
| 347 | ! clean |
---|
| 348 | DEALLOCATE( il_offset ) |
---|
| 349 | |
---|
[5832] | 350 | ALLOCATE( il_tmp(il_dim0(1),& |
---|
| 351 | & il_dim0(2),& |
---|
| 352 | & il_dim0(3)) ) |
---|
| 353 | il_tmp(:,:,:)=il_detect(:,:,:) |
---|
| 354 | ! select extra point depending on interpolation method |
---|
| 355 | ! compute point near grid point already inform |
---|
| 356 | DO jk0=1,il_dim0(3) |
---|
| 357 | DO jj0=1,il_dim0(2) |
---|
| 358 | DO ji0=1,il_dim0(1) |
---|
[4213] | 359 | |
---|
[5832] | 360 | IF( il_tmp(ji0,jj0,jk0) == 1 )THEN |
---|
| 361 | il_detect( & |
---|
| 362 | & MAX(1,ji0-il_ext(jp_I)):MIN(ji0+il_ext(jp_I),il_dim0(1)),& |
---|
| 363 | & MAX(1,jj0-il_ext(jp_J)):MIN(jj0+il_ext(jp_J),il_dim0(2)),& |
---|
| 364 | & MAX(1,jk0-il_ext(jp_K)):MIN(jk0+il_ext(jp_K),il_dim0(3)) & |
---|
| 365 | & ) = 1 |
---|
| 366 | ENDIF |
---|
[4213] | 367 | |
---|
[5832] | 368 | ENDDO |
---|
| 369 | ENDDO |
---|
| 370 | ENDDO |
---|
[4213] | 371 | |
---|
[5832] | 372 | ! clean |
---|
| 373 | DEALLOCATE( il_tmp ) |
---|
| 374 | |
---|
[4213] | 375 | ! do not compute grid point already inform |
---|
[5832] | 376 | DO jk0=1,td_var0%t_dim(3)%i_len |
---|
| 377 | DO jj0=1,td_var0%t_dim(2)%i_len |
---|
| 378 | DO ji0=1,td_var0%t_dim(1)%i_len |
---|
| 379 | IF( td_var0%d_value(ji0,jj0,jk0,1) /= td_var0%d_fill ) il_detect(ji0,jj0,jk0)=0 |
---|
| 380 | ENDDO |
---|
| 381 | ENDDO |
---|
| 382 | ENDDO |
---|
[4213] | 383 | |
---|
| 384 | ! save result |
---|
| 385 | extrap__detect(:,:,:)=il_detect(:,:,:) |
---|
| 386 | |
---|
| 387 | ! clean |
---|
| 388 | DEALLOCATE( il_dim0 ) |
---|
| 389 | DEALLOCATE( il_ext ) |
---|
| 390 | DEALLOCATE( il_detect ) |
---|
[5832] | 391 | DEALLOCATE( il_rho ) |
---|
[4213] | 392 | |
---|
| 393 | END FUNCTION extrap__detect |
---|
| 394 | !------------------------------------------------------------------- |
---|
| 395 | !> @brief |
---|
[5832] | 396 | !> This function sort variable to be extrapolated, depending on number of |
---|
| 397 | !> dimentsion, then detected point to be extrapolated. |
---|
[4213] | 398 | !> |
---|
[5832] | 399 | !> @author J.Paul |
---|
| 400 | !> - November, 2013- Initial Version |
---|
[4213] | 401 | !> |
---|
[5832] | 402 | !> @param[in] td_var coarse grid variable to extrapolate |
---|
| 403 | !> @param[in] td_level fine grid array of level |
---|
| 404 | !> @param[in] id_offset array of offset between fine and coarse grid |
---|
| 405 | !> @param[in] id_rho array of refinment factor |
---|
| 406 | !> @param[in] id_ext array of number of points to be extrapolated |
---|
| 407 | !> @return 3D array of point to be extrapolated |
---|
[4213] | 408 | !------------------------------------------------------------------- |
---|
| 409 | FUNCTION extrap__detect_wrapper( td_var, td_level, & |
---|
| 410 | & id_offset, id_rho, id_ext ) |
---|
| 411 | |
---|
| 412 | IMPLICIT NONE |
---|
| 413 | ! Argument |
---|
[5832] | 414 | TYPE(TVAR) , INTENT(IN ) :: td_var |
---|
| 415 | TYPE(TVAR) , DIMENSION(:) , INTENT(IN ), OPTIONAL :: td_level |
---|
[4213] | 416 | INTEGER(i4), DIMENSION(:,:), INTENT(IN ), OPTIONAL :: id_offset |
---|
| 417 | INTEGER(i4), DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_rho |
---|
| 418 | INTEGER(i4), DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_ext |
---|
| 419 | |
---|
| 420 | ! function |
---|
| 421 | INTEGER(i4), DIMENSION(td_var%t_dim(1)%i_len,& |
---|
| 422 | & td_var%t_dim(2)%i_len,& |
---|
| 423 | & td_var%t_dim(3)%i_len ) :: extrap__detect_wrapper |
---|
| 424 | |
---|
| 425 | ! local variable |
---|
| 426 | ! loop indices |
---|
| 427 | !---------------------------------------------------------------- |
---|
| 428 | ! init |
---|
| 429 | extrap__detect_wrapper(:,:,:)=0 |
---|
| 430 | |
---|
| 431 | IF( .NOT. ANY(td_var%t_dim(1:3)%l_use) )THEN |
---|
| 432 | ! no dimension I-J-K used |
---|
| 433 | CALL logger_debug(" EXTRAP DETECT: nothing done for variable"//& |
---|
| 434 | & TRIM(td_var%c_name) ) |
---|
| 435 | ELSE IF( ALL(td_var%t_dim(1:3)%l_use) )THEN |
---|
| 436 | |
---|
[5832] | 437 | ! detect point to be extrapolated on I-J-K |
---|
[4213] | 438 | CALL logger_debug(" EXTRAP DETECT: detect point "//& |
---|
| 439 | & " for variable "//TRIM(td_var%c_name) ) |
---|
| 440 | |
---|
| 441 | extrap__detect_wrapper(:,:,:)=extrap__detect( td_var, td_level, & |
---|
| 442 | & id_offset, & |
---|
| 443 | & id_rho, & |
---|
| 444 | & id_ext ) |
---|
| 445 | |
---|
| 446 | ELSE IF( ALL(td_var%t_dim(1:2)%l_use) )THEN |
---|
| 447 | |
---|
[5832] | 448 | ! detect point to be extrapolated on I-J |
---|
[4213] | 449 | CALL logger_debug(" EXTRAP DETECT: detect horizontal point "//& |
---|
| 450 | & " for variable "//TRIM(td_var%c_name) ) |
---|
| 451 | |
---|
| 452 | extrap__detect_wrapper(:,:,1:1)=extrap__detect( td_var , td_level,& |
---|
| 453 | & id_offset, & |
---|
| 454 | & id_rho, & |
---|
| 455 | & id_ext ) |
---|
| 456 | |
---|
| 457 | ELSE IF( td_var%t_dim(3)%l_use )THEN |
---|
| 458 | |
---|
[5832] | 459 | ! detect point to be extrapolated on K |
---|
[4213] | 460 | CALL logger_debug(" EXTRAP DETECT: detect vertical point "//& |
---|
| 461 | & " for variable "//TRIM(td_var%c_name) ) |
---|
| 462 | |
---|
| 463 | extrap__detect_wrapper(1:1,1:1,:)=extrap__detect( td_var , td_level, & |
---|
| 464 | & id_offset, & |
---|
| 465 | & id_rho, & |
---|
| 466 | & id_ext ) |
---|
| 467 | |
---|
| 468 | ENDIF |
---|
| 469 | |
---|
| 470 | CALL logger_debug(" EXTRAP DETECT: "//& |
---|
| 471 | & TRIM(fct_str(SUM(extrap__detect_wrapper(:,:,:))))//& |
---|
| 472 | & " points to be extrapolated" ) |
---|
| 473 | |
---|
| 474 | END FUNCTION extrap__detect_wrapper |
---|
| 475 | !------------------------------------------------------------------- |
---|
| 476 | !> @brief |
---|
[5832] | 477 | !> This subroutine select method to be used for extrapolation. |
---|
| 478 | !> If need be, increase number of points to be extrapolated. |
---|
| 479 | !> Finally launch extrap__fill_value. |
---|
[4213] | 480 | !> |
---|
[5832] | 481 | !> @details |
---|
| 482 | !> optionaly, you could specify :<br/> |
---|
| 483 | !> - refinment factor (default 1) |
---|
| 484 | !> - offset between fine and coarse grid (default compute from refinment factor |
---|
| 485 | !> as offset=(rho-1)/2) |
---|
| 486 | !> - number of point to be extrapolated in each direction (default im_minext) |
---|
| 487 | !> - radius of the halo used to compute extrapolation |
---|
| 488 | !> - maximum number of iteration |
---|
[4213] | 489 | !> |
---|
| 490 | !> @author J.Paul |
---|
| 491 | !> - Nov, 2013- Initial Version |
---|
| 492 | ! |
---|
[5832] | 493 | !> @param[inout] td_var variable structure |
---|
| 494 | !> @param[in] td_level fine grid array of level |
---|
| 495 | !> @param[in] id_offset array of offset between fine and coarse grid |
---|
| 496 | !> @param[in] id_rho array of refinment factor |
---|
| 497 | !> @param[in] id_iext number of points to be extrapolated in i-direction |
---|
| 498 | !> @param[in] id_jext number of points to be extrapolated in j-direction |
---|
| 499 | !> @param[in] id_kext number of points to be extrapolated in k-direction |
---|
| 500 | !> @param[in] id_radius radius of the halo used to compute extrapolation |
---|
| 501 | !> @param[in] id_maxiter maximum number of iteration |
---|
[4213] | 502 | !------------------------------------------------------------------- |
---|
| 503 | SUBROUTINE extrap__fill_value_wrapper( td_var, td_level, & |
---|
| 504 | & id_offset, & |
---|
| 505 | & id_rho, & |
---|
| 506 | & id_iext, id_jext, id_kext, & |
---|
| 507 | & id_radius, id_maxiter ) |
---|
| 508 | IMPLICIT NONE |
---|
| 509 | ! Argument |
---|
| 510 | TYPE(TVAR) , INTENT(INOUT) :: td_var |
---|
| 511 | TYPE(TVAR) , DIMENSION(:) , INTENT(IN ), OPTIONAL :: td_level |
---|
| 512 | INTEGER(i4), DIMENSION(:,:), INTENT(IN ), OPTIONAL :: id_offset |
---|
| 513 | INTEGER(i4), DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_rho |
---|
| 514 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext |
---|
| 515 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext |
---|
| 516 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_kext |
---|
| 517 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_radius |
---|
| 518 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_maxiter |
---|
| 519 | |
---|
| 520 | ! local variable |
---|
| 521 | INTEGER(i4) :: il_iext |
---|
| 522 | INTEGER(i4) :: il_jext |
---|
| 523 | INTEGER(i4) :: il_kext |
---|
| 524 | INTEGER(i4) :: il_radius |
---|
| 525 | INTEGER(i4) :: il_maxiter |
---|
| 526 | |
---|
| 527 | CHARACTER(LEN=lc) :: cl_method |
---|
| 528 | |
---|
| 529 | ! loop indices |
---|
| 530 | !---------------------------------------------------------------- |
---|
| 531 | IF( .NOT. ASSOCIATED(td_var%d_value) )THEN |
---|
[5832] | 532 | CALL logger_error("EXTRAP FILL VALUE: no value "//& |
---|
[4213] | 533 | & "associted to variable "//TRIM(td_var%c_name) ) |
---|
| 534 | ELSE |
---|
| 535 | |
---|
| 536 | SELECT CASE(TRIM(td_var%c_extrap(1))) |
---|
| 537 | CASE('min_error') |
---|
| 538 | cl_method='min_error' |
---|
| 539 | CASE DEFAULT |
---|
| 540 | cl_method='dist_weight' |
---|
| 541 | |
---|
| 542 | !update variable structure |
---|
| 543 | td_var%c_extrap(1)='dist_weight' |
---|
| 544 | END SELECT |
---|
| 545 | |
---|
| 546 | il_iext=im_minext |
---|
| 547 | IF( PRESENT(id_iext) ) il_iext=id_iext |
---|
| 548 | il_jext=im_minext |
---|
| 549 | IF( PRESENT(id_jext) ) il_jext=id_jext |
---|
| 550 | il_kext=0 |
---|
| 551 | IF( PRESENT(id_kext) ) il_kext=id_kext |
---|
| 552 | |
---|
| 553 | IF( TRIM(td_var%c_interp(1)) == 'cubic')THEN |
---|
| 554 | IF( il_iext > 0 .AND. il_iext < im_mincubic ) il_iext=im_mincubic |
---|
| 555 | IF( il_jext > 0 .AND. il_jext < im_mincubic ) il_jext=im_mincubic |
---|
| 556 | ENDIF |
---|
| 557 | |
---|
| 558 | IF( il_iext < 0 )THEN |
---|
| 559 | CALL logger_error("EXTRAP FILL VALUE: invalid "//& |
---|
| 560 | & " number of points to be extrapolated in i-direction "//& |
---|
| 561 | & "("//TRIM(fct_str(il_iext))//")") |
---|
| 562 | ENDIF |
---|
| 563 | |
---|
| 564 | IF( il_jext < 0 )THEN |
---|
| 565 | CALL logger_error("EXTRAP FILL VALUE: invalid "//& |
---|
| 566 | & " number of points to be extrapolated in j-direction "//& |
---|
| 567 | & "("//TRIM(fct_str(il_jext))//")") |
---|
| 568 | ENDIF |
---|
| 569 | |
---|
| 570 | IF( il_kext < 0 )THEN |
---|
| 571 | CALL logger_error("EXTRAP FILL VALUE: invalid "//& |
---|
| 572 | & " number of points to be extrapolated in k-direction "//& |
---|
| 573 | & "("//TRIM(fct_str(il_kext))//")") |
---|
| 574 | ENDIF |
---|
| 575 | |
---|
| 576 | IF( (il_iext /= 0 .AND. td_var%t_dim(1)%l_use) .OR. & |
---|
| 577 | & (il_jext /= 0 .AND. td_var%t_dim(2)%l_use) .OR. & |
---|
[5832] | 578 | & (il_kext /= 0 .AND. td_var%t_dim(3)%l_use) )THEN |
---|
[4213] | 579 | |
---|
| 580 | ! number of point use to compute box |
---|
| 581 | il_radius=1 |
---|
| 582 | IF( PRESENT(id_radius) ) il_radius=id_radius |
---|
| 583 | IF( il_radius < 0 )THEN |
---|
| 584 | CALL logger_error("EXTRAP FILL VALUE: invalid "//& |
---|
| 585 | & " radius of the box used to compute extrapolation "//& |
---|
| 586 | & "("//TRIM(fct_str(il_radius))//")") |
---|
| 587 | ENDIF |
---|
| 588 | |
---|
| 589 | ! maximum number of iteration |
---|
| 590 | il_maxiter=im_maxiter |
---|
| 591 | IF( PRESENT(id_maxiter) ) il_maxiter=id_maxiter |
---|
| 592 | IF( il_maxiter < 0 )THEN |
---|
| 593 | CALL logger_error("EXTRAP FILL VALUE: invalid "//& |
---|
| 594 | & " maximum nuber of iteration "//& |
---|
| 595 | & "("//TRIM(fct_str(il_maxiter))//")") |
---|
| 596 | ENDIF |
---|
| 597 | |
---|
| 598 | CALL logger_info("EXTRAP FILL: extrapolate "//TRIM(td_var%c_name)//& |
---|
| 599 | & " using "//TRIM(cl_method)//" method." ) |
---|
| 600 | |
---|
| 601 | CALL extrap__fill_value( td_var, cl_method, & |
---|
| 602 | & il_iext, il_jext, il_kext, & |
---|
| 603 | & il_radius, il_maxiter, & |
---|
| 604 | & td_level, & |
---|
| 605 | & id_offset, id_rho ) |
---|
| 606 | |
---|
| 607 | ENDIF |
---|
| 608 | |
---|
| 609 | ENDIF |
---|
| 610 | |
---|
| 611 | END SUBROUTINE extrap__fill_value_wrapper |
---|
| 612 | !------------------------------------------------------------------- |
---|
| 613 | !> @brief |
---|
[5832] | 614 | !> This subroutine compute point to be extrapolated, then extrapolate point. |
---|
[4213] | 615 | !> |
---|
| 616 | !> @details |
---|
[5832] | 617 | !> optionaly, you could specify :<br/> |
---|
| 618 | !> - refinment factor (default 1) |
---|
| 619 | !> - offset between fine and coarse grid (default compute from refinment factor |
---|
| 620 | !> as offset=(rho-1)/2) |
---|
[4213] | 621 | !> |
---|
| 622 | !> @author J.Paul |
---|
[5832] | 623 | !> - November, 2013- Initial Version |
---|
[4213] | 624 | ! |
---|
[5832] | 625 | !> @param[inout] td_var variable structure |
---|
| 626 | !> @param[in] cd_method extrapolation method |
---|
| 627 | !> @param[in] id_iext number of points to be extrapolated in i-direction |
---|
| 628 | !> @param[in] id_jext number of points to be extrapolated in j-direction |
---|
| 629 | !> @param[in] id_kext number of points to be extrapolated in k-direction |
---|
| 630 | !> @param[in] id_radius radius of the halo used to compute extrapolation |
---|
| 631 | !> @param[in] id_maxiter maximum number of iteration |
---|
| 632 | !> @param[in] td_level fine grid array of level |
---|
| 633 | !> @param[in] id_offset array of offset between fine and coarse grid |
---|
| 634 | !> @param[in] id_rho array of refinment factor |
---|
[4213] | 635 | !------------------------------------------------------------------- |
---|
| 636 | SUBROUTINE extrap__fill_value( td_var, cd_method, & |
---|
| 637 | & id_iext, id_jext, id_kext, & |
---|
| 638 | & id_radius, id_maxiter, & |
---|
| 639 | & td_level, & |
---|
| 640 | & id_offset, & |
---|
| 641 | & id_rho ) |
---|
| 642 | IMPLICIT NONE |
---|
| 643 | ! Argument |
---|
| 644 | TYPE(TVAR) , INTENT(INOUT) :: td_var |
---|
| 645 | CHARACTER(LEN=*), INTENT(IN ) :: cd_method |
---|
| 646 | INTEGER(i4) , INTENT(IN ) :: id_iext |
---|
| 647 | INTEGER(i4) , INTENT(IN ) :: id_jext |
---|
| 648 | INTEGER(i4) , INTENT(IN ) :: id_kext |
---|
| 649 | INTEGER(i4) , INTENT(IN ) :: id_radius |
---|
| 650 | INTEGER(i4) , INTENT(IN ) :: id_maxiter |
---|
| 651 | TYPE(TVAR) , DIMENSION(:) , INTENT(IN ), OPTIONAL :: td_level |
---|
| 652 | INTEGER(i4) , DIMENSION(:,:), INTENT(IN ), OPTIONAL :: id_offset |
---|
| 653 | INTEGER(i4) , DIMENSION(:) , INTENT(IN ), OPTIONAL :: id_rho |
---|
| 654 | |
---|
| 655 | ! local variable |
---|
| 656 | CHARACTER(LEN=lc) :: cl_extrap |
---|
| 657 | |
---|
| 658 | INTEGER(i4), DIMENSION(:,:,:) , ALLOCATABLE :: il_detect |
---|
| 659 | |
---|
| 660 | TYPE(TATT) :: tl_att |
---|
| 661 | |
---|
| 662 | ! loop indices |
---|
| 663 | !---------------------------------------------------------------- |
---|
| 664 | |
---|
| 665 | !1- detect point to be extrapolated |
---|
| 666 | ALLOCATE( il_detect( td_var%t_dim(1)%i_len, & |
---|
| 667 | & td_var%t_dim(2)%i_len, & |
---|
| 668 | & td_var%t_dim(3)%i_len) ) |
---|
| 669 | |
---|
| 670 | il_detect(:,:,:) = extrap_detect( td_var, td_level, & |
---|
| 671 | & id_offset, & |
---|
| 672 | & id_rho, & |
---|
| 673 | & id_ext=(/id_iext, id_jext, id_kext/) ) |
---|
| 674 | !2- add attribute to variable |
---|
| 675 | cl_extrap=fct_concat(td_var%c_extrap(:)) |
---|
| 676 | tl_att=att_init('extrapolation',cl_extrap) |
---|
| 677 | CALL var_move_att(td_var, tl_att) |
---|
| 678 | |
---|
[5832] | 679 | CALL att_clean(tl_att) |
---|
| 680 | |
---|
| 681 | CALL logger_info(" EXTRAP FILL: "//& |
---|
[4213] | 682 | & TRIM(fct_str(SUM(il_detect(:,:,:))))//& |
---|
| 683 | & " point(s) to extrapolate " ) |
---|
| 684 | |
---|
| 685 | !3- extrapolate |
---|
[5832] | 686 | CALL extrap__3D(td_var%d_value(:,:,:,:), td_var%d_fill, & |
---|
| 687 | & il_detect(:,:,:), & |
---|
| 688 | & cd_method, id_radius, id_maxiter ) |
---|
[4213] | 689 | |
---|
| 690 | DEALLOCATE(il_detect) |
---|
| 691 | |
---|
| 692 | END SUBROUTINE extrap__fill_value |
---|
| 693 | !------------------------------------------------------------------- |
---|
| 694 | !> @brief |
---|
[5832] | 695 | !> This subroutine compute point to be extrapolated in 3D array. |
---|
[4213] | 696 | !> |
---|
| 697 | !> @details |
---|
[5832] | 698 | !> in case of 'min_error' method:<br/> |
---|
| 699 | !> - compute derivative in i-, j- and k- direction |
---|
| 700 | !> - compute minimum error coefficient (distance to center of halo) |
---|
| 701 | !> - compute extrapolatd values by calculated minimum error using taylor expansion |
---|
| 702 | !> in case of 'dist_weight' method:<br/> |
---|
| 703 | !> - compute distance weight coefficient (inverse of distance to center of halo) |
---|
| 704 | !> - compute extrapolatd values using Inverse Distance Weighting |
---|
[4213] | 705 | !> |
---|
| 706 | !> @author J.Paul |
---|
| 707 | !> - Nov, 2013- Initial Version |
---|
| 708 | ! |
---|
[5832] | 709 | !> @param[inout] dd_value 3D array of variable to be extrapolated |
---|
| 710 | !> @param[in] dd_fill FillValue of variable |
---|
| 711 | !> @param[inout] id_detect array of point to extrapolate |
---|
| 712 | !> @param[in] cd_method extrapolation method |
---|
| 713 | !> @param[in] id_radius radius of the halo used to compute extrapolation |
---|
[4213] | 714 | !------------------------------------------------------------------- |
---|
| 715 | SUBROUTINE extrap__3D( dd_value, dd_fill, id_detect,& |
---|
[5832] | 716 | & cd_method, id_radius, id_maxiter ) |
---|
[4213] | 717 | IMPLICIT NONE |
---|
| 718 | ! Argument |
---|
[5832] | 719 | REAL(dp) , DIMENSION(:,:,:,:), INTENT(INOUT) :: dd_value |
---|
[4213] | 720 | REAL(dp) , INTENT(IN ) :: dd_fill |
---|
| 721 | INTEGER(i4), DIMENSION(:,:,:), INTENT(INOUT) :: id_detect |
---|
| 722 | CHARACTER(LEN=*), INTENT(IN ) :: cd_method |
---|
[5832] | 723 | INTEGER(i4), INTENT(IN ) :: id_radius |
---|
| 724 | INTEGER(i4), INTENT(IN ) :: id_maxiter |
---|
[4213] | 725 | |
---|
| 726 | ! local variable |
---|
| 727 | INTEGER(i4) :: il_imin |
---|
| 728 | INTEGER(i4) :: il_imax |
---|
| 729 | INTEGER(i4) :: il_jmin |
---|
| 730 | INTEGER(i4) :: il_jmax |
---|
| 731 | INTEGER(i4) :: il_kmin |
---|
| 732 | INTEGER(i4) :: il_kmax |
---|
[5832] | 733 | INTEGER(i4) :: il_iter |
---|
| 734 | INTEGER(i4) :: il_radius |
---|
[4213] | 735 | |
---|
[5832] | 736 | INTEGER(i4), DIMENSION(4) :: il_shape |
---|
[4213] | 737 | INTEGER(i4), DIMENSION(3) :: il_dim |
---|
| 738 | |
---|
[5832] | 739 | INTEGER(i4), DIMENSION(:,:,:), ALLOCATABLE :: il_detect |
---|
| 740 | |
---|
[4213] | 741 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_dfdx |
---|
| 742 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_dfdy |
---|
| 743 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_dfdz |
---|
| 744 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_coef |
---|
| 745 | |
---|
| 746 | ! loop indices |
---|
| 747 | INTEGER(i4) :: ji |
---|
| 748 | INTEGER(i4) :: jj |
---|
| 749 | INTEGER(i4) :: jk |
---|
[5832] | 750 | INTEGER(i4) :: jl |
---|
[4213] | 751 | !---------------------------------------------------------------- |
---|
| 752 | |
---|
| 753 | il_shape(:)=SHAPE(dd_value) |
---|
| 754 | |
---|
[5832] | 755 | ALLOCATE( il_detect( il_shape(1), il_shape(2), il_shape(3)) ) |
---|
| 756 | |
---|
[4213] | 757 | SELECT CASE(TRIM(cd_method)) |
---|
[5832] | 758 | CASE('min_error') |
---|
| 759 | DO jl=1,il_shape(4) |
---|
[4213] | 760 | |
---|
[5832] | 761 | ! intitialise table of poitn to be extrapolated |
---|
| 762 | il_detect(:,:,:)=id_detect(:,:,:) |
---|
[4213] | 763 | |
---|
[5832] | 764 | il_iter=1 |
---|
| 765 | DO WHILE( ANY(il_detect(:,:,:)==1) ) |
---|
| 766 | ! change extend value to minimize number of iteration |
---|
| 767 | il_radius=id_radius+(il_iter/id_maxiter) |
---|
[4213] | 768 | |
---|
[5832] | 769 | ALLOCATE( dl_dfdx(il_shape(1), il_shape(2), il_shape(3)) ) |
---|
| 770 | ALLOCATE( dl_dfdy(il_shape(1), il_shape(2), il_shape(3)) ) |
---|
| 771 | ALLOCATE( dl_dfdz(il_shape(1), il_shape(2), il_shape(3)) ) |
---|
[4213] | 772 | |
---|
[5832] | 773 | ! compute derivative in i-direction |
---|
| 774 | dl_dfdx(:,:,:)=dd_fill |
---|
| 775 | IF( il_shape(1) > 1 )THEN |
---|
| 776 | dl_dfdx(:,:,:)=extrap_deriv_3D( dd_value(:,:,:,jl), dd_fill, 'I' ) |
---|
| 777 | ENDIF |
---|
[4213] | 778 | |
---|
[5832] | 779 | ! compute derivative in j-direction |
---|
| 780 | dl_dfdy(:,:,:)=dd_fill |
---|
| 781 | IF( il_shape(2) > 1 )THEN |
---|
| 782 | dl_dfdy(:,:,:)=extrap_deriv_3D( dd_value(:,:,:,jl), dd_fill, 'J' ) |
---|
| 783 | ENDIF |
---|
| 784 | |
---|
| 785 | ! compute derivative in k-direction |
---|
| 786 | dl_dfdz(:,:,:)=dd_fill |
---|
| 787 | IF( il_shape(3) > 1 )THEN |
---|
| 788 | dl_dfdz(:,:,:)=extrap_deriv_3D( dd_value(:,:,:,jl), dd_fill, 'K' ) |
---|
| 789 | ENDIF |
---|
| 790 | |
---|
| 791 | il_dim(1)=2*il_radius+1 |
---|
| 792 | IF( il_shape(1) < 2*il_radius+1 ) il_dim(1)=1 |
---|
| 793 | il_dim(2)=2*il_radius+1 |
---|
| 794 | IF( il_shape(2) < 2*il_radius+1 ) il_dim(2)=1 |
---|
| 795 | il_dim(3)=2*il_radius+1 |
---|
| 796 | IF( il_shape(3) < 2*il_radius+1 ) il_dim(3)=1 |
---|
[4213] | 797 | |
---|
[5832] | 798 | ALLOCATE( dl_coef(il_dim(1), il_dim(2), il_dim(3)) ) |
---|
[4213] | 799 | |
---|
[5832] | 800 | dl_coef(:,:,:)=extrap__3D_min_error_coef(dd_value( 1:il_dim(1), & |
---|
| 801 | & 1:il_dim(2), & |
---|
| 802 | & 1:il_dim(3), & |
---|
| 803 | & jl )) |
---|
[4213] | 804 | |
---|
[5832] | 805 | DO jk=1,il_shape(3) |
---|
| 806 | IF( ALL(il_detect(:,:,jk) == 0) ) CYCLE |
---|
| 807 | DO jj=1,il_shape(2) |
---|
| 808 | IF( ALL(il_detect(:,jj,jk) == 0) ) CYCLE |
---|
| 809 | DO ji=1,il_shape(1) |
---|
[4213] | 810 | |
---|
[5832] | 811 | IF( il_detect(ji,jj,jk) == 1 )THEN |
---|
| 812 | |
---|
| 813 | il_imin=MAX(ji-il_radius,1) |
---|
| 814 | il_imax=MIN(ji+il_radius,il_shape(1)) |
---|
| 815 | IF( il_dim(1) == 1 )THEN |
---|
| 816 | il_imin=ji |
---|
| 817 | il_imax=ji |
---|
| 818 | ENDIF |
---|
[4213] | 819 | |
---|
[5832] | 820 | il_jmin=MAX(jj-il_radius,1) |
---|
| 821 | il_jmax=MIN(jj+il_radius,il_shape(2)) |
---|
| 822 | IF( il_dim(2) == 1 )THEN |
---|
| 823 | il_jmin=jj |
---|
| 824 | il_jmax=jj |
---|
| 825 | ENDIF |
---|
[4213] | 826 | |
---|
[5832] | 827 | il_kmin=MAX(jk-il_radius,1) |
---|
| 828 | il_kmax=MIN(jk+il_radius,il_shape(3)) |
---|
| 829 | IF( il_dim(3) == 1 )THEN |
---|
| 830 | il_kmin=jk |
---|
| 831 | il_kmax=jk |
---|
| 832 | ENDIF |
---|
[4213] | 833 | |
---|
[5832] | 834 | dd_value(ji,jj,jk,jl)=extrap__3D_min_error_fill( & |
---|
| 835 | & dd_value( il_imin:il_imax, & |
---|
| 836 | & il_jmin:il_jmax, & |
---|
| 837 | & il_kmin:il_kmax,jl ), dd_fill, il_radius, & |
---|
| 838 | & dl_dfdx( il_imin:il_imax, & |
---|
| 839 | & il_jmin:il_jmax, & |
---|
| 840 | & il_kmin:il_kmax ), & |
---|
| 841 | & dl_dfdy( il_imin:il_imax, & |
---|
| 842 | & il_jmin:il_jmax, & |
---|
| 843 | & il_kmin:il_kmax ), & |
---|
| 844 | & dl_dfdz( il_imin:il_imax, & |
---|
| 845 | & il_jmin:il_jmax, & |
---|
| 846 | & il_kmin:il_kmax ), & |
---|
| 847 | & dl_coef(:,:,:) ) |
---|
| 848 | |
---|
| 849 | IF( dd_value(ji,jj,jk,jl) /= dd_fill )THEN |
---|
| 850 | il_detect(ji,jj,jk)= 0 |
---|
| 851 | ENDIF |
---|
| 852 | |
---|
[4213] | 853 | ENDIF |
---|
| 854 | |
---|
[5832] | 855 | ENDDO |
---|
[4213] | 856 | ENDDO |
---|
| 857 | ENDDO |
---|
[5832] | 858 | |
---|
| 859 | DEALLOCATE( dl_dfdx ) |
---|
| 860 | DEALLOCATE( dl_dfdy ) |
---|
| 861 | DEALLOCATE( dl_dfdz ) |
---|
| 862 | DEALLOCATE( dl_coef ) |
---|
| 863 | |
---|
| 864 | il_iter=il_iter+1 |
---|
[4213] | 865 | ENDDO |
---|
[5832] | 866 | ENDDO |
---|
[4213] | 867 | |
---|
[5832] | 868 | CASE DEFAULT ! 'dist_weight' |
---|
| 869 | DO jl=1,il_shape(4) |
---|
[4213] | 870 | |
---|
[5832] | 871 | ! intitialise table of poitn to be extrapolated |
---|
| 872 | il_detect(:,:,:)=id_detect(:,:,:) |
---|
[4213] | 873 | |
---|
[5832] | 874 | il_iter=1 |
---|
| 875 | DO WHILE( ANY(il_detect(:,:,:)==1) ) |
---|
| 876 | ! change extend value to minimize number of iteration |
---|
| 877 | il_radius=id_radius+(il_iter/id_maxiter) |
---|
[4213] | 878 | |
---|
[5832] | 879 | il_dim(1)=2*il_radius+1 |
---|
| 880 | IF( il_shape(1) < 2*il_radius+1 ) il_dim(1)=1 |
---|
| 881 | il_dim(2)=2*il_radius+1 |
---|
| 882 | IF( il_shape(2) < 2*il_radius+1 ) il_dim(2)=1 |
---|
| 883 | il_dim(3)=2*il_radius+1 |
---|
| 884 | IF( il_shape(3) < 2*il_radius+1 ) il_dim(3)=1 |
---|
| 885 | |
---|
| 886 | ALLOCATE( dl_coef(il_dim(1), il_dim(2), il_dim(3)) ) |
---|
[4213] | 887 | |
---|
[5832] | 888 | dl_coef(:,:,:)=extrap__3D_dist_weight_coef(dd_value(1:il_dim(1), & |
---|
| 889 | & 1:il_dim(2), & |
---|
| 890 | & 1:il_dim(3), & |
---|
| 891 | & jl ) ) |
---|
[4213] | 892 | |
---|
[5832] | 893 | DO jk=1,il_shape(3) |
---|
| 894 | IF( ALL(il_detect(:,:,jk) == 0) ) CYCLE |
---|
| 895 | DO jj=1,il_shape(2) |
---|
| 896 | IF( ALL(il_detect(:,jj,jk) == 0) ) CYCLE |
---|
| 897 | DO ji=1,il_shape(1) |
---|
[4213] | 898 | |
---|
[5832] | 899 | IF( il_detect(ji,jj,jk) == 1 )THEN |
---|
| 900 | |
---|
| 901 | il_imin=MAX(ji-il_radius,1) |
---|
| 902 | il_imax=MIN(ji+il_radius,il_shape(1)) |
---|
| 903 | IF( il_dim(1) == 1 )THEN |
---|
| 904 | il_imin=ji |
---|
| 905 | il_imax=ji |
---|
| 906 | ENDIF |
---|
[4213] | 907 | |
---|
[5832] | 908 | il_jmin=MAX(jj-il_radius,1) |
---|
| 909 | il_jmax=MIN(jj+il_radius,il_shape(2)) |
---|
| 910 | IF( il_dim(2) == 1 )THEN |
---|
| 911 | il_jmin=jj |
---|
| 912 | il_jmax=jj |
---|
| 913 | ENDIF |
---|
[4213] | 914 | |
---|
[5832] | 915 | il_kmin=MAX(jk-il_radius,1) |
---|
| 916 | il_kmax=MIN(jk+il_radius,il_shape(3)) |
---|
| 917 | IF( il_dim(3) == 1 )THEN |
---|
| 918 | il_kmin=jk |
---|
| 919 | il_kmax=jk |
---|
| 920 | ENDIF |
---|
[4213] | 921 | |
---|
[5832] | 922 | dd_value(ji,jj,jk,jl)=extrap__3D_dist_weight_fill( & |
---|
| 923 | & dd_value( il_imin:il_imax, & |
---|
| 924 | & il_jmin:il_jmax, & |
---|
| 925 | & il_kmin:il_kmax, & |
---|
| 926 | & jl), dd_fill, il_radius, & |
---|
| 927 | & dl_coef(:,:,:) ) |
---|
| 928 | |
---|
| 929 | IF( dd_value(ji,jj,jk,jl) /= dd_fill )THEN |
---|
| 930 | il_detect(ji,jj,jk)= 0 |
---|
| 931 | ENDIF |
---|
| 932 | |
---|
[4213] | 933 | ENDIF |
---|
| 934 | |
---|
[5832] | 935 | ENDDO |
---|
[4213] | 936 | ENDDO |
---|
| 937 | ENDDO |
---|
[5832] | 938 | |
---|
| 939 | DEALLOCATE( dl_coef ) |
---|
| 940 | il_iter=il_iter+1 |
---|
[4213] | 941 | ENDDO |
---|
[5832] | 942 | ENDDO |
---|
[4213] | 943 | END SELECT |
---|
| 944 | |
---|
[5832] | 945 | DEALLOCATE( il_detect ) |
---|
| 946 | |
---|
[4213] | 947 | END SUBROUTINE extrap__3D |
---|
| 948 | !------------------------------------------------------------------- |
---|
| 949 | !> @brief |
---|
[5832] | 950 | !> This function compute derivative of 1D array. |
---|
[4213] | 951 | !> |
---|
| 952 | !> @details |
---|
[5832] | 953 | !> optionaly you could specify to take into account east west discontinuity |
---|
| 954 | !> (-180° 180° or 0° 360° for longitude variable) |
---|
[4213] | 955 | !> |
---|
| 956 | !> @author J.Paul |
---|
[5832] | 957 | !> - November, 2013- Initial Version |
---|
[4213] | 958 | ! |
---|
[5832] | 959 | !> @param[in] dd_value 1D array of variable to be extrapolated |
---|
| 960 | !> @param[in] dd_fill FillValue of variable |
---|
| 961 | !> @param[in] ld_discont logical to take into account east west discontinuity |
---|
[4213] | 962 | !------------------------------------------------------------------- |
---|
| 963 | PURE FUNCTION extrap_deriv_1D( dd_value, dd_fill, ld_discont ) |
---|
| 964 | |
---|
| 965 | IMPLICIT NONE |
---|
| 966 | ! Argument |
---|
| 967 | REAL(dp) , DIMENSION(:), INTENT(IN) :: dd_value |
---|
| 968 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
| 969 | LOGICAL , INTENT(IN), OPTIONAL :: ld_discont |
---|
| 970 | |
---|
| 971 | ! function |
---|
| 972 | REAL(dp), DIMENSION(SIZE(dd_value,DIM=1) ) :: extrap_deriv_1D |
---|
| 973 | |
---|
| 974 | ! local variable |
---|
| 975 | INTEGER(i4) :: il_imin |
---|
| 976 | INTEGER(i4) :: il_imax |
---|
| 977 | INTEGER(i4), DIMENSION(1) :: il_shape |
---|
| 978 | |
---|
| 979 | REAL(dp) :: dl_min |
---|
| 980 | REAL(dp) :: dl_max |
---|
| 981 | REAL(dp) , DIMENSION(:), ALLOCATABLE :: dl_value |
---|
| 982 | |
---|
| 983 | LOGICAL :: ll_discont |
---|
| 984 | |
---|
| 985 | ! loop indices |
---|
| 986 | INTEGER(i4) :: ji |
---|
| 987 | |
---|
| 988 | INTEGER(i4) :: i1 |
---|
| 989 | INTEGER(i4) :: i2 |
---|
| 990 | !---------------------------------------------------------------- |
---|
| 991 | ! init |
---|
| 992 | extrap_deriv_1D(:)=dd_fill |
---|
| 993 | |
---|
| 994 | ll_discont=.FALSE. |
---|
| 995 | IF( PRESENT(ld_discont) ) ll_discont=ld_discont |
---|
| 996 | |
---|
| 997 | il_shape(:)=SHAPE(dd_value(:)) |
---|
| 998 | |
---|
| 999 | ALLOCATE( dl_value(3)) |
---|
| 1000 | |
---|
| 1001 | ! compute derivative in i-direction |
---|
| 1002 | DO ji=1,il_shape(1) |
---|
| 1003 | |
---|
| 1004 | il_imin=MAX(ji-1,1) |
---|
| 1005 | il_imax=MIN(ji+1,il_shape(1)) |
---|
| 1006 | |
---|
| 1007 | IF( il_imin==ji-1 .AND. il_imax==ji+1 )THEN |
---|
| 1008 | i1=1 ; i2=3 |
---|
| 1009 | ELSEIF( il_imin==ji .AND. il_imax==ji+1 )THEN |
---|
| 1010 | i1=1 ; i2=2 |
---|
| 1011 | ELSEIF( il_imin==ji-1 .AND. il_imax==ji )THEN |
---|
| 1012 | i1=2 ; i2=3 |
---|
| 1013 | ENDIF |
---|
| 1014 | |
---|
| 1015 | dl_value(i1:i2)=dd_value(il_imin:il_imax) |
---|
| 1016 | IF( il_imin == 1 )THEN |
---|
| 1017 | dl_value(:)=EOSHIFT( dl_value(:), & |
---|
| 1018 | & DIM=1, & |
---|
| 1019 | & SHIFT=-1, & |
---|
| 1020 | & BOUNDARY=dl_value(1) ) |
---|
| 1021 | ENDIF |
---|
| 1022 | IF( il_imax == il_shape(1) )THEN |
---|
| 1023 | dl_value(:)=EOSHIFT( dl_value(:), & |
---|
| 1024 | & DIM=1, & |
---|
| 1025 | & SHIFT=1, & |
---|
| 1026 | & BOUNDARY=dl_value(3)) |
---|
| 1027 | ENDIF |
---|
| 1028 | |
---|
| 1029 | IF( ll_discont )THEN |
---|
| 1030 | dl_min=MINVAL( dl_value(:), dl_value(:)/=dd_fill ) |
---|
| 1031 | dl_max=MAXVAL( dl_value(:), dl_value(:)/=dd_fill ) |
---|
| 1032 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1033 | WHERE( dl_value(:) < 0._dp ) |
---|
| 1034 | dl_value(:) = dl_value(:)+360._dp |
---|
| 1035 | END WHERE |
---|
| 1036 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1037 | WHERE( dl_value(:) > 180._dp ) |
---|
| 1038 | dl_value(:) = dl_value(:)-180._dp |
---|
| 1039 | END WHERE |
---|
| 1040 | ENDIF |
---|
| 1041 | ENDIF |
---|
| 1042 | |
---|
| 1043 | IF( dl_value( 2) /= dd_fill .AND. & ! ji |
---|
| 1044 | & dl_value( 3) /= dd_fill .AND. & ! ji+1 |
---|
| 1045 | & dl_value( 1) /= dd_fill )THEN ! ji-1 |
---|
| 1046 | |
---|
| 1047 | extrap_deriv_1D(ji)=& |
---|
| 1048 | & ( dl_value(3) - dl_value(1) ) / & |
---|
| 1049 | & REAL( il_imax-il_imin ,dp) |
---|
| 1050 | |
---|
| 1051 | ENDIF |
---|
| 1052 | |
---|
| 1053 | ENDDO |
---|
| 1054 | |
---|
| 1055 | DEALLOCATE( dl_value ) |
---|
| 1056 | |
---|
| 1057 | END FUNCTION extrap_deriv_1D |
---|
| 1058 | !------------------------------------------------------------------- |
---|
| 1059 | !> @brief |
---|
[5832] | 1060 | !> This function compute derivative of 2D array. |
---|
| 1061 | !> you have to specify in which direction derivative have to be computed: |
---|
| 1062 | !> first (I) or second (J) dimension. |
---|
| 1063 | !> |
---|
[4213] | 1064 | !> @details |
---|
[5832] | 1065 | !> optionaly you could specify to take into account east west discontinuity |
---|
| 1066 | !> (-180° 180° or 0° 360° for longitude variable) |
---|
[4213] | 1067 | !> |
---|
| 1068 | !> @author J.Paul |
---|
[5832] | 1069 | !> - November, 2013- Initial Version |
---|
[4213] | 1070 | ! |
---|
[5832] | 1071 | !> @param[in] dd_value 2D array of variable to be extrapolated |
---|
| 1072 | !> @param[in] dd_fill FillValue of variable |
---|
| 1073 | !> @param[in] cd_dim compute derivative on first (I) or second (J) dimension |
---|
| 1074 | !> @param[in] ld_discont logical to take into account east west discontinuity |
---|
[4213] | 1075 | !------------------------------------------------------------------- |
---|
| 1076 | FUNCTION extrap_deriv_2D( dd_value, dd_fill, cd_dim, ld_discont ) |
---|
| 1077 | |
---|
| 1078 | IMPLICIT NONE |
---|
| 1079 | ! Argument |
---|
| 1080 | REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_value |
---|
| 1081 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
| 1082 | CHARACTER(LEN=*) , INTENT(IN) :: cd_dim |
---|
| 1083 | LOGICAL , INTENT(IN), OPTIONAL :: ld_discont |
---|
| 1084 | |
---|
| 1085 | ! function |
---|
| 1086 | REAL(dp), DIMENSION(SIZE(dd_value,DIM=1), & |
---|
| 1087 | & SIZE(dd_value,DIM=2) ) :: extrap_deriv_2D |
---|
| 1088 | |
---|
| 1089 | ! local variable |
---|
| 1090 | INTEGER(i4) :: il_imin |
---|
| 1091 | INTEGER(i4) :: il_imax |
---|
| 1092 | INTEGER(i4) :: il_jmin |
---|
| 1093 | INTEGER(i4) :: il_jmax |
---|
| 1094 | INTEGER(i4), DIMENSION(2) :: il_shape |
---|
| 1095 | |
---|
| 1096 | REAL(dp) :: dl_min |
---|
| 1097 | REAL(dp) :: dl_max |
---|
| 1098 | REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_value |
---|
| 1099 | |
---|
| 1100 | LOGICAL :: ll_discont |
---|
| 1101 | |
---|
| 1102 | ! loop indices |
---|
| 1103 | INTEGER(i4) :: ji |
---|
| 1104 | INTEGER(i4) :: jj |
---|
| 1105 | |
---|
| 1106 | INTEGER(i4) :: i1 |
---|
| 1107 | INTEGER(i4) :: i2 |
---|
| 1108 | |
---|
| 1109 | INTEGER(i4) :: j1 |
---|
| 1110 | INTEGER(i4) :: j2 |
---|
| 1111 | !---------------------------------------------------------------- |
---|
| 1112 | ! init |
---|
| 1113 | extrap_deriv_2D(:,:)=dd_fill |
---|
| 1114 | |
---|
| 1115 | ll_discont=.FALSE. |
---|
| 1116 | IF( PRESENT(ld_discont) ) ll_discont=ld_discont |
---|
| 1117 | |
---|
| 1118 | il_shape(:)=SHAPE(dd_value(:,:)) |
---|
| 1119 | |
---|
| 1120 | SELECT CASE(TRIM(fct_upper(cd_dim))) |
---|
| 1121 | |
---|
| 1122 | CASE('I') |
---|
| 1123 | |
---|
| 1124 | ALLOCATE( dl_value(3,il_shape(2)) ) |
---|
| 1125 | ! compute derivative in i-direction |
---|
| 1126 | DO ji=1,il_shape(1) |
---|
| 1127 | |
---|
| 1128 | ! init |
---|
| 1129 | dl_value(:,:)=dd_fill |
---|
| 1130 | |
---|
| 1131 | il_imin=MAX(ji-1,1) |
---|
| 1132 | il_imax=MIN(ji+1,il_shape(1)) |
---|
| 1133 | |
---|
| 1134 | IF( il_imin==ji-1 .AND. il_imax==ji+1 )THEN |
---|
| 1135 | i1=1 ; i2=3 |
---|
| 1136 | ELSEIF( il_imin==ji .AND. il_imax==ji+1 )THEN |
---|
| 1137 | i1=1 ; i2=2 |
---|
| 1138 | ELSEIF( il_imin==ji-1 .AND. il_imax==ji )THEN |
---|
| 1139 | i1=2 ; i2=3 |
---|
| 1140 | ENDIF |
---|
| 1141 | |
---|
| 1142 | dl_value(i1:i2,:)=dd_value(il_imin:il_imax,:) |
---|
| 1143 | IF( il_imin == 1 )THEN |
---|
| 1144 | dl_value(:,:)=EOSHIFT( dl_value(:,:), & |
---|
| 1145 | & DIM=1, & |
---|
| 1146 | & SHIFT=-1, & |
---|
| 1147 | & BOUNDARY=dl_value(1,:) ) |
---|
| 1148 | ENDIF |
---|
| 1149 | IF( il_imax == il_shape(1) )THEN |
---|
| 1150 | dl_value(:,:)=EOSHIFT( dl_value(:,:), & |
---|
| 1151 | & DIM=1, & |
---|
| 1152 | & SHIFT=1, & |
---|
| 1153 | & BOUNDARY=dl_value(3,:)) |
---|
| 1154 | ENDIF |
---|
| 1155 | |
---|
| 1156 | IF( ll_discont )THEN |
---|
| 1157 | dl_min=MINVAL( dl_value(:,:), dl_value(:,:)/=dd_fill ) |
---|
| 1158 | dl_max=MAXVAL( dl_value(:,:), dl_value(:,:)/=dd_fill ) |
---|
| 1159 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1160 | WHERE( dl_value(:,:) < 0_dp ) |
---|
| 1161 | dl_value(:,:) = dl_value(:,:)+360._dp |
---|
| 1162 | END WHERE |
---|
| 1163 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1164 | WHERE( dl_value(:,:) > 180 ) |
---|
| 1165 | dl_value(:,:) = dl_value(:,:)-180._dp |
---|
| 1166 | END WHERE |
---|
| 1167 | ENDIF |
---|
| 1168 | ENDIF |
---|
| 1169 | |
---|
| 1170 | WHERE( dl_value(2,:) /= dd_fill .AND. & ! ji |
---|
| 1171 | & dl_value(3,:) /= dd_fill .AND. & ! ji+1 |
---|
| 1172 | & dl_value(1,:) /= dd_fill ) ! ji-1 |
---|
| 1173 | |
---|
| 1174 | extrap_deriv_2D(ji,:)=& |
---|
| 1175 | & ( dl_value(3,:) - dl_value(1,:) ) / & |
---|
| 1176 | & REAL( il_imax-il_imin,dp) |
---|
| 1177 | |
---|
| 1178 | END WHERE |
---|
| 1179 | |
---|
[5832] | 1180 | ENDDO |
---|
[4213] | 1181 | |
---|
| 1182 | CASE('J') |
---|
| 1183 | |
---|
| 1184 | ALLOCATE( dl_value(il_shape(1),3) ) |
---|
| 1185 | ! compute derivative in j-direction |
---|
| 1186 | DO jj=1,il_shape(2) |
---|
| 1187 | |
---|
| 1188 | il_jmin=MAX(jj-1,1) |
---|
| 1189 | il_jmax=MIN(jj+1,il_shape(2)) |
---|
| 1190 | |
---|
| 1191 | IF( il_jmin==jj-1 .AND. il_jmax==jj+1 )THEN |
---|
| 1192 | j1=1 ; j2=3 |
---|
| 1193 | ELSEIF( il_jmin==jj .AND. il_jmax==jj+1 )THEN |
---|
| 1194 | j1=1 ; j2=2 |
---|
| 1195 | ELSEIF( il_jmin==jj-1 .AND. il_jmax==jj )THEN |
---|
| 1196 | j1=2 ; j2=3 |
---|
| 1197 | ENDIF |
---|
| 1198 | |
---|
| 1199 | dl_value(:,j1:j2)=dd_value(:,il_jmin:il_jmax) |
---|
| 1200 | IF( il_jmin == 1 )THEN |
---|
| 1201 | dl_value(:,:)=EOSHIFT( dl_value(:,:), & |
---|
| 1202 | & DIM=2, & |
---|
| 1203 | & SHIFT=-1, & |
---|
| 1204 | & BOUNDARY=dl_value(:,1)) |
---|
| 1205 | ENDIF |
---|
| 1206 | IF( il_jmax == il_shape(2) )THEN |
---|
| 1207 | dl_value(:,:)=EOSHIFT( dl_value(:,:), & |
---|
| 1208 | & DIM=2, & |
---|
| 1209 | & SHIFT=1, & |
---|
| 1210 | & BOUNDARY=dl_value(:,3)) |
---|
| 1211 | ENDIF |
---|
| 1212 | |
---|
| 1213 | IF( ll_discont )THEN |
---|
| 1214 | dl_min=MINVAL( dl_value(:,:), dl_value(:,:)/=dd_fill ) |
---|
| 1215 | dl_max=MAXVAL( dl_value(:,:), dl_value(:,:)/=dd_fill ) |
---|
| 1216 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1217 | WHERE( dl_value(:,:) < 0_dp ) |
---|
| 1218 | dl_value(:,:) = dl_value(:,:)+360._dp |
---|
| 1219 | END WHERE |
---|
| 1220 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1221 | WHERE( dl_value(:,:) > 180 ) |
---|
| 1222 | dl_value(:,:) = dl_value(:,:)-180._dp |
---|
| 1223 | END WHERE |
---|
| 1224 | ENDIF |
---|
| 1225 | ENDIF |
---|
| 1226 | |
---|
| 1227 | WHERE( dl_value(:, 2) /= dd_fill .AND. & ! jj |
---|
| 1228 | & dl_value(:, 3) /= dd_fill .AND. & ! jj+1 |
---|
| 1229 | & dl_value(:, 1) /= dd_fill ) ! jj-1 |
---|
| 1230 | |
---|
| 1231 | extrap_deriv_2D(:,jj)=& |
---|
| 1232 | & ( dl_value(:,3) - dl_value(:,1) ) / & |
---|
| 1233 | & REAL(il_jmax-il_jmin,dp) |
---|
| 1234 | |
---|
| 1235 | END WHERE |
---|
| 1236 | |
---|
| 1237 | ENDDO |
---|
| 1238 | |
---|
| 1239 | END SELECT |
---|
| 1240 | |
---|
| 1241 | DEALLOCATE( dl_value ) |
---|
| 1242 | |
---|
| 1243 | END FUNCTION extrap_deriv_2D |
---|
| 1244 | !------------------------------------------------------------------- |
---|
| 1245 | !> @brief |
---|
[5832] | 1246 | !> This function compute derivative of 3D array. |
---|
| 1247 | !> you have to specify in which direction derivative have to be computed: |
---|
| 1248 | !> first (I), second (J) or third (K) dimension. |
---|
[4213] | 1249 | !> |
---|
| 1250 | !> @details |
---|
[5832] | 1251 | !> optionaly you could specify to take into account east west discontinuity |
---|
| 1252 | !> (-180° 180° or 0° 360° for longitude variable) |
---|
[4213] | 1253 | !> |
---|
| 1254 | !> @author J.Paul |
---|
[5832] | 1255 | !> - November, 2013- Initial Version |
---|
[4213] | 1256 | ! |
---|
[5832] | 1257 | !> @param[inout] dd_value 3D array of variable to be extrapolated |
---|
| 1258 | !> @param[in] dd_fill FillValue of variable |
---|
| 1259 | !> @param[in] cd_dim compute derivative on first (I) second (J) or third (K) dimension |
---|
| 1260 | !> @param[in] ld_discont logical to take into account east west discontinuity |
---|
[4213] | 1261 | !------------------------------------------------------------------- |
---|
| 1262 | PURE FUNCTION extrap_deriv_3D( dd_value, dd_fill, cd_dim, ld_discont ) |
---|
| 1263 | |
---|
| 1264 | IMPLICIT NONE |
---|
| 1265 | ! Argument |
---|
| 1266 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_value |
---|
| 1267 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
| 1268 | CHARACTER(LEN=*) , INTENT(IN) :: cd_dim |
---|
| 1269 | LOGICAL , INTENT(IN), OPTIONAL :: ld_discont |
---|
| 1270 | |
---|
| 1271 | ! function |
---|
| 1272 | REAL(dp), DIMENSION(SIZE(dd_value,DIM=1), & |
---|
| 1273 | & SIZE(dd_value,DIM=2), & |
---|
| 1274 | & SIZE(dd_value,DIM=3)) :: extrap_deriv_3D |
---|
| 1275 | |
---|
| 1276 | ! local variable |
---|
| 1277 | INTEGER(i4) :: il_imin |
---|
| 1278 | INTEGER(i4) :: il_imax |
---|
| 1279 | INTEGER(i4) :: il_jmin |
---|
| 1280 | INTEGER(i4) :: il_jmax |
---|
| 1281 | INTEGER(i4) :: il_kmin |
---|
| 1282 | INTEGER(i4) :: il_kmax |
---|
| 1283 | INTEGER(i4), DIMENSION(3) :: il_shape |
---|
| 1284 | |
---|
| 1285 | REAL(dp) :: dl_min |
---|
| 1286 | REAL(dp) :: dl_max |
---|
| 1287 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_value |
---|
| 1288 | |
---|
| 1289 | LOGICAL :: ll_discont |
---|
| 1290 | |
---|
| 1291 | ! loop indices |
---|
| 1292 | INTEGER(i4) :: ji |
---|
| 1293 | INTEGER(i4) :: jj |
---|
| 1294 | INTEGER(i4) :: jk |
---|
| 1295 | |
---|
| 1296 | INTEGER(i4) :: i1 |
---|
| 1297 | INTEGER(i4) :: i2 |
---|
| 1298 | |
---|
| 1299 | INTEGER(i4) :: j1 |
---|
| 1300 | INTEGER(i4) :: j2 |
---|
| 1301 | |
---|
| 1302 | INTEGER(i4) :: k1 |
---|
| 1303 | INTEGER(i4) :: k2 |
---|
| 1304 | !---------------------------------------------------------------- |
---|
| 1305 | ! init |
---|
| 1306 | extrap_deriv_3D(:,:,:)=dd_fill |
---|
| 1307 | |
---|
| 1308 | ll_discont=.FALSE. |
---|
| 1309 | IF( PRESENT(ld_discont) ) ll_discont=ld_discont |
---|
| 1310 | |
---|
| 1311 | il_shape(:)=SHAPE(dd_value(:,:,:)) |
---|
| 1312 | |
---|
| 1313 | |
---|
| 1314 | SELECT CASE(TRIM(fct_upper(cd_dim))) |
---|
| 1315 | |
---|
| 1316 | CASE('I') |
---|
| 1317 | |
---|
| 1318 | ALLOCATE( dl_value(3,il_shape(2),il_shape(3)) ) |
---|
| 1319 | ! compute derivative in i-direction |
---|
| 1320 | DO ji=1,il_shape(1) |
---|
| 1321 | |
---|
| 1322 | il_imin=MAX(ji-1,1) |
---|
| 1323 | il_imax=MIN(ji+1,il_shape(1)) |
---|
| 1324 | |
---|
| 1325 | IF( il_imin==ji-1 .AND. il_imax==ji+1 )THEN |
---|
| 1326 | i1=1 ; i2=3 |
---|
| 1327 | ELSEIF( il_imin==ji .AND. il_imax==ji+1 )THEN |
---|
| 1328 | i1=1 ; i2=2 |
---|
| 1329 | ELSEIF( il_imin==ji-1 .AND. il_imax==ji )THEN |
---|
| 1330 | i1=2 ; i2=3 |
---|
| 1331 | ENDIF |
---|
| 1332 | |
---|
| 1333 | dl_value(i1:i2,:,:)=dd_value(il_imin:il_imax,:,:) |
---|
| 1334 | IF( il_imin == 1 )THEN |
---|
| 1335 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1336 | & DIM=1, & |
---|
| 1337 | & SHIFT=-1, & |
---|
| 1338 | & BOUNDARY=dl_value(1,:,:) ) |
---|
| 1339 | ENDIF |
---|
| 1340 | IF( il_imax == il_shape(1) )THEN |
---|
| 1341 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1342 | & DIM=1, & |
---|
| 1343 | & SHIFT=1, & |
---|
| 1344 | & BOUNDARY=dl_value(3,:,:)) |
---|
| 1345 | ENDIF |
---|
| 1346 | |
---|
| 1347 | IF( ll_discont )THEN |
---|
| 1348 | dl_min=MINVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1349 | dl_max=MAXVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1350 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1351 | WHERE( dl_value(:,:,:) < 0_dp ) |
---|
| 1352 | dl_value(:,:,:) = dl_value(:,:,:)+360._dp |
---|
| 1353 | END WHERE |
---|
| 1354 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1355 | WHERE( dl_value(:,:,:) > 180 ) |
---|
| 1356 | dl_value(:,:,:) = dl_value(:,:,:)-180._dp |
---|
| 1357 | END WHERE |
---|
| 1358 | ENDIF |
---|
| 1359 | ENDIF |
---|
| 1360 | |
---|
| 1361 | WHERE( dl_value(2,:,:) /= dd_fill .AND. & ! ji |
---|
| 1362 | & dl_value(3,:,:) /= dd_fill .AND. & !ji+1 |
---|
| 1363 | & dl_value(1,:,:) /= dd_fill ) !ji-1 |
---|
| 1364 | |
---|
| 1365 | extrap_deriv_3D(ji,:,:)= & |
---|
| 1366 | & ( dl_value(3,:,:) - dl_value(1,:,:) ) / & |
---|
| 1367 | & REAL( il_imax-il_imin ,dp) |
---|
| 1368 | |
---|
| 1369 | END WHERE |
---|
| 1370 | |
---|
| 1371 | ENDDO |
---|
| 1372 | |
---|
| 1373 | CASE('J') |
---|
| 1374 | |
---|
| 1375 | ALLOCATE( dl_value(il_shape(1),3,il_shape(3)) ) |
---|
| 1376 | ! compute derivative in j-direction |
---|
| 1377 | DO jj=1,il_shape(2) |
---|
| 1378 | |
---|
| 1379 | il_jmin=MAX(jj-1,1) |
---|
| 1380 | il_jmax=MIN(jj+1,il_shape(2)) |
---|
| 1381 | |
---|
| 1382 | IF( il_jmin==jj-1 .AND. il_jmax==jj+1 )THEN |
---|
| 1383 | j1=1 ; j2=3 |
---|
| 1384 | ELSEIF( il_jmin==jj .AND. il_jmax==jj+1 )THEN |
---|
| 1385 | j1=1 ; j2=2 |
---|
| 1386 | ELSEIF( il_jmin==jj-1 .AND. il_jmax==jj )THEN |
---|
| 1387 | j1=2 ; j2=3 |
---|
| 1388 | ENDIF |
---|
| 1389 | |
---|
| 1390 | dl_value(:,j1:j2,:)=dd_value(:,il_jmin:il_jmax,:) |
---|
| 1391 | IF( il_jmin == 1 )THEN |
---|
| 1392 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1393 | & DIM=2, & |
---|
| 1394 | & SHIFT=-1, & |
---|
| 1395 | & BOUNDARY=dl_value(:,1,:) ) |
---|
| 1396 | ENDIF |
---|
| 1397 | IF( il_jmax == il_shape(2) )THEN |
---|
| 1398 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1399 | & DIM=2, & |
---|
| 1400 | & SHIFT=1, & |
---|
| 1401 | & BOUNDARY=dl_value(:,3,:)) |
---|
| 1402 | ENDIF |
---|
| 1403 | |
---|
| 1404 | IF( ll_discont )THEN |
---|
| 1405 | dl_min=MINVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1406 | dl_max=MAXVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1407 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1408 | WHERE( dl_value(:,:,:) < 0_dp ) |
---|
| 1409 | dl_value(:,:,:) = dl_value(:,:,:)+360._dp |
---|
| 1410 | END WHERE |
---|
| 1411 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1412 | WHERE( dl_value(:,:,:) > 180 ) |
---|
| 1413 | dl_value(:,:,:) = dl_value(:,:,:)-180._dp |
---|
| 1414 | END WHERE |
---|
| 1415 | ENDIF |
---|
| 1416 | ENDIF |
---|
| 1417 | |
---|
| 1418 | WHERE( dl_value(:, 2,:) /= dd_fill .AND. & ! jj |
---|
| 1419 | & dl_value(:, 3,:) /= dd_fill .AND. & ! jj+1 |
---|
| 1420 | & dl_value(:, 1,:) /= dd_fill ) ! jj-1 |
---|
| 1421 | |
---|
| 1422 | extrap_deriv_3D(:,jj,:)=& |
---|
| 1423 | & ( dl_value(:,3,:) - dl_value(:,1,:) ) / & |
---|
| 1424 | & REAL( il_jmax - il_jmin ,dp) |
---|
| 1425 | |
---|
| 1426 | END WHERE |
---|
| 1427 | |
---|
| 1428 | ENDDO |
---|
| 1429 | |
---|
| 1430 | CASE('K') |
---|
| 1431 | ! compute derivative in k-direction |
---|
| 1432 | DO jk=1,il_shape(3) |
---|
| 1433 | |
---|
| 1434 | il_kmin=MAX(jk-1,1) |
---|
| 1435 | il_kmax=MIN(jk+1,il_shape(3)) |
---|
| 1436 | |
---|
| 1437 | IF( il_kmin==jk-1 .AND. il_kmax==jk+1 )THEN |
---|
| 1438 | k1=1 ; k2=3 |
---|
| 1439 | ELSEIF( il_kmin==jk .AND. il_kmax==jk+1 )THEN |
---|
| 1440 | k1=1 ; k2=2 |
---|
| 1441 | ELSEIF( il_kmin==jk-1 .AND. il_kmax==jk )THEN |
---|
| 1442 | k1=2 ; k2=3 |
---|
| 1443 | ENDIF |
---|
| 1444 | |
---|
| 1445 | dl_value(:,:,k1:k2)=dd_value(:,:,il_kmin:il_kmax) |
---|
| 1446 | IF( il_kmin == 1 )THEN |
---|
| 1447 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1448 | & DIM=3, & |
---|
| 1449 | & SHIFT=-1, & |
---|
| 1450 | & BOUNDARY=dl_value(:,:,1) ) |
---|
| 1451 | ENDIF |
---|
| 1452 | IF( il_kmax == il_shape(3) )THEN |
---|
| 1453 | dl_value(:,:,:)=EOSHIFT( dl_value(:,:,:), & |
---|
| 1454 | & DIM=3, & |
---|
| 1455 | & SHIFT=1, & |
---|
| 1456 | & BOUNDARY=dl_value(:,:,3)) |
---|
| 1457 | ENDIF |
---|
| 1458 | |
---|
| 1459 | IF( ll_discont )THEN |
---|
| 1460 | dl_min=MINVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1461 | dl_max=MAXVAL( dl_value(:,:,:), dl_value(:,:,:)/=dd_fill ) |
---|
| 1462 | IF( dl_min < -170_dp .AND. dl_max > 170_dp )THEN |
---|
| 1463 | WHERE( dl_value(:,:,:) < 0_dp ) |
---|
| 1464 | dl_value(:,:,:) = dl_value(:,:,:)+360._dp |
---|
| 1465 | END WHERE |
---|
| 1466 | ELSEIF( dl_min < 10_dp .AND. dl_max > 350_dp )THEN |
---|
| 1467 | WHERE( dl_value(:,:,:) > 180 ) |
---|
| 1468 | dl_value(:,:,:) = dl_value(:,:,:)-180._dp |
---|
| 1469 | END WHERE |
---|
| 1470 | ENDIF |
---|
| 1471 | ENDIF |
---|
| 1472 | |
---|
| 1473 | WHERE( dl_value(:,:, 2) /= dd_fill .AND. & ! jk |
---|
| 1474 | & dl_value(:,:, 3) /= dd_fill .AND. & ! jk+1 |
---|
| 1475 | & dl_value(:,:, 1) /= dd_fill ) ! jk-1 |
---|
| 1476 | |
---|
| 1477 | extrap_deriv_3D(:,:,jk)=& |
---|
| 1478 | & ( dl_value(:,:,3) - dl_value(:,:,1) ) / & |
---|
| 1479 | & REAL( il_kmax-il_kmin,dp) |
---|
| 1480 | |
---|
| 1481 | END WHERE |
---|
| 1482 | |
---|
| 1483 | ENDDO |
---|
| 1484 | |
---|
| 1485 | END SELECT |
---|
| 1486 | |
---|
| 1487 | DEALLOCATE( dl_value ) |
---|
| 1488 | |
---|
| 1489 | END FUNCTION extrap_deriv_3D |
---|
| 1490 | !------------------------------------------------------------------- |
---|
| 1491 | !> @brief |
---|
[5832] | 1492 | !> This function compute coefficient for min_error extrapolation. |
---|
[4213] | 1493 | !> |
---|
| 1494 | !> @details |
---|
[5832] | 1495 | !> coefficients are "grid distance" to the center of the box choosed to compute |
---|
| 1496 | !> extrapolation. |
---|
[4213] | 1497 | !> |
---|
| 1498 | !> @author J.Paul |
---|
[5832] | 1499 | !> - November, 2013- Initial Version |
---|
[4213] | 1500 | ! |
---|
[5832] | 1501 | !> @param[in] dd_value 3D array of variable to be extrapolated |
---|
| 1502 | !> @return 3D array of coefficient for minimum error extrapolation |
---|
[4213] | 1503 | !------------------------------------------------------------------- |
---|
| 1504 | PURE FUNCTION extrap__3D_min_error_coef( dd_value ) |
---|
| 1505 | |
---|
| 1506 | IMPLICIT NONE |
---|
| 1507 | ! Argument |
---|
| 1508 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_value |
---|
| 1509 | |
---|
| 1510 | ! function |
---|
| 1511 | REAL(dp), DIMENSION(SIZE(dd_value(:,:,:),DIM=1), & |
---|
| 1512 | & SIZE(dd_value(:,:,:),DIM=2), & |
---|
| 1513 | & SIZE(dd_value(:,:,:),DIM=3) ) :: extrap__3D_min_error_coef |
---|
| 1514 | |
---|
| 1515 | ! local variable |
---|
| 1516 | INTEGER(i4), DIMENSION(3) :: il_shape |
---|
| 1517 | |
---|
| 1518 | INTEGER(i4) :: il_imid |
---|
| 1519 | INTEGER(i4) :: il_jmid |
---|
| 1520 | INTEGER(i4) :: il_kmid |
---|
| 1521 | |
---|
| 1522 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_dist |
---|
| 1523 | |
---|
| 1524 | ! loop indices |
---|
| 1525 | INTEGER(i4) :: ji |
---|
| 1526 | INTEGER(i4) :: jj |
---|
| 1527 | INTEGER(i4) :: jk |
---|
| 1528 | !---------------------------------------------------------------- |
---|
| 1529 | |
---|
| 1530 | ! init |
---|
| 1531 | extrap__3D_min_error_coef(:,:,:)=0 |
---|
| 1532 | |
---|
| 1533 | il_shape(:)=SHAPE(dd_value(:,:,:)) |
---|
| 1534 | |
---|
| 1535 | il_imid=INT(REAL(il_shape(1),sp)*0.5+1) |
---|
| 1536 | il_jmid=INT(REAL(il_shape(2),sp)*0.5+1) |
---|
| 1537 | il_kmid=INT(REAL(il_shape(3),sp)*0.5+1) |
---|
| 1538 | |
---|
| 1539 | ALLOCATE( dl_dist(il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1540 | |
---|
| 1541 | DO jk=1,il_shape(3) |
---|
| 1542 | DO jj=1,il_shape(2) |
---|
| 1543 | DO ji=1,il_shape(1) |
---|
| 1544 | |
---|
| 1545 | ! compute distance |
---|
| 1546 | dl_dist(ji,jj,jk) = (ji-il_imid)**2 + & |
---|
| 1547 | & (jj-il_jmid)**2 + & |
---|
| 1548 | & (jk-il_kmid)**2 |
---|
| 1549 | |
---|
| 1550 | IF( dl_dist(ji,jj,jk) /= 0 )THEN |
---|
| 1551 | dl_dist(ji,jj,jk)=SQRT( dl_dist(ji,jj,jk) ) |
---|
| 1552 | ENDIF |
---|
| 1553 | |
---|
| 1554 | ENDDO |
---|
| 1555 | ENDDO |
---|
| 1556 | ENDDO |
---|
| 1557 | |
---|
| 1558 | WHERE( dl_dist(:,:,:) /= 0 ) |
---|
| 1559 | extrap__3D_min_error_coef(:,:,:)=dl_dist(:,:,:) |
---|
| 1560 | END WHERE |
---|
| 1561 | |
---|
| 1562 | DEALLOCATE( dl_dist ) |
---|
| 1563 | |
---|
| 1564 | END FUNCTION extrap__3D_min_error_coef |
---|
| 1565 | !------------------------------------------------------------------- |
---|
| 1566 | !> @brief |
---|
[5832] | 1567 | !> This function compute extrapolatd value by calculated minimum error using |
---|
[4213] | 1568 | !> taylor expansion |
---|
| 1569 | !> |
---|
[5832] | 1570 | !> @author J.Paul |
---|
| 1571 | !> - November, 2013- Initial Version |
---|
[4213] | 1572 | !> |
---|
[5832] | 1573 | !> @param[in] dd_value 3D array of variable to be extrapolated |
---|
| 1574 | !> @param[in] dd_fill FillValue of variable |
---|
| 1575 | !> @param[in] id_radius radius of the halo used to compute extrapolation |
---|
| 1576 | !> @param[in] dd_dfdx derivative of function in i-direction |
---|
| 1577 | !> @param[in] dd_dfdy derivative of function in j-direction |
---|
| 1578 | !> @param[in] dd_dfdz derivative of function in k-direction |
---|
| 1579 | !> @param[in] dd_coef array of coefficient for min_error extrapolation |
---|
| 1580 | !> @return extrapolatd value |
---|
[4213] | 1581 | !------------------------------------------------------------------- |
---|
[5832] | 1582 | PURE FUNCTION extrap__3D_min_error_fill( dd_value, dd_fill, id_radius, & |
---|
[4213] | 1583 | & dd_dfdx, dd_dfdy, dd_dfdz, & |
---|
| 1584 | & dd_coef ) |
---|
| 1585 | IMPLICIT NONE |
---|
| 1586 | ! Argument |
---|
| 1587 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_value |
---|
| 1588 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
[5832] | 1589 | INTEGER(i4), INTENT(IN) :: id_radius |
---|
[4213] | 1590 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_dfdx |
---|
| 1591 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_dfdy |
---|
| 1592 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_dfdz |
---|
| 1593 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_coef |
---|
| 1594 | |
---|
| 1595 | ! function |
---|
| 1596 | REAL(dp) :: extrap__3d_min_error_fill |
---|
| 1597 | |
---|
| 1598 | ! local variable |
---|
| 1599 | INTEGER(i4), DIMENSION(3) :: il_shape |
---|
| 1600 | INTEGER(i4), DIMENSION(3) :: il_ind |
---|
| 1601 | |
---|
| 1602 | INTEGER(i4), DIMENSION(:,:,:), ALLOCATABLE :: il_mask |
---|
| 1603 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_error |
---|
| 1604 | |
---|
| 1605 | INTEGER(i4) :: il_min |
---|
| 1606 | ! loop indices |
---|
| 1607 | |
---|
| 1608 | !---------------------------------------------------------------- |
---|
| 1609 | |
---|
| 1610 | ! init |
---|
| 1611 | extrap__3D_min_error_fill=dd_fill |
---|
| 1612 | |
---|
[5832] | 1613 | il_min=MAX(1,(SIZE(dd_value(:,:,:)))/(1+id_radius*2)) |
---|
[4213] | 1614 | |
---|
| 1615 | IF( COUNT(dd_value(:,:,:) /= dd_fill) >= il_min )THEN |
---|
| 1616 | |
---|
| 1617 | il_shape(:)=SHAPE(dd_value(:,:,:)) |
---|
| 1618 | ALLOCATE( il_mask( il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1619 | ALLOCATE( dl_error(il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1620 | |
---|
| 1621 | ! compute error |
---|
| 1622 | dl_error(:,:,:)=0. |
---|
| 1623 | il_mask(:,:,:)=0 |
---|
| 1624 | WHERE( dd_dfdx(:,:,:) /= dd_fill ) |
---|
| 1625 | dl_error(:,:,:)=dd_coef(:,:,:)*dd_dfdx(:,:,:) |
---|
| 1626 | il_mask(:,:,:)=1 |
---|
| 1627 | END WHERE |
---|
| 1628 | WHERE( dd_dfdy(:,:,:) /= dd_fill ) |
---|
| 1629 | dl_error(:,:,:)=(dl_error(:,:,:)+dd_coef(:,:,:)*dd_dfdy(:,:,:)) |
---|
| 1630 | il_mask(:,:,:)=1 |
---|
| 1631 | END WHERE |
---|
| 1632 | WHERE( dd_dfdz(:,:,:) /= dd_fill ) |
---|
| 1633 | dl_error(:,:,:)=(dl_error(:,:,:)+dd_coef(:,:,:)*dd_dfdz(:,:,:)) |
---|
| 1634 | il_mask(:,:,:)=1 |
---|
| 1635 | END WHERE |
---|
| 1636 | |
---|
| 1637 | ! get minimum error indices |
---|
| 1638 | il_ind(:)=MINLOC(dl_error(:,:,:),il_mask(:,:,:)==1) |
---|
| 1639 | |
---|
| 1640 | ! return value |
---|
| 1641 | IF( ALL(il_ind(:)/=0) )THEN |
---|
| 1642 | extrap__3D_min_error_fill=dd_value(il_ind(1),il_ind(2),il_ind(3)) |
---|
| 1643 | ENDIF |
---|
| 1644 | |
---|
| 1645 | DEALLOCATE( il_mask ) |
---|
| 1646 | DEALLOCATE( dl_error ) |
---|
| 1647 | |
---|
| 1648 | ENDIF |
---|
| 1649 | |
---|
| 1650 | END FUNCTION extrap__3D_min_error_fill |
---|
| 1651 | !------------------------------------------------------------------- |
---|
| 1652 | !> @brief |
---|
[5832] | 1653 | !> This function compute coefficient for inverse distance weighted method |
---|
[4213] | 1654 | !> |
---|
| 1655 | !> @details |
---|
[5832] | 1656 | !> coefficients are inverse "grid distance" to the center of the box choosed to compute |
---|
| 1657 | !> extrapolation. |
---|
[4213] | 1658 | !> |
---|
| 1659 | !> @author J.Paul |
---|
[5832] | 1660 | !> - November, 2013- Initial Version |
---|
[4213] | 1661 | ! |
---|
[5832] | 1662 | !> @param[in] dd_value 3D array of variable to be extrapolated |
---|
| 1663 | !> @return 3D array of coefficient for inverse distance weighted extrapolation |
---|
[4213] | 1664 | !------------------------------------------------------------------- |
---|
| 1665 | PURE FUNCTION extrap__3D_dist_weight_coef( dd_value ) |
---|
| 1666 | |
---|
| 1667 | IMPLICIT NONE |
---|
| 1668 | ! Argument |
---|
| 1669 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_value |
---|
| 1670 | |
---|
| 1671 | ! function |
---|
| 1672 | REAL(dp), DIMENSION(SIZE(dd_value(:,:,:),DIM=1), & |
---|
| 1673 | & SIZE(dd_value(:,:,:),DIM=2), & |
---|
| 1674 | & SIZE(dd_value(:,:,:),DIM=3) ) :: extrap__3D_dist_weight_coef |
---|
| 1675 | |
---|
| 1676 | ! local variable |
---|
| 1677 | INTEGER(i4), DIMENSION(3) :: il_shape |
---|
| 1678 | |
---|
| 1679 | INTEGER(i4) :: il_imid |
---|
| 1680 | INTEGER(i4) :: il_jmid |
---|
| 1681 | INTEGER(i4) :: il_kmid |
---|
| 1682 | |
---|
| 1683 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_dist |
---|
| 1684 | |
---|
| 1685 | ! loop indices |
---|
| 1686 | INTEGER(i4) :: ji |
---|
| 1687 | INTEGER(i4) :: jj |
---|
| 1688 | INTEGER(i4) :: jk |
---|
| 1689 | !---------------------------------------------------------------- |
---|
| 1690 | |
---|
| 1691 | ! init |
---|
| 1692 | extrap__3D_dist_weight_coef(:,:,:)=0 |
---|
| 1693 | |
---|
| 1694 | il_shape(:)=SHAPE(dd_value(:,:,:)) |
---|
| 1695 | |
---|
| 1696 | il_imid=INT(REAL(il_shape(1),sp)*0.5+1,i4) |
---|
| 1697 | il_jmid=INT(REAL(il_shape(2),sp)*0.5+1,i4) |
---|
| 1698 | il_kmid=INT(REAL(il_shape(3),sp)*0.5+1,i4) |
---|
| 1699 | |
---|
| 1700 | ALLOCATE( dl_dist(il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1701 | |
---|
| 1702 | DO jk=1,il_shape(3) |
---|
| 1703 | DO jj=1,il_shape(2) |
---|
| 1704 | DO ji=1,il_shape(1) |
---|
| 1705 | |
---|
| 1706 | ! compute distance |
---|
| 1707 | dl_dist(ji,jj,jk) = (ji-il_imid)**2 + & |
---|
| 1708 | & (jj-il_jmid)**2 + & |
---|
| 1709 | & (jk-il_kmid)**2 |
---|
| 1710 | |
---|
| 1711 | IF( dl_dist(ji,jj,jk) /= 0 )THEN |
---|
| 1712 | dl_dist(ji,jj,jk)=SQRT( dl_dist(ji,jj,jk) ) |
---|
| 1713 | ENDIF |
---|
| 1714 | |
---|
| 1715 | ENDDO |
---|
| 1716 | ENDDO |
---|
| 1717 | ENDDO |
---|
| 1718 | |
---|
| 1719 | WHERE( dl_dist(:,:,:) /= 0 ) |
---|
| 1720 | extrap__3D_dist_weight_coef(:,:,:)=1./dl_dist(:,:,:) |
---|
| 1721 | END WHERE |
---|
| 1722 | |
---|
| 1723 | DEALLOCATE( dl_dist ) |
---|
| 1724 | |
---|
| 1725 | END FUNCTION extrap__3D_dist_weight_coef |
---|
| 1726 | !------------------------------------------------------------------- |
---|
| 1727 | !> @brief |
---|
[5832] | 1728 | !> This function compute extrapolatd value using inverse distance weighted |
---|
| 1729 | !> method |
---|
[4213] | 1730 | !> |
---|
| 1731 | !> @details |
---|
| 1732 | !> |
---|
| 1733 | !> @author J.Paul |
---|
[5832] | 1734 | !> - November, 2013- Initial Version |
---|
[4213] | 1735 | ! |
---|
[5832] | 1736 | !> @param[in] dd_value 3D array of variable to be extrapolated |
---|
| 1737 | !> @param[in] dd_fill FillValue of variable |
---|
| 1738 | !> @param[in] id_radius radius of the halo used to compute extrapolation |
---|
| 1739 | !> @param[in] dd_coef 3D array of coefficient for inverse distance weighted extrapolation |
---|
| 1740 | !> @return extrapolatd value |
---|
[4213] | 1741 | !------------------------------------------------------------------- |
---|
[5832] | 1742 | FUNCTION extrap__3D_dist_weight_fill( dd_value, dd_fill, id_radius, & |
---|
| 1743 | & dd_coef ) |
---|
[4213] | 1744 | IMPLICIT NONE |
---|
| 1745 | ! Argument |
---|
| 1746 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_value |
---|
| 1747 | REAL(dp) , INTENT(IN) :: dd_fill |
---|
[5832] | 1748 | INTEGER(i4), INTENT(IN) :: id_radius |
---|
[4213] | 1749 | REAL(dp) , DIMENSION(:,:,:), INTENT(IN) :: dd_coef |
---|
| 1750 | |
---|
| 1751 | ! function |
---|
| 1752 | REAL(dp) :: extrap__3D_dist_weight_fill |
---|
| 1753 | |
---|
| 1754 | ! local variable |
---|
| 1755 | INTEGER(i4), DIMENSION(3) :: il_shape |
---|
| 1756 | |
---|
| 1757 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_value |
---|
| 1758 | REAL(dp) , DIMENSION(:,:,:), ALLOCATABLE :: dl_coef |
---|
| 1759 | |
---|
| 1760 | INTEGER(i4) :: il_min |
---|
| 1761 | ! loop indices |
---|
| 1762 | INTEGER(i4) :: ji |
---|
| 1763 | INTEGER(i4) :: jj |
---|
| 1764 | INTEGER(i4) :: jk |
---|
| 1765 | |
---|
| 1766 | !---------------------------------------------------------------- |
---|
| 1767 | |
---|
| 1768 | ! init |
---|
| 1769 | extrap__3D_dist_weight_fill=dd_fill |
---|
| 1770 | |
---|
[5832] | 1771 | il_min=MAX(1,(SIZE(dd_value(:,:,:)))/(1+id_radius*2)) |
---|
[4213] | 1772 | |
---|
| 1773 | IF( COUNT(dd_value(:,:,:)/= dd_fill) >= il_min )THEN |
---|
| 1774 | |
---|
| 1775 | il_shape(:)=SHAPE(dd_value(:,:,:)) |
---|
| 1776 | ALLOCATE( dl_value( il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1777 | ALLOCATE( dl_coef( il_shape(1),il_shape(2),il_shape(3)) ) |
---|
| 1778 | |
---|
| 1779 | dl_value(:,:,:)=0 |
---|
| 1780 | dl_coef(:,:,:)=0 |
---|
| 1781 | |
---|
[5832] | 1782 | DO jk=1,il_shape(3) |
---|
| 1783 | DO jj=1,il_shape(2) |
---|
| 1784 | DO ji=1,il_shape(1) |
---|
[4213] | 1785 | |
---|
[5832] | 1786 | IF( dd_value(ji,jj,jk) /= dd_fill )THEN |
---|
| 1787 | ! compute factor |
---|
| 1788 | dl_value(ji,jj,jk)=dd_coef(ji,jj,jk)*dd_value(ji,jj,jk) |
---|
| 1789 | dl_coef(ji,jj,jk)=dd_coef(ji,jj,jk) |
---|
| 1790 | ENDIF |
---|
[4213] | 1791 | |
---|
[5832] | 1792 | ENDDO |
---|
| 1793 | ENDDO |
---|
| 1794 | ENDDO |
---|
[4213] | 1795 | |
---|
| 1796 | ! return value |
---|
| 1797 | IF( SUM( dl_coef(:,:,:) ) /= 0 )THEN |
---|
[5832] | 1798 | extrap__3D_dist_weight_fill = & |
---|
| 1799 | & SUM( dl_value(:,:,:) )/SUM( dl_coef(:,:,:) ) |
---|
[4213] | 1800 | ENDIF |
---|
| 1801 | |
---|
| 1802 | DEALLOCATE( dl_value ) |
---|
| 1803 | DEALLOCATE( dl_coef ) |
---|
| 1804 | |
---|
| 1805 | ENDIF |
---|
| 1806 | |
---|
| 1807 | END FUNCTION extrap__3D_dist_weight_fill |
---|
| 1808 | !------------------------------------------------------------------- |
---|
| 1809 | !> @brief |
---|
| 1810 | !> This subroutine add to the variable (to be extrapolated) an |
---|
| 1811 | !> extraband of N points at north,south,east and west boundaries. |
---|
| 1812 | !> |
---|
[5832] | 1813 | !> @details |
---|
| 1814 | !> optionaly you could specify size of extra bands in i- and j-direction |
---|
| 1815 | !> |
---|
[4213] | 1816 | !> @author J.Paul |
---|
[5832] | 1817 | !> - November, 2013-Initial version |
---|
[4213] | 1818 | ! |
---|
[5832] | 1819 | !> @param[inout] td_var variable |
---|
| 1820 | !> @param[in] id_isize i-direction size of extra bands (default=im_minext) |
---|
| 1821 | !> @param[in] id_jsize j-direction size of extra bands (default=im_minext) |
---|
[4213] | 1822 | !> @todo |
---|
| 1823 | !> - invalid special case for grid with north fold |
---|
| 1824 | !------------------------------------------------------------------- |
---|
| 1825 | SUBROUTINE extrap_add_extrabands(td_var, id_isize, id_jsize ) |
---|
| 1826 | IMPLICIT NONE |
---|
| 1827 | ! Argument |
---|
| 1828 | TYPE(TVAR) , INTENT(INOUT) :: td_var |
---|
| 1829 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_isize |
---|
| 1830 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jsize |
---|
| 1831 | |
---|
| 1832 | ! local variable |
---|
| 1833 | REAL(dp), DIMENSION(:,:,:,:) , ALLOCATABLE :: dl_value |
---|
| 1834 | |
---|
| 1835 | INTEGER(i4) :: il_isize |
---|
| 1836 | INTEGER(i4) :: il_jsize |
---|
| 1837 | INTEGER(i4) :: il_tmp |
---|
| 1838 | |
---|
| 1839 | ! loop indices |
---|
| 1840 | INTEGER(i4) :: ji |
---|
| 1841 | INTEGER(i4) :: ij |
---|
| 1842 | !---------------------------------------------------------------- |
---|
| 1843 | il_isize=im_minext |
---|
| 1844 | IF(PRESENT(id_isize)) il_isize=id_isize |
---|
| 1845 | IF( il_isize < im_minext .AND. & |
---|
| 1846 | & TRIM(td_var%c_interp(1)) == 'cubic' )THEN |
---|
| 1847 | CALL logger_warn("EXTRAP ADD EXTRABANDS: size of extrabands "//& |
---|
| 1848 | & "should be at least "//TRIM(fct_str(im_minext))//" for "//& |
---|
| 1849 | & " cubic interpolation ") |
---|
| 1850 | ENDIF |
---|
| 1851 | |
---|
| 1852 | il_jsize=im_minext |
---|
| 1853 | IF(PRESENT(id_jsize)) il_jsize=id_jsize |
---|
| 1854 | IF( il_jsize < im_minext .AND. & |
---|
| 1855 | & TRIM(td_var%c_interp(1)) == 'cubic' )THEN |
---|
| 1856 | CALL logger_warn("EXTRAP ADD EXTRABANDS: size of extrabands "//& |
---|
| 1857 | & "should be at least "//TRIM(fct_str(im_minext))//" for "//& |
---|
| 1858 | & " cubic interpolation ") |
---|
| 1859 | ENDIF |
---|
| 1860 | |
---|
| 1861 | IF( .NOT. td_var%t_dim(1)%l_use ) il_isize=0 |
---|
| 1862 | IF( .NOT. td_var%t_dim(2)%l_use ) il_jsize=0 |
---|
| 1863 | |
---|
| 1864 | CALL logger_trace( "EXTRAP ADD EXTRABANDS: dimension change "//& |
---|
| 1865 | & "in variable "//TRIM(td_var%c_name) ) |
---|
| 1866 | |
---|
| 1867 | ! add extrabands in variable |
---|
| 1868 | ALLOCATE(dl_value( td_var%t_dim(1)%i_len, & |
---|
| 1869 | & td_var%t_dim(2)%i_len, & |
---|
| 1870 | & td_var%t_dim(3)%i_len, & |
---|
| 1871 | & td_var%t_dim(4)%i_len )) |
---|
| 1872 | |
---|
| 1873 | dl_value(:,:,:,:)=td_var%d_value(:,:,:,:) |
---|
| 1874 | |
---|
[5832] | 1875 | |
---|
[4213] | 1876 | td_var%t_dim(1)%i_len = td_var%t_dim(1)%i_len + 2*il_isize |
---|
| 1877 | td_var%t_dim(2)%i_len = td_var%t_dim(2)%i_len + 2*il_jsize |
---|
| 1878 | |
---|
| 1879 | DEALLOCATE(td_var%d_value) |
---|
| 1880 | ALLOCATE( td_var%d_value(td_var%t_dim(1)%i_len, & |
---|
| 1881 | & td_var%t_dim(2)%i_len, & |
---|
| 1882 | & td_var%t_dim(3)%i_len, & |
---|
| 1883 | & td_var%t_dim(4)%i_len ) ) |
---|
| 1884 | |
---|
| 1885 | ! intialise |
---|
| 1886 | td_var%d_value(:,:,:,:)=td_var%d_fill |
---|
| 1887 | |
---|
| 1888 | ! fill center |
---|
| 1889 | td_var%d_value( 1+il_isize:td_var%t_dim(1)%i_len-il_isize, & |
---|
| 1890 | & 1+il_jsize:td_var%t_dim(2)%i_len-il_jsize, & |
---|
| 1891 | & :,:) = dl_value(:,:,:,:) |
---|
| 1892 | |
---|
| 1893 | ! special case for overlap |
---|
| 1894 | IF( td_var%i_ew >= 0 .AND. il_isize /= 0 )THEN |
---|
| 1895 | DO ji=1,il_isize |
---|
| 1896 | ! from east to west |
---|
| 1897 | il_tmp=td_var%t_dim(1)%i_len-td_var%i_ew+ji-2*il_isize |
---|
| 1898 | td_var%d_value(ji,:,:,:) = td_var%d_value(il_tmp,:,:,:) |
---|
| 1899 | |
---|
| 1900 | ! from west to east |
---|
| 1901 | ij=td_var%t_dim(1)%i_len-ji+1 |
---|
| 1902 | il_tmp=td_var%i_ew-ji+2*il_isize+1 |
---|
| 1903 | td_var%d_value(ij,:,:,:) = td_var%d_value(il_tmp,:,:,:) |
---|
| 1904 | ENDDO |
---|
| 1905 | ENDIF |
---|
| 1906 | |
---|
| 1907 | DEALLOCATE( dl_value ) |
---|
| 1908 | |
---|
| 1909 | END SUBROUTINE extrap_add_extrabands |
---|
| 1910 | !------------------------------------------------------------------- |
---|
| 1911 | !> @brief |
---|
| 1912 | !> This subroutine remove of the variable an extraband |
---|
| 1913 | !> of N points at north,south,east and west boundaries. |
---|
| 1914 | !> |
---|
[5832] | 1915 | !> @details |
---|
| 1916 | !> optionaly you could specify size of extra bands in i- and j-direction |
---|
| 1917 | !> |
---|
[4213] | 1918 | !> @author J.Paul |
---|
[5832] | 1919 | !> - November, 2013-Initial version |
---|
| 1920 | !> |
---|
| 1921 | !> @param[inout] td_var variable |
---|
| 1922 | !> @param[in] id_isize i-direction size of extra bands (default=im_minext) |
---|
| 1923 | !> @param[in] id_jsize j-direction size of extra bands (default=im_minext) |
---|
[4213] | 1924 | !------------------------------------------------------------------- |
---|
| 1925 | SUBROUTINE extrap_del_extrabands(td_var, id_isize, id_jsize ) |
---|
| 1926 | IMPLICIT NONE |
---|
| 1927 | ! Argument |
---|
[5832] | 1928 | TYPE(TVAR) , INTENT(INOUT) :: td_var |
---|
[4213] | 1929 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_isize |
---|
| 1930 | INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jsize |
---|
| 1931 | |
---|
| 1932 | ! local variable |
---|
| 1933 | REAL(dp), DIMENSION(:,:,:,:) , ALLOCATABLE :: dl_value |
---|
| 1934 | |
---|
| 1935 | INTEGER(i4) :: il_isize |
---|
| 1936 | INTEGER(i4) :: il_jsize |
---|
[5832] | 1937 | |
---|
[4213] | 1938 | INTEGER(i4) :: il_imin |
---|
| 1939 | INTEGER(i4) :: il_imax |
---|
| 1940 | INTEGER(i4) :: il_jmin |
---|
| 1941 | INTEGER(i4) :: il_jmax |
---|
| 1942 | |
---|
| 1943 | ! loop indices |
---|
| 1944 | !---------------------------------------------------------------- |
---|
| 1945 | il_isize=im_minext |
---|
| 1946 | IF(PRESENT(id_isize)) il_isize=id_isize |
---|
| 1947 | |
---|
| 1948 | il_jsize=im_minext |
---|
| 1949 | IF(PRESENT(id_jsize)) il_jsize=id_jsize |
---|
| 1950 | |
---|
| 1951 | IF( .NOT. td_var%t_dim(1)%l_use ) il_isize=0 |
---|
| 1952 | IF( .NOT. td_var%t_dim(2)%l_use ) il_jsize=0 |
---|
| 1953 | |
---|
| 1954 | CALL logger_trace( "EXTRAP DEL EXTRABANDS: dimension change "//& |
---|
| 1955 | & "in variable "//TRIM(td_var%c_name) ) |
---|
| 1956 | |
---|
| 1957 | ! add extrabands in variable |
---|
| 1958 | ALLOCATE(dl_value( td_var%t_dim(1)%i_len, & |
---|
| 1959 | & td_var%t_dim(2)%i_len, & |
---|
| 1960 | & td_var%t_dim(3)%i_len, & |
---|
| 1961 | & td_var%t_dim(4)%i_len )) |
---|
| 1962 | |
---|
| 1963 | dl_value(:,:,:,:)=td_var%d_value(:,:,:,:) |
---|
| 1964 | |
---|
| 1965 | ! fill center |
---|
| 1966 | il_imin=1+il_isize |
---|
| 1967 | il_imax=td_var%t_dim(1)%i_len-il_isize |
---|
| 1968 | |
---|
| 1969 | il_jmin=1+il_jsize |
---|
| 1970 | il_jmax=td_var%t_dim(2)%i_len-il_jsize |
---|
| 1971 | |
---|
| 1972 | td_var%t_dim(1)%i_len = td_var%t_dim(1)%i_len - 2*il_isize |
---|
| 1973 | td_var%t_dim(2)%i_len = td_var%t_dim(2)%i_len - 2*il_jsize |
---|
| 1974 | |
---|
| 1975 | DEALLOCATE(td_var%d_value) |
---|
| 1976 | ALLOCATE( td_var%d_value(td_var%t_dim(1)%i_len, & |
---|
| 1977 | & td_var%t_dim(2)%i_len, & |
---|
| 1978 | & td_var%t_dim(3)%i_len, & |
---|
| 1979 | & td_var%t_dim(4)%i_len ) ) |
---|
| 1980 | |
---|
| 1981 | ! intialise |
---|
| 1982 | td_var%d_value(:,:,:,:)=td_var%d_fill |
---|
| 1983 | |
---|
| 1984 | td_var%d_value(:,:,:,:)=dl_value(il_imin:il_imax,& |
---|
| 1985 | & il_jmin:il_jmax,& |
---|
| 1986 | & :,:) |
---|
| 1987 | |
---|
| 1988 | DEALLOCATE( dl_value ) |
---|
| 1989 | |
---|
| 1990 | END SUBROUTINE extrap_del_extrabands |
---|
| 1991 | END MODULE extrap |
---|