Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:37:02 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
4416d8eda1 Gael*0001 #include "CPP_OPTIONS.h"
                0002 
                0003 C--  File rotate_uv2en.F: Routines to handle a vector coordinate system rotation.
                0004 C--   Contents
                0005 C--   o ROTATE_UV2EN_RL
                0006 C--   o ROTATE_UV2EN_RS
                0007 
                0008       subroutine rotate_uv2en_rl(
                0009      U          uFldX, vFldY,
                0010      U          uFldE, vFldN,
                0011      I          xy2en, switchGrid, applyMask, kSize, mythid
                0012      &                     )
                0013 
                0014 c     ==================================================================
                0015 c     SUBROUTINE rotate_uv2en_rl
                0016 c     ==================================================================
                0017 c
                0018 c     o uFldX/vFldY are in the model grid directions. 
                0019 c     o uFldE/vFldN are eastward/northward.
                0020 c     o This routine goes from uFldX/vFldY to uFldE/vFldN (for xy2en=.TRUE.) 
                0021 c         or vice versa (for xy2en=.FALSE.).
                0022 c     o uFldX/vFldY may be at the C grid velocity points, or at the A grid
                0023 c         velocity points (i.e. the C grid cell center). uFldE/vFldN are always 
                0024 c         at the cell center. If switchGrid=.TRUE. we go from C grid uFldX/vFldY
                0025 c         to A grid uFldE/vFldN, or vice versa. If switchGrid=.FALSE. we go
                0026 c         from A grid uFldX/vFldY to A grid uFldE/vFldN, or vice versa.
                0027 c     o If applyMask=.TRUE. then masks are applied to the output.
                0028 c         If kSize=1 (resp. nr) we then use the surface (resp. 3D) masks.
                0029 c     o In any case it is assumed that exchanges are done on the outside.
                0030 c
                0031 c     ==================================================================
                0032 c     SUBROUTINE rotate_uv2en_rl
                0033 c     ==================================================================
                0034 
                0035       implicit none
                0036 
                0037 c     == global variables ==
                0038 
                0039 #include "EEPARAMS.h"
                0040 #include "SIZE.h"
                0041 #include "PARAMS.h"
                0042 #include "GRID.h"
                0043 
                0044 c     == routine arguments ==
                0045 
                0046       integer kSize 
                0047       logical xy2en, switchGrid, applyMask
                0048       _RL     uFldX(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0049       _RL     vFldY(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0050       _RL     uFldE(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0051       _RL     vFldN(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0052 
                0053       integer mythid
                0054 
                0055 c     == local variables ==
                0056 
                0057       integer bi,bj
                0058       integer i,j,k,kk
                0059       _RL     tmpU(1-olx:snx+olx,1-oly:sny+oly)
                0060       _RL     tmpV(1-olx:snx+olx,1-oly:sny+oly)
                0061       CHARACTER*(MAX_LEN_MBUF) msgBuf
                0062 
                0063 c     == end of interface ==
                0064 
                0065       if ( (kSize.NE.1).AND.(kSize.NE.nr)
                0066      &   .AND.(applyMask) ) then
                0067         WRITE(msgBuf,'(2A,I4,A)') ' ROTATE_UV2EN: ',
                0068      &       'no mask has ',kSize,' levels'
                0069         CALL PRINT_ERROR(msgBuf, myThid)      
                0070         STOP 'ABNROMAL END: S/R ROTATE_UV2EN' 
                0071       endif
                0072 
                0073       do bj = mybylo(mythid),mybyhi(mythid)
                0074        do bi = mybxlo(mythid),mybxhi(mythid)
                0075         do k = 1,kSize
                0076 
                0077         if ( (kSize.EQ.1).AND.(usingPCoords) ) then
                0078           kk=nr
                0079         else
                0080           kk=k
                0081         endif
                0082 
                0083         if ( xy2en ) then
                0084 c go from uFldX/vFldY to uFldE/vFldN
                0085         if ( switchGrid ) then
                0086 C 1a) go from C grid velocity points to A grid velocity points
                0087         do i = 1-olx,snx+olx
                0088           tmpU(i,sny+Oly)=0.
                0089           tmpV(i,sny+Oly)=0.
                0090         enddo
                0091         do j = 1-oly,sny+oly-1
                0092           tmpU(snx+Olx,j)=0.
                0093           tmpV(snx+Olx,j)=0.
                0094         do i = 1-olx,snx+olx-1
                0095           tmpU(i,j) = 0.5 _d 0
                0096      &          *( uFldX(i+1,j,k,bi,bj) + uFldX(i,j,k,bi,bj) )
                0097           tmpV(i,j) = 0.5 _d 0
                0098      &          *( vFldY(i,j+1,k,bi,bj) + vFldY(i,j,k,bi,bj) )
                0099           if (applyMask) then
                0100             tmpU(i,j) = tmpU(i,j) * maskC(i,j,kk,bi,bj)
                0101             tmpV(i,j) = tmpV(i,j) * maskC(i,j,kk,bi,bj)
                0102           endif
                0103           enddo
                0104           enddo
                0105         else
                0106 C 1b) stay at A grid velocity points (i.e. at the C grid cell center)
                0107           do j = 1-oly,sny+oly
                0108           do i = 1-olx,snx+olx
                0109             tmpU(i,j) = uFldX(i,j,k,bi,bj)
                0110             tmpV(i,j) = vFldY(i,j,k,bi,bj)
                0111             if (applyMask) then
                0112               tmpU(i,j) = tmpU(i,j) * maskC(i,j,kk,bi,bj)
                0113               tmpV(i,j) = tmpV(i,j) * maskC(i,j,kk,bi,bj)
                0114             endif
                0115           enddo
                0116           enddo
                0117         endif!if ( switchGrid ) then
                0118 
                0119 C 2) rotation
                0120           do j = 1-oly,sny+oly
                0121           do i = 1-olx,snx+olx
                0122             uFldE(i,j,k,bi,bj) = 
                0123      &         angleCosC(i,j,bi,bj)*tmpU(i,j)
                0124      &        -angleSinC(i,j,bi,bj)*tmpV(i,j)
                0125             vFldN(i,j,k,bi,bj) = 
                0126      &         angleSinC(i,j,bi,bj)*tmpU(i,j)
                0127      &        +angleCosC(i,j,bi,bj)*tmpV(i,j)
                0128           enddo
                0129           enddo
                0130 
                0131       else!if (xy2en) then
                0132 c go from uFldE/vFldN to uFldX/vFldY
                0133 
                0134 C 1) rotation
                0135           do j = 1-oly,sny+oly
                0136           do i = 1-olx,snx+olx
                0137             tmpU(i,j) = 
                0138      &         angleCosC(i,j,bi,bj)*uFldE(i,j,k,bi,bj)
                0139      &        +angleSinC(i,j,bi,bj)*vFldN(i,j,k,bi,bj)
                0140             tmpV(i,j) =
                0141      &        -angleSinC(i,j,bi,bj)*uFldE(i,j,k,bi,bj)
                0142      &        +angleCosC(i,j,bi,bj)*vFldN(i,j,k,bi,bj)
                0143           enddo
                0144           enddo
                0145 
                0146         if ( switchGrid ) then
                0147 C 2a) go from A grid velocity points to C grid velocity points
                0148           do i = 1-olx,snx+olx
                0149             uFldX(i,1,k,bi,bj)=0.
                0150             vFldY(i,1,k,bi,bj)=0.
                0151           enddo
                0152           do j = 1-oly+1,sny+oly
                0153              uFldX(1,j,k,bi,bj)=0.
                0154              vFldY(1,j,k,bi,bj)=0.
                0155           do i = 1-olx+1,snx+olx
                0156             uFldX(i,j,k,bi,bj) = 0.5 _d 0
                0157      &         *( tmpU(i-1,j) + tmpU(i,j) )
                0158             vFldY(i,j,k,bi,bj) = 0.5 _d 0
                0159      &         *( tmpV(i,j-1) + tmpV(i,j) )
                0160             if (applyMask) then
                0161               uFldX(i,j,k,bi,bj)=uFldX(i,j,k,bi,bj)*maskW(i,j,kk,bi,bj)
                0162               vFldY(i,j,k,bi,bj)=vFldY(i,j,k,bi,bj)*maskS(i,j,kk,bi,bj)
                0163             endif
                0164           enddo
                0165           enddo
                0166         else
                0167 C 2b) stay at A grid velocity points (i.e. at the C grid cell center)
                0168           do j = 1-oly,sny+oly
                0169           do i = 1-olx,snx+olx
                0170             uFldX(i,j,k,bi,bj) = tmpU(i,j)
                0171             vFldY(i,j,k,bi,bj) = tmpV(i,j)
                0172             if (applyMask) then
                0173               uFldX(i,j,k,bi,bj)=uFldX(i,j,k,bi,bj)*maskC(i,j,kk,bi,bj)
                0174               vFldY(i,j,k,bi,bj)=vFldY(i,j,k,bi,bj)*maskC(i,j,kk,bi,bj)
                0175             endif
                0176           enddo
                0177           enddo
                0178         endif!if ( switchGrid ) then
                0179 
                0180         endif!if (xy2en) then
                0181 
                0182         enddo
                0183        enddo
                0184       enddo
                0185 
                0186       return
                0187       end
                0188 
                0189       subroutine rotate_uv2en_rs(
                0190      U          uFldX, vFldY,
                0191      U          uFldE, vFldN,
                0192      I          xy2en, switchGrid, applyMask, kSize, mythid
                0193      &                     )
                0194 
                0195 c     ==================================================================
                0196 c     SUBROUTINE rotate_uv2en_rs
                0197 c     ==================================================================
                0198 c
                0199 c     o uFldX/vFldY are in the model grid directions. 
                0200 c     o uFldE/vFldN are eastward/northward.
                0201 c     o This routine goes from uFldX/vFldY to uFldE/vFldN (for xy2en=.TRUE.) 
                0202 c         or vice versa (for xy2en=.FALSE.).
                0203 c     o uFldX/vFldY may be at the C grid velocity points, or at the A grid
                0204 c         velocity points (i.e. the C grid cell center). uFldE/vFldN are always 
                0205 c         at the cell center. If switchGrid=.TRUE. we go from C grid uFldX/vFldY
                0206 c         to A grid uFldE/vFldN, or vice versa. If switchGrid=.FALSE. we go
                0207 c         from A grid uFldX/vFldY to A grid uFldE/vFldN, or vice versa.
                0208 c     o If applyMask=.TRUE. then masks are applied to the output.
                0209 c         If kSize=1 (resp. nr) we then use the surface (resp. 3D) masks.
                0210 c     o In any case it is assumed that exchanges are done on the outside.
                0211 c
                0212 c     ==================================================================
                0213 c     SUBROUTINE rotate_uv2en_rs
                0214 c     ==================================================================
                0215 
                0216       implicit none
                0217 
                0218 c     == global variables ==
                0219 
                0220 #include "EEPARAMS.h"
                0221 #include "SIZE.h"
                0222 #include "PARAMS.h"
                0223 #include "GRID.h"
                0224 
                0225 c     == routine arguments ==
                0226 
                0227       integer kSize 
                0228       logical xy2en, switchGrid, applyMask
                0229       _RS     uFldX(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0230       _RS     vFldY(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0231       _RS     uFldE(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0232       _RS     vFldN(1-olx:snx+olx,1-oly:sny+oly,kSize,nsx,nsy)
                0233 
                0234       integer mythid
                0235 
                0236 c     == local variables ==
                0237 
                0238       integer bi,bj
                0239       integer i,j,k,kk
                0240       _RS     tmpU(1-olx:snx+olx,1-oly:sny+oly)
                0241       _RS     tmpV(1-olx:snx+olx,1-oly:sny+oly)
                0242       CHARACTER*(MAX_LEN_MBUF) msgBuf
                0243 
                0244 c     == end of interface ==
                0245 
                0246       if ( (kSize.NE.1).AND.(kSize.NE.nr)
                0247      &   .AND.(applyMask) ) then
                0248         WRITE(msgBuf,'(2A,I4,A)') ' ROTATE_UV2EN: ',
                0249      &       'no mask has ',kSize,' levels'
                0250         CALL PRINT_ERROR(msgBuf, myThid)      
                0251         STOP 'ABNROMAL END: S/R ROTATE_UV2EN' 
                0252       endif
                0253 
                0254       do bj = mybylo(mythid),mybyhi(mythid)
                0255        do bi = mybxlo(mythid),mybxhi(mythid)
                0256         do k = 1,kSize
                0257 
                0258         if ( (kSize.EQ.1).AND.(usingPCoords) ) then
                0259           kk=nr
                0260         else
                0261           kk=k
                0262         endif
                0263 
                0264         if ( xy2en ) then
                0265 c go from uFldX/vFldY to uFldE/vFldN
                0266         if ( switchGrid ) then
                0267 C 1a) go from C grid velocity points to A grid velocity points
                0268         do i = 1-olx,snx+olx
                0269           tmpU(i,sny+Oly)=0.
                0270           tmpV(i,sny+Oly)=0.
                0271         enddo
                0272         do j = 1-oly,sny+oly-1
                0273           tmpU(snx+Olx,j)=0.
                0274           tmpV(snx+Olx,j)=0.
                0275         do i = 1-olx,snx+olx-1
                0276           tmpU(i,j) = 0.5 _d 0
                0277      &          *( uFldX(i+1,j,k,bi,bj) + uFldX(i,j,k,bi,bj) )
                0278           tmpV(i,j) = 0.5 _d 0
                0279      &          *( vFldY(i,j+1,k,bi,bj) + vFldY(i,j,k,bi,bj) )
                0280           if (applyMask) then
                0281             tmpU(i,j) = tmpU(i,j) * maskC(i,j,kk,bi,bj)
                0282             tmpV(i,j) = tmpV(i,j) * maskC(i,j,kk,bi,bj)
                0283           endif
                0284           enddo
                0285           enddo
                0286         else
                0287 C 1b) stay at A grid velocity points (i.e. at the C grid cell center)
                0288           do j = 1-oly,sny+oly
                0289           do i = 1-olx,snx+olx
                0290             tmpU(i,j) = uFldX(i,j,k,bi,bj)
                0291             tmpV(i,j) = vFldY(i,j,k,bi,bj)
                0292             if (applyMask) then
                0293               tmpU(i,j) = tmpU(i,j) * maskC(i,j,kk,bi,bj)
                0294               tmpV(i,j) = tmpV(i,j) * maskC(i,j,kk,bi,bj)
                0295             endif
                0296           enddo
                0297           enddo
                0298         endif!if ( switchGrid ) then
                0299 
                0300 C 2) rotation
                0301           do j = 1-oly,sny+oly
                0302           do i = 1-olx,snx+olx
                0303             uFldE(i,j,k,bi,bj) = 
                0304      &         angleCosC(i,j,bi,bj)*tmpU(i,j)
                0305      &        -angleSinC(i,j,bi,bj)*tmpV(i,j)
                0306             vFldN(i,j,k,bi,bj) = 
                0307      &         angleSinC(i,j,bi,bj)*tmpU(i,j)
                0308      &        +angleCosC(i,j,bi,bj)*tmpV(i,j)
                0309           enddo
                0310           enddo
                0311 
                0312       else!if (xy2en) then
                0313 c go from uFldE/vFldN to uFldX/vFldY
                0314 
                0315 C 1) rotation
                0316           do j = 1-oly,sny+oly
                0317           do i = 1-olx,snx+olx
                0318             tmpU(i,j) = 
                0319      &         angleCosC(i,j,bi,bj)*uFldE(i,j,k,bi,bj)
                0320      &        +angleSinC(i,j,bi,bj)*vFldN(i,j,k,bi,bj)
                0321             tmpV(i,j) =
                0322      &        -angleSinC(i,j,bi,bj)*uFldE(i,j,k,bi,bj)
                0323      &        +angleCosC(i,j,bi,bj)*vFldN(i,j,k,bi,bj)
                0324           enddo
                0325           enddo
                0326 
                0327         if ( switchGrid ) then
                0328 C 2a) go from A grid velocity points to C grid velocity points
                0329           do i = 1-olx,snx+olx
                0330             uFldX(i,1,k,bi,bj)=0.
                0331             vFldY(i,1,k,bi,bj)=0.
                0332           enddo
                0333           do j = 1-oly+1,sny+oly
                0334              uFldX(1,j,k,bi,bj)=0.
                0335              vFldY(1,j,k,bi,bj)=0.
                0336           do i = 1-olx+1,snx+olx
                0337             uFldX(i,j,k,bi,bj) = 0.5 _d 0
                0338      &         *( tmpU(i-1,j) + tmpU(i,j) )
                0339             vFldY(i,j,k,bi,bj) = 0.5 _d 0
                0340      &         *( tmpV(i,j-1) + tmpV(i,j) )
                0341             if (applyMask) then
                0342               uFldX(i,j,k,bi,bj)=uFldX(i,j,k,bi,bj)*maskW(i,j,kk,bi,bj)
                0343               vFldY(i,j,k,bi,bj)=vFldY(i,j,k,bi,bj)*maskS(i,j,kk,bi,bj)
                0344             endif
                0345           enddo
                0346           enddo
                0347         else
                0348 C 2b) stay at A grid velocity points (i.e. at the C grid cell center)
                0349           do j = 1-oly,sny+oly
                0350           do i = 1-olx,snx+olx
                0351             uFldX(i,j,k,bi,bj) = tmpU(i,j)
                0352             vFldY(i,j,k,bi,bj) = tmpV(i,j)
                0353             if (applyMask) then
                0354               uFldX(i,j,k,bi,bj)=uFldX(i,j,k,bi,bj)*maskC(i,j,kk,bi,bj)
                0355               vFldY(i,j,k,bi,bj)=vFldY(i,j,k,bi,bj)*maskC(i,j,kk,bi,bj)
                0356             endif
                0357           enddo
                0358           enddo
                0359         endif!if ( switchGrid ) then
                0360 
                0361         endif!if (xy2en) then
                0362 
                0363         enddo
                0364        enddo
                0365       enddo
                0366 
                0367       return
                0368       end
                0369