/[MITgcm]/MITgcm_contrib/ifenty/Curvi/code/ini_curvilinear_grid.F
ViewVC logotype

Annotation of /MITgcm_contrib/ifenty/Curvi/code/ini_curvilinear_grid.F

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


Revision 1.1 - (hide annotations) (download)
Wed Jul 5 16:36:02 2006 UTC (19 years ago) by ifenty
Branch: MAIN
CVS Tags: HEAD
Initial checkin of curvilinear configuration of Labrador Sea region

1 ifenty 1.1 C $Header: /u/gcmpack/MITgcm/model/src/ini_curvilinear_grid.F,v 1.24 2005/09/28 01:43:55 edhill 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     #ifndef ALLOW_EXCH2
39     C- note: default is to use "new" grid files (OLD_GRID_IO undef) with EXCH2
40     C but can still use (on 1 cpu) OLD_GRID_IO and EXCH2 independently
41     #ifdef ALLOW_MDSIO
42     #define OLD_GRID_IO
43     #endif
44     #endif /* ALLOW_EXCH2 */
45    
46     C !INPUT/OUTPUT PARAMETERS:
47     C == Routine arguments ==
48     C myThid - Number of this instance of INI_CURVILINEAR_GRID
49     INTEGER myThid
50    
51     C !LOCAL VARIABLES:
52     C == Local variables ==
53     INTEGER bi,bj, myIter
54     INTEGER I,J
55     CHARACTER*(MAX_LEN_FNAM) fName
56     #ifdef ALLOW_MNC
57     CHARACTER*(80) mncFn
58     #endif
59     #ifdef ALLOW_EXCH2
60     _RL buf(sNx*nSx*nPx+1)
61     INTEGER myTile
62     #else
63     _RL buf(sNx+1,sNy+1)
64     #endif
65     INTEGER iG, iL, iLen
66     CHARACTER*(MAX_LEN_MBUF) msgBuf, tmpBuf
67     INTEGER ILNBLNK
68     EXTERNAL ILNBLNK
69     CEOP
70    
71     C-- Set everything to zero everywhere
72     DO bj = myByLo(myThid), myByHi(myThid)
73     DO bi = myBxLo(myThid), myBxHi(myThid)
74    
75     DO J=1-Oly,sNy+Oly
76     DO I=1-Olx,sNx+Olx
77     XC(i,j,bi,bj)=0.
78     YC(i,j,bi,bj)=0.
79     XG(i,j,bi,bj)=0.
80     YG(i,j,bi,bj)=0.
81     DXC(i,j,bi,bj)=0.
82     DYC(i,j,bi,bj)=0.
83     DXG(i,j,bi,bj)=0.
84     DYG(i,j,bi,bj)=0.
85     DXF(i,j,bi,bj)=0.
86     DYF(i,j,bi,bj)=0.
87     DXV(i,j,bi,bj)=0.
88     DYU(i,j,bi,bj)=0.
89     RA(i,j,bi,bj)=0.
90     RAZ(i,j,bi,bj)=0.
91     RAW(i,j,bi,bj)=0.
92     RAS(i,j,bi,bj)=0.
93     tanPhiAtU(i,j,bi,bj)=0.
94     tanPhiAtV(i,j,bi,bj)=0.
95     angleCosC(i,j,bi,bj)=1.
96     angleSinC(i,j,bi,bj)=0.
97     cosFacU(J,bi,bj)=1.
98     cosFacV(J,bi,bj)=1.
99     sqcosFacU(J,bi,bj)=1.
100     sqcosFacV(J,bi,bj)=1.
101     ENDDO
102     ENDDO
103    
104     ENDDO
105     ENDDO
106    
107    
108     #ifdef ALLOW_MNC
109     IF (useMNC .AND. readgrid_mnc) THEN
110    
111     _BEGIN_MASTER(myThid)
112     DO i = 1,80
113     mncFn(i:i) = ' '
114     ENDDO
115     write(mncFn,'(a)') 'mitgrid'
116     DO i = 1,MAX_LEN_MBUF
117     msgBuf(i:i) = ' '
118     ENDDO
119     WRITE(msgBuf,'(2A)') msgBuf,' ; Reading grid info using MNC'
120     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
121     & SQUEEZE_RIGHT , myThid)
122     CALL MNC_FILE_CLOSE_ALL_MATCHING(mncFn, myThid)
123     CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
124     CALL MNC_CW_SET_CITER(mncFn, 2, -1, -1, -1, myThid)
125     CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
126     CALL MNC_CW_RS_R('D',mncFn,0,0,'XC', XC, myThid)
127     CALL MNC_CW_RS_R('D',mncFn,0,0,'XG', XG, myThid)
128     CALL MNC_CW_RS_R('D',mncFn,0,0,'YC', YC, myThid)
129     CALL MNC_CW_RS_R('D',mncFn,0,0,'YG', YG, myThid)
130     CALL MNC_CW_RS_R('D',mncFn,0,0,'dxC',DXC, myThid)
131     CALL MNC_CW_RS_R('D',mncFn,0,0,'dyC',DYC, myThid)
132     CALL MNC_CW_RS_R('D',mncFn,0,0,'dxF',DXF, myThid)
133     CALL MNC_CW_RS_R('D',mncFn,0,0,'dyF',DYF, myThid)
134     CALL MNC_CW_RS_R('D',mncFn,0,0,'dxG',DXG, myThid)
135     CALL MNC_CW_RS_R('D',mncFn,0,0,'dyG',DYG, myThid)
136     CALL MNC_CW_RS_R('D',mncFn,0,0,'dxV',DXV, myThid)
137     CALL MNC_CW_RS_R('D',mncFn,0,0,'dyU',DYU, myThid)
138     CALL MNC_CW_RS_R('D',mncFn,0,0,'rA', RA, myThid)
139     CALL MNC_CW_RS_R('D',mncFn,0,0,'rAz',RAZ, myThid)
140     CALL MNC_CW_RS_R('D',mncFn,0,0,'rAw',RAW, myThid)
141     CALL MNC_CW_RS_R('D',mncFn,0,0,'rAs',RAS, myThid)
142     CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleCS',angleCosC,myThid)
143     CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleSN',angleSinC,myThid)
144    
145     _END_MASTER(myThid)
146    
147     CALL EXCH_XY_RS(XC,myThid)
148     CALL EXCH_XY_RS(YC,myThid)
149     #ifdef HRCUBE
150     CALL EXCH_XY_RS(DXF,myThid)
151     CALL EXCH_XY_RS(DYF,myThid)
152     #endif
153     CALL EXCH_XY_RS(RA,myThid )
154     CALL EXCH_Z_XY_RS(XG,myThid)
155     CALL EXCH_Z_XY_RS(YG,myThid)
156     CALL EXCH_Z_XY_RS(RAZ,myThid)
157     CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
158     CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
159     CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
160     CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
161    
162     ELSE
163     #endif
164    
165     C Here we make no assumptions about grid symmetry and simply
166     C read the raw grid data from files
167    
168     #ifdef MONKEY_BUTT_OLD_GRID_IO
169    
170     C- Cell centered quantities
171     CALL MDSREADFIELD('LONC.bin',readBinaryPrec,'RS',1,XC, 1,myThid)
172     CALL MDSREADFIELD('LATC.bin',readBinaryPrec,'RS',1,YC, 1,myThid)
173     _EXCH_XY_R4(XC,myThid)
174     _EXCH_XY_R4(YC,myThid)
175    
176     CALL MDSREADFIELD('DXF.bin',readBinaryPrec,'RS',1,DXF, 1,myThid)
177     CALL MDSREADFIELD('DYF.bin',readBinaryPrec,'RS',1,DYF, 1,myThid)
178     C !!! _EXCH_OUV_XY_R4(DXF, DYF, unSigned, myThid )
179     cs! this is not correct! <= need paired exchange for DXF,DYF
180     _EXCH_XY_R4(DXF,myThid)
181     _EXCH_XY_R4(DYF,myThid)
182     IF (useCubedSphereExchange) THEN
183     cs! fix overlaps:
184     DO bj = myByLo(myThid), myByHi(myThid)
185     DO bi = myBxLo(myThid), myBxHi(myThid)
186     DO j=1,sNy
187     DO i=1,Olx
188     DXF(1-i,j,bi,bj)=DXF(i,j,bi,bj)
189     DXF(sNx+i,j,bi,bj)=DXF(sNx+1-i,j,bi,bj)
190     DYF(1-i,j,bi,bj)=DYF(i,j,bi,bj)
191     DYF(sNx+i,j,bi,bj)=DYF(sNx+1-i,j,bi,bj)
192     ENDDO
193     ENDDO
194     DO j=1,Oly
195     DO i=1,sNx
196     DXF(i,1-j,bi,bj)=DXF(i,j,bi,bj)
197     DXF(i,sNy+j,bi,bj)=DXF(i,sNy+1-j,bi,bj)
198     DYF(i,1-j,bi,bj)=DYF(i,j,bi,bj)
199     DYF(i,sNy+j,bi,bj)=DYF(i,sNy+1-j,bi,bj)
200     ENDDO
201     ENDDO
202     ENDDO
203     ENDDO
204     ENDIF
205     cs
206    
207     CALL MDSREADFIELD('RA.bin',readBinaryPrec,'RS',1,RA, 1,myThid)
208     _EXCH_XY_R4(RA,myThid )
209    
210     C- Corner quantities
211     C *********** this are not degbugged ************
212     CALL MDSREADFIELD('LONG.bin',readBinaryPrec,'RS',1,XG, 1,myThid)
213     CALL MDSREADFIELD('LATG.bin',readBinaryPrec,'RS',1,YG, 1,myThid)
214     IF (useCubedSphereExchange) THEN
215     cs- this block needed by cubed sphere until we write more useful I/O routines
216     bi=3
217     bj=1
218     YG(1,sNy+1,bj,1)=YG(1,1,bi,1)
219     bj=bj+2
220     YG(1,sNy+1,bj,1)=YG(1,1,bi,1)
221     bj=bj+2
222     YG(1,sNy+1,bj,1)=YG(1,1,bi,1)
223     bi=6
224     bj=2
225     YG(sNx+1,1,bj,1)=YG(1,1,bi,1)
226     bj=bj+2
227     YG(sNx+1,1,bj,1)=YG(1,1,bi,1)
228     bj=bj+2
229     YG(sNx+1,1,bj,1)=YG(1,1,bi,1)
230     cs- end block
231     ENDIF
232     CALL EXCH_Z_XY_RS(XG,myThid)
233     CALL EXCH_Z_XY_RS(YG,myThid)
234    
235     CALL MDSREADFIELD('DXV.bin',readBinaryPrec,'RS',1,DXV, 1,myThid)
236     CALL MDSREADFIELD('DYU.bin',readBinaryPrec,'RS',1,DYU, 1,myThid)
237     cs- this block needed by cubed sphere until we write more useful I/O routines
238     C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)
239     cs! this is not correct <= need paired exchange for dxv,dyu
240     IF (.NOT.useCubedSphereExchange) THEN
241     CALL EXCH_Z_XY_RS(DXV,myThid)
242     CALL EXCH_Z_XY_RS(DYU,myThid)
243     ELSE
244     DO bj = myByLo(myThid), myByHi(myThid)
245     DO bi = myBxLo(myThid), myBxHi(myThid)
246     cs! fix overlaps:
247     DO j=1,sNy
248     DO i=1,Olx
249     DXV(1-i,j,bi,bj)=DXV(1+i,j,bi,bj)
250     DXV(sNx+i,j,bi,bj)=DXV(i,j,bi,bj)
251     DYU(1-i,j,bi,bj)=DYU(1+i,j,bi,bj)
252     DYU(sNx+i,j,bi,bj)=DYU(i,j,bi,bj)
253     ENDDO
254     ENDDO
255     DO j=1,Oly
256     DO i=1-Olx,sNx+Olx
257     DXV(i,1-j,bi,bj)=DXV(i,1+j,bi,bj)
258     DXV(i,sNy+j,bi,bj)=DXV(i,j,bi,bj)
259     DYU(i,1-j,bi,bj)=DYU(i,1+j,bi,bj)
260     DYU(i,sNy+j,bi,bj)=DYU(i,j,bi,bj)
261     ENDDO
262     ENDDO
263     ENDDO
264     ENDDO
265     cs- end block
266     ENDIF
267    
268     CALL MDSREADFIELD('RAZ.bin',readBinaryPrec,'RS',1,RAZ, 1,myThid)
269     IF (useCubedSphereExchange) THEN
270     cs- this block needed by cubed sphere until we write more useful I/O routines
271     CALL EXCH_Z_XY_RS(RAZ , myThid )
272     DO bj = myByLo(myThid), myByHi(myThid)
273     DO bi = myBxLo(myThid), myBxHi(myThid)
274     RAZ(sNx+1,1,bi,bj)=RAZ(1,1,bi,bj)
275     RAZ(1,sNy+1,bi,bj)=RAZ(1,1,bi,bj)
276     ENDDO
277     ENDDO
278     cs- end block
279     ENDIF
280     CALL EXCH_Z_XY_RS(RAZ,myThid)
281    
282     C- Staggered (u,v pairs) quantities
283     CALL MDSREADFIELD('DXC.bin',readBinaryPrec,'RS',1,DXC, 1,myThid)
284     CALL MDSREADFIELD('DYC.bin',readBinaryPrec,'RS',1,DYC, 1,myThid)
285     CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
286    
287     CALL MDSREADFIELD('RAW.bin',readBinaryPrec,'RS',1,RAW, 1,myThid)
288     CALL MDSREADFIELD('RAS.bin',readBinaryPrec,'RS',1,RAS, 1,myThid)
289     IF (useCubedSphereExchange) THEN
290     cs- this block needed by cubed sphere until we write more useful I/O routines
291     DO bj = myByLo(myThid), myByHi(myThid)
292     DO bi = myBxLo(myThid), myBxHi(myThid)
293     DO J = 1,sNy
294     c RAW(sNx+1,J,bi,bj)=RAW(1,J,bi,bj)
295     c RAS(J,sNy+1,bi,bj)=RAS(J,1,bi,bj)
296     ENDDO
297     ENDDO
298     ENDDO
299     cs- end block
300     ENDIF
301     CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
302    
303     CALL MDSREADFIELD('DXG.bin',readBinaryPrec,'RS',1,DXG, 1,myThid)
304     CALL MDSREADFIELD('DYG.bin',readBinaryPrec,'RS',1,DYG, 1,myThid)
305     IF (useCubedSphereExchange) THEN
306     cs- this block needed by cubed sphere until we write more useful I/O routines
307     DO bj = myByLo(myThid), myByHi(myThid)
308     DO bi = myBxLo(myThid), myBxHi(myThid)
309     DO J = 1,sNy
310     c DYG(sNx+1,J,bi,bj)=DYG(1,J,bi,bj)
311     c DXG(J,sNy+1,bi,bj)=DXG(J,1,bi,bj)
312     ENDDO
313     ENDDO
314     ENDDO
315     cs- end block
316     ENDIF
317     CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
318     CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
319    
320     c write(10) XC
321     c write(10) YC
322     c write(10) DXF
323     c write(10) DYF
324     c write(10) RA
325     c write(10) XG
326     c write(10) YG
327     c write(10) DXV
328     c write(10) DYU
329     c write(10) RAZ
330     c write(10) DXC
331     c write(10) DYC
332     c write(10) RAW
333     c write(10) RAS
334     c write(10) DXG
335     c write(10) DYG
336    
337     #else /* ifndef OLD_GRID_IO */
338    
339     C-- Only do I/O if I am the master thread
340     _BEGIN_MASTER(myThid)
341    
342     write(msgbuf,'(3a)') 'Open Bondary File'
343     call print_message( msgbuf, standardmessageunit,
344     & SQUEEZE_RIGHT , mythid)
345    
346     DO bj = 1,nSy
347     DO bi = 1,nSx
348     iG=bi+(myXGlobalLo-1)/sNx
349     WRITE(tmpBuf,'(A,I4)') 'tile:',iG
350     #ifdef ALLOW_EXCH2
351     myTile = W2_myTileList(bi)
352     WRITE(tmpBuf,'(A,I4)') 'tile:',myTile
353     iG = exch2_myface(myTile)
354     #endif
355     iLen = ILNBLNK(horizGridFile)
356     IF ( iLen .EQ. 0 ) THEN
357     WRITE(fName,'("tile",I3.3,".mitgrid")') iG
358     ELSE
359     WRITE(fName,'(2A,I3.3,A)') horizGridFile(1:iLen),
360     & '.face',iG,'.bin'
361     ENDIF
362     iLen = ILNBLNK(fName)
363     iL = ILNBLNK(tmpBuf)
364     WRITE(msgBuf,'(3A)') tmpBuf(1:iL),
365     & ' ; Read from file ',fName(1:iLen)
366     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
367     & SQUEEZE_RIGHT , myThid)
368     WRITE(msgBuf,'(A)') ' =>'
369     CALL READSYMTILE_RS(fName,1,XC,bi,bj,buf,myThid)
370     iL = ILNBLNK(msgBuf)
371     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'XC'
372     CALL READSYMTILE_RS(fName,2,YC,bi,bj,buf,myThid)
373     iL = ILNBLNK(tmpBuf)
374     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'YC'
375     CALL READSYMTILE_RS(fName,3,DXF,bi,bj,buf,myThid)
376     iL = ILNBLNK(msgBuf)
377     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXF'
378     CALL READSYMTILE_RS(fName,4,DYF,bi,bj,buf,myThid)
379     iL = ILNBLNK(tmpBuf)
380     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYF'
381     CALL READSYMTILE_RS(fName,5,RA,bi,bj,buf,myThid)
382     iL = ILNBLNK(msgBuf)
383     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'RA'
384     CALL READSYMTILE_RS(fName,6,XG,bi,bj,buf,myThid)
385     iL = ILNBLNK(tmpBuf)
386     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'XG'
387     CALL READSYMTILE_RS(fName,7,YG,bi,bj,buf,myThid)
388     iL = ILNBLNK(msgBuf)
389     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'YG'
390     CALL READSYMTILE_RS(fName,8,DXV,bi,bj,buf,myThid)
391     iL = ILNBLNK(tmpBuf)
392     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DXV'
393     CALL READSYMTILE_RS(fName,9,DYU,bi,bj,buf,myThid)
394     iL = ILNBLNK(msgBuf)
395     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DYU'
396     CALL READSYMTILE_RS(fName,10,RAZ,bi,bj,buf,myThid)
397     iL = ILNBLNK(tmpBuf)
398     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'RAZ'
399     CALL READSYMTILE_RS(fName,11,DXC,bi,bj,buf,myThid)
400     iL = ILNBLNK(msgBuf)
401     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXC'
402     CALL READSYMTILE_RS(fName,12,DYC,bi,bj,buf,myThid)
403     iL = ILNBLNK(tmpBuf)
404     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYC'
405     CALL READSYMTILE_RS(fName,13,RAW,bi,bj,buf,myThid)
406     iL = ILNBLNK(msgBuf)
407     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'RAW'
408     CALL READSYMTILE_RS(fName,14,RAS,bi,bj,buf,myThid)
409     iL = ILNBLNK(tmpBuf)
410     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'RAS'
411     CALL READSYMTILE_RS(fName,15,DXG,bi,bj,buf,myThid)
412     iL = ILNBLNK(msgBuf)
413     WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXG'
414     CALL READSYMTILE_RS(fName,16,DYG,bi,bj,buf,myThid)
415     iL = ILNBLNK(tmpBuf)
416     WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYG'
417     c iLen = ILNBLNK(horizGridFile)
418     c IF ( iLen.GT.0 ) THEN
419     c CALL READSYMTILE_RS(fName,17,angleCosC,bi,bj,buf,myThid)
420     c iL = ILNBLNK(msgBuf)
421     c WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'AngleCS'
422     c CALL READSYMTILE_RS(fName,18,angleSinC,bi,bj,buf,myThid)
423     c iL = ILNBLNK(tmpBuf)
424     c WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'AngleSN'
425     c ENDIF
426    
427     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
428     & SQUEEZE_RIGHT , myThid)
429    
430     ENDDO
431     ENDDO
432    
433     _END_MASTER(myThid)
434    
435     CALL EXCH_XY_RS(XC,myThid)
436     CALL EXCH_XY_RS(YC,myThid)
437     C !!! _EXCH_OUV_XY_R4(DXF, DYF, unSigned, myThid )
438     #ifdef HRCUBE
439     CALL EXCH_XY_RS(DXF,myThid)
440     CALL EXCH_XY_RS(DYF,myThid)
441     #endif
442     CALL EXCH_XY_RS(RA,myThid )
443     CALL EXCH_Z_XY_RS(XG,myThid)
444     CALL EXCH_Z_XY_RS(YG,myThid)
445     C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)
446     c CALL EXCH_Z_XY_RS(DXV,myThid)
447     c CALL EXCH_Z_XY_RS(DYU,myThid)
448     CALL EXCH_Z_XY_RS(RAZ,myThid)
449     CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
450     CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
451     CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
452     CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
453    
454     #endif /* OLD_GRID_IO */
455    
456     #ifdef ALLOW_MNC
457     ENDIF
458     #endif /* ALLOW_MNC */
459    
460     c CALL WRITE_FULLARRAY_RL('DXV',DXV,1,0,0,0,myThid)
461     c CALL WRITE_FULLARRAY_RL('DYU',DYU,1,0,0,0,myThid)
462     c CALL WRITE_FULLARRAY_RL('RAZ',RAZ,1,0,0,0,myThid)
463     c CALL WRITE_FULLARRAY_RL('XG',XG,1,0,0,0,myThid)
464     c CALL WRITE_FULLARRAY_RL('YG',YG,1,0,0,0,myThid)
465    
466     C-- Require that 0 <= longitude < 360 if using exf package
467     #ifdef ALLOW_EXF
468     DO bj = 1,nSy
469     DO bi = 1,nSx
470     DO J=1-Oly,sNy+Oly
471     DO I=1-Olx,sNx+Olx
472     IF (XC(i,j,bi,bj).lt.0.) XC(i,j,bi,bj)=XC(i,j,bi,bj)+360.
473     IF (XG(i,j,bi,bj).lt.0.) XG(i,j,bi,bj)=XG(i,j,bi,bj)+360.
474     ENDDO
475     ENDDO
476     ENDDO
477     ENDDO
478     #endif /* ALLOW_EXF */
479    
480     C-- Now let's look at all these beasts
481     IF ( debugLevel .GE. debLevB ) THEN
482     myIter = 1
483     CALL PLOT_FIELD_XYRL( XC , 'Current XC ' ,
484     & myIter, myThid )
485     CALL PLOT_FIELD_XYRL( YC , 'Current YC ' ,
486     & myIter, myThid )
487     CALL PLOT_FIELD_XYRL( DXF , 'Current DXF ' ,
488     & myIter, myThid )
489     CALL PLOT_FIELD_XYRL( XC , 'Current XC ' ,
490     & myIter, myThid )
491     CALL PLOT_FIELD_XYRL( DYF , 'Current DYF ' ,
492     & myIter, myThid )
493     CALL PLOT_FIELD_XYRL( RA , 'Current RA ' ,
494     & myIter, myThid )
495     CALL PLOT_FIELD_XYRL( XG , 'Current XG ' ,
496     & myIter, myThid )
497     CALL PLOT_FIELD_XYRL( YG , 'Current YG ' ,
498     & myIter, myThid )
499     CALL PLOT_FIELD_XYRL( DXV , 'Current DXV ' ,
500     & myIter, myThid )
501     CALL PLOT_FIELD_XYRL( DYU , 'Current DYU ' ,
502     & myIter, myThid )
503     CALL PLOT_FIELD_XYRL( RAZ , 'Current RAZ ' ,
504     & myIter, myThid )
505     CALL PLOT_FIELD_XYRL( DXC , 'Current DXC ' ,
506     & myIter, myThid )
507     CALL PLOT_FIELD_XYRL( DYC , 'Current DYC ' ,
508     & myIter, myThid )
509     CALL PLOT_FIELD_XYRL( RAW , 'Current RAW ' ,
510     & myIter, myThid )
511     CALL PLOT_FIELD_XYRL( RAS , 'Current RAS ' ,
512     & myIter, myThid )
513     CALL PLOT_FIELD_XYRL( DXG , 'Current DXG ' ,
514     & myIter, myThid )
515     CALL PLOT_FIELD_XYRL( DYG , 'Current DYG ' ,
516     & myIter, myThid )
517     CALL PLOT_FIELD_XYRL(angleCosC, 'Current AngleCS ' ,
518     & myIter, myThid )
519     CALL PLOT_FIELD_XYRL(angleSinC, 'Current AngleSN ' ,
520     & myIter, myThid )
521     ENDIF
522    
523     RETURN
524     END
525    
526     C --------------------------------------------------------------------------
527    
528     SUBROUTINE READSYMTILE_RS(fName,irec,array,bi,bj,buf,myThid)
529     C /==========================================================\
530     C | SUBROUTINE READSYMTILE_RS |
531     C |==========================================================|
532     C \==========================================================/
533     IMPLICIT NONE
534    
535     C === Global variables ===
536     #include "SIZE.h"
537     #include "EEPARAMS.h"
538     #ifdef ALLOW_EXCH2
539     #include "W2_EXCH2_TOPOLOGY.h"
540     #include "W2_EXCH2_PARAMS.h"
541     #endif /* ALLOW_EXCH2 */
542    
543     C == Routine arguments ==
544     CHARACTER*(*) fName
545     INTEGER irec
546     _RS array(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
547     INTEGER bi,bj,myThid
548     #ifdef ALLOW_EXCH2
549     _RL buf(1:sNx*nSx*nPx+1)
550     #else
551     _RL buf(1:sNx+1,1:sNy+1)
552     #endif /* ALLOW_EXCH2 */
553    
554     C == Local variables ==
555     INTEGER I,J,dUnit, iLen
556     INTEGER length_of_rec
557     INTEGER MDS_RECLEN
558     #ifdef ALLOW_EXCH2
559     INTEGER TN, dNx, dNy, TBX, TBY, TNX, TNY, II, iBase
560     #endif
561     INTEGER ILNBLNK
562     EXTERNAL ILNBLNK
563    
564     iLen = ILNBLNK(fName)
565     #ifdef ALLOW_EXCH2
566     C Figure out offset of tile within face
567     TN = W2_myTileList(bi)
568     dNx = exch2_mydnx(TN)
569     dNy = exch2_mydny(TN)
570     TBX = exch2_tbasex(TN)
571     TBY = exch2_tbasey(TN)
572     TNX = exch2_tnx(TN)
573     TNY = exch2_tny(TN)
574    
575     CALL MDSFINDUNIT( dUnit, myThid )
576     length_of_rec=MDS_RECLEN( 64, (dNx+1), myThid )
577     OPEN( dUnit, file=fName(1:iLen), status='old',
578     & access='direct', recl=length_of_rec )
579     J=0
580     iBase=(irec-1)*(dny+1)
581     DO I=1+TBY,sNy+1+TBY
582     READ(dUnit,rec=I+iBase)(buf(ii),ii=1,dNx+1)
583     #ifdef _BYTESWAPIO
584     #ifdef REAL4_IS_SLOW
585     CALL MDS_BYTESWAPR8((dNx+1), buf)
586     #else
587     CALL MDS_BYTESWAPR4((dNx+1), buf)
588     #endif
589     #endif
590     J=J+1
591     DO II=1,sNx+1
592     array(II,J,bi,bj)=buf(II+TBX)
593     ENDDO
594     ENDDO
595     CLOSE( dUnit )
596    
597     #else /* ALLOW_EXCH2 */
598    
599     CALL MDSFINDUNIT( dUnit, myThid )
600     length_of_rec=MDS_RECLEN( 64, (sNx+1)*(sNy+1), myThid )
601     OPEN( dUnit, file=fName(1:iLen), status='old',
602     & access='direct', recl=length_of_rec )
603     READ(dUnit,rec=irec) buf
604     CLOSE( dUnit )
605    
606     #ifdef _BYTESWAPIO
607     #ifdef REAL4_IS_SLOW
608     CALL MDS_BYTESWAPR8((sNx+1)*(sNy+1), buf)
609     #else
610     CALL MDS_BYTESWAPR4((sNx+1)*(sNy+1), buf)
611     #endif
612     #endif
613    
614     DO J=1,sNy+1
615     DO I=1,sNx+1
616     array(I,J,bi,bj)=buf(I,J)
617     ENDDO
618     ENDDO
619     c write(0,*) irec,buf(1,1),array(1,1,1,1)
620    
621     #endif /* ALLOW_EXCH2 */
622    
623     RETURN
624     END

  ViewVC Help
Powered by ViewVC 1.1.22