/[MITgcm]/MITgcm/pkg/exf/exf_interp.F
ViewVC logotype

Diff of /MITgcm/pkg/exf/exf_interp.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.4 by dimitri, Sun Dec 7 07:41:45 2003 UTC revision 1.5 by dimitri, Wed Dec 10 19:37:25 2003 UTC
# Line 76  C subroutine variables Line 76  C subroutine variables
76  C local variables  C local variables
77        integer  ierr        integer  ierr
78        real*8   ne_fac,nw_fac,se_fac,sw_fac        real*8   ne_fac,nw_fac,se_fac,sw_fac
79        integer  e_ind(snx),w_ind(snx)        integer  e_ind(snx,sny),w_ind(snx,sny)
80        integer  n_ind(sny),s_ind(sny)        integer  n_ind(snx,sny),s_ind(snx,sny)
81        real*8   px_ind(4), py_ind(4), ew_val(4)        real*8   px_ind(4), py_ind(4), ew_val(4)
82        external lagran        external lagran
83        real*8   lagran        real*8   lagran
84        real*4   arrayin(-1:nx_in+2 ,      -1:ny_in+2)        real*4   arrayin(-1:nx_in+2 ,      -1:ny_in+2)
85        real*8   x_in   (-1:nx_in+2), y_in(-1:ny_in+2)        real*8   x_in   (-1:nx_in+2), y_in(-1:ny_in+2)
86        integer  i, j, k, l, js, bi, bj, interp_unit, sp        integer  i, j, k, l, js, bi, bj, sp, interp_unit
87        real*4   global(nx_in,ny_in)        real*4   global(nx_in,ny_in)
88    
89        _BEGIN_MASTER( myThid )        _BEGIN_MASTER( myThid )
# Line 177  C check validity of input/output coordin Line 177  C check validity of input/output coordin
177    
178  C compute interpolation indices  C compute interpolation indices
179          do i=1,snx          do i=1,snx
180           if (xG(i,1,bi,bj)-x_in(1) .ge. 0.) then           do j=1,sny
181              w_ind(i) = int((xG(i,1,bi,bj)-x_in(1))/lon_inc) + 1            if (xG(i,j,bi,bj)-x_in(1) .ge. 0.) then
182           else             w_ind(i,j) = int((xG(i,j,bi,bj)-x_in(1))/lon_inc) + 1
183              w_ind(i) = int((xG(i,1,bi,bj)-x_in(1))/lon_inc)            else
184           endif             w_ind(i,j) = int((xG(i,j,bi,bj)-x_in(1))/lon_inc)
185           e_ind(i) = w_ind(i) + 1            endif
186          enddo            e_ind(i,j) = w_ind(i,j) + 1
187          js = ny_in/2            js = ny_in/2
188          do j=1,sny            do while (yG(i,j,bi,bj) .lt. y_in(js))
189           do while (yG(1,j,bi,bj) .lt. y_in(js))             js = (js + 1)/2
190              js = (js + 1)/2            enddo
191           enddo            do while (yG(i,j,bi,bj) .ge. y_in(js+1))
192           do while (yG(1,j,bi,bj) .ge. y_in(js+1))             js = js + 1
193              js = js + 1            enddo
194              s_ind(i,j) = js
195              n_ind(i,j) = js + 1
196           enddo           enddo
          s_ind(j) = js  
          n_ind(j) = js + 1  
197          enddo          enddo
198    
199          if (method .eq. 1) then          if (method .eq. 1) then
# Line 204  C bilinear interpolation Line 204  C bilinear interpolation
204            do i=1,snx            do i=1,snx
205             arrayout(i,j,bi,bj) = 0.             arrayout(i,j,bi,bj) = 0.
206             do l=0,1             do l=0,1
207              px_ind(l+1) = x_in(w_ind(i)+l)              px_ind(l+1) = x_in(w_ind(i,j)+l)
208              py_ind(l+1) = y_in(s_ind(j)+l)              py_ind(l+1) = y_in(s_ind(i,j)+l)
209             enddo             enddo
210             do k=1,2             do k=1,2
211              ew_val(k) = arrayin(w_ind(i),s_ind(j)+k-1)              ew_val(k) = arrayin(w_ind(i,j),s_ind(i,j)+k-1)
212       &             *lagran(1,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(1,xG(i,j,bi,bj),px_ind,sp)
213       &             +arrayin(e_ind(i),s_ind(j)+k-1)       &             +arrayin(e_ind(i,j),s_ind(i,j)+k-1)
214       &             *lagran(2,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(2,xG(i,j,bi,bj),px_ind,sp)
215              arrayout(i,j,bi,bj)=arrayout(i,j,bi,bj)              arrayout(i,j,bi,bj)=arrayout(i,j,bi,bj)
216       &             +ew_val(k)*lagran(k,yG(1,j,bi,bj),py_ind,sp)       &             +ew_val(k)*lagran(k,yG(i,j,bi,bj),py_ind,sp)
217             enddo             enddo
218            enddo            enddo
219           enddo           enddo
# Line 225  C bicubic interpolation Line 225  C bicubic interpolation
225            do i=1,snx            do i=1,snx
226             arrayout(i,j,bi,bj) = 0.             arrayout(i,j,bi,bj) = 0.
227             do l=-1,2             do l=-1,2
228              px_ind(l+2) = x_in(w_ind(i)+l)              px_ind(l+2) = x_in(w_ind(i,j)+l)
229              py_ind(l+2) = y_in(s_ind(j)+l)              py_ind(l+2) = y_in(s_ind(i,j)+l)
230             enddo             enddo
231             do k=1,4             do k=1,4
232              ew_val(k) =              ew_val(k) =
233       &             arrayin(w_ind(i)-1,s_ind(j)+k-2)       &             arrayin(w_ind(i,j)-1,s_ind(i,j)+k-2)
234       &             *lagran(1,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(1,xG(i,j,bi,bj),px_ind,sp)
235       &             +arrayin(w_ind(i)  ,s_ind(j)+k-2)       &             +arrayin(w_ind(i,j)  ,s_ind(i,j)+k-2)
236       &             *lagran(2,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(2,xG(i,j,bi,bj),px_ind,sp)
237       &             +arrayin(e_ind(i)  ,s_ind(j)+k-2)       &             +arrayin(e_ind(i,j)  ,s_ind(i,j)+k-2)
238       &             *lagran(3,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(3,xG(i,j,bi,bj),px_ind,sp)
239       &             +arrayin(e_ind(i)+1,s_ind(j)+k-2)       &             +arrayin(e_ind(i,j)+1,s_ind(i,j)+k-2)
240       &             *lagran(4,xG(i,1,bi,bj),px_ind,sp)       &             *lagran(4,xG(i,j,bi,bj),px_ind,sp)
241              arrayout(i,j,bi,bj)=arrayout(i,j,bi,bj)              arrayout(i,j,bi,bj)=arrayout(i,j,bi,bj)
242       &             +ew_val(k)*lagran(k,yG(1,j,bi,bj),py_ind,sp)       &             +ew_val(k)*lagran(k,yG(i,j,bi,bj),py_ind,sp)
243             enddo             enddo
244            enddo            enddo
245           enddo           enddo
Legend:
Removed from v.1.4  
changed lines
  Added in v.1.5
  ViewVC Help
Powered by ViewVC 1.1.22