      PROGRAM VPINTV                                                                                                                
C     PROGRAM VPINTV (XIN,       XOUT,       PSURF,       INPUT,                C2                                                  
C    1                                                    OUTPUT,       )       C2                                                  
C    2          TAPE1=XIN, TAPE2=XOUT, TAPE3=PSURF, TAPE5=INPUT,                                                                    
C    3                                              TAPE6=OUTPUT)                                                                   
C     -----------------------------------------------------------               C2                                                  
C                                                                               C2                                                  
C     JUL 03/03 - A.IRWIN,F.MAJAESS (REVISED FOR F90 AND LINUX PLATFORMS)       C2                                                  
C     JAN 14/93 - E. CHAN  (DECODE LEVELS IN 8-WORD LABEL)                                                                          
C     JAN 29/92 - E. CHAN   (CONVERT HOLLERITH LITERALS TO ASCII)                                                                   
C     JUN 09/88 - F.MAJAESS (INCREASE THE DIMENSION SO THAT IT CAN                                                                  
C                            HANDLE 20 LEVELS OF T-30 RUNS)                     C2                                                  
C     NOV 29/83 - B.DUGAS                                                                                                           
C                                                                               C2                                                  
CVPINTV  - VERTICAL INTEGRAL WITH VARIABLE BOUNDARIES                   2  1 C GC1                                                  
C                                                                               C3                                                  
CAUTHOR  - B.DUGAS                                                              C3                                                  
C                                                                               C3                                                  
CPURPOSE - DOES VERTICAL INTEGRAL OF XIN WHILE OPTIONNALY TAKING                C3                                                  
C          INTO ACCOUNT VARIABLE SURFACE PRESSURE AND/OR SPECIFIC               C3                                                  
C          INTEGRATION LIMITS.                                                  C3                                                  
C          NOTE - PSURF FILE MAY NOT BE USED DEPENDING ON THE VALUE READ        C3                                                  
C                 FOR MODE PARAMETER, (SEE THE DESCRIPTION OF MODE PARA-        C3                                                  
C                 METER BELOW). AND IF PSURF IS TO BE USED, IT HAS TO BE        C3                                                  
C                 OF THE SAME KIND AS XIN, (EITHER ZONL OR GRID).               C3                                                  
C                 MINIMUM NUMBER OF LEVELS IS 2, MAX IS $PL$.                   C3                                                  
C                                                                               C3                                                  
CINPUT FILES...                                                                 C3                                                  
C                                                                               C3                                                  
C      XIN   = INPUT SERIES OF PRESSURE LEVEL SETS TO BE INTEGRATED VERICALLY.  C3                                                  
C      PSURF = LOG OF SURFACE PRESSURE GRIDS OR ZONAL AVERAGES,                 C3                                                  
C              (OPTIONAL, CHECK MODE PARAMETER DESCRIPTION BELOW).              C3                                                  
C                                                                               C3                                                  
COUTPUT FILE...                                                                 C3                                                  
C                                                                               C3                                                  
C      XOUT  = OUTPUT ONE LEVEL VERTICAL INTEGRAL OF XIN.                       C3                                                  
C                                                                                                                                   
CINPUT PARAMETERS...                                                                                                                
C                                                                               C5                                                  
C      PB   = 0, DEFULTS TO 1013.3 MB,                                          C5                                                  
C             OTHERWISE, PB IS THE BOTTOM INTEGRATION LIMIT IN MB.              C5                                                  
C      PT   = TOP INTEGRATION LIMIT IN MB IF PB.NE.0                            C5                                                  
C           = DEFAULTS TO LEV(1)/2        IF PB.EQ.0                            C5                                                  
C      MODE = 0, VPINTV DOES NOT TAKE INTO ACCOUNT VARYING SURFACE              C5                                                  
C                PRESSURES AND XIN CAN BE ANY KIND OF FIELD. AND                C5                                                  
C                IN THIS CASE, VPINTV DOES NOT EVEN READ PSURF.                 C5                                                  
C             OTHERWISE, VPINTV TAKES INTO ACCOUNT VARYING SURFACE              C5                                                  
C                PRESSURES AND XIN HAS TO BE SAME KIND AS PSURF                 C5                                                  
C                (EITHER ZONL OR GRID).                                         C5                                                  
C                                                                               C5                                                  
CEXAMPLE OF INPUT CARD...                                                       C5                                                  
C                                                                               C5                                                  
C* VPINTV.      850.      100.         0                                        C5                                                  
C                                                                               C5                                                  
C  IN THIS CASE VPINTV WILL INTEGRATE FROM 850 TO 100 MB WITHOUT                C5                                                  
C  CONSIDERING LNSP.                                                            C5                                                  
C-----------------------------------------------------------------------------                                                      
C                                                                                                                                   
      IMPLICIT REAL (A-H,O-Z),                                                                                                      
     +INTEGER (I-N)                                                                                                                 
      COMMON/GG/GG(1852800)                                                                                                         
      COMMON/PS/PS(18528)                                                                                                           
      COMMON/VOUT/VOUT(18528)                                                                                                       
                                                                                                                                    
      LOGICAL OK,OLD                                                                                                                
      INTEGER LEV(100)                                                                                                              
      REAL PR(100),PRH(100)                                                                                                         
                                                                                                                                    
      COMMON/ICOM/IBUF(8),IDAT(18528)                                                                                               
      COMMON/JCOM/JBUF(8),JDAT(18528)                                                                                               
      DATA MAXX/18528/, MAXL/100/                                                                                                   
C---------------------------------------------------------------------                                                              
      NFF=5                                                                                                                         
      CALL JCLPNT(NFF,1,2,3,5,6)                                                                                                    
      GINV=1.E0/9.80616E0                                                                                                           
      REWIND 1                                                                                                                      
      REWIND 2                                                                                                                      
                                                                                                                                    
C     * READ MODE, PB, PT.                                                                                                          
                                                                                                                                    
      READ(5,5000,END=900)PB,PT,MODE                                            C4                                                  
      IF( INT(PB).EQ.0 ) THEN                                                                                                       
        OLD=.TRUE.                                                                                                                  
      ELSE                                                                                                                          
        OLD=.FALSE.                                                                                                                 
      ENDIF                                                                                                                         
      IF( .NOT.OLD.AND.PT.GT.PB ) CALL             XIT('VPINTV',-1)                                                                 
                                                                                                                                    
C     * READ THE NEXT SETS.                                                                                                         
                                                                                                                                    
      NSETS=0                                                                                                                       
  100 CALL GETSET2(1,GG,LEV,NLEV,IBUF,MAXX,OK)                                                                                      
      IF(.NOT.OK)THEN                                                                                                               
        IF(NSETS.EQ.0)THEN                                                                                                          
          WRITE(6,6010)                                                                                                             
          CALL                                     XIT('VPINTV',-2)                                                                 
        ENDIF                                                                                                                       
        WRITE(6,6020) NSETS                                                                                                         
        CALL                                       XIT('VPINTV',0)                                                                  
      ENDIF                                                                                                                         
      IF(NSETS.EQ.0)THEN                                                                                                            
        WRITE(6,6025) IBUF                                                                                                          
        WRITE(6,6005)IBUF(3),NLEV,(LEV(L),L=1,NLEV)                                                                                 
      ENDIF                                                                                                                         
      IF(NLEV.LT.2 .OR. NLEV.GT.MAXL) CALL         XIT('VPINTV',-3)                                                                 
      KIND=IBUF(1)                                                                                                                  
      IF(MODE.EQ.0)                         GO TO 110                                                                               
      IF(KIND.NE.NC4TO8("GRID").AND.                                                                                                
     +   KIND.NE.NC4TO8("ZONL")     ) CALL         XIT('VPINTV',-4)                                                                 
  110 NWDS=IBUF(5)*IBUF(6)                                                                                                          
      IF(KIND.EQ.NC4TO8("SPEC").OR.KIND.EQ.NC4TO8("FOUR"))  NWDS=NWDS*2                                                             
                                                                                                                                    
C     * INITIALISE PS AND DETERMINE PTOP, PBOT.                                                                                     
                                                                                                                                    
      IF(OLD)THEN                                                                                                                   
        CALL LVDCODE(PTOP2,LEV(1),1)                                                                                                
        PTOP=PTOP2/2.E0                                                                                                             
        PBOT=1013.3E0                                                                                                               
      ELSE                                                                                                                          
        PTOP=PT                                                                                                                     
        PBOT=PB                                                                                                                     
      ENDIF                                                                                                                         
      IF(NSETS.EQ.0) WRITE(6,6030) PTOP,PBOT                                                                                        
      DO 120 N=1,NWDS                                                                                                               
        PS(N)=PBOT                                                                                                                  
  120 CONTINUE                                                                                                                      
                                                                                                                                    
C     * READ THE SURFACE PRESSURE IF MODE .NE. 0 .                                                                                  
                                                                                                                                    
      IF( MODE.EQ.0 ) GO TO 140                                                                                                     
      IF ( NSETS.EQ.0 ) REWIND 3                                                                                                    
      CALL GETFLD2 (3,PS,-1,IBUF(2),NC4TO8("LNSP"),-1,JBUF,MAXX,OK)                                                                 
      IF(NSETS.EQ.0)      WRITE(6,6025) JBUF                                                                                        
      IF(KIND.NE.JBUF(1)) CALL                     XIT('VPINTV',-5)                                                                 
      IF( .NOT.OK )  CALL                          XIT('VPINTV',-6)                                                                 
      DO 130 N=1,NWDS                                                                                                               
        PS(N)   =EXP(PS(N))                                                                                                         
  130 CONTINUE                                                                                                                      
                                                                                                                                    
C     * COMPUTE PRESSURE HALF LEVELS FOR INTEGRAL IN P,                                                                             
C     * TAKING PT AND PB INTO ACCOUNT.                                                                                              
                                                                                                                                    
  140 PTOP=PTOP*100.E0                                                                                                              
      PBOT=PBOT*100.E0                                                                                                              
      CALL LVDCODE(PR,LEV,NLEV)                                                                                                     
      DO 150 L=1,NLEV                                                                                                               
        PR(L)=PR(L)*100.E0                                                                                                          
  150 CONTINUE                                                                                                                      
                                                                                                                                    
      PRH(1)=PR(1)*0.5E0                                                                                                            
      IF(PRH(1).LT.PTOP) PRH(1)=PTOP                                                                                                
      DO 160 L=2,NLEV                                                                                                               
        P=0.5E0*(PR(L-1)+PR(L))                                                                                                     
        IF(P.LE.PTOP) P=PTOP                                                                                                        
        IF(P.GE.PBOT) P=PBOT                                                                                                        
        PRH(L)=P                                                                                                                    
  160 CONTINUE                                                                                                                      
      PRH(NLEV+1)=PBOT                                                                                                              
                                                                                                                                    
      DO 170 N=1,NWDS                                                                                                               
  170 PS(N)=PS(N)*100.E0                                                                                                            
                                                                                                                                    
C     * INTEGRATE THE GRIDS IN THE VERTICAL AND DIVIDE BY G.                                                                        
                                                                                                                                    
      IF(MODE.NE.0) GO TO 200                                                                                                       
      DO 180 N=1,NWDS                                                                                                               
        VOUT(N)=0.E0                                                                                                                
  180 CONTINUE                                                                                                                      
      LCOUNT=0                                                                                                                      
      DO 190 L=1,NLEV                                                                                                               
        DP=GINV*(PRH(L+1)-PRH(L))                                                                                                   
        DO 190 N=1,NWDS                                                                                                             
          LCOUNT=LCOUNT+1                                                                                                           
          VOUT(N)=VOUT(N)+GG(LCOUNT)*DP                                                                                             
  190 CONTINUE                                                                                                                      
      GO TO 250                                                                                                                     
                                                                                                                                    
  200 DO 240 I=1,NWDS                                                                                                               
      VAL=0.E0                                                                                                                      
      DO 220 L=1,NLEV                                                                                                               
      N=(L-1)*NWDS+I                                                                                                                
      GGN=GG(N)                                                                                                                     
      PH=PRH(L+1)                                                                                                                   
      PL=PRH(L)                                                                                                                     
      IF(PS(I).LE.PRH(L))   GO TO 230                                                                                               
      IF(PS(I).GE.PRH(L+1)) GO TO 210                                                                                               
      PH=PS(I)                                                                                                                      
  210 VAL=VAL+GGN*(PH-PL)                                                                                                           
  220 CONTINUE                                                                                                                      
  230 VOUT(I)=VAL*GINV                                                                                                              
  240 CONTINUE                                                                                                                      
                                                                                                                                    
C     * PUT THE RESULT ONTO FILE 2.                                                                                                 
                                                                                                                                    
  250 IBUF(4)=1                                                                                                                     
      CALL PUTFLD2(2,VOUT,IBUF,MAXX)                                                                                                
      IF(NSETS.EQ.0) WRITE(6,6025) IBUF                                                                                             
      NSETS=NSETS+1                                                                                                                 
      GO TO 100                                                                                                                     
                                                                                                                                    
C     * E.O.F. ON INPUT.                                                                                                            
                                                                                                                                    
  900 CALL                                         XIT('VPINTV',-7)                                                                 
C---------------------------------------------------------------------                                                              
 5000 FORMAT(10X,2E10.0,5X,I5)                                                  C4                                                  
 6005 FORMAT('0 VPINTV ON ',A4,I5,' LEVELS =',20I5/(30X,20I5)/)                                                                     
 6010 FORMAT('0.. VPINTV INPUT GRID FILE IS EMPTY')                                                                                 
 6020 FORMAT('0',I6,' SETS WERE PROCESSED')                                                                                         
 6025 FORMAT(' ',A4,I10,2X,A4,I10,4I6)                                                                                              
 6030 FORMAT('0 INTEGRATION LIMITS =',E10.4,' ,',E10.4)                                                                             
      END                                                                                                                           
