      PROGRAM GRIDO                                                                                                                 
C     PROGRAM GRIDO(GRIDI,        BETA,          DA,           DX,              D2                                                  
C    1                 DY,          DZ,         OUTPUT,                )        D2                                                  
C    2       TAPE11=GRIDI, TAPE12=BETA,   TAPE13=DA,    TAPE14=DX,                                                                  
C    3          TAPE15=DY,   TAPE16=DZ, , TAPE6=OUTPUT)                                                                             
C     ------------------------------------------------------------              D2                                                  
C                                                                               D2                                                  
C     JUL 03/03 - A.IRWIN,F.MAJAESS (REVISED FOR F90 AND LINUX PLATFORMS)       D2                                                  
C     APR 04/2002 - BM  (MODIFIED FOR NEW DATA DESCRIPTION FORMAT)              D2                                                  
C     OCT 12/2001 - WGL (CODE RE-WRITTEN TO SEPARATE OUT SECTION THAT READS IN                                                      
C                        THE OCEAN GRIDINFO FILE)                                                                                   
C     JUN 06/2001 - BM - ORIGINAL VERSION                                                                                           
C                                                                               D2                                                  
CGRIDO   - CALCULATE GRID INFO FOR OGCM MODELS                          1  6    D1                                                  
C                                                                               D3                                                  
CAUTHOR  - B. MIVILLE                                                           D3                                                  
C                                                                               D3                                                  
CPURPOSE - COMPUTE BETA, DA,DX,DY AND DZ FOR FULL INDEX GRID CELL (TRACER CELL) D3                                                  
C          FOR OGCM MODELS. SEE INTERNAL WEB PAGE FOR FURTHER INFORMATION       D3                                                  
C                                                                               D3                                                  
CINPUT FILE...                                                                  D3                                                  
C                                                                               D3                                                  
C  GRIDI = GRID INFORMATION FILE CONTAINING:                                    D3                                                  
C                                                                               D3                                                  
C           DSRC = DATA SOURCE  (CURRENTLY MODEL TYPE OGCM1,2,3)                D3                                                  
C           GTYP = GRID TYPE (CURRENTLY ONLY LL)                                D3                                                  
C           RADO = RADIUS OF THE EARTH (DEPENDS ON MODEL)                       D3                                                  
C           LAT  = LATITUDE OF FULL INDEX (J) POINTS (CURRENTLY T, S, UI)       D3                                                  
C           LATH = LATITUDE OF HALF INDEX (J+1/2) POINTS (CURRENTLY U, V, VI)   D3                                                  
C           LON  = LONGITUDE OF FULL INDEX (I) POINTS (CURRENTLY T, S, VI)      D3                                                  
C           LONH = LONGITUDE OF HALF INDEX (I+1/2) POINTS (U, V, UI)            D3                                                  
C           Z    = DEPTH OF FULL LEVEL (K) POINTS (CURRENTLY T,S,U,V,UI,VI)     D3                                                  
C           ZH   = DEPTH OF HELP LEVEL (K+1/2) POINTS (CURRENTLY  W)            D3                                                  
C           ZBOT = DEPTH OF BOTTOM OF OCEAN                                     D3                                                  
C           BASN = BASIN NUMBER ASSIGNED TO EACH OCEAN GRID                     D3                                                  
C                                                                               D3                                                  
COUTPUT FILES...                                                                D3                                                  
C                                                                               D3                                                  
C  BETA = 3-D BETA (0 FOR LAND, 1 FOR OCEAN, FRACTION FOR BOTTOM SHAVED CELL)   D3                                                  
C  DA   = AREA OF FULL INDEX (J) GRID CELL, 2-D                                 D3                                                  
C  DX   = EAST-WEST EXTENT OF FULL INDEX (J) GRID CELL, 2-D                     D3                                                  
C  DY   = NORTH-SOUTH EXTENT OF FULL INDEX (J) GRID CELL, 2-D                   D3                                                  
C  DZ   = VERTICAL EXTENT OF FULL INDEX (J) GRID CELL, 3-D                      D3                                                  
C------------------------------------------------------------------------------                                                     
C                                                                                                                                   
C     * MAXIMUM 2-D GRID SIZE (I+2)*(J+2)=IM*JM --> IJM                                                                             
C     * FOR MAX GRID SIZE OF I*J PLUS ONE MORE ON EACH SIDE.                                                                        
C                                                                                                                                   
C     * MAXIMUM NUMBER OF OCEAN DEPTH IS "KM".                                                                                      
C                                                                                                                                   
      IMPLICIT REAL (A-H,O-Z),                                                                                                      
     +INTEGER (I-N)                                                                                                                 
      PARAMETER (IM=771,JM=386,IJM=IM*JM,IJMV=IJM*1)                                                                                
      PARAMETER (KM=200)                                                                                                            
C                                                                                                                                   
      LOGICAL OK,SPEC                                                                                                               
C                                                                                                                                   
      REAL LAT(IJM), LATH(IJM), Z(KM), ZH(KM), ZB(IJM)                                                                              
      REAL BETA(IJM), DA(IJM), DX(IJM), DY(IJM), DZ(IJM)                                                                            
      REAL LON(IJM), LONH(IJM), BASN(IJM), DZ3D(IJM)                                                                                
      INTEGER IBASIN(IJM)                                                                                                           
      CHARACTER*8 DSOURCE,GTYPE                                                                                                     
C                                                                                                                                   
      INTEGER IBUF,IDAT,MAXX                                                                                                        
      COMMON/ICOM/IBUF(8),IDAT(IJMV)                                                                                                
      DATA MAXX/IJMV/                                                                                                               
C---------------------------------------------------------------------                                                              
C                                                                                                                                   
      NF = 7                                                                                                                        
C                                                                                                                                   
      CALL JCLPNT(NF,11,12,13,14,15,16,6)                                                                                           
c                                                                                                                                   
      REWIND 11                                                                                                                     
      REWIND 12                                                                                                                     
      REWIND 13                                                                                                                     
      REWIND 14                                                                                                                     
      REWIND 15                                                                                                                     
      REWIND 16                                                                                                                     
C                                                                                                                                   
C     * GRID INFO INPUT FILE READ UNIT NUMBER                                                                                       
C                                                                                                                                   
      NUNIT = 11                                                                                                                    
C                                                                                                                                   
C     * READ GRID INFO FILE GRIDI                                                                                                   
C                                                                                                                                   
      CALL READOGI(NUNIT,IM,JM,IJM,KM,IBUF,MAXX,DSOURCE,GTYPE,RADIUS,                                                               
     1             LAT,LATH,LON,LONH,Z,ZH,ZB,BASN,IBASIN,NLON,NLAT,NLEV)                                                            
C                                                                                                                                   
C     * WRITE MODEL VERSION, GRID SIZE INFORMATION TO STANDARD OUTPUT                                                               
C                                                                                                                                   
      WRITE(6,6001) DSOURCE, NLON, NLAT, NLEV                                                                                       
C                                                                                                                                   
      NWRDS = NLON * NLAT                                                                                                           
C                                                                                                                                   
C     * CALCULATE GRID GEOMETRIC INFORMATION                                                                                        
C                                                                                                                                   
      CALL GRIDINF(RADIUS,NLON,NLAT,LAT,LATH,LONH,DA,DX,DY)                                                                         
C                                                                                                                                   
C     * CALCULATE THE DZ                                                                                                            
C                                                                                                                                   
      DZ(1)=ZH(1)                                                                                                                   
      DO 160 K = 2,NLEV                                                                                                             
         DZ(K) = ZH(K) - ZH(K-1)                                                                                                    
 160  CONTINUE                                                                                                                      
C                                                                                                                                   
C     * WRITE DATA TO INDIVIDUAL FILE                                                                                               
C                                                                                                                                   
      CALL SETLAB(IBUF,NC4TO8("GRID"),0,-1,0,NLON,NLAT,0,1)                                                                         
C                                                                                                                                   
      IBUF(3) = NC4TO8("  DA")                                                                                                      
      CALL  PUTFLD2(13,DA,IBUF,MAXX)                                                                                                
C                                                                                                                                   
      IBUF(3) = NC4TO8("  DX")                                                                                                      
      CALL  PUTFLD2(14,DX,IBUF,MAXX)                                                                                                
C                                                                                                                                   
      IBUF(3) = NC4TO8("  DY")                                                                                                      
      CALL  PUTFLD2(15,DY,IBUF,MAXX)                                                                                                
C                                                                                                                                   
      IBUF(3) = NC4TO8("  DZ")                                                                                                      
      DO 180 K = 1,NLEV                                                                                                             
         IBUF(4) = INT(Z(K)*10.0E0+0.5E0)                                                                                           
         DO 190 IJ = 1,NWRDS                                                                                                        
            DZ3D(IJ)=DZ(K)                                                                                                          
 190     CONTINUE                                                                                                                   
         CALL  PUTFLD2(16,DZ3D,IBUF,MAXX)                                                                                           
 180  CONTINUE                                                                                                                      
C                                                                                                                                   
C     * 3-D BETA                                                                                                                    
C                                                                                                                                   
      IBUF(5) = NLON                                                                                                                
      IBUF(6) = NLAT                                                                                                                
      ZHM     = 0.0E0                                                                                                               
C                                                                                                                                   
      IBUF(3) = NC4TO8("BETA")                                                                                                      
C                                                                                                                                   
      DO 200 K = 1,NLEV                                                                                                             
         ZHP = ZH(K)                                                                                                                
C                                                                                                                                   
C        * THIS IS TO MAKE SURE THAT THERE IS NO DIVIDE BY ZERO                                                                     
C        * IN BETAINF                                                                                                               
C                                                                                                                                   
         IF(ZHP.EQ.ZHM) THEN                                                                                                        
            WRITE(6,'(A,E12.5,A,E12.5)')'PREVIOUS LEVEL: '                                                                          
     1           ,ZHM,' EQUALS CURRENT LEVEL: ',ZHP                                                                                 
            CALL                                   XIT('GRIDO',-1)                                                                  
         ENDIF                                                                                                                      
C                                                                                                                                   
         CALL BETAINF(NLON,NLAT,ZHP,ZHM,ZB,BETA)                                                                                    
C                                                                                                                                   
C        * THE 0.5 IS TO MAKE SURE THAT THERE IS NO ROUNDOFF ERROR                                                                  
C                                                                                                                                   
         IBUF(4) = INT(Z(K)*10.0E0+0.5E0)                                                                                           
         CALL PUTFLD2(12,BETA,IBUF,MAXX)                                                                                            
         ZHM = ZHP                                                                                                                  
 200  CONTINUE                                                                                                                      
C                                                                                                                                   
      CALL                                         XIT('GRIDO',0)                                                                   
C---------------------------------------------------------------------                                                              
 6001 FORMAT(' DATA SOURCE: ',A8,'  NLON: ',I5,'  NLAT: ',I5,                                                                       
     1       '  NLEV: ',I5)                                                                                                         
C                                                                                                                                   
        END                                                                                                                         
