/[MITgcm]/MITgcm/pkg/flt/flt_mapping.F
ViewVC logotype

Contents of /MITgcm/pkg/flt/flt_mapping.F

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

Revision 1.4 - (show annotations) (download)
Wed Dec 22 21:21:24 2010 UTC (13 years, 4 months ago) by jahn
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64, checkpoint65, checkpoint63, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x
Changes since 1.3: +18 -1 lines 
add mapping for curvilinear grid
1 C $Header: /u/gcmpack/MITgcm/pkg/flt/flt_mapping.F,v 1.3 2009/02/13 04:22:22 jmc Exp $
2 C $Name: $
3
4 #include "FLT_OPTIONS.h"
5
6 C-- Contents
7 C-- o FLT_MAP_XY2IJLOCAL
8 C-- o FLT_MAP_IJLOCAL2XY
9 C-- o FLT_MAP_R2K (Function)
10 C-- o FLT_MAP_K2R (Function)
11
12 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
13
14 SUBROUTINE FLT_MAP_XY2IJLOCAL(
15 O ix, jy,
16 I xx, yy, bi, bj, myThid )
17
18 C ==================================================================
19 C SUBROUTINE FLT_MAP_XY2IJLOCAL
20 C ==================================================================
21 C o Converts global x,y-coordinates (grid) to corresponding
22 C local fractional horizontal indices for specific tile
23 C Range: [1/2 , sNx+1/2] , [1/2 , sNy+1/2]
24 C Center (Tracer Pt) <-> integer , integer
25 C U-velocity Pt <-> half integer , integer
26 C V-velocity Pt <-> integer , half integer
27 C Vorticity Pt <-> half integer , half integer
28 C ==================================================================
29
30 C !USES:
31 IMPLICIT NONE
32
33 C == global variables ==
34 #include "SIZE.h"
35 #include "EEPARAMS.h"
36 #include "GRID.h"
37 #include "PARAMS.h"
38
39 C == routine arguments ==
40 _RL ix, jy
41 _RL xx, yy
42 INTEGER bi, bj, myThid
43
44 C == local variables ==
45 _RL fm, dist
46 INTEGER i, j
47
48 C == end of interface ==
49
50 IF ( usingCartesianGrid .OR.
51 & usingSphericalPolarGrid .AND. .NOT.rotateGrid
52 & ) THEN
53
54 ix = -1. _d 0
55 jy = -1. _d 0
56
57 j = 1
58 DO i=0,sNx+1
59 IF ( ix.EQ.-1. _d 0 ) THEN
60 IF ( xG(i,j,bi,bj).LE.xx .AND. xx.LT.xG(i+1,j,bi,bj) ) THEN
61 dist = xG(i+1,j,bi,bj) - xG(i,j,bi,bj)
62 fm = ( xx - xG(i,j,bi,bj) ) / dist
63 ix = DFLOAT(i)+fm-0.5 _d 0
64 ENDIF
65 ENDIF
66 ENDDO
67
68 i = 1
69 DO j=0,sNy+1
70 IF ( jy.EQ.-1. _d 0 ) THEN
71 IF ( yG(i,j,bi,bj).LE.yy .AND. yy.LT.yG(i,j+1,bi,bj) ) THEN
72 dist = yG(i,j+1,bi,bj) - yG(i,j,bi,bj)
73 fm = ( yy - yG(i,j,bi,bj) ) / dist
74 jy = DFLOAT(j)+fm-0.5 _d 0
75 ENDIF
76 ENDIF
77 ENDDO
78
79 ELSE
80 STOP 'FLT_MAP_XY2IJLOCAL: not yet coded for this grid'
81 ENDIF
82
83 RETURN
84 END
85
86 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
87
88 SUBROUTINE FLT_MAP_IJLOCAL2XY(
89 O xx, yy,
90 I ix, jy, bi, bj, myThid )
91
92 C ==================================================================
93 C SUBROUTINE FLT_MAP_IJLOCAL2XY
94 C ==================================================================
95 C o Converts local fractional horizontal indices for specific tile
96 C to corresponding global x,y-coordinates (grid)
97 C Range: [1/2 , sNx+1/2] , [1/2 , sNy+1/2]
98 C Center (Tracer Pt) <-> integer , integer
99 C U-velocity Pt <-> half integer , integer
100 C V-velocity Pt <-> integer , half integer
101 C Vorticity Pt <-> half integer , half integer
102 C ==================================================================
103
104 C !USES:
105 IMPLICIT NONE
106
107 C == global variables ==
108 #include "SIZE.h"
109 #include "EEPARAMS.h"
110 #include "GRID.h"
111 #include "PARAMS.h"
112
113 C == routine arguments ==
114 _RL xx, yy
115 _RL ix, jy
116 INTEGER bi, bj, myThid
117
118 C == local variables ==
119 _RL ddx, ddy
120 INTEGER i, j
121
122 C == end of interface ==
123
124 IF ( usingCartesianGrid .OR.
125 & usingSphericalPolarGrid .AND. .NOT.rotateGrid
126 & ) THEN
127
128 i = NINT(ix)
129 j = NINT(jy)
130 ddx = 0.5 _d 0 + ix - DFLOAT(i)
131 ddy = 0.5 _d 0 + jy - DFLOAT(j)
132
133 xx = xG(i,j,bi,bj) + ddx*( xG(i+1,j,bi,bj) - xG(i,j,bi,bj) )
134 yy = yG(i,j,bi,bj) + ddy*( yG(i,j+1,bi,bj) - yG(i,j,bi,bj) )
135
136 ELSEIF ( usingCurvilinearGrid ) THEN
137
138 i = NINT(ix)
139 j = NINT(jy)
140 ddx = 0.5 _d 0 + ix - DFLOAT(i)
141 ddy = 0.5 _d 0 + jy - DFLOAT(j)
142
143 C bilinear interpolation within grid cell (should use arcs instead?)
144 xx = xG(i,j,bi,bj) + ddx*( xG(i+1,j,bi,bj) - xG(i,j,bi,bj) )
145 & + ddy*( xG(i,j+1,bi,bj) - xG(i,j,bi,bj) )
146 & + ddx*ddy*( xG(i+1,j+1,bi,bj) - xG(i+1,j,bi,bj)
147 & - xG(i,j+1,bi,bj) + xG(i,j,bi,bj) )
148 yy = yG(i,j,bi,bj) + ddx*( yG(i+1,j,bi,bj) - yG(i,j,bi,bj) )
149 & + ddy*( yG(i,j+1,bi,bj) - yG(i,j,bi,bj) )
150 & + ddx*ddy*( yG(i+1,j+1,bi,bj) - yG(i+1,j,bi,bj)
151 & - yG(i,j+1,bi,bj) + yG(i,j,bi,bj) )
152
153 ELSE
154 STOP 'FLT_MAP_IJLOCAL2XY: not yet coded for this grid'
155 ENDIF
156
157 RETURN
158 END
159
160 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
161
162 _RL FUNCTION FLT_MAP_R2K(
163 I rr, bi, bj, myThid )
164
165 C ==================================================================
166 C FUNCTION FLT_MAP_R2K
167 C ==================================================================
168 C o Converts r-coordinate (grid) to corresponding
169 C fractional vertical index
170 C Range: [1/2 , Nr+1/2],
171 C Center (Tracer Pt) <-> integer
172 C Interface (wVel Pt) <-> half integer
173 C ==================================================================
174
175 C !USES:
176 IMPLICIT NONE
177
178 C == global variables ==
179 #include "SIZE.h"
180 #include "EEPARAMS.h"
181 #include "GRID.h"
182
183 C == routine arguments ==
184 _RL rr
185 INTEGER bi, bj, myThid
186
187 C == local variables ==
188 _RL fm
189 INTEGER k
190
191 C == end of interface ==
192
193 FLT_MAP_R2K = 0. _d 0
194 DO k=1,Nr
195 IF ( FLT_MAP_R2K .EQ. 0. _d 0 ) THEN
196 C- r decreases when k increases (rkSign < 0):
197 IF ( rF(k) .GE. rr .AND. rr.GT.rF(k+1) ) THEN
198 fm = ( rr - rF(k) ) * recip_drF(k)*rkSign
199 FLT_MAP_R2K = DFLOAT(k)+fm-0.5 _d 0
200 ENDIF
201 ENDIF
202 ENDDO
203
204 RETURN
205 END
206
207 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
208
209 _RL FUNCTION FLT_MAP_K2R(
210 I kr, bi, bj, myThid )
211
212 C ==================================================================
213 C FUNCTION FLT_MAP_K2R
214 C ==================================================================
215 C o Converts fractional vertical index to corresponding
216 C r-coordinate (grid)
217 C Range: [1/2 , Nr+1/2],
218 C Center (Tracer Pt) <-> integer
219 C Interface (wVel Pt) <-> half integer
220 C ==================================================================
221
222 C !USES:
223 IMPLICIT NONE
224
225 C == global variables ==
226 #include "SIZE.h"
227 #include "EEPARAMS.h"
228 #include "GRID.h"
229
230 C == routine arguments ==
231 _RL kr
232 INTEGER bi, bj, myThid
233
234 C == local variables ==
235 _RL ddz
236 INTEGER k
237
238 C == end of interface ==
239
240 k = NINT(kr)
241 IF ( k.LT.1 ) THEN
242 FLT_MAP_K2R = rF(1)
243 ELSEIF ( k.GT.Nr ) THEN
244 FLT_MAP_K2R = rF(Nr+1)
245 ELSE
246 ddz = 0.5 _d 0 + kr - DFLOAT(k)
247 FLT_MAP_K2R = rF(k) + ddz*drF(k)*rkSign
248 ENDIF
249
250 RETURN
251 END
  ViewVC Help
Powered by ViewVC 1.1.22