/[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.14 by jmc, Sat Jan 14 20:47:43 2006 UTC revision 1.15 by dimitri, Fri Jun 30 12:05:42 2006 UTC
# Line 52  C     lon_0, lat_0 = lon and lat of sw c Line 52  C     lon_0, lat_0 = lon and lat of sw c
52  C     lon_inc      = scalar x-grid increment  C     lon_inc      = scalar x-grid increment
53  C     lat_inc      = vector y-grid increments  C     lat_inc      = vector y-grid increments
54  C     nx_in, ny_in = input x-grid and y-grid size  C     nx_in, ny_in = input x-grid and y-grid size
55  C     method       = 1 for bilinear 2 for bicubic  C     method       = 1,11,21 for bilinear; 2,12,22 for bicubic
56    C                    1,2 for tracer; 11,12 for U; 21,22 for V
57  C     mythid       = thread id  C     mythid       = thread id
58  C  C
59    
# Line 71  C subroutine variables Line 72  C subroutine variables
72        integer       method, mythid        integer       method, mythid
73    
74  C local variables  C local variables
 c     real*8   ne_fac,nw_fac,se_fac,sw_fac  
75        integer  e_ind(snx,sny),w_ind(snx,sny)        integer  e_ind(snx,sny),w_ind(snx,sny)
76        integer  n_ind(snx,sny),s_ind(snx,sny)        integer  n_ind(snx,sny),s_ind(snx,sny)
77        real*8   px_ind(4), py_ind(4), ew_val(4)        real*8   px_ind(4), py_ind(4), ew_val(4)
# Line 81  c     real*8   ne_fac,nw_fac,se_fac,sw_f Line 81  c     real*8   ne_fac,nw_fac,se_fac,sw_f
81        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)
82        integer  i, j, k, l, js, bi, bj, sp, interp_unit        integer  i, j, k, l, js, bi, bj, sp, interp_unit
83        _RS      xG(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)        _RS      xG(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
84        _RS      threeSixtyRS        _RS      threeSixtyRS, NorthValue
85        PARAMETER ( threeSixtyRS = 360. )        PARAMETER ( threeSixtyRS = 360. )
86    
87  C     put xG in interval [ lon_0 , lon_0+360 [  C     put xG in interval [ lon_0 , lon_0+360 [
# Line 98  C     put xG in interval [ lon_0 , lon_0 Line 98  C     put xG in interval [ lon_0 , lon_0
98        enddo        enddo
99    
100         call exf_interp_read(         call exf_interp_read(
101       I   infile,       I   infile, filePrec,
      I   filePrec,  
102       O   arrayin,       O   arrayin,
103       I   irecord, xG, yG,       I   irecord, nx_in, ny_in, mythid)
      I   lon_0, lon_inc,  
      I   lat_0, lat_inc,  
      I   nx_in, ny_in, method, mythid)  
104        _BARRIER        _BARRIER
105    
106  C     _BEGIN_MASTER( myThid )  C     _BEGIN_MASTER( myThid )
# Line 120  C setup input grid Line 116  C setup input grid
116         do j=2,ny_in         do j=2,ny_in
117          y_in(j) = y_in(j-1) + lat_inc(j-1)          y_in(j) = y_in(j-1) + lat_inc(j-1)
118         enddo         enddo
119         y_in(ny_in+1) = y_in(ny_in) + lat_inc(ny_in-1)  c       y_in(ny_in+1) = y_in(ny_in) + lat_inc(ny_in-1)
120         y_in(ny_in+2) = y_in(ny_in) + 2.*lat_inc(ny_in-1)  c       y_in(ny_in+2) = y_in(ny_in) + 2.*lat_inc(ny_in-1)
121           y_in(ny_in+1) = min( y_in(ny_in) + lat_inc(ny_in-1), 90 )
122           y_in(ny_in+2) = min( y_in(ny_in) + 2.*lat_inc(ny_in-1), 90)
123    
124  C enlarge boundary  C enlarge boundary
125         do j=1,ny_in         do j=1,ny_in
# Line 137  C enlarge boundary Line 135  C enlarge boundary
135          arrayin(i,ny_in+2) = arrayin(i,ny_in)          arrayin(i,ny_in+2) = arrayin(i,ny_in)
136         enddo         enddo
137    
138  C     _END_MASTER( myThid )  C     For tracer (method=1,2) set to northernmost zonal-mean value
139    C     at 90N to avoid sharp zonal gradients near the Pole.
140    C     For U (method=11,12) set to zero at 90N to minimize velocity
141    C     gradient at North Pole
142    C     For V (method=11,12) set to northernmost zonal value at 90N,
143    C     as is already done above in order to allow cross-PoleArctic flow
144           if (y_in(ny_in+1).eq.90) then
145            if (method.eq.1 .or. method.eq.2) then
146             NorthValue = 0
147             do i=1,nx_in
148              NorthValue = NorthValue + arrayin(i,ny_in)
149             enddo
150             NorthValue = NorthValue / nx_in
151             do i=-1,nx_in+2
152              arrayin(i,ny_in+1) = NorthValue
153             enddo
154            elseif (method.eq.11 .or. method.eq.12) then
155             do i=-1,nx_in+2
156              arrayin(i,ny_in+1) = 0
157             enddo
158            endif
159           endif
160           if (y_in(ny_in+2).eq.90) then
161            if (method.eq.1 .or. method.eq.2) then
162             NorthValue = 0
163             do i=1,nx_in
164              NorthValue = NorthValue + arrayin(i,ny_in)
165             enddo
166             NorthValue = NorthValue / nx_in
167             do i=-1,nx_in+2
168              arrayin(i,ny_in+2) = NorthValue
169             enddo
170            elseif (method.eq.11 .or. method.eq.12) then
171             do i=-1,nx_in+2
172              arrayin(i,ny_in+2) = 0
173             enddo
174            endif
175           endif
176    
177    C     _END_MASTER( myThid )
178          
179        do bj = mybylo(mythid), mybyhi(mythid)        do bj = mybylo(mythid), mybyhi(mythid)
180         do bi = mybxlo(mythid), mybxhi(mythid)         do bi = mybxlo(mythid), mybxhi(mythid)
181    
# Line 186  C compute interpolation indices Line 223  C compute interpolation indices
223           enddo           enddo
224          enddo          enddo
225    
226          if (method .eq. 1) then          if (method.eq.1 .or. method.eq.11 .or. method.eq.21) then
227    
228  C bilinear interpolation  C bilinear interpolation
229           sp = 2           sp = 2
# Line 207  C bilinear interpolation Line 244  C bilinear interpolation
244             enddo             enddo
245            enddo            enddo
246           enddo           enddo
247          elseif (method .eq. 2) then          elseif (method .eq. 2 .or. method.eq.12 .or. method.eq.22) then
248    
249  C bicubic interpolation  C bicubic interpolation
250           sp = 4           sp = 4
Legend:
Removed from v.1.14  
changed lines
  Added in v.1.15
  ViewVC Help
Powered by ViewVC 1.1.22