      PROGRAM RZONAVG                                                                                                               
C     PROGRAM RZONAVG (X,       D,       RZX,       XS,       RZXS2,            D2                                                  
C    1                                                        OUTPUT,   )       D2                                                  
C    2          TAPE11=X,TAPE12=D,TAPE13=RZX,TAPE14=XS,TAPE15=RZXS2,            D2                                                  
C    3                                                 TAPE6 =OUTPUT)           D2                                                  
C     ---------------------------------------------------------------           D2                                                  
C                                                                               D2                                                  
C     JUL 03/03 - A.IRWIN,F.MAJAESS (REVISED FOR F90 AND LINUX PLATFORMS)       D2                                                  
C     APR 13/00 - S.KHARIN (OPTIONALLY CALCULATE DEVIATIONS AND ZONAL           D2                                                  
C                           AVERAGES OF SQUARED DEVIATIONS.                     D2                                                  
C                           USE MAX PACKING DENSITY = 2)                        D2                                                  
C     JAN 29/92 - E. CHAN  (CONVERT HOLLERITH LITERALS TO ASCII)                                                                    
C     MAY 13/83 - R.LAPRISE.                                                                                                        
C                                                                               D2                                                  
CRZONAVG - COMPUTES THE REPRESENTATIVE ZONAL AVERAGE OF A FIELD         2  3    D1                                                  
C                                                                               D3                                                  
CAUTHOR  - R.LAPRISE                                                            D3                                                  
C                                                                               D3                                                  
CPURPOSE - COMPUTE THE REPRESENTATIVE ZONAL AVERAGE (RZX) OF A SET              D3                                                  
C          OF GRIDS (X) USING (D) AS A MASK.                                    D3                                                  
C          OPTIONALLY, ZONAL DEVIATIONS AND THE REPRESENTATIVE ZONAL            D3                                                  
C          AVERAGE OF  SQUARED DEVIATIONS ARE COMPUTED.                         D3                                                  
C          NOTE - THERE MUST BE AS MANY TIME STEPS OF (D) AS OF (X),            D3                                                  
C                 AND THEIR VERTICAL LEVELS MUST MATCH.                         D3                                                  
C                                                                               D3                                                  
CINPUT FILES...                                                                 D3                                                  
C                                                                               D3                                                  
C      X    = GRIDS TO BE ZONALLY AVERAGED.                                     D3                                                  
C      D    = GRIDS OF MASK (RANGING FROM 0 BELOW GROUND TO 1 ABOVE)            D3                                                  
C                                                                               D3                                                  
COUTPUT FILE...                                                                 D3                                                  
C                                                                               D3                                                  
C      RZX  = CONTAINS RZX VECTORS OF REPRESENTATIVE ZONAL AVERAGES,            D3                                                  
C             WHERE,                                                            D3                                                  
C                RZX = Z(D*X) / Z(D),  IF Z(D).NE.0.,                           D3                                                  
C                    = Z(X)         ,  IF Z(D).EQ.0.                            D3                                                  
C      XS   = (OPTIONAL) GRIDS OF DEVIATIONS FROM THE ZONAL MEANS.              D3                                                  
C      RZXS2= (OPTIONAL) REPRESENTATIVE ZONAL AVERAGES OF SQUARED DEVIATIONS    D3                                                  
C             FROM THE ZONAL MEANS. IF ONLY RZXS2' IS REQUIRED AND 'XS' IS NOT  D3                                                  
C             NEEDED, USE THE UNDERSCORE '_' INSTEAD OF 'XS'.                   D3                                                  
CEXAMPLES:                                                                      D3                                                  
C    rzonavg x d rzx          # CALCULATE REPRESENTATIVE ZONAL AVERAGE.         D3                                                  
C    rzonavg x d rzx xs       # CALCULATE REPR. ZONAL AVERAGE AND DEVIATIONS.   D3                                                  
C    rzonavg x d rzx xs rzxs2 # CALCULATE REPR. ZONAL AVERAGE, DEVIATIONS AND   D3                                                  
C                             # VARIANCE OF STATIONARY EDDIES.                  D3                                                  
C    rzonavg x d rzx _  rzxs2 # CALCULATE REPR. ZONAL AVERAGE AND STAT. EDDIES. D3                                                  
C-------------------------------------------------------------------------                                                          
C                                                                                                                                   
      IMPLICIT REAL (A-H,O-Z),                                                                                                      
     +INTEGER (I-N)                                                                                                                 
      COMMON/BLANCK/ X(18528),D(18528),RZX(96),RZXS2(96)                                                                            
C                                                                                                                                   
      LOGICAL OK,ZONDEV,ZONVAR                                                                                                      
C                                                                                                                                   
      COMMON /ICOM/ IBUF(8),IDAT(18528)                                                                                             
      COMMON /JCOM/ JBUF(8),JDAT(18528)                                                                                             
C                                                                                                                                   
C     * IUA/IUST  - ARRAY KEEPING TRACK OF UNITS ASSIGNED/STATUS.                                                                   
C                                                                                                                                   
      INTEGER IUA(100),IUST(100)                                                                                                    
      COMMON /JCLPNTU/ IUA,IUST                                                                                                     
C                                                                                                                                   
      DATA MAXX/18528/                                                                                                              
C---------------------------------------------------------------------                                                              
      NF=6                                                                                                                          
      CALL JCLPNT(NF,11,12,13,14,15,6)                                                                                              
      IF (NF.LT.4) CALL                            XIT('RZONAVG',-1)                                                                
      REWIND 11                                                                                                                     
      REWIND 12                                                                                                                     
      REWIND 13                                                                                                                     
C                                                                                                                                   
C     * CHECK IF ZONAL DEVIATIONS ARE REQUESTED.                                                                                    
C                                                                                                                                   
      ZONDEV=.FALSE.                                                                                                                
      IF (IUST(14).EQ.1) THEN                                                                                                       
         WRITE(6,6005)                                                                                                              
         ZONDEV=.TRUE.                                                                                                              
         REWIND 14                                                                                                                  
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * CHECK IF VARIANCE OF STATIONARY EDDIES IS REQUESTED.                                                                        
C                                                                                                                                   
      ZONVAR=.FALSE.                                                                                                                
      IF (IUST(15).EQ.1) THEN                                                                                                       
         WRITE(6,6006)                                                                                                              
         ZONVAR=.TRUE.                                                                                                              
         REWIND 15                                                                                                                  
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * READ THE NEXT PAIR OF X AND D FIELDS. STOP AT EOF.                                                                          
C                                                                                                                                   
      NR=0                                                                                                                          
  100 CALL GETFLD2(11,X,NC4TO8("GRID"),-1,-1,-1,IBUF,MAXX,OK)                                                                       
      IF(.NOT.OK)THEN                                                                                                               
        IF(NR.EQ.0) CALL                           XIT('RZONAVG',-2)                                                                
        WRITE(6,6010) NR                                                                                                            
        CALL                                       XIT('RZONAVG',0)                                                                 
      ENDIF                                                                                                                         
      IF(NR.EQ.0) THEN                                                                                                              
        CALL PRTLAB (IBUF)                                                                                                          
        NPACK=MIN(2,IBUF(8))                                                                                                        
      ENDIF                                                                                                                         
C                                                                                                                                   
      CALL GETFLD2(12,D,NC4TO8("GRID"),-1,-1,-1,JBUF,MAXX,OK)                                                                       
      IF(.NOT.OK) CALL                             XIT('RZONAVG',-3)                                                                
      IF(NR.EQ.0) CALL PRTLAB (JBUF)                                                                                                
C                                                                                                                                   
C     * CHECK SIZES.                                                                                                                
C                                                                                                                                   
      CALL CMPLBL (0,IBUF,0,JBUF,OK)                                                                                                
      IF(.NOT.OK .OR. IBUF(4).NE.JBUF(4)) THEN                                                                                      
        CALL PRTLAB (IBUF)                                                                                                          
        CALL PRTLAB (JBUF)                                                                                                          
        CALL                                       XIT('RZONAVG',-4)                                                                
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * COMPUTE THE REPRESENTATIVE ZONAL AVERAGE FOR THIS GRID.                                                                     
C     * REMEMBER THAT THE LAST VALUE IN EACH ROW IS A COPY OF THE FIRST.                                                            
C                                                                                                                                   
      NLG=IBUF(5)                                                                                                                   
      NLAT=IBUF(6)                                                                                                                  
      NLGM=NLG-1                                                                                                                    
C                                                                                                                                   
      DO 300 J=1,NLAT                                                                                                               
      ZDJ =0.E0                                                                                                                     
      ZXJ =0.E0                                                                                                                     
      ZDXJ=0.E0                                                                                                                     
      N=(J-1)*NLG                                                                                                                   
      DO 200 I=1,NLGM                                                                                                               
      ZDJ  =ZDJ +D(N+I)                                                                                                             
      ZXJ  =ZXJ +X(N+I)                                                                                                             
      ZDXJ =ZDXJ+X(N+I)*D(N+I)                                                                                                      
  200 CONTINUE                                                                                                                      
      ZDJ =ZDJ /FLOAT(NLGM)                                                                                                         
      ZXJ =ZXJ /FLOAT(NLGM)                                                                                                         
      ZDXJ=ZDXJ/FLOAT(NLGM)                                                                                                         
      IF(ZDJ.NE.0.E0)THEN                                                                                                           
         RZX(J)=ZDXJ/ZDJ                                                                                                            
      ELSE                                                                                                                          
         RZX(J)=ZXJ                                                                                                                 
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * COMPUTE THE ZONAL DEVIATIONS FOR X.                                                                                         
C                                                                                                                                   
      IF (ZONDEV.OR.ZONVAR) THEN                                                                                                    
        DO 240 I=1,NLG                                                                                                              
  240   X(N+I)=X(N+I)-RZX(J)                                                                                                        
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * COMPUTE THE REPRESENTATIVE ZONAL AVERAGE OF SQUARED DEVIATIONS.                                                             
C                                                                                                                                   
      IF (ZONVAR) THEN                                                                                                              
         ZXJ =0.E0                                                                                                                  
         ZDXJ=0.E0                                                                                                                  
         DO 250 I=1,NLGM                                                                                                            
         ZXJ =ZXJ +X(N+I)*X(N+I)                                                                                                    
         ZDXJ=ZDXJ+X(N+I)*X(N+I)*D(N+I)                                                                                             
  250    CONTINUE                                                                                                                   
         ZXJ =ZXJ /FLOAT(NLGM)                                                                                                      
         ZDXJ=ZDXJ/FLOAT(NLGM)                                                                                                      
         IF(ZDJ.NE.0.E0)THEN                                                                                                        
            RZXS2(J)=ZDXJ/ZDJ                                                                                                       
         ELSE                                                                                                                       
            RZXS2(J)=ZXJ                                                                                                            
         ENDIF                                                                                                                      
      ENDIF                                                                                                                         
  300 CONTINUE                                                                                                                      
C                                                                                                                                   
C     * WRITE OUT THE REPRESENTATIVE ZONAL AVERAGE VECTOR.                                                                          
C                                                                                                                                   
      CALL SETLAB(IBUF,NC4TO8("ZONL"),-1,-1,-1,NLAT,1,-1,NPACK)                                                                     
      CALL PUTFLD2(13,RZX,IBUF,MAXX)                                                                                                
      IF(NR.EQ.0) CALL PRTLAB (IBUF)                                                                                                
C                                                                                                                                   
C     * WRITE OUT THE GRID OF ZONAL DEVIATIONS.                                                                                     
C                                                                                                                                   
      IF (ZONDEV) THEN                                                                                                              
        CALL SETLAB(IBUF,NC4TO8("GRID"),-1,-1,-1,NLG,NLAT,-1,NPACK)                                                                 
        CALL PUTFLD2(14,X,IBUF,MAXX)                                                                                                
        IF(NR.EQ.0) CALL PRTLAB (IBUF)                                                                                              
      ENDIF                                                                                                                         
C                                                                                                                                   
C     * WRITE OUT THE REPRESENTATIVE ZONAL AVERAGES OF SQUARED DEVIATIONS.                                                          
C                                                                                                                                   
      IF (ZONVAR) THEN                                                                                                              
        CALL SETLAB(IBUF,NC4TO8("ZONL"),-1,-1,-1,NLAT,1,-1,NPACK)                                                                   
        CALL PUTFLD2(15,RZXS2,IBUF,MAXX)                                                                                            
        IF(NR.EQ.0) CALL PRTLAB (IBUF)                                                                                              
      ENDIF                                                                                                                         
      NR=NR+1                                                                                                                       
      GO TO 100                                                                                                                     
C---------------------------------------------------------------------                                                              
 6005 FORMAT('0 CALCULATE ALSO ZONAL DEVIATIONS.')                                                                                  
 6006 FORMAT('0 CALCULATE ALSO ZONAL AVERAGES OF SQUARED DEVIATIONS.')                                                              
 6010 FORMAT('0',I6,'  PAIRS OF RECORDS PROCESSED')                                                                                 
      END                                                                                                                           
