Oasis3 4.0.2
|
00001 C**** 00002 C ***************************** 00003 C * OASIS MODULE - LEVEL ? * 00004 C * ------------- ------- * 00005 C ***************************** 00006 C 00007 C**** remap_bilinear_reduced - calculate reduced grid bilinear remapping 00008 C 00009 C Purpose: 00010 C ------- 00011 C Adaptation of SCRIP 1.4 remap_bilinear module to calculate 00012 C bilinear remapping for reduced grids. 00013 C 00014 C** Interface: 00015 C --------- 00016 C *CALL* *remap_bilin_reduced* 00017 C 00018 C Input: 00019 C ----- 00020 C 00021 C Output: 00022 C ------ 00023 C 00024 C History: 00025 C ------- 00026 C Version Programmer Date Description 00027 C ------- ---------- ---- ----------- 00028 C 2.5 D. Declat 2002/07 created 00029 C 2.5 S. Valcke 2002/09 modified 00030 C 00031 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00032 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00033 ! 00034 ! this module contains necessary routines for performing an 00035 ! bilinear interpolation on reduced grids. 00036 ! 00037 !----------------------------------------------------------------------- 00038 ! 00039 ! CVS:$Id: remap_bilinear_reduced.f 3500 2012-03-22 14:58:30Z coquart $ 00040 ! 00041 ! Copyright (c) 1997, 1998 the Regents of the University of 00042 ! California. 00043 ! 00044 ! This software and ancillary information (herein called software) 00045 ! called SCRIP is made available under the terms described here. 00046 ! The software has been approved for release with associated 00047 ! LA-CC Number 98-45. 00048 ! 00049 ! Unless otherwise indicated, this software has been authored 00050 ! by an employee or employees of the University of California, 00051 ! operator of the Los Alamos National Laboratory under Contract 00052 ! No. W-7405-ENG-36 with the U.S. Department of Energy. The U.S. 00053 ! Government has rights to use, reproduce, and distribute this 00054 ! software. The public may copy and use this software without 00055 ! charge, provided that this Notice and any statement of authorship 00056 ! are reproduced on all copies. Neither the Government nor the 00057 ! University makes any warranty, express or implied, or assumes 00058 ! any liability or responsibility for the use of this software. 00059 ! 00060 ! If software is modified to produce derivative works, such modified 00061 ! software should be clearly marked, so as not to confuse it with 00062 ! the version available from Los Alamos National Laboratory. 00063 ! 00064 !*********************************************************************** 00065 00066 module remap_bilinear_reduced 00067 00068 !----------------------------------------------------------------------- 00069 00070 use kinds_mod ! defines common data types 00071 use constants ! defines common constants 00072 use grids ! module containing grid info 00073 use remap_vars ! module containing remap info 00074 USE mod_unit 00075 USE mod_printing 00076 00077 implicit none 00078 00079 !----------------------------------------------------------------------- 00080 00081 integer (kind=int_kind), parameter :: 00082 & max_iter = 100 ! max iteration count for i,j iteration 00083 00084 real (kind=dbl_kind), parameter :: 00085 & converge = 1.e-10_dbl_kind ! convergence criterion 00086 00087 !*********************************************************************** 00088 00089 contains 00090 00091 !*********************************************************************** 00092 00093 subroutine remap_bilin_reduced (lextrapdone) 00094 00095 !----------------------------------------------------------------------- 00096 ! 00097 ! this routine computes the weights for a bilinear interpolation. 00098 ! 00099 !----------------------------------------------------------------------- 00100 00101 LOGICAL :: 00102 & lextrapdone, ! logical, true if EXTRAP done on field 00103 & ll_nnei ! true if extra search is done 00104 !----------------------------------------------------------------------- 00105 ! 00106 ! local variables 00107 ! 00108 !----------------------------------------------------------------------- 00109 00110 integer (kind=int_kind) :: n, icount, i, 00111 & dst_add, ! destination address 00112 & iter, ! iteration counter 00113 & nmap ! index of current map being computed 00114 00115 integer (kind=int_kind), dimension(4) :: 00116 & src_add ! address for the four source points 00117 00118 real (kind=dbl_kind), dimension(4) :: 00119 & src_lats, ! latitudes of four bilinear corners 00120 & src_lons, ! longitudes of four bilinear corners 00121 & wgts ! bilinear weights for four corners 00122 00123 real (kind=dbl_kind) :: 00124 & plat, plon, ! lat/lon coords of destination point 00125 & iguess, jguess, ! current guess for bilinear coordinate 00126 & deli, delj, ! corrections to i,j 00127 & dth1, dth2, dth3, ! some latitude differences 00128 & dph1, dph2, dph3, ! some longitude differences 00129 & dthp, dphp, ! difference between point and sw corner 00130 & mat1, mat2, mat3, mat4, ! matrix elements 00131 & determinant, ! matrix determinant 00132 & sum_wgts ! sum of weights for normalization 00133 00134 real (kind=dbl_kind) :: 00135 & coslat_dst, sinlat_dst, coslon_dst, sinlon_dst, 00136 & dist_min, distance, ! for computing dist-weighted avg 00137 & src_latsnn, arg 00138 00139 integer (kind=int_kind) :: min_add, max_add, srch_add, src_addnn 00140 00141 !----------------------------------------------------------------------- 00142 ! 00143 ! compute mappings from grid1 to grid2 00144 ! 00145 !----------------------------------------------------------------------- 00146 ! 00147 IF (nlogprt .GE. 2) THEN 00148 WRITE (UNIT = nulou,FMT = *)' ' 00149 WRITE (UNIT = nulou,FMT = *) 00150 & 'Entering routine remap_bilin_reduced' 00151 CALL FLUSH(nulou) 00152 ENDIF 00153 ! 00154 ll_nnei = .true. 00155 nmap = 1 00156 if (grid1_rank /= 2) then 00157 stop 'Can not do bilinear interpolation when grid_rank /= 2' 00158 endif 00159 00160 !*** 00161 !*** loop over destination grid 00162 !*** 00163 00164 grid_loop1: do dst_add = 1, grid2_size 00165 00166 if (.not. grid2_mask(dst_add)) cycle grid_loop1 00167 00168 plat = grid2_center_lat(dst_add) 00169 plon = grid2_center_lon(dst_add) 00170 00171 !*** 00172 !*** find nearest square of grid points on source grid 00173 !*** 00174 00175 call grid_search_bilin_rd(src_add, src_lats, src_lons, 00176 & min_add, max_add, 00177 & plat, plon, grid1_dims, 00178 & grid1_center_lat, grid1_center_lon, 00179 & lextrapdone) 00180 if (src_add(1) > 0) THEN 00181 00182 !*** 00183 !*** if the 4 surrounding points have been found and are 00184 !*** non-masked, do bilinear interpolation 00185 !*** 00186 00187 grid2_frac(dst_add) = one 00188 00189 !*** 00190 !*** iterate to find i,j for bilinear approximation 00191 !*** 00192 00193 dth1 = src_lats(2) - src_lats(1) 00194 dth2 = src_lats(4) - src_lats(1) 00195 dth3 = src_lats(3) - src_lats(2) - dth2 00196 00197 dph1 = src_lons(2) - src_lons(1) 00198 dph2 = src_lons(4) - src_lons(1) 00199 dph3 = src_lons(3) - src_lons(2) 00200 00201 if (dph1 > three*pih) dph1 = dph1 - pi2 00202 if (dph2 > three*pih) dph2 = dph2 - pi2 00203 if (dph3 > three*pih) dph3 = dph3 - pi2 00204 if (dph1 < -three*pih) dph1 = dph1 + pi2 00205 if (dph2 < -three*pih) dph2 = dph2 + pi2 00206 if (dph3 < -three*pih) dph3 = dph3 + pi2 00207 00208 dph3 = dph3 - dph2 00209 00210 iguess = half 00211 jguess = half 00212 00213 iter_loop1: do iter=1,max_iter 00214 00215 dthp = plat - src_lats(1) - dth1*iguess - 00216 & dth2*jguess - dth3*iguess*jguess 00217 dphp = plon - src_lons(1) 00218 00219 if (dphp > three*pih) dphp = dphp - pi2 00220 if (dphp < -three*pih) dphp = dphp + pi2 00221 00222 dphp = dphp - dph1*iguess - dph2*jguess - 00223 & dph3*iguess*jguess 00224 00225 mat1 = dth1 + dth3*jguess 00226 mat2 = dth2 + dth3*iguess 00227 mat3 = dph1 + dph3*jguess 00228 mat4 = dph2 + dph3*iguess 00229 00230 determinant = mat1*mat4 - mat2*mat3 00231 00232 deli = (dthp*mat4 - mat2*dphp)/determinant 00233 delj = (mat1*dphp - dthp*mat3)/determinant 00234 00235 if (abs(deli) < converge .and. 00236 & abs(delj) < converge) exit iter_loop1 00237 00238 iguess = iguess + deli 00239 jguess = jguess + delj 00240 00241 end do iter_loop1 00242 00243 if (iter <= max_iter) then 00244 00245 !*** 00246 !*** successfully found i,j - compute weights 00247 !*** 00248 00249 wgts(1) = (one-iguess)*(one-jguess) 00250 wgts(2) = iguess*(one-jguess) 00251 wgts(3) = iguess*jguess 00252 wgts(4) = (one-iguess)*jguess 00253 call store_link_bilin(dst_add, src_add, wgts, nmap) 00254 else 00255 write(nulou,*) 'Point coords: ',plat,plon 00256 write(nulou,*) 'Dest grid lats: ',src_lats 00257 write(nulou,*) 'Dest grid lons: ',src_lons 00258 write(nulou,*) 'Dest grid addresses: ',src_add 00259 write(nulou,*) 'Current i,j : ',iguess, jguess 00260 write(nulou,*) 00261 & 'Iteration for i,j exceed max iteration count' 00262 stop 00263 endif 00264 00265 else if (src_add(1) < 0) THEN 00266 00267 !*** 00268 !*** We are in the first or last bin or at least one of the 4 00269 !*** neighbours was masked. Do distance-weighted average using 00270 !*** the non-masked points among the 4 closest ones. 00271 00272 IF (nlogprt .eq. 2) then 00273 WRITE(nulou,*) ' ' 00274 WRITE(nulou,*) 00275 & 'WARNING: Cannot make bilinear interpolation for target point' 00276 & ,dst_add 00277 WRITE(nulou,*) 00278 & 'Using non-masked points among 4 nearest neighbors.' 00279 CALL FLUSH(nulou) 00280 ENDIF 00281 00282 !*** 00283 !*** Find the 4 closest points 00284 !*** 00285 00286 coslat_dst = cos(plat) 00287 sinlat_dst = sin(plat) 00288 coslon_dst = cos(plon) 00289 sinlon_dst = sin(plon) 00290 src_add = 0 00291 dist_min = bignum 00292 src_lats = bignum 00293 IF (min_add == 0) min_add = 1 00294 IF (max_add > grid1_size) max_add = grid1_size 00295 do srch_add = min_add,max_add 00296 arg = coslat_dst*cos(grid1_center_lat(srch_add))* 00297 & (coslon_dst*cos(grid1_center_lon(srch_add)) + 00298 & sinlon_dst*sin(grid1_center_lon(srch_add)))+ 00299 & sinlat_dst*sin(grid1_center_lat(srch_add)) 00300 IF (arg < -1.0d0) THEN 00301 arg = -1.0d0 00302 ELSE IF (arg > 1.0d0) THEN 00303 arg = 1.0d0 00304 END IF 00305 distance=acos(arg) 00306 00307 if (distance < dist_min) then 00308 sort_loop: do n=1,4 00309 if (distance < src_lats(n)) then 00310 do i=4,n+1,-1 00311 src_add (i) = src_add (i-1) 00312 src_lats(i) = src_lats(i-1) 00313 end do 00314 src_add (n) = srch_add 00315 src_lats(n) = distance 00316 dist_min = src_lats(4) 00317 exit sort_loop 00318 endif 00319 end do sort_loop 00320 endif 00321 end do 00322 00323 src_lons = one/(src_lats + tiny) 00324 distance = sum(src_lons) 00325 src_lats = src_lons/distance 00326 00327 !*** 00328 !*** Among 4 closest points, keep only the non-masked ones 00329 !*** 00330 00331 icount = 0 00332 do n=1,4 00333 if (grid1_mask(src_add(n)) .or. 00334 & (.not. grid1_mask(src_add(n)) .and. lextrapdone)) then 00335 icount = icount + 1 00336 else 00337 src_lats(n) = zero 00338 endif 00339 end do 00340 00341 if (icount > 0) then 00342 !*** renormalize weights 00343 sum_wgts = sum(src_lats) 00344 wgts(1) = src_lats(1)/sum_wgts 00345 wgts(2) = src_lats(2)/sum_wgts 00346 wgts(3) = src_lats(3)/sum_wgts 00347 wgts(4) = src_lats(4)/sum_wgts 00348 00349 grid2_frac(dst_add) = one 00350 call store_link_bilin(dst_add, src_add, wgts, nmap) 00351 ELSE 00352 IF (ll_nnei .eqv. .true. ) then 00353 IF (nlogprt .ge. 2) THEN 00354 WRITE(nulou,*) ' ' 00355 WRITE(nulou,*) 00356 & 'All 4 surrounding source grid points are masked' 00357 WRITE(nulou,*) 'for target point ',dst_add 00358 WRITE(nulou,*) 'with longitude and latitude', 00359 & plon, plat 00360 WRITE(nulou,*) 00361 & 'Using the nearest non-masked neighbour.' 00362 CALL FLUSH(nulou) 00363 ENDIF 00364 src_latsnn = bignum 00365 !cdir novector 00366 do srch_add = min_add,max_add 00367 if (grid1_mask(srch_add) .or. 00368 & (.not. grid1_mask(srch_add) .and. lextrapdone)) THEN 00369 arg = coslat_dst*cos(grid1_center_lat(srch_add))* 00370 & (coslon_dst*cos(grid1_center_lon(srch_add)) + 00371 & sinlon_dst*sin(grid1_center_lon(srch_add)))+ 00372 & sinlat_dst*sin(grid1_center_lat(srch_add)) 00373 IF (arg < -1.0d0) THEN 00374 arg = -1.0d0 00375 ELSE IF (arg > 1.0d0) THEN 00376 arg = 1.0d0 00377 END IF 00378 distance=acos(arg) 00379 00380 if (distance < src_latsnn) then 00381 src_addnn = srch_add 00382 src_latsnn = distance 00383 endif 00384 endif 00385 end DO 00386 IF (nlogprt .ge. 2) THEN 00387 WRITE(nulou,*) 00388 & 'Nearest non masked neighbour is source point ' 00389 & ,src_addnn 00390 WRITE(nulou,*) 'with longitude and latitude', 00391 & grid1_center_lon(src_addnn), 00392 & grid1_center_lat(src_addnn) 00393 WRITE(nulou,*) ' ' 00394 CALL FLUSH(nulou) 00395 ENDIF 00396 wgts(1) = 1. 00397 wgts(2) = 0. 00398 wgts(3) = 0. 00399 wgts(4) = 0. 00400 src_add(1) = src_addnn 00401 src_add(2) = 0 00402 src_add(3) = 0 00403 src_add(4) = 0 00404 00405 grid2_frac(dst_add) = one 00406 call store_link_bilin(dst_add, src_add, wgts, nmap) 00407 00408 ENDIF 00409 endif 00410 ENDIF 00411 end do grid_loop1 00412 ! 00413 !----------------------------------------------------------------------- 00414 00415 end subroutine remap_bilin_reduced 00416 00417 !*********************************************************************** 00418 00419 subroutine grid_search_bilin_rd(src_add, src_lats, src_lons, 00420 & min_add, max_add, 00421 & plat, plon, src_grid_dims, 00422 & src_center_lat, src_center_lon, 00423 & lextrapdone) 00424 00425 !----------------------------------------------------------------------- 00426 ! 00427 ! this routine finds the location of the search point plat, plon 00428 ! in the source grid and returns the corners needed for a bilinear 00429 ! interpolation. The target grid is a reduced grid. 00430 ! 00431 !----------------------------------------------------------------------- 00432 00433 !----------------------------------------------------------------------- 00434 ! 00435 ! output variables 00436 ! 00437 !----------------------------------------------------------------------- 00438 00439 integer (kind=int_kind), dimension(4), intent(out) :: 00440 & src_add ! address of each corner point enclosing P 00441 00442 real (kind=dbl_kind), dimension(4), intent(out) :: 00443 & src_lats, ! latitudes of the four corner points 00444 & src_lons ! longitudes of the four corner points 00445 00446 integer (kind=int_kind) :: min_add, max_add 00447 00448 !----------------------------------------------------------------------- 00449 ! 00450 ! input variables 00451 ! 00452 !----------------------------------------------------------------------- 00453 00454 real (kind=dbl_kind), intent(in) :: 00455 & plat, ! latitude of the search point 00456 & plon ! longitude of the search point 00457 00458 integer (kind=int_kind), dimension(2), intent(in) :: 00459 & src_grid_dims ! size of each src grid dimension 00460 00461 real (kind=dbl_kind), dimension(:), intent(in) :: 00462 & src_center_lat, ! latitude of each src grid center 00463 & src_center_lon ! longitude of each src grid center 00464 00465 LOGICAL :: 00466 & lextrapdone ! logical, true if EXTRAP done on field 00467 00468 !----------------------------------------------------------------------- 00469 ! 00470 ! local variables 00471 ! 00472 !----------------------------------------------------------------------- 00473 00474 integer (kind=int_kind) :: n, next_n, srch_add, ni, ! dummy indices 00475 & nx, ny, ntotmask, ! dimensions of src grid 00476 & inter_add ! add for restricting search 00477 ! 00478 integer (kind=int_kind), DIMENSION(4) :: src_bid 00479 00480 !----------------------------------------------------------------------- 00481 ! 00482 ! restrict search first using bins 00483 ! 00484 !----------------------------------------------------------------------- 00485 00486 nx = src_grid_dims(1) 00487 inter_add = 0 00488 00489 src_add = 0 00490 00491 min_add = size(src_center_lat) + 1 00492 max_add = 1 00493 if (plat >= bin_lats_r(1,1)) then 00494 min_add = 0 00495 max_add = bin_addr1_r(4,1) 00496 inter_add = bin_addr1_r(3,1) 00497 else if (plat <= bin_lats_r(1,num_srch_red)) then 00498 max_add = nx + 1 00499 min_add = bin_addr1_r(1,num_srch_red) 00500 inter_add = bin_addr1_r(3,num_srch_red) 00501 else 00502 srch_loopred: do n=1,num_srch_red 00503 if (plat <= bin_lats_r(1,n) 00504 & .and. plat >= bin_lats_r(2,n)) then 00505 min_add = bin_addr1_r(1,n) 00506 max_add = bin_addr1_r(4,n) 00507 inter_add = bin_addr1_r(3,n) 00508 exit srch_loopred 00509 endif 00510 end DO srch_loopred 00511 ENDIF 00512 00513 !----------------------------------------------------------------------- 00514 ! 00515 ! now perform a more detailed search 00516 ! 00517 !----------------------------------------------------------------------- 00518 if (min_add .ne. 0 .and. max_add .ne. nx+1) THEN 00519 !*** The concerned bins are not the top north or south ones. 00520 !*** We should be able to find the four corners 00521 !*** for the bilinear interpolation. 00522 00523 IF ( plon <= src_center_lon(min_add) ) THEN 00524 src_add(1) = inter_add-1 00525 src_add(2) = min_add 00526 ELSE IF ( plon > src_center_lon(inter_add-1) ) then 00527 src_add(1) = inter_add-1 00528 src_add(2) = min_add 00529 else 00530 srch_loop2: do srch_add = min_add, inter_add-2 00531 if ( (plon > src_center_lon(srch_add)) .and. 00532 & (plon <= src_center_lon(srch_add+1)) )then 00533 src_add(1) = srch_add 00534 src_add(2) = srch_add+1 00535 exit srch_loop2 00536 endif 00537 end do srch_loop2 00538 ENDIF 00539 IF ( plon <= src_center_lon(inter_add) ) THEN 00540 src_add(3) = max_add 00541 src_add(4) = inter_add 00542 ELSE IF ( plon >= src_center_lon(max_add) ) then 00543 src_add(3) = max_add 00544 src_add(4) = inter_add 00545 else 00546 srch_loop3: do srch_add = inter_add, max_add 00547 if ( (plon >= src_center_lon(srch_add)) .and. 00548 & (plon <= src_center_lon(srch_add+1)) )then 00549 src_add(3) = srch_add 00550 src_add(4) = srch_add+1 00551 exit srch_loop3 00552 endif 00553 enddo srch_loop3 00554 ENDIF 00555 src_lats(1) = src_center_lat(src_add(3)) 00556 src_lats(2) = src_center_lat(src_add(4)) 00557 src_lats(3) = src_center_lat(src_add(2)) 00558 src_lats(4) = src_center_lat(src_add(1)) 00559 00560 src_lons(1) = src_center_lon(src_add(3)) 00561 src_lons(2) = src_center_lon(src_add(4)) 00562 src_lons(3) = src_center_lon(src_add(2)) 00563 src_lons(4) = src_center_lon(src_add(1)) 00564 00565 src_bid=src_add 00566 00567 src_add(1) = src_bid(3) 00568 src_add(2) = src_bid(4) 00569 src_add(3) = src_bid(2) 00570 src_add(4) = src_bid(1) 00571 00572 ! Check if one point is masked; IF so, nearest-neighbour 00573 ! interpolation will be used 00574 00575 ntotmask = 0 00576 do ni=1,4 00577 if (.not. grid1_mask(src_add(ni)).and. 00578 & .not. lextrapdone) 00579 & ntotmask = ntotmask + 1 00580 end DO 00581 IF (ntotmask .gt. 0) src_add(1) = -src_add(1) 00582 00583 ELSE 00584 00585 !*** We are in the first or last bin. Put src_add = -1 so that 00586 !*** distance-weighted average of the 4 non-masked closest points 00587 !*** is done in calling routine. 00588 00589 src_add = -1 00590 00591 ENDIF 00592 00593 !----------------------------------------------------------------------- 00594 00595 end subroutine grid_search_bilin_rd 00596 00597 !*********************************************************************** 00598 00599 subroutine store_link_bilin(dst_add, src_add, weights, nmap) 00600 00601 !----------------------------------------------------------------------- 00602 ! 00603 ! this routine stores the address and weight for four links 00604 ! associated with one destination point in the appropriate address 00605 ! and weight arrays and resizes those arrays if necessary. 00606 ! 00607 !----------------------------------------------------------------------- 00608 00609 !----------------------------------------------------------------------- 00610 ! 00611 ! input variables 00612 ! 00613 !----------------------------------------------------------------------- 00614 00615 integer (kind=int_kind), intent(in) :: 00616 & dst_add, ! address on destination grid 00617 & nmap ! identifies which direction for mapping 00618 00619 integer (kind=int_kind), dimension(4), intent(in) :: 00620 & src_add ! addresses on source grid 00621 00622 real (kind=dbl_kind), dimension(4), intent(in) :: 00623 & weights ! array of remapping weights for these links 00624 00625 !----------------------------------------------------------------------- 00626 ! 00627 ! local variables 00628 ! 00629 !----------------------------------------------------------------------- 00630 00631 integer (kind=int_kind) :: n, ! dummy index 00632 & num_links_old ! placeholder for old link number 00633 00634 !----------------------------------------------------------------------- 00635 ! 00636 ! increment number of links and check to see if remap arrays need 00637 ! to be increased to accomodate the new link. then store the 00638 ! link. 00639 ! 00640 !----------------------------------------------------------------------- 00641 00642 select case (nmap) 00643 case(1) 00644 00645 num_links_old = num_links_map1 00646 num_links_map1 = num_links_old + 4 00647 00648 if (num_links_map1 > max_links_map1) 00649 & call resize_remap_vars(1,resize_increment) 00650 00651 do n=1,4 00652 grid1_add_map1(num_links_old+n) = src_add(n) 00653 grid2_add_map1(num_links_old+n) = dst_add 00654 wts_map1 (1,num_links_old+n) = weights(n) 00655 end do 00656 00657 case(2) 00658 00659 num_links_old = num_links_map2 00660 num_links_map2 = num_links_old + 4 00661 00662 if (num_links_map2 > max_links_map2) 00663 & call resize_remap_vars(2,resize_increment) 00664 00665 do n=1,4 00666 grid1_add_map2(num_links_old+n) = dst_add 00667 grid2_add_map2(num_links_old+n) = src_add(n) 00668 wts_map2 (1,num_links_old+n) = weights(n) 00669 end do 00670 00671 end select 00672 00673 !----------------------------------------------------------------------- 00674 00675 end subroutine store_link_bilin 00676 00677 !*********************************************************************** 00678 00679 end module remap_bilinear_reduced 00680 00681 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!