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

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

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

revision 1.2 by adcroft, Tue May 29 14:01:37 2001 UTC revision 1.23 by edhill, Sat Sep 17 03:17:05 2005 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4    #include "PACKAGES_CONFIG.h"
5  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
6    
7    CBOP
8    C     !ROUTINE: INI_CURVILINEAR_GRID
9    C     !INTERFACE:
10        SUBROUTINE INI_CURVILINEAR_GRID( myThid )        SUBROUTINE INI_CURVILINEAR_GRID( myThid )
11  C     /==========================================================\  C     !DESCRIPTION: \bv
12  C     | SUBROUTINE INI_CURVILINEAR_GRID                          |  C     *==========================================================*
13  C     | o Initialise curvilinear coordinate system               |  C     | SUBROUTINE INI_CURVILINEAR_GRID                          
14  C     |==========================================================|  C     | o Initialise curvilinear coordinate system                
15  C     \==========================================================/  C     *==========================================================*
16        IMPLICIT NONE  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 ===  C     === Global variables ===
26  #include "SIZE.h"  #include "SIZE.h"
27  #include "EEPARAMS.h"  #include "EEPARAMS.h"
28  #include "PARAMS.h"  #include "PARAMS.h"
29  #include "GRID.h"  #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 ==  C     == Routine arguments ==
48  C     myThid -  Number of this instance of INI_CARTESIAN_GRID  C     myThid -  Number of this instance of INI_CURVILINEAR_GRID
49        INTEGER myThid        INTEGER myThid
50    
51    C     !LOCAL VARIABLES:
52  C     == Local variables ==  C     == Local variables ==
53        INTEGER bi,bj        INTEGER bi,bj, myIter
54        INTEGER I,J        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  C--   Set everything to zero everywhere
72        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
# Line 49  C--   Set everything to zero everywhere Line 92  C--   Set everything to zero everywhere
92            RAS(i,j,bi,bj)=0.            RAS(i,j,bi,bj)=0.
93            tanPhiAtU(i,j,bi,bj)=0.            tanPhiAtU(i,j,bi,bj)=0.
94            tanPhiAtV(i,j,bi,bj)=0.            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.            cosFacU(J,bi,bj)=1.
98            cosFacV(J,bi,bj)=1.            cosFacV(J,bi,bj)=1.
99            sqcosFacU(J,bi,bj)=1.            sqcosFacU(J,bi,bj)=1.
# Line 56  C--   Set everything to zero everywhere Line 101  C--   Set everything to zero everywhere
101           ENDDO           ENDDO
102          ENDDO          ENDDO
103    
104         ENDDO ! bi         ENDDO
105        ENDDO ! bj        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('R',mncFn,0,0,'XC', XC,  myThid)
127            CALL MNC_CW_RS_R('R',mncFn,0,0,'XG', XG,  myThid)
128            CALL MNC_CW_RS_R('R',mncFn,0,0,'YC', YC,  myThid)
129            CALL MNC_CW_RS_R('R',mncFn,0,0,'YG', YG,  myThid)
130            CALL MNC_CW_RS_R('R',mncFn,0,0,'dxC',DXC, myThid)
131            CALL MNC_CW_RS_R('R',mncFn,0,0,'dyC',DYC, myThid)
132            CALL MNC_CW_RS_R('R',mncFn,0,0,'dxF',DXF, myThid)
133            CALL MNC_CW_RS_R('R',mncFn,0,0,'dyF',DYF, myThid)
134            CALL MNC_CW_RS_R('R',mncFn,0,0,'dxG',DXG, myThid)
135            CALL MNC_CW_RS_R('R',mncFn,0,0,'dyG',DYG, myThid)
136            CALL MNC_CW_RS_R('R',mncFn,0,0,'dxV',DXV, myThid)
137            CALL MNC_CW_RS_R('R',mncFn,0,0,'dyU',DYU, myThid)
138            CALL MNC_CW_RS_R('R',mncFn,0,0,'rA', RA,  myThid)
139            CALL MNC_CW_RS_R('R',mncFn,0,0,'rAz',RAZ, myThid)
140            CALL MNC_CW_RS_R('R',mncFn,0,0,'rAw',RAW, myThid)
141            CALL MNC_CW_RS_R('R',mncFn,0,0,'rAs',RAS, myThid)
142    
143            _END_MASTER(myThid)
144    
145            CALL EXCH_XY_RS(XC,myThid)
146            CALL EXCH_XY_RS(YC,myThid)
147    #ifdef HRCUBE
148            CALL EXCH_XY_RS(DXF,myThid)
149            CALL EXCH_XY_RS(DYF,myThid)
150    #endif
151            CALL EXCH_XY_RS(RA,myThid )
152            CALL EXCH_Z_XY_RS(XG,myThid)
153            CALL EXCH_Z_XY_RS(YG,myThid)
154            CALL EXCH_Z_XY_RS(RAZ,myThid)
155            CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
156            CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
157            CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
158            CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
159    
160          ELSE
161    #endif
162    
163  C     Here we make no assumptions about grid symmetry and simply  C     Here we make no assumptions about grid symmetry and simply
164  C     read the raw grid data from files  C     read the raw grid data from files
165    
166    #ifdef OLD_GRID_IO
167    
168  C-    Cell centered quantities  C-    Cell centered quantities
169        CALL MDSREADFIELD('LONC.bin',readBinaryPrec,'RL',1,XC,  1,myThid)        CALL MDSREADFIELD('LONC.bin',readBinaryPrec,'RS',1,XC,  1,myThid)
170        CALL MDSREADFIELD('LATC.bin',readBinaryPrec,'RL',1,YC,  1,myThid)        CALL MDSREADFIELD('LATC.bin',readBinaryPrec,'RS',1,YC,  1,myThid)
171        _EXCH_XY_R4(XC,myThid)        _EXCH_XY_R4(XC,myThid)
172        _EXCH_XY_R4(YC,myThid)        _EXCH_XY_R4(YC,myThid)
173    
174        CALL MDSREADFIELD('DXF.bin',readBinaryPrec,'RL',1,DXF,  1,myThid)        CALL MDSREADFIELD('DXF.bin',readBinaryPrec,'RS',1,DXF,  1,myThid)
175        CALL MDSREADFIELD('DYF.bin',readBinaryPrec,'RL',1,DYF,  1,myThid)        CALL MDSREADFIELD('DYF.bin',readBinaryPrec,'RS',1,DYF,  1,myThid)
176  C !!! _EXCH_OUV_XY_R4(DXF, DYF, unSigned, myThid )  C !!! _EXCH_OUV_XY_R4(DXF, DYF, unSigned, myThid )
177  cs!   this is not correct! <= need paired exchange for DXF,DYF  cs!   this is not correct! <= need paired exchange for DXF,DYF
178        _EXCH_XY_R4(DXF,myThid)        _EXCH_XY_R4(DXF,myThid)
179        _EXCH_XY_R4(DYF,myThid)        _EXCH_XY_R4(DYF,myThid)
180          IF (useCubedSphereExchange) THEN
181    cs! fix overlaps:
182          DO bj = myByLo(myThid), myByHi(myThid)
183           DO bi = myBxLo(myThid), myBxHi(myThid)
184            DO j=1,sNy
185             DO i=1,Olx
186              DXF(1-i,j,bi,bj)=DXF(i,j,bi,bj)
187              DXF(sNx+i,j,bi,bj)=DXF(sNx+1-i,j,bi,bj)
188              DYF(1-i,j,bi,bj)=DYF(i,j,bi,bj)
189              DYF(sNx+i,j,bi,bj)=DYF(sNx+1-i,j,bi,bj)
190             ENDDO
191            ENDDO
192            DO j=1,Oly
193             DO i=1,sNx
194              DXF(i,1-j,bi,bj)=DXF(i,j,bi,bj)
195              DXF(i,sNy+j,bi,bj)=DXF(i,sNy+1-j,bi,bj)
196              DYF(i,1-j,bi,bj)=DYF(i,j,bi,bj)
197              DYF(i,sNy+j,bi,bj)=DYF(i,sNy+1-j,bi,bj)
198             ENDDO
199            ENDDO
200           ENDDO
201          ENDDO
202          ENDIF
203    cs
204    
205        CALL MDSREADFIELD('RA.bin',readBinaryPrec,'RL',1,RA,  1,myThid)        CALL MDSREADFIELD('RA.bin',readBinaryPrec,'RS',1,RA,  1,myThid)
206        _EXCH_XY_R4(RA,myThid )        _EXCH_XY_R4(RA,myThid )
207    
208  C-    Corner quantities  C-    Corner quantities
209  C       *********** this are not degbugged ************  C       *********** this are not degbugged ************
210        CALL MDSREADFIELD('LONG.bin',readBinaryPrec,'RL',1,XG,  1,myThid)        CALL MDSREADFIELD('LONG.bin',readBinaryPrec,'RS',1,XG,  1,myThid)
211        CALL MDSREADFIELD('LATG.bin',readBinaryPrec,'RL',1,YG,  1,myThid)        CALL MDSREADFIELD('LATG.bin',readBinaryPrec,'RS',1,YG,  1,myThid)
212          IF (useCubedSphereExchange) THEN
213  cs-   this block needed by cubed sphere until we write more useful I/O routines  cs-   this block needed by cubed sphere until we write more useful I/O routines
214        bi=3        bi=3
215        bj=1        bj=1
# Line 98  cs-   this block needed by cubed sphere Line 226  cs-   this block needed by cubed sphere
226        bj=bj+2        bj=bj+2
227        YG(sNx+1,1,bj,1)=YG(1,1,bi,1)        YG(sNx+1,1,bj,1)=YG(1,1,bi,1)
228  cs-   end block  cs-   end block
229          ENDIF
230        CALL EXCH_Z_XY_RS(XG,myThid)        CALL EXCH_Z_XY_RS(XG,myThid)
231        CALL EXCH_Z_XY_RS(YG,myThid)        CALL EXCH_Z_XY_RS(YG,myThid)
232    
233        CALL MDSREADFIELD('DXV.bin',readBinaryPrec,'RL',1,DXV,  1,myThid)        CALL MDSREADFIELD('DXV.bin',readBinaryPrec,'RS',1,DXV,  1,myThid)
234        CALL MDSREADFIELD('DYU.bin',readBinaryPrec,'RL',1,DYU,  1,myThid)        CALL MDSREADFIELD('DYU.bin',readBinaryPrec,'RS',1,DYU,  1,myThid)
235  cs-   this block needed by cubed sphere until we write more useful I/O routines  cs-   this block needed by cubed sphere until we write more useful I/O routines
236  C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)  C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)
237  cs!   this is not correct <= need paired exchange for dxv,dyu  cs!   this is not correct <= need paired exchange for dxv,dyu
238          IF (.NOT.useCubedSphereExchange) THEN
239        CALL EXCH_Z_XY_RS(DXV,myThid)        CALL EXCH_Z_XY_RS(DXV,myThid)
240        CALL EXCH_Z_XY_RS(DYU,myThid)        CALL EXCH_Z_XY_RS(DYU,myThid)
241          ELSE
242        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
243         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
244          DXV(sNx+1,1,bi,bj)=DXV(1,1,bi,bj)  cs! fix overlaps:
245          DXV(1,sNy+1,bi,bj)=DXV(1,1,bi,bj)          DO j=1,sNy
246          DYU(sNx+1,1,bi,bj)=DYU(1,1,bi,bj)           DO i=1,Olx
247          DYU(1,sNy+1,bi,bj)=DYU(1,1,bi,bj)            DXV(1-i,j,bi,bj)=DXV(1+i,j,bi,bj)
248              DXV(sNx+i,j,bi,bj)=DXV(i,j,bi,bj)
249              DYU(1-i,j,bi,bj)=DYU(1+i,j,bi,bj)
250              DYU(sNx+i,j,bi,bj)=DYU(i,j,bi,bj)
251             ENDDO
252            ENDDO
253            DO j=1,Oly
254             DO i=1-Olx,sNx+Olx
255              DXV(i,1-j,bi,bj)=DXV(i,1+j,bi,bj)
256              DXV(i,sNy+j,bi,bj)=DXV(i,j,bi,bj)
257              DYU(i,1-j,bi,bj)=DYU(i,1+j,bi,bj)
258              DYU(i,sNy+j,bi,bj)=DYU(i,j,bi,bj)
259             ENDDO
260            ENDDO
261         ENDDO         ENDDO
262        ENDDO        ENDDO
263  cs-   end block  cs-   end block
264  C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)        ENDIF
 cs!   this is not correct <= need paired exchange for dxv,dyu  
       CALL EXCH_Z_XY_RS(DXV,myThid)  
       CALL EXCH_Z_XY_RS(DYU,myThid)  
265    
266        CALL MDSREADFIELD('RAZ.bin',readBinaryPrec,'RL',1,RAZ,  1,myThid)        CALL MDSREADFIELD('RAZ.bin',readBinaryPrec,'RS',1,RAZ,  1,myThid)
267          IF (useCubedSphereExchange) THEN
268  cs-   this block needed by cubed sphere until we write more useful I/O routines  cs-   this block needed by cubed sphere until we write more useful I/O routines
269        CALL EXCH_Z_XY_RS(RAZ , myThid )        CALL EXCH_Z_XY_RS(RAZ , myThid )
270        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
# Line 132  cs-   this block needed by cubed sphere Line 274  cs-   this block needed by cubed sphere
274         ENDDO         ENDDO
275        ENDDO        ENDDO
276  cs-   end block  cs-   end block
277          ENDIF
278        CALL EXCH_Z_XY_RS(RAZ,myThid)        CALL EXCH_Z_XY_RS(RAZ,myThid)
279    
280  C-    Staggered (u,v pairs) quantities  C-    Staggered (u,v pairs) quantities
281        CALL MDSREADFIELD('DXC.bin',readBinaryPrec,'RL',1,DXC,  1,myThid)        CALL MDSREADFIELD('DXC.bin',readBinaryPrec,'RS',1,DXC,  1,myThid)
282        CALL MDSREADFIELD('DYC.bin',readBinaryPrec,'RL',1,DYC,  1,myThid)        CALL MDSREADFIELD('DYC.bin',readBinaryPrec,'RS',1,DYC,  1,myThid)
283        CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)        CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
284    
285        CALL MDSREADFIELD('RAW.bin',readBinaryPrec,'RL',1,RAW,  1,myThid)        CALL MDSREADFIELD('RAW.bin',readBinaryPrec,'RS',1,RAW,  1,myThid)
286        CALL MDSREADFIELD('RAS.bin',readBinaryPrec,'RL',1,RAS,  1,myThid)        CALL MDSREADFIELD('RAS.bin',readBinaryPrec,'RS',1,RAS,  1,myThid)
287          IF (useCubedSphereExchange) THEN
288  cs-   this block needed by cubed sphere until we write more useful I/O routines  cs-   this block needed by cubed sphere until we write more useful I/O routines
289        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
290         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
# Line 151  c        RAS(J,sNy+1,bi,bj)=RAS(J,1,bi,b Line 295  c        RAS(J,sNy+1,bi,bj)=RAS(J,1,bi,b
295         ENDDO         ENDDO
296        ENDDO        ENDDO
297  cs-   end block  cs-   end block
298          ENDIF
299        CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)        CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
300    
301        CALL MDSREADFIELD('DXG.bin',readBinaryPrec,'RL',1,DXG,  1,myThid)        CALL MDSREADFIELD('DXG.bin',readBinaryPrec,'RS',1,DXG,  1,myThid)
302        CALL MDSREADFIELD('DYG.bin',readBinaryPrec,'RL',1,DYG,  1,myThid)        CALL MDSREADFIELD('DYG.bin',readBinaryPrec,'RS',1,DYG,  1,myThid)
303          IF (useCubedSphereExchange) THEN
304  cs-   this block needed by cubed sphere until we write more useful I/O routines  cs-   this block needed by cubed sphere until we write more useful I/O routines
305        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
306         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
# Line 165  c        DXG(J,sNy+1,bi,bj)=DXG(J,1,bi,b Line 311  c        DXG(J,sNy+1,bi,bj)=DXG(J,1,bi,b
311         ENDDO         ENDDO
312        ENDDO        ENDDO
313  cs-   end block  cs-   end block
314          ENDIF
315        CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)        CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
316          CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
317    
318  c     write(10) XC  c     write(10) XC
319  c     write(10) YC  c     write(10) YC
# Line 184  c     write(10) RAS Line 332  c     write(10) RAS
332  c     write(10) DXG  c     write(10) DXG
333  c     write(10) DYG  c     write(10) DYG
334    
335    #else /* ifndef OLD_GRID_IO */
336    
337    C--   Only do I/O if I am the master thread
338          _BEGIN_MASTER(myThid)
339    
340          DO bj = 1,nSy
341           DO bi = 1,nSx
342            iG=bi+(myXGlobalLo-1)/sNx
343            WRITE(tmpBuf,'(A,I4)') 'tile:',iG
344    #ifdef ALLOW_EXCH2
345            myTile = W2_myTileList(bi)
346            WRITE(tmpBuf,'(A,I4)') 'tile:',myTile
347            iG = exch2_myface(myTile)
348    #endif
349            iLen = ILNBLNK(horizGridFile)
350            IF ( iLen .EQ. 0 ) THEN
351              WRITE(fName,'("tile",I3.3,".mitgrid")') iG
352            ELSE
353              WRITE(fName,'(2A,I3.3,A)') horizGridFile(1:iLen),
354         &                              '.face',iG,'.bin'
355            ENDIF
356            iLen = ILNBLNK(fName)
357            iL = ILNBLNK(tmpBuf)
358            WRITE(msgBuf,'(3A)') tmpBuf(1:iL),
359         &                   ' ; Read from file ',fName(1:iLen)
360            CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
361         &                      SQUEEZE_RIGHT , myThid)
362            WRITE(msgBuf,'(A)') '  =>'
363    
364            CALL READSYMTILE_RS(fName,1,XC,bi,bj,buf,myThid)
365            iL = ILNBLNK(msgBuf)
366            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'XC'
367            CALL READSYMTILE_RS(fName,2,YC,bi,bj,buf,myThid)
368            iL = ILNBLNK(tmpBuf)
369            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'YC'
370            CALL READSYMTILE_RS(fName,3,DXF,bi,bj,buf,myThid)
371            iL = ILNBLNK(msgBuf)
372            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXF'
373            CALL READSYMTILE_RS(fName,4,DYF,bi,bj,buf,myThid)
374            iL = ILNBLNK(tmpBuf)
375            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYF'
376            CALL READSYMTILE_RS(fName,5,RA,bi,bj,buf,myThid)
377            iL = ILNBLNK(msgBuf)
378            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'RA'
379            CALL READSYMTILE_RS(fName,6,XG,bi,bj,buf,myThid)
380            iL = ILNBLNK(tmpBuf)
381            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'XG'
382            CALL READSYMTILE_RS(fName,7,YG,bi,bj,buf,myThid)
383            iL = ILNBLNK(msgBuf)
384            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'YG'
385            CALL READSYMTILE_RS(fName,8,DXV,bi,bj,buf,myThid)
386            iL = ILNBLNK(tmpBuf)
387            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DXV'
388            CALL READSYMTILE_RS(fName,9,DYU,bi,bj,buf,myThid)
389            iL = ILNBLNK(msgBuf)
390            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DYU'
391            CALL READSYMTILE_RS(fName,10,RAZ,bi,bj,buf,myThid)
392            iL = ILNBLNK(tmpBuf)
393            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'RAZ'
394            CALL READSYMTILE_RS(fName,11,DXC,bi,bj,buf,myThid)
395            iL = ILNBLNK(msgBuf)
396            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXC'
397            CALL READSYMTILE_RS(fName,12,DYC,bi,bj,buf,myThid)
398            iL = ILNBLNK(tmpBuf)
399            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYC'
400            CALL READSYMTILE_RS(fName,13,RAW,bi,bj,buf,myThid)
401            iL = ILNBLNK(msgBuf)
402            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'RAW'
403            CALL READSYMTILE_RS(fName,14,RAS,bi,bj,buf,myThid)
404            iL = ILNBLNK(tmpBuf)
405            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'RAS'
406            CALL READSYMTILE_RS(fName,15,DXG,bi,bj,buf,myThid)
407            iL = ILNBLNK(msgBuf)
408            WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'DXG'
409            CALL READSYMTILE_RS(fName,16,DYG,bi,bj,buf,myThid)
410            iL = ILNBLNK(tmpBuf)
411            WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'DYG'
412    
413            iLen = ILNBLNK(horizGridFile)
414            IF ( iLen.GT.0 ) THEN
415             CALL READSYMTILE_RS(fName,17,angleCosC,bi,bj,buf,myThid)
416             iL = ILNBLNK(msgBuf)
417             WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'AngleCS'
418             CALL READSYMTILE_RS(fName,18,angleSinC,bi,bj,buf,myThid)
419             iL = ILNBLNK(tmpBuf)
420             WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'AngleSN'
421            ENDIF
422    
423            CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
424         &                      SQUEEZE_RIGHT , myThid)
425    
426           ENDDO
427          ENDDO
428    
429          _END_MASTER(myThid)
430    
431          CALL EXCH_XY_RS(XC,myThid)
432          CALL EXCH_XY_RS(YC,myThid)
433    C !!! _EXCH_OUV_XY_R4(DXF, DYF, unSigned, myThid )
434    #ifdef HRCUBE
435          CALL EXCH_XY_RS(DXF,myThid)
436          CALL EXCH_XY_RS(DYF,myThid)
437    #endif
438          CALL EXCH_XY_RS(RA,myThid )
439          CALL EXCH_Z_XY_RS(XG,myThid)
440          CALL EXCH_Z_XY_RS(YG,myThid)
441    C !!! _EXCH_ZUV_XY_R4(DXV, DYU, unSigned, myThid)
442    c     CALL EXCH_Z_XY_RS(DXV,myThid)
443    c     CALL EXCH_Z_XY_RS(DYU,myThid)
444          CALL EXCH_Z_XY_RS(RAZ,myThid)
445          CALL EXCH_UV_XY_RS(DXC,DYC,.FALSE.,myThid)
446          CALL EXCH_UV_XY_RS(RAW,RAS,.FALSE.,myThid)
447          CALL EXCH_UV_XY_RS(DYG,DXG,.FALSE.,myThid)
448          CALL EXCH_UV_AGRID_XY_RS(angleSinC,angleCosC,.TRUE.,myThid)
449    
450    #endif /* OLD_GRID_IO */
451    
452    #ifdef ALLOW_MNC
453          ENDIF
454    #endif /* ALLOW_MNC */
455    
456    c     CALL WRITE_FULLARRAY_RL('DXV',DXV,1,0,0,0,myThid)
457    c     CALL WRITE_FULLARRAY_RL('DYU',DYU,1,0,0,0,myThid)
458    c     CALL WRITE_FULLARRAY_RL('RAZ',RAZ,1,0,0,0,myThid)
459    c     CALL WRITE_FULLARRAY_RL('XG',XG,1,0,0,0,myThid)
460    c     CALL WRITE_FULLARRAY_RL('YG',YG,1,0,0,0,myThid)
461    
462    C--   Require that 0 <= longitude < 360 if using exf package
463    #ifdef ALLOW_EXF
464          DO bj = 1,nSy
465           DO bi = 1,nSx
466            DO J=1-Oly,sNy+Oly
467             DO I=1-Olx,sNx+Olx
468              IF (XC(i,j,bi,bj).lt.0.) XC(i,j,bi,bj)=XC(i,j,bi,bj)+360.
469              IF (XG(i,j,bi,bj).lt.0.) XG(i,j,bi,bj)=XG(i,j,bi,bj)+360.
470             ENDDO
471            ENDDO
472           ENDDO
473          ENDDO
474    #endif /* ALLOW_EXF */
475    
476    C--   Now let's look at all these beasts
477          IF ( debugLevel .GE. debLevB ) THEN
478             myIter = 1
479             CALL PLOT_FIELD_XYRL( XC      , 'Current XC      ' ,
480         &        myIter, myThid )
481             CALL PLOT_FIELD_XYRL( YC      , 'Current YC      ' ,
482         &        myIter, myThid )
483             CALL PLOT_FIELD_XYRL( DXF     , 'Current DXF     ' ,
484         &        myIter, myThid )
485             CALL PLOT_FIELD_XYRL( XC      , 'Current XC      ' ,
486         &        myIter, myThid )
487             CALL PLOT_FIELD_XYRL( DYF     , 'Current DYF     ' ,
488         &        myIter, myThid )
489             CALL PLOT_FIELD_XYRL( RA      , 'Current RA      ' ,
490         &        myIter, myThid )
491             CALL PLOT_FIELD_XYRL( XG      , 'Current XG      ' ,
492         &        myIter, myThid )
493             CALL PLOT_FIELD_XYRL( YG      , 'Current YG      ' ,
494         &        myIter, myThid )
495             CALL PLOT_FIELD_XYRL( DXV     , 'Current DXV     ' ,
496         &        myIter, myThid )
497             CALL PLOT_FIELD_XYRL( DYU     , 'Current DYU     ' ,
498         &        myIter, myThid )
499             CALL PLOT_FIELD_XYRL( RAZ     , 'Current RAZ     ' ,
500         &        myIter, myThid )
501             CALL PLOT_FIELD_XYRL( DXC     , 'Current DXC     ' ,
502         &        myIter, myThid )
503             CALL PLOT_FIELD_XYRL( DYC     , 'Current DYC     ' ,
504         &        myIter, myThid )
505             CALL PLOT_FIELD_XYRL( RAW     , 'Current RAW     ' ,
506         &        myIter, myThid )
507             CALL PLOT_FIELD_XYRL( RAS     , 'Current RAS     ' ,
508         &        myIter, myThid )
509             CALL PLOT_FIELD_XYRL( DXG     , 'Current DXG     ' ,
510         &        myIter, myThid )
511             CALL PLOT_FIELD_XYRL( DYG     , 'Current DYG     ' ,
512         &        myIter, myThid )
513             CALL PLOT_FIELD_XYRL(angleCosC, 'Current AngleCS ' ,
514         &        myIter, myThid )
515             CALL PLOT_FIELD_XYRL(angleSinC, 'Current AngleSN ' ,
516         &        myIter, myThid )
517          ENDIF
518    
519          RETURN
520          END
521    
522    C --------------------------------------------------------------------------
523    
524          SUBROUTINE READSYMTILE_RS(fName,irec,array,bi,bj,buf,myThid)
525    C     /==========================================================\
526    C     | SUBROUTINE READSYMTILE_RS                                |
527    C     |==========================================================|
528    C     \==========================================================/
529          IMPLICIT NONE
530    
531    C     === Global variables ===
532    #include "SIZE.h"
533    #include "EEPARAMS.h"
534    #ifdef ALLOW_EXCH2
535    #include "W2_EXCH2_TOPOLOGY.h"
536    #include "W2_EXCH2_PARAMS.h"
537    #endif /* ALLOW_EXCH2 */
538    
539    C     == Routine arguments ==
540          CHARACTER*(*) fName
541          INTEGER irec
542          _RS array(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
543          INTEGER bi,bj,myThid
544    #ifdef ALLOW_EXCH2
545          _RL buf(1:sNx*nSx*nPx+1)
546    #else
547          _RL buf(1:sNx+1,1:sNy+1)
548    #endif /* ALLOW_EXCH2 */
549    
550    C     == Local variables ==
551          INTEGER I,J,dUnit, iLen
552          INTEGER length_of_rec
553          INTEGER MDS_RECLEN
554    #ifdef ALLOW_EXCH2
555          INTEGER TN, dNx, dNy, TBX, TBY, TNX, TNY, II, iBase
556    #endif
557          INTEGER  ILNBLNK
558          EXTERNAL ILNBLNK
559    
560          iLen = ILNBLNK(fName)
561    #ifdef ALLOW_EXCH2
562    C     Figure out offset of tile within face
563          TN  = W2_myTileList(bi)
564          dNx = exch2_mydnx(TN)
565          dNy = exch2_mydny(TN)
566          TBX = exch2_tbasex(TN)
567          TBY = exch2_tbasey(TN)
568          TNX = exch2_tnx(TN)
569          TNY = exch2_tny(TN)
570    
571          CALL MDSFINDUNIT( dUnit, myThid )
572          length_of_rec=MDS_RECLEN( 64, (dNx+1), myThid )
573          OPEN( dUnit, file=fName(1:iLen), status='old',
574         &             access='direct', recl=length_of_rec )
575          J=0
576          iBase=(irec-1)*(dny+1)
577          DO I=1+TBY,sNy+1+TBY
578           READ(dUnit,rec=I+iBase)(buf(ii),ii=1,dNx+1)
579    #ifdef _BYTESWAPIO
580    #ifdef REAL4_IS_SLOW
581           CALL MDS_BYTESWAPR8((dNx+1), buf)
582    #else
583           CALL MDS_BYTESWAPR4((dNx+1), buf)
584    #endif
585    #endif
586           J=J+1
587           DO II=1,sNx+1
588            array(II,J,bi,bj)=buf(II+TBX)
589           ENDDO
590          ENDDO
591          CLOSE( dUnit )
592          
593    #else /* ALLOW_EXCH2 */
594    
595          CALL MDSFINDUNIT( dUnit, myThid )
596          length_of_rec=MDS_RECLEN( 64, (sNx+1)*(sNy+1), myThid )
597          OPEN( dUnit, file=fName(1:iLen), status='old',
598         &             access='direct', recl=length_of_rec )
599          READ(dUnit,rec=irec) buf
600          CLOSE( dUnit )
601    
602    #ifdef _BYTESWAPIO
603    #ifdef REAL4_IS_SLOW
604          CALL MDS_BYTESWAPR8((sNx+1)*(sNy+1), buf)
605    #else
606          CALL MDS_BYTESWAPR4((sNx+1)*(sNy+1), buf)
607    #endif
608    #endif
609    
610          DO J=1,sNy+1
611           DO I=1,sNx+1
612            array(I,J,bi,bj)=buf(I,J)
613           ENDDO
614          ENDDO
615    c       write(0,*) irec,buf(1,1),array(1,1,1,1)
616    
617    #endif /* ALLOW_EXCH2 */
618    
619        RETURN        RETURN
620        END        END
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.23
  ViewVC Help
Powered by ViewVC 1.1.22