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