      SUBROUTINE VFAST(F,nlevf,nlat,ilh,
     1                  S,nlevs,la,
     2                  WLP,iram,
     3                  ilev,LSR,LM)
C
C     * OCT 25/93 - M.LAZARE. MOVE DEFINITION OF "na3" IN "NON-LEFTOVER"
C     *                       CODE AT END TO INSIDE LOOP 270, SO THAT
C     *                       IT GETS RE-INITIALIZAED TO "NA" FOR EACH
C     *                       LOOP ITTERATION OVER "N". THIS WAS A BUG. 
C     * OCT 04/93 - M.LAZARE. DECLARE "WLP" TO BE REAL*8, AS SHOULD BE.
C     * OCT 30/92 - J. STACEY ("FAST" EXTENSIVELY REWRITTEN FOR PERFORMANCE.
C     *                        MULTIPLE LATITUDES, UNIT-STRIDE INPUT DATA)
C     * JUL 14/92 - E. CHAN. PREVIOUS VERSION "FAST".
C     *
C     * MULTIPLE LEGENDRE TRANSFORMS FROM FOURIER TO SPECTRAL.
C     * ACCUMULATE THE CONTRIBUTIONS OF FOURIER COEFFICIENTS
C     * ON ONE LATITUDE CIRCLE INTO THE LATITUDE INTEGRAL 
C     * TO OBTAIN SPHERICAL HARMONIC COEFFICIENTS.
C     *
C     * S   = ACCUMULATED SPHERICAL HARMONIC COEFFICIENTS,
C     * F   = FOURIER COEFFICIENTS AT ONE LATITUDE, 
C     * WLP = LEGENDRE POLYNOMIALS NORMALIZED BY GAUSSIAN WEIGHT, 
C     * LM  = NUMBER OF E-W WAVE NUMBERS (M=1,LM),
C     * ILH = FIRST DIMENSION OF COMPLEX F, 
C     * LA  = FIRST DIMENSION OF COMPLEX S, 
C     * ILEV=NUMBER OF LEVELS FOR WHICH TRANSFORM IS PERFORMED. 
C 
      IMPLICIT REAL (A-H,O-Z),
     +INTEGER (I-N)

C-------------------------- ORIGINAL CODE ------------------------------
C
C     REAL S(2,LA,ILEV),F(2,ILH,ILEV) 
C     REAL*8 WLP(1) 
C     INTEGER LSR(2,1)
C
C     DO 500 M=1,LM 
C     NA=LSR(2,M) 
C     NL=LSR(1,M) 
C     NR=LSR(1,M+1)-1 
C 
C     IF(ILEV.GT.NR-NL+1) THEN
C 
C        DO 200 N=NL,NR 
C        WP=WLP(NA) 
C           DO 100 L=1,ILEV 
C           S(1,N,L)=S(1,N,L)+WP*F(1,M,L) 
C           S(2,N,L)=S(2,N,L)+WP*F(2,M,L) 
C 100       CONTINUE
C        NA=NA+1
C 200    CONTINUE 
C 
C     ELSE
C 
C        DO 400 L=1,ILEV
C        NB=NA
C        F1=F(1,M,L)
C        F2=F(2,M,L)
C           DO 300 N=NL,NR
C           S(1,N,L)=S(1,N,L) +WLP(NB)*F1 
C           S(2,N,L)=S(2,N,L) +WLP(NB)*F2 
C           NB=NB+1 
C 300    CONTINUE 
C 400    CONTINUE 
C 
C     ENDIF 
C 
C 500 CONTINUE
C     RETURN
C-----------------------------------------------------------------------
C 
      DIMENSION S(2*nlevs,LA),F(2*nlevf,nlat,ILH)
      REAL*8 WLP(iram,nlat)
      DIMENSION LSR(2,LM)
C-----------------------------------------------------------------------
C     * OPTIMIZED CODE TO DO THE LM MATRIX MULTIPLIES.
C
      DO 400 M=1,LM
         NA = LSR(2,M) 
         NL = LSR(1,M) 
         NR = LSR(1,M+1)-1 
         ntot  = nlat
C
C        * First pass, initialization.
C  
         nleft = mod(nlat,8)
         ntot2 = ntot - nleft
         nlat2 = nleft + 1
         na2   = na

         if (nleft .eq. 1) then
           do 110 n=nl,nr
             do 105 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
  105        continue
             na2 = na2 + 1
  110      continue

         else if (nleft .eq. 2) then
           do 120 n=nl,nr
             do 115 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
  115        continue
             na2 = na2 + 1
  120      continue

         else if (nleft .eq. 3) then
           do 130 n=nl,nr
             do 125 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
     2                         + wlp(na2,3)*f(l,3,m)
  125        continue
             na2 = na2 + 1
  130      continue

         else if (nleft .eq. 4) then
           do 140 n=nl,nr
             do 135 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
     2                         + wlp(na2,3)*f(l,3,m)
     3                         + wlp(na2,4)*f(l,4,m)
  135        continue
             na2 = na2 + 1
  140      continue

         else if (nleft .eq. 5) then
           do 150 n=nl,nr
             do 145 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
     2                         + wlp(na2,3)*f(l,3,m)
     3                         + wlp(na2,4)*f(l,4,m)
     4                         + wlp(na2,5)*f(l,5,m)
  145        continue
             na2 = na2 + 1
  150      continue

         else if (nleft .eq. 6) then
           do 160 n=nl,nr
             do 155 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
     2                         + wlp(na2,3)*f(l,3,m)
     3                         + wlp(na2,4)*f(l,4,m)
     4                         + wlp(na2,5)*f(l,5,m)
     5                         + wlp(na2,6)*f(l,6,m)
  155        continue
             na2 = na2 + 1
  160      continue

         else if (nleft .eq. 7) then
           do 170 n=nl,nr
             do 165 l=1,2*ilev
               s(l,n) = s(l,n) + wlp(na2,1)*f(l,1,m)
     1                         + wlp(na2,2)*f(l,2,m)
     2                         + wlp(na2,3)*f(l,3,m)
     3                         + wlp(na2,4)*f(l,4,m)
     4                         + wlp(na2,5)*f(l,5,m)
     5                         + wlp(na2,6)*f(l,6,m)
     6                         + wlp(na2,7)*f(l,7,m)
  165        continue
             na2 = na2 + 1
  170      continue

         end if
         nlat3 = nlat2
C 
C        * Subsequent passes, use combined computational loop.
C        * This loop is designed to increase the computational 
C        * intensity for multipipe machines (SX-3, CRAY Y-MP/EL).
C
         ntim  = ntot2/8
         if (ntim .ge. 1) then
           do 270 lat=nlat3,nlat,8
             na3 = na
             do 260 n=nl,nr
               do 250 l=1,2*ilev
                 s(l,n) = s(l,n) + wlp(na3,lat  )*f(l,lat  ,m)
     1                           + wlp(na3,lat+1)*f(l,lat+1,m)
     2                           + wlp(na3,lat+2)*f(l,lat+2,m)
     3                           + wlp(na3,lat+3)*f(l,lat+3,m)
     4                           + wlp(na3,lat+4)*f(l,lat+4,m)
     5                           + wlp(na3,lat+5)*f(l,lat+5,m)
     6                           + wlp(na3,lat+6)*f(l,lat+6,m)
     7                           + wlp(na3,lat+7)*f(l,lat+7,m)
  250          continue
               na3 = na3 + 1
  260        continue
  270      continue
         end if

  400 CONTINUE
C
      RETURN 
      END 
