model/src/ini_curvilinear_grid.F

aa335bcda3 Ed H*0001 #include "PACKAGES_CONFIG.h"
aea29c8517 Alis*0002 #include "CPP_OPTIONS.h"
                0003
9366854e02 Chri*0004 CBOP
                0005 C     !ROUTINE: INI_CURVILINEAR_GRID
                0006 C     !INTERFACE:
aea29c8517 Alis*0007       SUBROUTINE INI_CURVILINEAR_GRID( myThid )
9744e36521 Jean*0008
9366854e02 Chri*0009 C     !DESCRIPTION: \bv
                0010 C     *==========================================================*
a7f6b46be7 Jean*0011 C     | SUBROUTINE INI_CURVILINEAR_GRID
                0012 C     | o Initialise curvilinear coordinate system
9366854e02 Chri*0013 C     *==========================================================*
                0014 C     | Curvilinear grid settings are read from a file rather
                0015 C     | than coded in-line as for cartesian and spherical polar.
                0016 C     | This is more general but you have to create the grid
                0017 C     | yourself.
                0018 C     *==========================================================*
                0019 C     \ev
aea29c8517 Alis*0020
9366854e02 Chri*0021 C     !USES:
                0022       IMPLICIT NONE
aea29c8517 Alis*0023 C     === Global variables ===
                0024 #include "SIZE.h"
                0025 #include "EEPARAMS.h"
                0026 #include "PARAMS.h"
                0027 #include "GRID.h"
4749c74143 Alis*0028 #ifdef ALLOW_EXCH2
1cc6effca6 Jean*0029 #include "W2_EXCH2_SIZE.h"
bdac319189 Jean*0030 #include "W2_EXCH2_TOPOLOGY.h"
                0031 #endif
e808c5b132 Ed H*0032 #ifdef ALLOW_MNC
                0033 #include "MNC_PARAMS.h"
                0034 #endif
bdac319189 Jean*0035
9366854e02 Chri*0036 C     !INPUT/OUTPUT PARAMETERS:
aea29c8517 Alis*0037 C     == Routine arguments ==
530d89e069 Jean*0038 C     myThid -  Number of this instance of INI_CURVILINEAR_GRID
aea29c8517 Alis*0039       INTEGER myThid
                0040
9366854e02 Chri*0041 C     !LOCAL VARIABLES:
54a8fe0a22 Jean*0042 C     == Shared Local variables ==
                0043       LOGICAL anglesAreSet
                0044       COMMON /LOCAL_INI_CURVILINEAR_GRID/ anglesAreSet
aea29c8517 Alis*0045 C     == Local variables ==
ae73826a86 Jean*0046       INTEGER bi,bj
7d173b8d56 Jean*0047       INTEGER i,j
8996cf5a3c Jean*0048       CHARACTER*(MAX_LEN_MBUF) msgBuf
30d424dc37 Jean*0049       INTEGER fp
9744e36521 Jean*0050       _RL tmpFac, tmpFac2
e808c5b132 Ed H*0051 #ifdef ALLOW_MNC
de57a2ec4b Mart*0052       CHARACTER*(MAX_LEN_FNAM) mncFn
e808c5b132 Ed H*0053 #endif
8996cf5a3c Jean*0054 #ifndef OLD_GRID_IO
623695a379 Jean*0055       INTEGER iG, jG, iL, iLen
8996cf5a3c Jean*0056       CHARACTER*(MAX_LEN_FNAM) fName
                0057       CHARACTER*(MAX_LEN_MBUF) tmpBuf
56cf4bd603 Jean*0058       INTEGER  ILNBLNK
                0059       EXTERNAL ILNBLNK
8996cf5a3c Jean*0060 #endif
b9dadda204 Mart*0061 #ifdef ALLOW_EXCH2
                0062       CHARACTER*(6) fmtStr
                0063       INTEGER iTmp
                0064 #endif
9366854e02 Chri*0065 CEOP
aea29c8517 Alis*0066
                0067 C--   Set everything to zero everywhere
be9aaf9124 Jean*0068 C     Note: this is now done earlier in main S/R INI_GRID
aea29c8517 Alis*0069
                0070 C     Here we make no assumptions about grid symmetry and simply
                0071 C     read the raw grid data from files
                0072
bdac319189 Jean*0073 #ifdef OLD_GRID_IO
30d424dc37 Jean*0074 C--   File Precision is different from "new grid IO" (always 64-bits precision)
                0075 C     which should probably be changed to the standard file-prec (= readBinaryPrec)
                0076       fp = readBinaryPrec
bdac319189 Jean*0077
4d2c31cbda Jean*0078 # ifdef ALLOW_MDSIO
aea29c8517 Alis*0079 C-    Cell centered quantities
30d424dc37 Jean*0080       CALL READ_REC_3D_RS( 'LONC.bin', fp, 1, xC, 1, 0, myThid )
                0081       CALL READ_REC_3D_RS( 'LATC.bin', fp, 1, yC, 1, 0, myThid )
12c8b75709 Jean*0082       _EXCH_XY_RS(xC,myThid)
                0083       _EXCH_XY_RS(yC,myThid)
aea29c8517 Alis*0084
30d424dc37 Jean*0085       CALL READ_REC_3D_RS( 'DXF.bin', fp, 1, dxF, 1, 0, myThid )
                0086       CALL READ_REC_3D_RS( 'DYF.bin', fp, 1, dyF, 1, 0, myThid )
d8a7c4feca Jean*0087       CALL EXCH_UV_AGRID_3D_RS( dxF, dyF, .FALSE., 1, myThid )
a7f6b46be7 Jean*0088
30d424dc37 Jean*0089       CALL READ_REC_3D_RS( 'RA.bin' , fp, 1, rA,  1, 0, myThid )
12c8b75709 Jean*0090       _EXCH_XY_RS(rA,myThid )
aea29c8517 Alis*0091
54a8fe0a22 Jean*0092       _BEGIN_MASTER(myThid)
                0093       anglesAreSet = .FALSE.
                0094       _END_MASTER(myThid)
                0095
aea29c8517 Alis*0096 C-    Corner quantities
30d424dc37 Jean*0097       CALL READ_REC_3D_RS( 'LONG.bin', fp, 1, xG, 1, 0, myThid )
                0098       CALL READ_REC_3D_RS( 'LATG.bin', fp, 1, yG, 1, 0, myThid )
f424074f4f Dimi*0099       IF (useCubedSphereExchange) THEN
aea29c8517 Alis*0100 cs-   this block needed by cubed sphere until we write more useful I/O routines
30d424dc37 Jean*0101        IF ( nPx*nPy*nSy.EQ.1 .AND. nSx.EQ.6 ) THEN
                0102         _BARRIER
                0103         _BEGIN_MASTER(myThid)
                0104          bi=3
                0105          bj=1
                0106          yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
                0107          bj=bj+2
                0108          yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
                0109          bj=bj+2
                0110          yG(1,sNy+1,bj,1)=yG(1,1,bi,1)
                0111          bi=6
                0112          bj=2
                0113          yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
                0114          bj=bj+2
                0115          yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
                0116          bj=bj+2
                0117          yG(sNx+1,1,bj,1)=yG(1,1,bi,1)
                0118         _END_MASTER(myThid)
                0119         _BARRIER
                0120        ELSE
                0121         WRITE(msgBuf,'(2A)') 'INI_CURVILINEAR_GRID:',
                0122      &        ' OLD_GRID_IO only works for 6 tiles on 1 proc'
                0123         CALL PRINT_ERROR( msgBuf, myThid )
                0124         STOP 'ABNORMAL END: S/R INI_CURVILINEAR_GRID'
                0125        ENDIF
aea29c8517 Alis*0126 cs-   end block
f424074f4f Dimi*0127       ENDIF
d8a7c4feca Jean*0128       CALL EXCH_Z_3D_RS( xG, 1, myThid )
                0129       CALL EXCH_Z_3D_RS( yG, 1, myThid )
aea29c8517 Alis*0130
30d424dc37 Jean*0131       CALL READ_REC_3D_RS( 'DXV.bin', fp, 1, dxV, 1, 0, myThid )
                0132       CALL READ_REC_3D_RS( 'DYU.bin', fp, 1, dyU, 1, 0, myThid )
aea29c8517 Alis*0133 cs-   this block needed by cubed sphere until we write more useful I/O routines
30d424dc37 Jean*0134       IF ( useCubedSphereExchange ) THEN
                0135        IF ( nPx*nPy*nSx*nSy.EQ.6 .AND. sNx.EQ.sNy ) THEN
                0136         DO bj = myByLo(myThid), myByHi(myThid)
                0137          DO bi = myBxLo(myThid), myBxHi(myThid)
                0138           dxV(sNx+1,1,bi,bj)=dxV(1,1,bi,bj)
                0139           dxV(1,sNy+1,bi,bj)=dxV(1,1,bi,bj)
                0140           dyU(sNx+1,1,bi,bj)=dyU(1,1,bi,bj)
                0141           dyU(1,sNy+1,bi,bj)=dyU(1,1,bi,bj)
f9915c1d33 Alis*0142          ENDDO
                0143         ENDDO
30d424dc37 Jean*0144        ELSE
                0145         WRITE(msgBuf,'(2A)') 'INI_CURVILINEAR_GRID:',
                0146      &        ' OLD_GRID_IO only works with 1 tile per face'
                0147         CALL PRINT_ERROR( msgBuf, myThid )
                0148         STOP 'ABNORMAL END: S/R INI_CURVILINEAR_GRID'
                0149        ENDIF
aea29c8517 Alis*0150 cs-   end block
f424074f4f Dimi*0151       ENDIF
30d424dc37 Jean*0152       CALL EXCH_UV_BGRID_3D_RS( dxV, dyU, .FALSE., 1, myThid )
aea29c8517 Alis*0153
30d424dc37 Jean*0154       CALL READ_REC_3D_RS( 'RAZ.bin', fp, 1, rAz, 1, 0, myThid )
f424074f4f Dimi*0155       IF (useCubedSphereExchange) THEN
aea29c8517 Alis*0156 cs-   this block needed by cubed sphere until we write more useful I/O routines
30d424dc37 Jean*0157        DO bj = myByLo(myThid), myByHi(myThid)
                0158         DO bi = myBxLo(myThid), myBxHi(myThid)
                0159          rAz(sNx+1,1,bi,bj)=rAz(1,1,bi,bj)
                0160          rAz(1,sNy+1,bi,bj)=rAz(1,1,bi,bj)
                0161         ENDDO
aea29c8517 Alis*0162        ENDDO
                0163 cs-   end block
f424074f4f Dimi*0164       ENDIF
d8a7c4feca Jean*0165       CALL EXCH_Z_3D_RS( rAz, 1, myThid )
aea29c8517 Alis*0166
                0167 C-    Staggered (u,v pairs) quantities
30d424dc37 Jean*0168       CALL READ_REC_3D_RS( 'DXC.bin', fp, 1, dxC, 1, 0, myThid )
                0169       CALL READ_REC_3D_RS( 'DYC.bin', fp, 1, dyC, 1, 0, myThid )
a7f6b46be7 Jean*0170       CALL EXCH_UV_XY_RS(dxC,dyC,.FALSE.,myThid)
aea29c8517 Alis*0171
30d424dc37 Jean*0172       CALL READ_REC_3D_RS( 'RAW.bin', fp, 1, rAw, 1, 0, myThid )
                0173       CALL READ_REC_3D_RS( 'RAS.bin', fp, 1, rAs, 1, 0, myThid )
a7f6b46be7 Jean*0174       CALL EXCH_UV_XY_RS(rAw,rAs,.FALSE.,myThid)
aea29c8517 Alis*0175
30d424dc37 Jean*0176       CALL READ_REC_3D_RS( 'DXG.bin', fp, 1, dxG, 1, 0, myThid )
                0177       CALL READ_REC_3D_RS( 'DYG.bin', fp, 1, dyG, 1, 0, myThid )
a7f6b46be7 Jean*0178       CALL EXCH_UV_XY_RS(dyG,dxG,.FALSE.,myThid)
4d2c31cbda Jean*0179 # else /* ALLOW_MDSIO */
                0180       WRITE(msgBuf,'(2A)')
                0181      &  'INI_CURVILINEAR_GRID: In order to use OLD_GRID_IO code,'
                0182       CALL PRINT_ERROR( msgBuf, myThid )
                0183       WRITE(msgBuf,'(2A)')
                0184      &  'INI_CURVILINEAR_GRID: needs to compile MDSIO pkg'
                0185       CALL PRINT_ERROR( msgBuf, myThid )
                0186       STOP 'ABNORMAL END: S/R INI_CURVILINEAR_GRID'
                0187 # endif /* ALLOW_MDSIO */
aea29c8517 Alis*0188
a7f6b46be7 Jean*0189 c     write(10) xC
                0190 c     write(10) yC
                0191 c     write(10) dxF
                0192 c     write(10) dyF
                0193 c     write(10) rA
                0194 c     write(10) xG
                0195 c     write(10) yG
                0196 c     write(10) dxV
                0197 c     write(10) dyU
                0198 c     write(10) rAz
                0199 c     write(10) dxC
                0200 c     write(10) dyC
                0201 c     write(10) rAw
                0202 c     write(10) rAs
                0203 c     write(10) dxG
                0204 c     write(10) dyG
aea29c8517 Alis*0205
f2ad5f0777 Dimi*0206 #else /* ifndef OLD_GRID_IO */
bdac319189 Jean*0207
285a3c0cbd Jean*0208 #ifdef ALLOW_MNC
                0209       IF (useMNC .AND. readgrid_mnc) THEN
                0210 C--   read NetCDF files:
                0211
                0212         DO i = 1,80
                0213           mncFn(i:i) = ' '
                0214         ENDDO
                0215         write(mncFn,'(a)') 'mitgrid'
                0216         DO i = 1,MAX_LEN_MBUF
                0217           msgBuf(i:i) = ' '
                0218         ENDDO
                0219         WRITE(msgBuf,'(2A)') msgBuf,' ; Reading grid info using MNC'
                0220         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0221      &       SQUEEZE_RIGHT , myThid)
                0222         CALL MNC_FILE_CLOSE_ALL_MATCHING(mncFn, myThid)
                0223         CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
                0224         CALL MNC_CW_SET_CITER(mncFn, 2, -1, -1, -1, myThid)
                0225         CALL MNC_CW_SET_UDIM(mncFn, 1, myThid)
                0226         CALL MNC_CW_RS_R('D',mncFn,0,0,'XC', xC,  myThid)
                0227         CALL MNC_CW_RS_R('D',mncFn,0,0,'XG', xG,  myThid)
                0228         CALL MNC_CW_RS_R('D',mncFn,0,0,'YC', yC,  myThid)
                0229         CALL MNC_CW_RS_R('D',mncFn,0,0,'YG', yG,  myThid)
                0230         CALL MNC_CW_RS_R('D',mncFn,0,0,'dxC',dxC, myThid)
                0231         CALL MNC_CW_RS_R('D',mncFn,0,0,'dyC',dyC, myThid)
                0232         CALL MNC_CW_RS_R('D',mncFn,0,0,'dxF',dxF, myThid)
                0233         CALL MNC_CW_RS_R('D',mncFn,0,0,'dyF',dyF, myThid)
                0234         CALL MNC_CW_RS_R('D',mncFn,0,0,'dxG',dxG, myThid)
                0235         CALL MNC_CW_RS_R('D',mncFn,0,0,'dyG',dyG, myThid)
                0236         CALL MNC_CW_RS_R('D',mncFn,0,0,'dxV',dxV, myThid)
                0237         CALL MNC_CW_RS_R('D',mncFn,0,0,'dyU',dyU, myThid)
                0238         CALL MNC_CW_RS_R('D',mncFn,0,0,'rA', rA,  myThid)
                0239         CALL MNC_CW_RS_R('D',mncFn,0,0,'rAz',rAz, myThid)
                0240         CALL MNC_CW_RS_R('D',mncFn,0,0,'rAw',rAw, myThid)
                0241         CALL MNC_CW_RS_R('D',mncFn,0,0,'rAs',rAs, myThid)
                0242         CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleCS',angleCosC,myThid)
                0243         CALL MNC_CW_RS_R('D',mncFn,0,0,'AngleSN',angleSinC,myThid)
                0244         anglesAreSet = .TRUE.
                0245
                0246       ELSE
                0247 C--   read Binary files:
                0248 #endif /* ALLOW_MNC */
                0249
7d173b8d56 Jean*0250 C--   File Precision: keep 64-bits precision (as it used to be)
                0251 C     but should probably change it to the standard file-prec (= readBinaryPrec)
                0252       fp = precFloat64
                0253 c     fp = readBinaryPrec
                0254
30d424dc37 Jean*0255 C--   Everyone must wait for the initialisation to be done
                0256       _BARRIER
                0257
                0258 C--   Only do I/O if I am the master thread
                0259       _BEGIN_MASTER(myThid)
                0260
56cf4bd603 Jean*0261       DO bj = 1,nSy
                0262        DO bi = 1,nSx
4749c74143 Alis*0263 #ifdef ALLOW_EXCH2
623695a379 Jean*0264 C-    Use face number:
8adbfea2f8 Jean*0265         jG = W2_myTileList(bi,bj)
623695a379 Jean*0266         iG = exch2_myface(jG)
b9dadda204 Mart*0267         iTmp = MAX(4,1 + INT(LOG10(DFLOAT(W2_maxNbTiles))))
                0268         WRITE(fmtStr,'(A,I1,A)') '(A,I',iTmp,')'
                0269         WRITE(tmpBuf,fmtStr) 'tile:',jG
623695a379 Jean*0270 #else
                0271 C-    Tile Id number = Bi + (Bj-1)*(nSx*nPx)  with tile global-indices Bi,Bj
                0272         iG = bi+(myXGlobalLo-1)/sNx
                0273         jG = bj+(myYGlobalLo-1)/sNy
                0274         WRITE(tmpBuf,'(2(A,I3))') 'tile:',iG,' ,',jG
                0275         iG = iG + (jG-1)*(nSx*nPx)
bdac319189 Jean*0276 #endif
7d173b8d56 Jean*0277
f23e97a1aa Jean*0278         iLen = ILNBLNK(horizGridFile)
                0279         IF ( iLen .EQ. 0 ) THEN
                0280           WRITE(fName,'("tile",I3.3,".mitgrid")') iG
                0281         ELSE
                0282           WRITE(fName,'(2A,I3.3,A)') horizGridFile(1:iLen),
                0283      &                              '.face',iG,'.bin'
                0284         ENDIF
                0285         iLen = ILNBLNK(fName)
                0286         iL = ILNBLNK(tmpBuf)
                0287         WRITE(msgBuf,'(3A)') tmpBuf(1:iL),
                0288      &                   ' ; Read from file ',fName(1:iLen)
56cf4bd603 Jean*0289         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0290      &                      SQUEEZE_RIGHT , myThid)
                0291         WRITE(msgBuf,'(A)') '  =>'
                0292
7d173b8d56 Jean*0293 #ifdef ALLOW_MDSIO
                0294         CALL MDS_FACEF_READ_RS( fName, fp, 1,  xC, bi, bj, myThid )
56cf4bd603 Jean*0295         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0296         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'xC'
7d173b8d56 Jean*0297         CALL MDS_FACEF_READ_RS( fName, fp, 2,  yC, bi, bj, myThid )
f23e97a1aa Jean*0298         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0299         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'yC'
7d173b8d56 Jean*0300         CALL MDS_FACEF_READ_RS( fName, fp, 3, dxF, bi, bj, myThid )
56cf4bd603 Jean*0301         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0302         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxF'
7d173b8d56 Jean*0303         CALL MDS_FACEF_READ_RS( fName, fp, 4, dyF, bi, bj, myThid )
f23e97a1aa Jean*0304         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0305         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyF'
7d173b8d56 Jean*0306         CALL MDS_FACEF_READ_RS( fName, fp, 5,  rA, bi, bj, myThid )
56cf4bd603 Jean*0307         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0308         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'rA'
7d173b8d56 Jean*0309         CALL MDS_FACEF_READ_RS( fName, fp, 6,  xG, bi, bj, myThid )
f23e97a1aa Jean*0310         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0311         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'xG'
7d173b8d56 Jean*0312         CALL MDS_FACEF_READ_RS( fName, fp, 7,  yG, bi, bj, myThid )
56cf4bd603 Jean*0313         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0314         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'yG'
7d173b8d56 Jean*0315         CALL MDS_FACEF_READ_RS( fName, fp, 8, dxV, bi, bj, myThid )
f23e97a1aa Jean*0316         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0317         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dxV'
7d173b8d56 Jean*0318         CALL MDS_FACEF_READ_RS( fName, fp, 9, dyU, bi, bj, myThid )
56cf4bd603 Jean*0319         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0320         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dyU'
7d173b8d56 Jean*0321         CALL MDS_FACEF_READ_RS( fName, fp,10, rAz, bi, bj, myThid )
f23e97a1aa Jean*0322         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0323         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'rAz'
7d173b8d56 Jean*0324         CALL MDS_FACEF_READ_RS( fName, fp,11, dxC, bi, bj, myThid )
56cf4bd603 Jean*0325         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0326         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxC'
7d173b8d56 Jean*0327         CALL MDS_FACEF_READ_RS( fName, fp,12, dyC, bi, bj, myThid )
f23e97a1aa Jean*0328         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0329         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyC'
7d173b8d56 Jean*0330         CALL MDS_FACEF_READ_RS( fName, fp,13, rAw, bi, bj, myThid )
56cf4bd603 Jean*0331         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0332         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'rAw'
7d173b8d56 Jean*0333         CALL MDS_FACEF_READ_RS( fName, fp,14, rAs, bi, bj, myThid )
f23e97a1aa Jean*0334         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0335         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'rAs'
7d173b8d56 Jean*0336         CALL MDS_FACEF_READ_RS( fName, fp,15, dxG, bi, bj, myThid )
56cf4bd603 Jean*0337         iL = ILNBLNK(msgBuf)
a7f6b46be7 Jean*0338         WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'dxG'
7d173b8d56 Jean*0339         CALL MDS_FACEF_READ_RS( fName, fp,16, dyG, bi, bj, myThid )
f23e97a1aa Jean*0340         iL = ILNBLNK(tmpBuf)
a7f6b46be7 Jean*0341         WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'dyG'
f23e97a1aa Jean*0342
                0343         iLen = ILNBLNK(horizGridFile)
                0344         IF ( iLen.GT.0 ) THEN
7d173b8d56 Jean*0345          CALL MDS_FACEF_READ_RS(fName,fp,17,angleCosC,bi,bj,myThid)
f23e97a1aa Jean*0346          iL = ILNBLNK(msgBuf)
                0347          WRITE(tmpBuf,'(A,1X,A)') msgBuf(1:iL),'AngleCS'
7d173b8d56 Jean*0348          CALL MDS_FACEF_READ_RS(fName,fp,18,angleSinC,bi,bj,myThid)
f23e97a1aa Jean*0349          iL = ILNBLNK(tmpBuf)
                0350          WRITE(msgBuf,'(A,1X,A)') tmpBuf(1:iL),'AngleSN'
b3e3d8ed80 Jean*0351          anglesAreSet = .TRUE.
                0352         ELSE
                0353          anglesAreSet = .FALSE.
f23e97a1aa Jean*0354         ENDIF
7d173b8d56 Jean*0355 #else /* ALLOW_MDSIO */
                0356         WRITE(msgBuf,'(2A)')
                0357      &   'INI_CURVILINEAR_GRID: Needs to compile MDSIO pkg'
                0358         CALL PRINT_ERROR( msgBuf, myThid )
                0359         STOP 'ABNORMAL END: S/R INI_CURVILINEAR_GRID'
                0360 #endif /* ALLOW_MDSIO */
56cf4bd603 Jean*0361
                0362         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0363      &                      SQUEEZE_RIGHT , myThid)
bdac319189 Jean*0364
                0365        ENDDO
                0366       ENDDO
e808c5b132 Ed H*0367
30d424dc37 Jean*0368       _END_MASTER(myThid)
                0369
285a3c0cbd Jean*0370 #ifdef ALLOW_MNC
                0371       ENDIF
                0372 #endif /* ALLOW_MNC */
                0373
a7f6b46be7 Jean*0374       CALL EXCH_XY_RS(xC,myThid)
                0375       CALL EXCH_XY_RS(yC,myThid)
d8a7c4feca Jean*0376       CALL EXCH_UV_AGRID_3D_RS( dxF, dyF, .FALSE., 1, myThid )
a7f6b46be7 Jean*0377       CALL EXCH_XY_RS(rA,myThid )
d8a7c4feca Jean*0378       CALL EXCH_Z_3D_RS( xG, 1, myThid )
                0379       CALL EXCH_Z_3D_RS( yG, 1, myThid )
d04bef6537 Jean*0380       CALL EXCH_UV_BGRID_3D_RS( dxV, dyU, .FALSE., 1, myThid)
d8a7c4feca Jean*0381       CALL EXCH_Z_3D_RS( rAz, 1, myThid )
a7f6b46be7 Jean*0382       CALL EXCH_UV_XY_RS(dxC,dyC,.FALSE.,myThid)
                0383       CALL EXCH_UV_XY_RS(rAw,rAs,.FALSE.,myThid)
                0384       CALL EXCH_UV_XY_RS(dyG,dxG,.FALSE.,myThid)
bdac319189 Jean*0385
f2ad5f0777 Dimi*0386 #endif /* OLD_GRID_IO */
bdac319189 Jean*0387
9744e36521 Jean*0388 C--   Scale all grid-factor when original grid-file corresponds to
                0389 C     a different planet radius (radius_fromHorizGrid <> rSphere)
                0390       IF ( rSphere.NE.radius_fromHorizGrid ) THEN
                0391         tmpFac = rSphere / radius_fromHorizGrid
                0392         tmpFac2 = tmpFac*tmpFac
                0393         DO bj = myByLo(myThid), myByHi(myThid)
                0394          DO bi = myBxLo(myThid), myBxHi(myThid)
c679ff8f37 Jean*0395           DO j=1-OLy,sNy+OLy
                0396            DO i=1-OLx,sNx+OLx
9744e36521 Jean*0397             dxC(i,j,bi,bj) = dxC(i,j,bi,bj)*tmpFac
                0398             dyC(i,j,bi,bj) = dyC(i,j,bi,bj)*tmpFac
                0399             dxG(i,j,bi,bj) = dxG(i,j,bi,bj)*tmpFac
                0400             dyG(i,j,bi,bj) = dyG(i,j,bi,bj)*tmpFac
                0401             dxF(i,j,bi,bj) = dxF(i,j,bi,bj)*tmpFac
                0402             dyF(i,j,bi,bj) = dyF(i,j,bi,bj)*tmpFac
                0403             dxV(i,j,bi,bj) = dxV(i,j,bi,bj)*tmpFac
                0404             dyU(i,j,bi,bj) = dyU(i,j,bi,bj)*tmpFac
                0405             rA (i,j,bi,bj) = rA (i,j,bi,bj)*tmpFac2
                0406             rAz(i,j,bi,bj) = rAz(i,j,bi,bj)*tmpFac2
                0407             rAw(i,j,bi,bj) = rAw(i,j,bi,bj)*tmpFac2
                0408             rAs(i,j,bi,bj) = rAs(i,j,bi,bj)*tmpFac2
                0409            ENDDO
                0410           ENDDO
                0411          ENDDO
                0412         ENDDO
                0413       ENDIF
                0414
30d424dc37 Jean*0415 C--   Calculate (sines and cosines of) angles of grid north with
be9aaf9124 Jean*0416 C     geographical north when they have not been read from a file
                0417       CALL CALC_GRID_ANGLES( anglesAreSet, myThid )
                0418
54a8fe0a22 Jean*0419 C--   Exchange Angle (either loaded from file or computed)
                0420       CALL EXCH_UV_AGRID_3D_RS(angleSinC,angleCosC,.TRUE., 1, myThid)
                0421
aa6b2555c8 Jean*0422 c     CALL WRITE_FULLARRAY_RS('dxV',dxV,1,0,0,1,0,myThid)
                0423 c     CALL WRITE_FULLARRAY_RS('dyU',dyU,1,0,0,1,0,myThid)
                0424 c     CALL WRITE_FULLARRAY_RS('rAz',rAz,1,0,0,1,0,myThid)
                0425 c     CALL WRITE_FULLARRAY_RS('xG' ,xG ,1,0,0,1,0,myThid)
                0426 c     CALL WRITE_FULLARRAY_RS('yG' ,yG ,1,0,0,1,0,myThid)
bdac319189 Jean*0427
ff02675122 Jean*0428 C--   Now let us look at all these beasts
522c728681 Jean*0429       IF ( plotLevel.GE.debLevC ) THEN
42b7ef933b Jean*0430         CALL PLOT_FIELD_XYRS( xC      , 'Current xC      ', 0, myThid )
                0431         CALL PLOT_FIELD_XYRS( yC      , 'Current yC      ', 0, myThid )
                0432         CALL PLOT_FIELD_XYRS( dxF     , 'Current dxF     ', 0, myThid )
                0433         CALL PLOT_FIELD_XYRS( dyF     , 'Current dyF     ', 0, myThid )
                0434         CALL PLOT_FIELD_XYRS( rA      , 'Current rA      ', 0, myThid )
                0435         CALL PLOT_FIELD_XYRS( xG      , 'Current xG      ', 0, myThid )
                0436         CALL PLOT_FIELD_XYRS( yG      , 'Current yG      ', 0, myThid )
                0437         CALL PLOT_FIELD_XYRS( dxV     , 'Current dxV     ', 0, myThid )
                0438         CALL PLOT_FIELD_XYRS( dyU     , 'Current dyU     ', 0, myThid )
                0439         CALL PLOT_FIELD_XYRS( rAz     , 'Current rAz     ', 0, myThid )
                0440         CALL PLOT_FIELD_XYRS( dxC     , 'Current dxC     ', 0, myThid )
                0441         CALL PLOT_FIELD_XYRS( dyC     , 'Current dyC     ', 0, myThid )
                0442         CALL PLOT_FIELD_XYRS( rAw     , 'Current rAw     ', 0, myThid )
                0443         CALL PLOT_FIELD_XYRS( rAs     , 'Current rAs     ', 0, myThid )
                0444         CALL PLOT_FIELD_XYRS( dxG     , 'Current dxG     ', 0, myThid )
                0445         CALL PLOT_FIELD_XYRS( dyG     , 'Current dyG     ', 0, myThid )
                0446         CALL PLOT_FIELD_XYRS(angleCosC, 'Current AngleCS ', 0, myThid )
                0447         CALL PLOT_FIELD_XYRS(angleSinC, 'Current AngleSN ', 0, myThid )
3e0af136db Dimi*0448       ENDIF
bdac319189 Jean*0449
                0450       RETURN
                0451       END