/[MITgcm]/MITgcm/model/src/ini_curvilinear_grid.F
ViewVC logotype

Contents of /MITgcm/model/src/ini_curvilinear_grid.F

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


Revision 1.35 - (show annotations) (download)
Mon Aug 20 22:08:15 2007 UTC (16 years, 10 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint59q, checkpoint59p, checkpoint59g, checkpoint59f, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59h, checkpoint59k, checkpoint59j
Changes since 1.34: +7 -6 lines
call the new UV_BGRID exchange for (dxV,dyU), for now, only with exch2 and
hiden from TAF.

1 C $Header: /u/gcmpack/MITgcm/model/src/ini_curvilinear_grid.F,v 1.34 2007/01/09 23:27:43 dimitri Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 CBOP
8 C !ROUTINE: INI_CURVILINEAR_GRID
9 C !INTERFACE:
10 SUBROUTINE INI_CURVILINEAR_GRID( myThid )
11 C !DESCRIPTION: \bv
12 C *==========================================================*
13 C | SUBROUTINE INI_CURVILINEAR_GRID
14 C | o Initialise curvilinear coordinate system
15 C *==========================================================*
16 C | Curvilinear grid settings are read from a file rather
17 C | than coded in-line as for cartesian and spherical polar.
18 C | This is more general but you have to create the grid
19 C | yourself.
20 C *==========================================================*
21 C \ev
22
23 C !USES:
24 IMPLICIT NONE
25 C === Global variables ===
26 #include "SIZE.h"
27 #include "EEPARAMS.h"
28 #include "PARAMS.h"
29 #include "GRID.h"
30 #ifdef ALLOW_EXCH2
31 #include "W2_EXCH2_TOPOLOGY.h"
32 #include "W2_EXCH2_PARAMS.h"
33 #endif
34 #ifdef ALLOW_MNC
35 #include "MNC_PARAMS.h"
36 #endif
37
38 C !INPUT/OUTPUT PARAMETERS:
39 C == Routine arguments ==
40 C myThid - Number of this instance of INI_CURVILINEAR_GRID
41 INTEGER myThid
42
43 C !LOCAL VARIABLES:
44 C == Local variables ==
45 INTEGER bi,bj
46 INTEGER I,J
47 CHARACTER*(MAX_LEN_MBUF) msgBuf
48 LOGICAL anglesAreSet
49 #ifdef ALLOW_MNC
50 CHARACTER*(80) mncFn
51 #endif
52 #ifndef OLD_GRID_IO
53 # ifdef ALLOW_EXCH2
54 _RL buf(sNx*nSx*nPx+1)
55 INTEGER myTile
56 # else
57 _RL buf(sNx+1,sNy+1)
58 # endif
59 INTEGER iG, iL, iLen
60 CHARACTER*(MAX_LEN_FNAM) fName
61 CHARACTER*(MAX_LEN_MBUF) tmpBuf
62 INTEGER ILNBLNK
63 EXTERNAL ILNBLNK
64 #endif
65 CEOP
66
67 C-- Set everything to zero everywhere
68 DO bj = myByLo(myThid), myByHi(myThid)
69 DO bi = myBxLo(myThid), myBxHi(myThid)
70
71 DO J=1-Oly,sNy+Oly
72 DO I=1-Olx,sNx+Olx
73 xC(i,j,bi,bj)=0.
74 yC(i,j,bi,bj)=0.
75 xG(i,j,bi,bj)=0.
76 yG(i,j,bi,bj)=0.
77 dxC(i,j,bi,bj)=0.
78 dyC(i,j,bi,bj)=0.
79 dxG(i,j,bi,bj)=0.
80 dyG(i,j,bi,bj)=0.
81 dxF(i,j,bi,bj)=0.
82 dyF(i,j,bi,bj)=0.
83 dxV(i,j,bi,bj)=0.
84 dyU(i,j,bi,bj)=0.
85 rA(i,j,bi,bj)=0.
86 rAz(i,j,bi,bj)=0.
87 rAw(i,j,bi,bj)=0.
88 rAs(i,j,bi,bj)=0.
89 tanPhiAtU(i,j,bi,bj)=0.
90 tanPhiAtV(i,j,bi,bj)=0.
91 angleCosC(i,j,bi,bj)=1.
92 angleSinC(i,j,bi,bj)=0.
93 cosFacU(J,bi,bj)=1.
94 cosFacV(J,bi,bj)=1.
95 sqCosFacU(J,bi,bj)=1.
96 sqCosFacV(J,bi,bj)=1.
97 ENDDO
98 ENDDO
99
100 ENDDO
101 ENDDO
102
103 C-- Everyone must wait for the initialisation to be done
104 _BARRIER
105
106 C Here we make no assumptions about grid symmetry and simply
107 C read the raw grid data from files
108
109 #ifdef OLD_GRID_IO
110
111 C- Cell centered quantities
112 CALL MDSREADFIELD('LONC.bin',readBinaryPrec,'RS',1,xC, 1,myThid)
113 CALL MDSREADFIELD('LATC.bin',readBinaryPrec,'RS',1,yC, 1,myThid)
114 _EXCH_XY_R4(xC,myThid)
115 _EXCH_XY_R4(yC,myThid)
116
117 CALL MDSREADFIELD('DXF.bin',readBinaryPrec,'RS',1,dxF, 1,myThid)
118 CALL MDSREADFIELD('DYF.bin',readBinaryPrec,'RS',1,dyF, 1,myThid)
119 CALL EXCH_UV_AGRID_3D_RS( dxF, dyF, .FALSE., 1, myThid )
120
121 CALL MDSREADFIELD('RA.bin',readBinaryPrec,'RS',1,rA, 1,myThid)
122 _EXCH_XY_R4(rA,myThid )
123
124 C- Corner quantities
125 CALL MDSREADFIELD('LONG.bin',readBinaryPrec,'RS',1,xG, 1,myThid)
126 CALL MDSREADFIELD('LATG.bin',readBinaryPrec,'RS',1,yG, 1,myThid)
127 IF (useCubedSphereExchange) THEN
128 cs- this block needed by cubed sphere until we write more useful I/O routines
129 bi=3
130 bj=1
131 yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
132 bj=bj+2
133 yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
134 bj=bj+2
135 yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
136 bi=6
137 bj=2
138 yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
139 bj=bj+2
140 yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
141 bj=bj+2
142 yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
143 cs- end block
144 ENDIF
145 CALL EXCH_Z_3D_RS( xG, 1, myThid )
146 CALL EXCH_Z_3D_RS( yG, 1, myThid )
147
148 CALL MDSREADFIELD('DXV.bin',readBinaryPrec,'RS',1,dxV, 1,myThid)
149 CALL MDSREADFIELD('DYU.bin',readBinaryPrec,'RS',1,dyU, 1,myThid)
150 c CALL EXCH_UV_BGRID_3D_RS( dxV, dyU, .FALSE., 1, myThid)
151 cs- this block needed by cubed sphere until we write more useful I/O routines
152 IF (.NOT.useCubedSphereExchange) THEN
153 CALL EXCH_Z_3D_RS( dxV, 1, myThid )
154 CALL EXCH_Z_3D_RS( dyU, 1, myThid )
155 ELSE
156 DO bj = myByLo(myThid), myByHi(myThid)
157 DO bi = myBxLo(myThid), myBxHi(myThid)
158 cs! fix overlaps:
159 DO j=1,sNy
160 DO i=1,Olx
161 dxV(1-i,j,bi,bj)=dxV(1+i,j,bi,bj)
162 dxV(sNx+i,j,bi,bj)=dxV(i,j,bi,bj)
163 dyU(1-i,j,bi,bj)=dyU(1+i,j,bi,bj)
164 dyU(sNx+i,j,bi,bj)=dyU(i,j,bi,bj)
165 ENDDO
166 ENDDO
167 DO j=1,Oly
168 DO i=1-Olx,sNx+Olx
169 dxV(i,1-j,bi,bj)=dxV(i,1+j,bi,bj)
170 dxV(i,sNy+j,bi,bj)=dxV(i,j,bi,bj)
171 dyU(i,1-j,bi,bj)=dyU(i,1+j,bi,bj)
172 dyU(i,sNy+j,bi,bj)=dyU(i,j,bi,bj)
173 ENDDO
174 ENDDO
175 ENDDO
176 ENDDO
177 cs- end block
178 ENDIF
179
180 CALL MDSREADFIELD('RAZ.bin',readBinaryPrec,'RS',1,rAz, 1,myThid)
181 IF (useCubedSphereExchange) THEN
182 cs- this block needed by cubed sphere until we write more useful I/O routines
183 CALL EXCH_Z_3D_RS( rAz, 1, myThid )
184 DO bj = myByLo(myThid), myByHi(myThid)
185 DO bi = myBxLo(myThid), myBxHi(myThid)
186 rAz(sNx+1,1,bi,bj)=rAz(1,1,bi,bj)
187 rAz(1,sNy+1,bi,bj)=rAz(1,1,bi,bj)
188 ENDDO
189 ENDDO
190 cs- end block
191 ENDIF
192 CALL EXCH_Z_3D_RS( rAz, 1, myThid )
193
194 C- Staggered (u,v pairs) quantities
195 CALL MDSREADFIELD('DXC.bin',readBinaryPrec,'RS',1,dxC, 1,myThid)
196 CALL MDSREADFIELD('DYC.bin',readBinaryPrec,'RS',1,dyC, 1,myThid)
197 CALL EXCH_UV_XY_RS(dxC,dyC,.FALSE.,myThid)
198
199 CALL MDSREADFIELD('RAW.bin',readBinaryPrec,'RS',1,rAw, 1,myThid)
200 CALL MDSREADFIELD('RAS.bin',readBinaryPrec,'RS',1,rAs, 1,myThid)
201 CALL EXCH_UV_XY_RS(rAw,rAs,.FALSE.,myThid)
202
203 CALL MDSREADFIELD('DXG.bin',readBinaryPrec,'RS',1,dxG, 1,myThid)
204 CALL MDSREADFIELD('DYG.bin',readBinaryPrec,'RS',1,dyG, 1,myThid)
205 CALL EXCH_UV_XY_RS(dyG,dxG,.FALSE.,myThid)
206 CALL EXCH_UV_AGRID_3D_RS(angleSinC,angleCosC,.TRUE., 1, myThid)
207 anglesAreSet = .FALSE.
208
209 c write(10) xC
210 c write(10) yC
211 c write(10) dxF
212 c write(10) dyF
213 c write(10) rA
214 c write(10) xG
215 c write(10) yG
216 c write(10) dxV
217 c write(10) dyU
218 c write(10) rAz
219 c write(10) dxC
220 c write(10) dyC
221 c write(10) rAw
222 c write(10) rAs
223 c write(10) dxG
224 c write(10) dyG
225
226 #else /* ifndef OLD_GRID_IO */
227
228 C-- Only do I/O if I am the master thread
229 _BEGIN_MASTER(myThid)
230
231 #ifdef ALLOW_MNC
232 IF (useMNC .AND. readgrid_mnc) THEN
233 C-- read NetCDF files:
234
235 DO i = 1,80
236 mncFn(i:i) = ' '
237 ENDDO
238 write(mncFn,'(a)') 'mitgrid'
239 DO i = 1,MAX_LEN_MBUF
240 msgBuf(i:i) = ' '
241 ENDDO
242 WRITE(msgBuf,'(2A)') msgBuf,' ; Reading grid info using MNC'
243 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
244 & SQUEEZE_RIGHT , myThid)
245 CALL MNC_FILE_CLOSE_ALL_MATCHING(mncFn, myThid)
246 CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
247 CALL MNC_CW_SET_CITER(mncFn, 2, -1, -1, -1, myThid)
248 CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
249 CALL MNC_CW_RS_R('D',mncFn,0,0,'XC', xC, myThid)
250 CALL MNC_CW_RS_R('D',mncFn,0,0,'XG', xG, myThid)
251 CALL MNC_CW_RS_R('D',mncFn,0,0,'YC', yC, myThid)
252 CALL MNC_CW_RS_R('D',mncFn,0,0,'YG', yG, myThid)
253 CALL MNC_CW_RS_R('D',mncFn,0,0,'dxC',dxC, myThid)
254 CALL MNC_CW_RS_R('D',mncFn,0,0,'dyC',dyC, myThid)
255 CALL MNC_CW_RS_R('D',mncFn,0,0,'dxF',dxF, myThid)
256 CALL MNC_CW_RS_R('D',mncFn,0,0,'dyF',dyF, myThid)
257 CALL MNC_CW_RS_R('D',mncFn,0,0,'dxG',dxG, myThid)
258 CALL MNC_CW_RS_R('D',mncFn,0,0,'dyG',dyG, myThid)
259 CALL MNC_CW_RS_R('D',mncFn,0,0,'dxV',dxV, myThid)
260 CALL MNC_CW_RS_R('D',mncFn,0,0,'dyU',dyU, myThid)
261 CALL MNC_CW_RS_R('D',mncFn,0,0,'rA', rA, myThid)
262 CALL MNC_CW_RS_R('D',mncFn,0,0,'rAz',rAz, myThid)
263 CALL MNC_CW_RS_R('D',mncFn,0,0,'rAw',rAw, myThid)
264 CALL MNC_CW_RS_R('D',mncFn,0,0,'rAs',rAs, myThid)
265 CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleCS',angleCosC,myThid)
266 CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleSN',angleSinC,myThid)
267 anglesAreSet = .TRUE.
268
269 ELSE
270 C-- read Binary files:
271 #endif /* ALLOW_MNC */
272
273 DO bj = 1,nSy
274 DO bi = 1,nSx
275 iG=bi+(myXGlobalLo-1)/sNx
276 WRITE(tmpBuf,'(A,I4)') 'tile:',iG
277 #ifdef ALLOW_EXCH2
278 myTile = W2_myTileList(bi)
279 WRITE(tmpBuf,'(A,I4)') 'tile:',myTile
280 iG = exch2_myface(myTile)
281 #endif
282 iLen = ILNBLNK(horizGridFile)
283 IF ( iLen .EQ. 0 ) THEN
284 WRITE(fName,'("tile",I3.3,".mitgrid")') iG
285 ELSE
286 WRITE(fName,'(2A,I3.3,A)') horizGridFile(1:iLen),
287 & '.face',iG,'.bin'
288 ENDIF
289 iLen = ILNBLNK(fName)
290 iL = ILNBLNK(tmpBuf)
291 WRITE(msgBuf,'(3A)') tmpBuf(1:iL),
292 & ' ; Read from file ',fName(1:iLen)
293 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
294 & SQUEEZE_RIGHT , myThid)
295 WRITE(msgBuf,'(A)') ' =>'
296
297 CALL READSYMTILE_RS(fName,1,xC,buf,bi,bj,myThid)
298 iL = ILNBLNK(msgBuf)
299 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'xC'
300 CALL READSYMTILE_RS(fName,2,yC,buf,bi,bj,myThid)
301 iL = ILNBLNK(tmpBuf)
302 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'yC'
303 CALL READSYMTILE_RS(fName,3,dxF,buf,bi,bj,myThid)
304 iL = ILNBLNK(msgBuf)
305 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxF'
306 CALL READSYMTILE_RS(fName,4,dyF,buf,bi,bj,myThid)
307 iL = ILNBLNK(tmpBuf)
308 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyF'
309 CALL READSYMTILE_RS(fName,5,rA,buf,bi,bj,myThid)
310 iL = ILNBLNK(msgBuf)
311 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'rA'
312 CALL READSYMTILE_RS(fName,6,xG,buf,bi,bj,myThid)
313 iL = ILNBLNK(tmpBuf)
314 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'xG'
315 CALL READSYMTILE_RS(fName,7,yG,buf,bi,bj,myThid)
316 iL = ILNBLNK(msgBuf)
317 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'yG'
318 CALL READSYMTILE_RS(fName,8,dxV,buf,bi,bj,myThid)
319 iL = ILNBLNK(tmpBuf)
320 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dxV'
321 CALL READSYMTILE_RS(fName,9,dyU,buf,bi,bj,myThid)
322 iL = ILNBLNK(msgBuf)
323 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dyU'
324 CALL READSYMTILE_RS(fName,10,rAz,buf,bi,bj,myThid)
325 iL = ILNBLNK(tmpBuf)
326 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'rAz'
327 CALL READSYMTILE_RS(fName,11,dxC,buf,bi,bj,myThid)
328 iL = ILNBLNK(msgBuf)
329 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxC'
330 CALL READSYMTILE_RS(fName,12,dyC,buf,bi,bj,myThid)
331 iL = ILNBLNK(tmpBuf)
332 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyC'
333 CALL READSYMTILE_RS(fName,13,rAw,buf,bi,bj,myThid)
334 iL = ILNBLNK(msgBuf)
335 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'rAw'
336 CALL READSYMTILE_RS(fName,14,rAs,buf,bi,bj,myThid)
337 iL = ILNBLNK(tmpBuf)
338 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'rAs'
339 CALL READSYMTILE_RS(fName,15,dxG,buf,bi,bj,myThid)
340 iL = ILNBLNK(msgBuf)
341 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxG'
342 CALL READSYMTILE_RS(fName,16,dyG,buf,bi,bj,myThid)
343 iL = ILNBLNK(tmpBuf)
344 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyG'
345
346 iLen = ILNBLNK(horizGridFile)
347 IF ( iLen.GT.0 ) THEN
348 CALL READSYMTILE_RS(fName,17,angleCosC,buf,bi,bj,myThid)
349 iL = ILNBLNK(msgBuf)
350 WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'AngleCS'
351 CALL READSYMTILE_RS(fName,18,angleSinC,buf,bi,bj,myThid)
352 iL = ILNBLNK(tmpBuf)
353 WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'AngleSN'
354 anglesAreSet = .TRUE.
355 ELSE
356 anglesAreSet = .FALSE.
357 ENDIF
358
359 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
360 & SQUEEZE_RIGHT , myThid)
361
362 ENDDO
363 ENDDO
364
365 #ifdef ALLOW_MNC
366 ENDIF
367 #endif /* ALLOW_MNC */
368
369 _END_MASTER(myThid)
370
371 CALL EXCH_XY_RS(xC,myThid)
372 CALL EXCH_XY_RS(yC,myThid)
373 CALL EXCH_UV_AGRID_3D_RS( dxF, dyF, .FALSE., 1, myThid )
374 CALL EXCH_XY_RS(rA,myThid )
375 CALL EXCH_Z_3D_RS( xG, 1, myThid )
376 CALL EXCH_Z_3D_RS( yG, 1, myThid )
377 #ifdef ALLOW_EXCH2
378 # ifndef ALLOW_AUTODIFF_TAMC
379 CALL EXCH2_UV_BGRID_3D_RS( dxV, dyU, .FALSE., 1, myThid)
380 # endif
381 #endif
382 CALL EXCH_Z_3D_RS( rAz, 1, myThid )
383 CALL EXCH_UV_XY_RS(dxC,dyC,.FALSE.,myThid)
384 CALL EXCH_UV_XY_RS(rAw,rAs,.FALSE.,myThid)
385 CALL EXCH_UV_XY_RS(dyG,dxG,.FALSE.,myThid)
386 CALL EXCH_UV_AGRID_3D_RS(angleSinC,angleCosC,.TRUE., 1, myThid)
387
388 #endif /* OLD_GRID_IO */
389
390 C-- Stop if Angle have not been loaded but are needed :
391 _BEGIN_MASTER(myThid)
392 IF ( .NOT.anglesAreSet .AND. use3dCoriolis ) THEN
393 WRITE(msgBuf,'(2A)')
394 & 'INI_CURVILINEAR_GRID: Angle of CurvilinearGrid not set',
395 & ' but needed for 3-D Coriolis'
396 CALL PRINT_ERROR( msgBuf , myThid)
397 STOP 'ABNORMAL END: S/R INI_CURVILINEAR_GRID'
398 ENDIF
399 _END_MASTER(myThid)
400
401 c CALL WRITE_FULLARRAY_RL('dxV',dxV,1,0,0,0,myThid)
402 c CALL WRITE_FULLARRAY_RL('dyU',dyU,1,0,0,0,myThid)
403 c CALL WRITE_FULLARRAY_RL('rAz',rAz,1,0,0,0,myThid)
404 c CALL WRITE_FULLARRAY_RL('xG',xG,1,0,0,0,myThid)
405 c CALL WRITE_FULLARRAY_RL('yG',yG,1,0,0,0,myThid)
406
407 C-- Now let's look at all these beasts
408 IF ( debugLevel .GE. debLevB ) THEN
409 CALL PLOT_FIELD_XYRS( xC , 'Current xC ', 0, myThid )
410 CALL PLOT_FIELD_XYRS( yC , 'Current yC ', 0, myThid )
411 CALL PLOT_FIELD_XYRS( dxF , 'Current dxF ', 0, myThid )
412 CALL PLOT_FIELD_XYRS( dyF , 'Current dyF ', 0, myThid )
413 CALL PLOT_FIELD_XYRS( rA , 'Current rA ', 0, myThid )
414 CALL PLOT_FIELD_XYRS( xG , 'Current xG ', 0, myThid )
415 CALL PLOT_FIELD_XYRS( yG , 'Current yG ', 0, myThid )
416 CALL PLOT_FIELD_XYRS( dxV , 'Current dxV ', 0, myThid )
417 CALL PLOT_FIELD_XYRS( dyU , 'Current dyU ', 0, myThid )
418 CALL PLOT_FIELD_XYRS( rAz , 'Current rAz ', 0, myThid )
419 CALL PLOT_FIELD_XYRS( dxC , 'Current dxC ', 0, myThid )
420 CALL PLOT_FIELD_XYRS( dyC , 'Current dyC ', 0, myThid )
421 CALL PLOT_FIELD_XYRS( rAw , 'Current rAw ', 0, myThid )
422 CALL PLOT_FIELD_XYRS( rAs , 'Current rAs ', 0, myThid )
423 CALL PLOT_FIELD_XYRS( dxG , 'Current dxG ', 0, myThid )
424 CALL PLOT_FIELD_XYRS( dyG , 'Current dyG ', 0, myThid )
425 CALL PLOT_FIELD_XYRS(angleCosC, 'Current AngleCS ', 0, myThid )
426 CALL PLOT_FIELD_XYRS(angleSinC, 'Current AngleSN ', 0, myThid )
427 ENDIF
428
429 RETURN
430 END
431
432 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
433
434 CBOP
435 C !ROUTINE: READSYMTILE_RS
436 C !INTERFACE:
437 SUBROUTINE READSYMTILE_RS(
438 I fName, irec,
439 U array, buf,
440 I bi,bj, myThid )
441 C !DESCRIPTION: \bv
442 C *==========================================================*
443 C | SUBROUTINE READSYMTILE_RS
444 C *==========================================================*
445 C *==========================================================*
446 C \ev
447
448 C !USES:
449 IMPLICIT NONE
450 C === Global variables ===
451 #include "SIZE.h"
452 #include "EEPARAMS.h"
453 #ifdef ALLOW_EXCH2
454 #include "W2_EXCH2_TOPOLOGY.h"
455 #include "W2_EXCH2_PARAMS.h"
456 #endif /* ALLOW_EXCH2 */
457
458 C !INPUT/OUTPUT PARAMETERS:
459 C == Routine arguments ==
460 CHARACTER*(*) fName
461 INTEGER irec
462 _RS array(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
463 #ifdef ALLOW_EXCH2
464 _RL buf(1:sNx*nSx*nPx+1)
465 #else
466 _RL buf(1:sNx+1,1:sNy+1)
467 #endif /* ALLOW_EXCH2 */
468 INTEGER bi,bj, myThid
469 CEOP
470
471 C !LOCAL VARIABLES:
472 C == Local variables ==
473 INTEGER I,J,dUnit, iLen
474 INTEGER length_of_rec
475 INTEGER MDS_RECLEN
476 #ifdef ALLOW_EXCH2
477 INTEGER TN, dNx, dNy, TBX, TBY, TNX, TNY, II, iBase
478 #endif
479 INTEGER ILNBLNK
480 EXTERNAL ILNBLNK
481
482 iLen = ILNBLNK(fName)
483 #ifdef ALLOW_EXCH2
484 C Figure out offset of tile within face
485 TN = W2_myTileList(bi)
486 dNx = exch2_mydnx(TN)
487 dNy = exch2_mydny(TN)
488 TBX = exch2_tbasex(TN)
489 TBY = exch2_tbasey(TN)
490 TNX = exch2_tnx(TN)
491 TNY = exch2_tny(TN)
492
493 CALL MDSFINDUNIT( dUnit, myThid )
494 length_of_rec=MDS_RECLEN( 64, (dNx+1), myThid )
495 OPEN( dUnit, file=fName(1:iLen), status='old',
496 & access='direct', recl=length_of_rec )
497 J=0
498 iBase=(irec-1)*(dny+1)
499 DO I=1+TBY,sNy+1+TBY
500 READ(dUnit,rec=I+iBase)(buf(ii),ii=1,dNx+1)
501 #ifdef _BYTESWAPIO
502 #ifdef REAL4_IS_SLOW
503 CALL MDS_BYTESWAPR8((dNx+1), buf)
504 #else
505 CALL MDS_BYTESWAPR4((dNx+1), buf)
506 #endif
507 #endif
508 J=J+1
509 DO II=1,sNx+1
510 array(II,J,bi,bj)=buf(II+TBX)
511 ENDDO
512 ENDDO
513 CLOSE( dUnit )
514
515 #else /* ALLOW_EXCH2 */
516
517 CALL MDSFINDUNIT( dUnit, myThid )
518 length_of_rec=MDS_RECLEN( 64, (sNx+1)*(sNy+1), myThid )
519 OPEN( dUnit, file=fName(1:iLen), status='old',
520 & access='direct', recl=length_of_rec )
521 READ(dUnit,rec=irec) buf
522 CLOSE( dUnit )
523
524 #ifdef _BYTESWAPIO
525 #ifdef REAL4_IS_SLOW
526 CALL MDS_BYTESWAPR8((sNx+1)*(sNy+1), buf)
527 #else
528 CALL MDS_BYTESWAPR4((sNx+1)*(sNy+1), buf)
529 #endif
530 #endif
531
532 DO J=1,sNy+1
533 DO I=1,sNx+1
534 array(I,J,bi,bj)=buf(I,J)
535 ENDDO
536 ENDDO
537 c write(0,*) irec,buf(1,1),array(1,1,1,1)
538
539 #endif /* ALLOW_EXCH2 */
540
541 RETURN
542 END

  ViewVC Help
Powered by ViewVC 1.1.22