Index: Source/Fortran/changelog.txt
===================================================================
--- Source/Fortran/changelog.txt	(revision dfa7dbd2f70a8aa071dbb9f8b4928984af159da3)
+++ Source/Fortran/changelog.txt	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -23,3 +23,12 @@
      makefile_local_gfortran -> makefile_fast
    
-   
+2022-07-21 AT
+   solved issue with emoslib!
+   see ticket #312 on flexpart.eu
+
+   identified the source files from emoslib which calc_etadot needs
+   converted them from f77 to f90 code
+   adapted makefile to get rid of emoslib and make use of the fft source files
+
+   new files: set99.f90, fft99.f90, qpassm.f90, rpassm.f90, jsppole.f90
+      
Index: Source/Fortran/fft99.f90
===================================================================
--- Source/Fortran/fft99.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
+++ Source/Fortran/fft99.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -0,0 +1,218 @@
+!c Copyright 1981-2016 ECMWF.
+!c
+!c This software is licensed under the terms of the Apache Licence 
+!c Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+!c
+!c In applying this licence, ECMWF does not waive the privileges and immunities 
+!c granted to it by virtue of its status as an intergovernmental organisation 
+!c nor does it submit to any jurisdiction.
+!c
+
+      SUBROUTINE FFT99(A,WORK,TRIGS,IFAX,INC,JUMP,N,LOT,ISIGN)
+      DIMENSION A(N),WORK(N),TRIGS(N),IFAX(*)
+!C
+!C     SUBROUTINE 'FFT99' - MULTIPLE FAST REAL PERIODIC TRANSFORM
+!C     SUPERSEDES PREVIOUS ROUTINE 'FFT99'
+!C
+!C     REAL TRANSFORM OF LENGTH N PERFORMED BY REMOVING REDUNDANT
+!C     OPERATIONS FROM COMPLEX TRANSFORM OF LENGTH N
+!C
+!C     A IS THE ARRAY CONTAINING INPUT & OUTPUT DATA
+!C     WORK IS AN AREA OF SIZE (N+1)*MIN(LOT,64)
+!C     TRIGS IS A PREVIOUSLY PREPARED LIST OF TRIG FUNCTION VALUES
+!C     IFAX IS A PREVIOUSLY PREPARED LIST OF FACTORS OF N
+!C     INC IS THE INCREMENT WITHIN EACH DATA 'VECTOR'
+!C         (E.G. INC=1 FOR CONSECUTIVELY STORED DATA)
+!C     JUMP IS THE INCREMENT BETWEEN THE START OF EACH DATA VECTOR
+!C     N IS THE LENGTH OF THE DATA VECTORS
+!C     LOT IS THE NUMBER OF DATA VECTORS
+!C     ISIGN = +1 FOR TRANSFORM FROM SPECTRAL TO GRIDPOINT
+!C           = -1 FOR TRANSFORM FROM GRIDPOINT TO SPECTRAL
+!C
+!C     ORDERING OF COEFFICIENTS:
+!C         A(0),B(0),A(1),B(1),A(2),B(2),...,A(N/2),B(N/2)
+!C         WHERE B(0)=B(N/2)=0; (N+2) LOCATIONS REQUIRED
+!C
+!C     ORDERING OF DATA:
+!C         X(N-1),X(0),X(1),X(2),...,X(N-1),X(0)
+!C         I.E. EXPLICIT CYCLIC CONTINUITY; (N+2) LOCATIONS REQUIRED
+!C
+!C     VECTORIZATION IS ACHIEVED ON CRAY BY DOING THE TRANSFORMS
+!C     IN PARALLEL
+!C
+!C     N MUST BE COMPOSED OF FACTORS 2,3 & 5 BUT DOES NOT HAVE TO BE EVEN
+!C
+!C     DEFINITION OF TRANSFORMS:
+!C     -------------------------
+!C
+!C     ISIGN=+1: X(J)=SUM(K=0,...,N-1)(C(K)*EXP(2*I*J*K*PI/N))
+!C         WHERE C(K)=A(K)+I*B(K) AND C(N-K)=A(K)-I*B(K)
+!C
+!C     ISIGN=-1: A(K)=(1/N)*SUM(J=0,...,N-1)(X(J)*COS(2*J*K*PI/N))
+!C               B(K)=-(1/N)*SUM(J=0,...,N-1)(X(J)*SIN(2*J*K*PI/N))
+!C
+      IF(IFAX(10).NE.N) CALL SET99(TRIGS,IFAX,N)
+      NFAX=IFAX(1)
+      NX=N+1
+      IF (MOD(N,2).EQ.1) NX=N
+      NBLOX=1+(LOT-1)/64
+      NVEX=LOT-(NBLOX-1)*64
+      IF (ISIGN.EQ.-1) GO TO 300
+!C
+!C     ISIGN=+1, SPECTRAL TO GRIDPOINT TRANSFORM
+!C     -----------------------------------------
+  100 CONTINUE
+      ISTART=1
+      DO 220 NB=1,NBLOX
+      IA=ISTART
+      I=ISTART
+      DO 110 J=1,NVEX
+      A(I+INC)=0.5*A(I)
+      I=I+JUMP
+  110 CONTINUE
+      IF (MOD(N,2).EQ.1) GO TO 130
+      I=ISTART+N*INC
+      DO 120 J=1,NVEX
+      A(I)=0.5*A(I)
+      I=I+JUMP
+  120 CONTINUE
+  130 CONTINUE
+      IA=ISTART+INC
+      LA=1
+      IGO=+1
+!C
+      DO 160 K=1,NFAX
+      IFAC=IFAX(K+1)
+      IERR=-1
+      IF (IGO.EQ.-1) GO TO 140
+      CALL RPASSM(A(IA),A(IA+LA*INC),WORK(1),WORK(IFAC*LA+1),TRIGS,&
+         INC,1,JUMP,NX,NVEX,N,IFAC,LA,IERR)
+      GO TO 150
+  140 CONTINUE
+      CALL RPASSM(WORK(1),WORK(LA+1),A(IA),A(IA+IFAC*LA*INC),TRIGS,&
+         1,INC,NX,JUMP,NVEX,N,IFAC,LA,IERR)
+  150 CONTINUE
+      IF (IERR.NE.0) GO TO 500
+      LA=IFAC*LA
+      IGO=-IGO
+      IA=ISTART+INC
+  160 CONTINUE
+!C
+!C     IF NECESSARY, COPY RESULTS BACK TO A
+!C     ------------------------------------
+      IF (MOD(NFAX,2).EQ.0) GO TO 190
+      IBASE=1
+      JBASE=IA
+      DO 180 JJ=1,NVEX
+      I=IBASE
+      J=JBASE
+      DO 170 II=1,N
+      A(J)=WORK(I)
+      I=I+1
+      J=J+INC
+  170 CONTINUE
+      IBASE=IBASE+NX
+      JBASE=JBASE+JUMP
+  180 CONTINUE
+  190 CONTINUE
+!C
+!C     FILL IN CYCLIC BOUNDARY VALUES
+!C     ------------------------------
+      IX=ISTART
+      IZ=ISTART+N*INC
+!DIR$ IVDEP
+      DO 200 J=1,NVEX
+      A(IX)=A(IZ)
+      A(IZ+INC)=A(IX+INC)
+      IX=IX+JUMP
+      IZ=IZ+JUMP
+  200 CONTINUE
+!C
+      ISTART=ISTART+NVEX*JUMP
+      NVEX=64
+  220 CONTINUE
+      RETURN
+!C
+!C     ISIGN=-1, GRIDPOINT TO SPECTRAL TRANSFORM
+!C     -----------------------------------------
+  300 CONTINUE
+      ISTART=1
+      DO 410 NB=1,NBLOX
+      IA=ISTART+INC
+      LA=N
+      IGO=+1
+!C
+      DO 340 K=1,NFAX
+      IFAC=IFAX(NFAX+2-K)
+      LA=LA/IFAC
+      IERR=-1
+      IF (IGO.EQ.-1) GO TO 320
+      CALL QPASSM(A(IA),A(IA+IFAC*LA*INC),WORK(1),WORK(LA+1),TRIGS,&
+         INC,1,JUMP,NX,NVEX,N,IFAC,LA,IERR)
+      GO TO 330
+  320 CONTINUE
+      CALL QPASSM(WORK(1),WORK(IFAC*LA+1),A(IA),A(IA+LA*INC),TRIGS,&
+         1,INC,NX,JUMP,NVEX,N,IFAC,LA,IERR)
+  330 CONTINUE
+      IF (IERR.NE.0) GO TO 500
+      IGO=-IGO
+      IA=ISTART+INC
+  340 CONTINUE
+!C
+!C     IF NECESSARY, COPY RESULTS BACK TO A
+!C     ------------------------------------
+      IF (MOD(NFAX,2).EQ.0) GO TO 370
+      IBASE=1
+      JBASE=IA
+      DO 360 JJ=1,NVEX
+      I=IBASE
+      J=JBASE
+      DO 350 II=1,N
+      A(J)=WORK(I)
+      I=I+1
+      J=J+INC
+  350 CONTINUE
+      IBASE=IBASE+NX
+      JBASE=JBASE+JUMP
+  360 CONTINUE
+  370 CONTINUE
+!C
+!C     SHIFT A(0) & FILL IN ZERO IMAG PARTS
+!C     ------------------------------------
+      IX=ISTART
+      DO 380 J=1,NVEX
+      A(IX)=A(IX+INC)
+      A(IX+INC)=0.0
+      IX=IX+JUMP
+  380 CONTINUE
+      IF (MOD(N,2).EQ.1) GO TO 400
+      IZ=ISTART+(N+1)*INC
+      DO 390 J=1,NVEX
+      A(IZ)=0.0
+      IZ=IZ+JUMP
+  390 CONTINUE
+  400 CONTINUE
+!C
+      ISTART=ISTART+NVEX*JUMP
+      NVEX=64
+  410 CONTINUE
+      RETURN
+!C
+!C     ERROR MESSAGES
+!C     --------------
+  500 CONTINUE
+      GO TO (510,530,550) IERR
+  510 CONTINUE
+      WRITE(6,520) NVEX
+  520 FORMAT('1VECTOR LENGTH =',I4,', GREATER THAN 64')
+      GO TO 570
+  530 CONTINUE
+      WRITE(6,540) IFAC
+  540 FORMAT( '1FACTOR =',I3,', NOT CATERED FOR')
+      GO TO 570
+  550 CONTINUE
+      WRITE(6,560) IFAC
+  560 FORMAT('1FACTOR =',I3,', ONLY CATERED FOR IF LA*IFAC=N')
+  570 CONTINUE
+      RETURN
+      END
Index: Source/Fortran/jsppole.f90
===================================================================
--- Source/Fortran/jsppole.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
+++ Source/Fortran/jsppole.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -0,0 +1,182 @@
+! Copyright 1981-2016 ECMWF.
+!
+! This software is licensed under the terms of the Apache Licence 
+! Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+!
+! In applying this licence, ECMWF does not waive the privileges and immunities 
+! granted to it by virtue of its status as an intergovernmental organisation 
+! nor does it submit to any jurisdiction.
+!
+
+      SUBROUTINE JSPPOLE(PSHUP,KNUMB,KTRUNC,OMARS,PXF)
+      IMPLICIT NONE
+!
+!---->
+!**** *JSPPOLE* - Calculates fourier coefficient for U or V at pole
+!
+!     Purpose
+!     -------
+!
+!     Calculates fourier coefficient for first harmonic only
+!     for U and V wind component at the pole.
+!
+!     Interface
+!     ---------
+!
+!     CALL JSPPOLE(PSHUP,KNUMB,KTRUNC,OMARS,PXF)
+!
+!     Input parameters
+!     ----------------
+!
+!     PSHUP    - Unpacked harmonics field, unpacked
+!     KNUMB    - 1 for North Pole, otherwise South Pole
+!     KTRUNC   - Number (value) of the trucation 
+!     OMARS    - .TRUE. if data is from MARS 
+!     PXF      - Fourier coefficients (zero on input)
+!
+!
+!     Output parameters
+!     -----------------
+!
+!     PXF(2)   - Single fourier coefficient calculated
+!
+!
+!     Common block usage
+!     -----------------
+!
+!     None.
+!
+!
+!     Externals
+!     ---------
+!
+!     None.
+!
+!
+!     Author
+!     ------
+!
+!     J.D.Chambers     *ECMWF*      Oct 1993
+!
+!
+!     Modifications
+!     -------------
+!
+!     None.
+!
+!
+!     Comments
+!     --------
+!
+!     Created from SPPOLE. 
+!     Changed to provide all parameters in the call, i.e. no common
+!     blocks are used.
+!
+!
+!     Method
+!     ------
+!
+!     None.
+!
+!
+!     Reference
+!     _________
+!
+!     None.
+!
+!----<
+!     _______________________________________________________
+!
+!*    Section 0. Definition of variables.
+!     _______________________________________________________
+!
+!*    Prefix conventions for variable names
+!
+!     Logical      L (but not LP), global or common.
+!                  O, dummy argument
+!                  G, local variable
+!                  LP, parameter.
+!     Character    C, global or common.
+!                  H, dummy argument
+!                  Y (but not YP), local variable
+!                  YP, parameter.
+!     Integer      M and N, global or common.
+!                  K, dummy argument
+!                  I, local variable
+!                  J (but not JP), loop control
+!                  JP, parameter.
+!     REAL         A to F and Q to X, global or common.
+!                  P (but not PP), dummy argument
+!                  Z, local variable
+!                  PP, parameter.
+!
+!     Dummy arguments
+!
+      COMPLEX   PSHUP
+      INTEGER   KNUMB
+      INTEGER   KTRUNC
+      LOGICAL   OMARS
+      COMPLEX   PXF
+      DIMENSION PSHUP(*)
+      DIMENSION PXF(*)
+!
+!     Local variables
+!
+      INTEGER   I1, ITIN1, ITOUT1, JN
+      REAL      Z1, Z2, ZNORM, ZP1, ZP2, ZPOL
+! 
+!     -----------------------------------------------------------
+!
+!*    1.    Set initial values
+!           ------------------
+!
+ 100  CONTINUE
+!
+      ITIN1  = KTRUNC + 1
+      ITOUT1 = KTRUNC
+!
+      ZPOL = 1.
+      IF (KNUMB .NE. 1) ZPOL = -1.0
+!
+      ZP1  = -1.0
+      ZP2  = -3.0 * ZPOL
+      I1   = ITIN1 + 1
+!
+!*    2.    Change normalisation (if flagged as necessary)
+!           --------------------
+!
+ 200  CONTINUE
+!
+      IF (OMARS) THEN
+         ZNORM = -SQRT(2.0)
+      ELSE
+         ZNORM = 1
+      ENDIF
+!
+!
+!*    3.    Calculation
+!           -----------
+!
+ 300  CONTINUE
+      PXF(2) = (0.0,0.0)
+!
+!     Calculate the fourier coefficient for the first harmonic only.
+      DO 310 JN = 1,ITOUT1,2
+!
+        Z1 = SQRT( (2.0*JN + 1.0)/(2.0*JN*(JN + 1.0)) )
+        Z2 = SQRT( (2.0*(JN + 1.0) +1.0)/(2.0*(JN +1.0)*(JN +2.0)) )
+!
+        IF (JN .EQ. ITOUT1) Z2 = 0.0
+!
+        PXF(2) = PXF(2) +(Z1*ZP1*PSHUP(I1) +Z2*ZP2*PSHUP(I1+1))*ZNORM
+        ZP1   = ZP1 - 2.0*(JN + 1.0) - 1.0
+        ZP2   = ZP2 - (2.0*(JN + 2.0) + 1.0)*ZPOL
+        I1    = I1 + 2
+!
+ 310  CONTINUE
+!
+!     -------------------------------------------------------------
+!
+      RETURN
+!
+      END
Index: Source/Fortran/makefile_cray
===================================================================
--- Source/Fortran/makefile_cray	(revision 4d68c4a1835c900e4801f53e90e813553087073e)
+++ Source/Fortran/makefile_cray	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -15,5 +15,5 @@
 EXE      =  calc_etadot
 
-LIB     =  $(ECCODES_LIB) $(EMOSLIB)   
+LIB     =  $(ECCODES_LIB)   
 INC = -I. -I$(ECCODES_INCLUDE_DIR)
 
Index: Source/Fortran/makefile_debug
===================================================================
--- Source/Fortran/makefile_debug	(revision d2d6cf918895f938ac0288f15d1b1032d4323b03)
+++ Source/Fortran/makefile_debug	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -16,6 +16,5 @@
 
 ECCODES_LIB=  -Bstatic -leccodes_f90 -leccodes -Bdynamic -lm 
-EMOSLIB=-lemosR64
-LIB =  $(ECCODES_LIB) $(EMOSLIB)
+LIB =  $(ECCODES_LIB) 
 
 ECCODES_INCLUDE_DIR=/usr/lib/x86_64-linux-gnu/fortran/gfortran-mod-15
Index: Source/Fortran/makefile_ecgate
===================================================================
--- Source/Fortran/makefile_ecgate	(revision 1610f7325d2d44586c844cd85efb2b15a502abd3)
+++ Source/Fortran/makefile_ecgate	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -15,5 +15,5 @@
 EXE      =  calc_etadot_fast.out
 
-LIB     =  $(ECCODES_LIB) $(EMOSLIB)   
+LIB     =  $(ECCODES_LIB)    
 INC = -I. -I$(ECCODES_INCLUDE_DIR)
 
Index: Source/Fortran/makefile_fast
===================================================================
--- Source/Fortran/makefile_fast	(revision 14d96186895955cf5ee8d4076dafba8ab936bff6)
+++ Source/Fortran/makefile_fast	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -16,6 +16,5 @@
 
 ECCODES_LIB =  -Bstatic -leccodes_f90 -leccodes -Bdynamic -lm 
-EMOSLIB=-lemosR64
-LIB =  $(ECCODES_LIB) $(EMOSLIB)
+LIB =  $(ECCODES_LIB) 
 
 ECCODES_INCLUDE_DIR=/usr/lib/x86_64-linux-gnu/fortran/gfortran-mod-15
Index: Source/Fortran/qpassm.f90
===================================================================
--- Source/Fortran/qpassm.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
+++ Source/Fortran/qpassm.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -0,0 +1,778 @@
+!c Copyright 1981-2016 ECMWF.
+!c
+!c This software is licensed under the terms of the Apache Licence 
+!c Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+!c
+!c In applying this licence, ECMWF does not waive the privileges and immunities 
+!c granted to it by virtue of its status as an intergovernmental organisation 
+!c nor does it submit to any jurisdiction.
+!c
+
+      SUBROUTINE QPASSM(A,B,C,D,TRIGS,INC1,INC2,INC3,INC4,ILOT,N,IFAC,ILA,IERR)
+      REAL A(*),B(*),C(*),D(*),TRIGS(*)
+!C
+!C     SUBROUTINE 'QPASSM' - PERFORMS ONE PASS THROUGH DATA AS PART
+!C     OF MULTIPLE REAL FFT (FOURIER ANALYSIS) ROUTINE
+!C
+!C     A IS FIRST REAL INPUT VECTOR
+!C         EQUIVALENCE B(1) WITH A(IFAC*ILA*INC1+1)
+!C     C IS FIRST REAL OUTPUT VECTOR
+!C         EQUIVALENCE D(1) WITH C(ILA*INC2+1)
+!C     TRIGS IS A PRECALCULATED LIST OF SINES & COSINES
+!C     INC1 IS THE ADDRESSING INCREMENT FOR A
+!C     INC2 IS THE ADDRESSING INCREMENT FOR C
+!C     INC3 IS THE INCREMENT BETWEEN INPUT VECTORS A
+!C     INC4 IS THE INCREMENT BETWEEN OUTPUT VECTORS C
+!C     ILOT IS THE NUMBER OF VECTORS
+!C     N IS THE LENGTH OF THE VECTORS
+!C     IFAC IS THE CURRENT FACTOR OF N
+!C     ILA = N/(PRODUCT OF FACTORS USED SO FAR)
+!C     IERR IS AN ERROR INDICATOR:
+!C              0 - PASS COMPLETED WITHOUT ERROR
+!C              1 - ILOT GREATER THAN 64
+!C              2 - IFAC NOT CATERED FOR
+!C              3 - IFAC ONLY CATERED FOR IF ILA=N/IFAC
+!C
+!C-----------------------------------------------------------------------
+!C
+      SAVE SIN36, SIN72, QRT5, SIN60
+
+      DATA SIN36/0.587785252292473/,SIN72/0.951056516295154/,&
+           QRT5/0.559016994374947/,SIN60/0.866025403784437/
+
+      M=N/IFAC
+      IINK=ILA*INC1
+      JINK=ILA*INC2
+      IJUMP=(IFAC-1)*IINK
+      KSTOP=(N-IFAC)/(2*IFAC)
+
+      IBAD=1
+      IF (ILOT.GT.512) GO TO 910
+      IBASE=0
+      JBASE=0
+      IGO=IFAC-1
+      IF (IGO.EQ.7) IGO=6
+      IBAD=2
+      IF (IGO.LT.1.OR.IGO.GT.6) GO TO 910
+      GO TO (200,300,400,500,600,800),IGO
+
+!     CODING FOR FACTOR 2
+!     -------------------
+  200 CONTINUE
+      IA=1
+      IB=IA+IINK
+      JA=1
+      JB=JA+(2*M-ILA)*INC2
+
+      IF (ILA.EQ.M) GO TO 290
+
+      DO 220 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 210 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      C(JB+J)=A(IA+I)-A(IB+I)
+      I=I+INC3
+      J=J+INC4
+  210 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  220 CONTINUE
+      JA=JA+JINK
+      JINK=2*JINK
+      JB=JB-JINK
+      IBASE=IBASE+IJUMP
+      IJUMP=2*IJUMP+IINK
+      IF (JA.EQ.JB) GO TO 260
+      DO 250 K=ILA,KSTOP,ILA
+      KB=K+K
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      JBASE=0
+      DO 240 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 230 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+(C1*A(IB+I)+S1*B(IB+I))
+      C(JB+J)=A(IA+I)-(C1*A(IB+I)+S1*B(IB+I))
+      D(JA+J)=(C1*B(IB+I)-S1*A(IB+I))+B(IA+I)
+      D(JB+J)=(C1*B(IB+I)-S1*A(IB+I))-B(IA+I)
+      I=I+INC3
+      J=J+INC4
+  230 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  240 CONTINUE
+      IBASE=IBASE+IJUMP
+      JA=JA+JINK
+      JB=JB-JINK
+  250 CONTINUE
+      IF (JA.GT.JB) GO TO 900
+  260 CONTINUE
+      JBASE=0
+      DO 280 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 270 IJK=1,ILOT
+      C(JA+J)=A(IA+I)
+      D(JA+J)=-A(IB+I)
+      I=I+INC3
+      J=J+INC4
+  270 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  280 CONTINUE
+      GO TO 900
+
+  290 CONTINUE
+      Z=1.0/FLOAT(N)
+      DO 294 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 292 IJK=1,ILOT
+      C(JA+J)=Z*(A(IA+I)+A(IB+I))
+      C(JB+J)=Z*(A(IA+I)-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  292 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  294 CONTINUE
+      GO TO 900
+
+!     CODING FOR FACTOR 3
+!     -------------------
+  300 CONTINUE
+      IA=1
+      IB=IA+IINK
+      IC=IB+IINK
+      JA=1
+      JB=JA+(2*M-ILA)*INC2
+      JC=JB
+
+      IF (ILA.EQ.M) GO TO 390
+
+      DO 320 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 310 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+(A(IB+I)+A(IC+I))
+      C(JB+J)=A(IA+I)-0.5*(A(IB+I)+A(IC+I))
+      D(JB+J)=SIN60*(A(IC+I)-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  310 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  320 CONTINUE
+      JA=JA+JINK
+      JINK=2*JINK
+      JB=JB+JINK
+      JC=JC-JINK
+      IBASE=IBASE+IJUMP
+      IJUMP=2*IJUMP+IINK
+      IF (JA.EQ.JC) GO TO 360
+      DO 350 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      JBASE=0
+      DO 340 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 330 IJK=1,ILOT
+      A1=(C1*A(IB+I)+S1*B(IB+I))+(C2*A(IC+I)+S2*B(IC+I))
+      B1=(C1*B(IB+I)-S1*A(IB+I))+(C2*B(IC+I)-S2*A(IC+I))
+      A2=A(IA+I)-0.5*A1
+      B2=B(IA+I)-0.5*B1
+      A3=SIN60*((C1*A(IB+I)+S1*B(IB+I))-(C2*A(IC+I)+S2*B(IC+I)))
+      B3=SIN60*((C1*B(IB+I)-S1*A(IB+I))-(C2*B(IC+I)-S2*A(IC+I)))
+      C(JA+J)=A(IA+I)+A1
+      D(JA+J)=B(IA+I)+B1
+      C(JB+J)=A2+B3
+      D(JB+J)=B2-A3
+      C(JC+J)=A2-B3
+      D(JC+J)=-(B2+A3)
+      I=I+INC3
+      J=J+INC4
+  330 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  340 CONTINUE
+      IBASE=IBASE+IJUMP
+      JA=JA+JINK
+      JB=JB+JINK
+      JC=JC-JINK
+  350 CONTINUE
+      IF (JA.GT.JC) GO TO 900
+  360 CONTINUE
+      JBASE=0
+      DO 380 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 370 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+0.5*(A(IB+I)-A(IC+I))
+      D(JA+J)=-SIN60*(A(IB+I)+A(IC+I))
+      C(JB+J)=A(IA+I)-(A(IB+I)-A(IC+I))
+      I=I+INC3
+      J=J+INC4
+  370 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  380 CONTINUE
+      GO TO 900
+
+  390 CONTINUE
+      Z=1.0/FLOAT(N)
+      ZSIN60=Z*SIN60
+      DO 394 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 392 IJK=1,ILOT
+      C(JA+J)=Z*(A(IA+I)+(A(IB+I)+A(IC+I)))
+      C(JB+J)=Z*(A(IA+I)-0.5*(A(IB+I)+A(IC+I)))
+      D(JB+J)=ZSIN60*(A(IC+I)-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  392 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  394 CONTINUE
+      GO TO 900
+
+!     CODING FOR FACTOR 4
+!     -------------------
+  400 CONTINUE
+      IA=1
+      IB=IA+IINK
+      IC=IB+IINK
+      ID=IC+IINK
+      JA=1
+      JB=JA+(2*M-ILA)*INC2
+      JC=JB+2*M*INC2
+      JD=JB
+
+      IF (ILA.EQ.M) GO TO 490
+
+      DO 420 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 410 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+A(IC+I))+(A(IB+I)+A(ID+I))
+      C(JC+J)=(A(IA+I)+A(IC+I))-(A(IB+I)+A(ID+I))
+      C(JB+J)=A(IA+I)-A(IC+I)
+      D(JB+J)=A(ID+I)-A(IB+I)
+      I=I+INC3
+      J=J+INC4
+  410 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  420 CONTINUE
+      JA=JA+JINK
+      JINK=2*JINK
+      JB=JB+JINK
+      JC=JC-JINK
+      JD=JD-JINK
+      IBASE=IBASE+IJUMP
+      IJUMP=2*IJUMP+IINK
+      IF (JB.EQ.JC) GO TO 460
+      DO 450 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      JBASE=0
+      DO 440 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 430 IJK=1,ILOT
+      A0=A(IA+I)+(C2*A(IC+I)+S2*B(IC+I))
+      A2=A(IA+I)-(C2*A(IC+I)+S2*B(IC+I))
+      A1=(C1*A(IB+I)+S1*B(IB+I))+(C3*A(ID+I)+S3*B(ID+I))
+      A3=(C1*A(IB+I)+S1*B(IB+I))-(C3*A(ID+I)+S3*B(ID+I))
+      B0=B(IA+I)+(C2*B(IC+I)-S2*A(IC+I))
+      B2=B(IA+I)-(C2*B(IC+I)-S2*A(IC+I))
+      B1=(C1*B(IB+I)-S1*A(IB+I))+(C3*B(ID+I)-S3*A(ID+I))
+      B3=(C1*B(IB+I)-S1*A(IB+I))-(C3*B(ID+I)-S3*A(ID+I))
+      C(JA+J)=A0+A1
+      C(JC+J)=A0-A1
+      D(JA+J)=B0+B1
+      D(JC+J)=B1-B0
+      C(JB+J)=A2+B3
+      C(JD+J)=A2-B3
+      D(JB+J)=B2-A3
+      D(JD+J)=-(B2+A3)
+      I=I+INC3
+      J=J+INC4
+  430 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  440 CONTINUE
+      IBASE=IBASE+IJUMP
+      JA=JA+JINK
+      JB=JB+JINK
+      JC=JC-JINK
+      JD=JD-JINK
+  450 CONTINUE
+      IF (JB.GT.JC) GO TO 900
+  460 CONTINUE
+      SIN45=SQRT(0.5)
+      JBASE=0
+      DO 480 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 470 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+SIN45*(A(IB+I)-A(ID+I))
+      C(JB+J)=A(IA+I)-SIN45*(A(IB+I)-A(ID+I))
+      D(JA+J)=-A(IC+I)-SIN45*(A(IB+I)+A(ID+I))
+      D(JB+J)=A(IC+I)-SIN45*(A(IB+I)+A(ID+I))
+      I=I+INC3
+      J=J+INC4
+  470 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  480 CONTINUE
+      GO TO 900
+
+  490 CONTINUE
+      Z=1.0/FLOAT(N)
+      DO 494 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 492 IJK=1,ILOT
+      C(JA+J)=Z*((A(IA+I)+A(IC+I))+(A(IB+I)+A(ID+I)))
+      C(JC+J)=Z*((A(IA+I)+A(IC+I))-(A(IB+I)+A(ID+I)))
+      C(JB+J)=Z*(A(IA+I)-A(IC+I))
+      D(JB+J)=Z*(A(ID+I)-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  492 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  494 CONTINUE
+      GO TO 900
+
+!     CODING FOR FACTOR 5
+!     -------------------
+  500 CONTINUE
+      IA=1
+      IB=IA+IINK
+      IC=IB+IINK
+      ID=IC+IINK
+      IE=ID+IINK
+      JA=1
+      JB=JA+(2*M-ILA)*INC2
+      JC=JB+2*M*INC2
+      JD=JC
+      JE=JB
+
+      IF (ILA.EQ.M) GO TO 590
+
+      DO 520 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 510 IJK=1,ILOT
+      A1=A(IB+I)+A(IE+I)
+      A3=A(IB+I)-A(IE+I)
+      A2=A(IC+I)+A(ID+I)
+      A4=A(IC+I)-A(ID+I)
+      A5=A(IA+I)-0.25*(A1+A2)
+      A6=QRT5*(A1-A2)
+      C(JA+J)=A(IA+I)+(A1+A2)
+      C(JB+J)=A5+A6
+      C(JC+J)=A5-A6
+      D(JB+J)=-SIN72*A3-SIN36*A4
+      D(JC+J)=-SIN36*A3+SIN72*A4
+      I=I+INC3
+      J=J+INC4
+  510 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  520 CONTINUE
+      JA=JA+JINK
+      JINK=2*JINK
+      JB=JB+JINK
+      JC=JC+JINK
+      JD=JD-JINK
+      JE=JE-JINK
+      IBASE=IBASE+IJUMP
+      IJUMP=2*IJUMP+IINK
+      IF (JB.EQ.JD) GO TO 560
+      DO 550 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      KE=KD+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      C4=TRIGS(KE+1)
+      S4=TRIGS(KE+2)
+      JBASE=0
+      DO 540 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 530 IJK=1,ILOT
+      A1=(C1*A(IB+I)+S1*B(IB+I))+(C4*A(IE+I)+S4*B(IE+I))
+      A3=(C1*A(IB+I)+S1*B(IB+I))-(C4*A(IE+I)+S4*B(IE+I))
+      A2=(C2*A(IC+I)+S2*B(IC+I))+(C3*A(ID+I)+S3*B(ID+I))
+      A4=(C2*A(IC+I)+S2*B(IC+I))-(C3*A(ID+I)+S3*B(ID+I))
+      B1=(C1*B(IB+I)-S1*A(IB+I))+(C4*B(IE+I)-S4*A(IE+I))
+      B3=(C1*B(IB+I)-S1*A(IB+I))-(C4*B(IE+I)-S4*A(IE+I))
+      B2=(C2*B(IC+I)-S2*A(IC+I))+(C3*B(ID+I)-S3*A(ID+I))
+      B4=(C2*B(IC+I)-S2*A(IC+I))-(C3*B(ID+I)-S3*A(ID+I))
+      A5=A(IA+I)-0.25*(A1+A2)
+      A6=QRT5*(A1-A2)
+      B5=B(IA+I)-0.25*(B1+B2)
+      B6=QRT5*(B1-B2)
+      A10=A5+A6
+      A20=A5-A6
+      B10=B5+B6
+      B20=B5-B6
+      A11=SIN72*B3+SIN36*B4
+      A21=SIN36*B3-SIN72*B4
+      B11=SIN72*A3+SIN36*A4
+      B21=SIN36*A3-SIN72*A4
+      C(JA+J)=A(IA+I)+(A1+A2)
+      C(JB+J)=A10+A11
+      C(JE+J)=A10-A11
+      C(JC+J)=A20+A21
+      C(JD+J)=A20-A21
+      D(JA+J)=B(IA+I)+(B1+B2)
+      D(JB+J)=B10-B11
+      D(JE+J)=-(B10+B11)
+      D(JC+J)=B20-B21
+      D(JD+J)=-(B20+B21)
+      I=I+INC3
+      J=J+INC4
+  530 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  540 CONTINUE
+      IBASE=IBASE+IJUMP
+      JA=JA+JINK
+      JB=JB+JINK
+      JC=JC+JINK
+      JD=JD-JINK
+      JE=JE-JINK
+  550 CONTINUE
+      IF (JB.GT.JD) GO TO 900
+  560 CONTINUE
+      JBASE=0
+      DO 580 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 570 IJK=1,ILOT
+      A1=A(IB+I)+A(IE+I)
+      A3=A(IB+I)-A(IE+I)
+      A2=A(IC+I)+A(ID+I)
+      A4=A(IC+I)-A(ID+I)
+      A5=A(IA+I)+0.25*(A3-A4)
+      A6=QRT5*(A3+A4)
+      C(JA+J)=A5+A6
+      C(JB+J)=A5-A6
+      C(JC+J)=A(IA+I)-(A3-A4)
+      D(JA+J)=-SIN36*A1-SIN72*A2
+      D(JB+J)=-SIN72*A1+SIN36*A2
+      I=I+INC3
+      J=J+INC4
+  570 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  580 CONTINUE
+      GO TO 900
+
+  590 CONTINUE
+      Z=1.0/FLOAT(N)
+      ZQRT5=Z*QRT5
+      ZSIN36=Z*SIN36
+      ZSIN72=Z*SIN72
+      DO 594 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 592 IJK=1,ILOT
+      A1=A(IB+I)+A(IE+I)
+      A3=A(IB+I)-A(IE+I)
+      A2=A(IC+I)+A(ID+I)
+      A4=A(IC+I)-A(ID+I)
+      A5=Z*(A(IA+I)-0.25*(A1+A2))
+      A6=ZQRT5*(A1-A2)
+      C(JA+J)=Z*(A(IA+I)+(A1+A2))
+      C(JB+J)=A5+A6
+      C(JC+J)=A5-A6
+      D(JB+J)=-ZSIN72*A3-ZSIN36*A4
+      D(JC+J)=-ZSIN36*A3+ZSIN72*A4
+      I=I+INC3
+      J=J+INC4
+  592 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  594 CONTINUE
+      GO TO 900
+
+!     CODING FOR FACTOR 6
+!     -------------------
+  600 CONTINUE
+      IA=1
+      IB=IA+IINK
+      IC=IB+IINK
+      ID=IC+IINK
+      IE=ID+IINK
+      IF=IE+IINK
+      JA=1
+      JB=JA+(2*M-ILA)*INC2
+      JC=JB+2*M*INC2
+      JD=JC+2*M*INC2
+      JE=JC
+      JF=JB
+
+      IF (ILA.EQ.M) GO TO 690
+
+      DO 620 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 610 IJK=1,ILOT
+      A11=(A(IC+I)+A(IF+I))+(A(IB+I)+A(IE+I))
+      C(JA+J)=(A(IA+I)+A(ID+I))+A11
+      C(JC+J)=(A(IA+I)+A(ID+I)-0.5*A11)
+      D(JC+J)=SIN60*((A(IC+I)+A(IF+I))-(A(IB+I)+A(IE+I)))
+      A11=(A(IC+I)-A(IF+I))+(A(IE+I)-A(IB+I))
+      C(JB+J)=(A(IA+I)-A(ID+I))-0.5*A11
+      D(JB+J)=SIN60*((A(IE+I)-A(IB+I))-(A(IC+I)-A(IF+I)))
+      C(JD+J)=(A(IA+I)-A(ID+I))+A11
+      I=I+INC3
+      J=J+INC4
+  610 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  620 CONTINUE
+      JA=JA+JINK
+      JINK=2*JINK
+      JB=JB+JINK
+      JC=JC+JINK
+      JD=JD-JINK
+      JE=JE-JINK
+      JF=JF-JINK
+      IBASE=IBASE+IJUMP
+      IJUMP=2*IJUMP+IINK
+      IF (JC.EQ.JD) GO TO 660
+      DO 650 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      KE=KD+KB
+      KF=KE+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      C4=TRIGS(KE+1)
+      S4=TRIGS(KE+2)
+      C5=TRIGS(KF+1)
+      S5=TRIGS(KF+2)
+      JBASE=0
+      DO 640 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 630 IJK=1,ILOT
+      A1=C1*A(IB+I)+S1*B(IB+I)
+      B1=C1*B(IB+I)-S1*A(IB+I)
+      A2=C2*A(IC+I)+S2*B(IC+I)
+      B2=C2*B(IC+I)-S2*A(IC+I)
+      A3=C3*A(ID+I)+S3*B(ID+I)
+      B3=C3*B(ID+I)-S3*A(ID+I)
+      A4=C4*A(IE+I)+S4*B(IE+I)
+      B4=C4*B(IE+I)-S4*A(IE+I)
+      A5=C5*A(IF+I)+S5*B(IF+I)
+      B5=C5*B(IF+I)-S5*A(IF+I)
+      A11=(A2+A5)+(A1+A4)
+      A20=(A(IA+I)+A3)-0.5*A11
+      A21=SIN60*((A2+A5)-(A1+A4))
+      B11=(B2+B5)+(B1+B4)
+      B20=(B(IA+I)+B3)-0.5*B11
+      B21=SIN60*((B2+B5)-(B1+B4))
+      C(JA+J)=(A(IA+I)+A3)+A11
+      D(JA+J)=(B(IA+I)+B3)+B11
+      C(JC+J)=A20-B21
+      D(JC+J)=A21+B20
+      C(JE+J)=A20+B21
+      D(JE+J)=A21-B20
+      A11=(A2-A5)+(A4-A1)
+      A20=(A(IA+I)-A3)-0.5*A11
+      A21=SIN60*((A4-A1)-(A2-A5))
+      B11=(B5-B2)-(B4-B1)
+      B20=(B3-B(IA+I))-0.5*B11
+      B21=SIN60*((B5-B2)+(B4-B1))
+      C(JB+J)=A20-B21
+      D(JB+J)=A21-B20
+      C(JD+J)=A11+(A(IA+I)-A3)
+      D(JD+J)=B11+(B3-B(IA+I))
+      C(JF+J)=A20+B21
+      D(JF+J)=A21+B20
+      I=I+INC3
+      J=J+INC4
+  630 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  640 CONTINUE
+      IBASE=IBASE+IJUMP
+      JA=JA+JINK
+      JB=JB+JINK
+      JC=JC+JINK
+      JD=JD-JINK
+      JE=JE-JINK
+      JF=JF-JINK
+  650 CONTINUE
+      IF (JC.GT.JD) GO TO 900
+  660 CONTINUE
+      JBASE=0
+      DO 680 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 670 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+0.5*(A(IC+I)-A(IE+I)))+ SIN60*(A(IB+I)-A(IF+I))
+      D(JA+J)=-(A(ID+I)+0.5*(A(IB+I)+A(IF+I)))-SIN60*(A(IC+I)+A(IE+I))
+      C(JB+J)=A(IA+I)-(A(IC+I)-A(IE+I))
+      D(JB+J)=A(ID+I)-(A(IB+I)+A(IF+I))
+      C(JC+J)=(A(IA+I)+0.5*(A(IC+I)-A(IE+I)))-SIN60*(A(IB+I)-A(IF+I))
+      D(JC+J)=-(A(ID+I)+0.5*(A(IB+I)+A(IF+I)))+SIN60*(A(IC+I)+A(IE+I))
+      I=I+INC3
+      J=J+INC4
+  670 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  680 CONTINUE
+      GO TO 900
+
+  690 CONTINUE
+      Z=1.0/FLOAT(N)
+      ZSIN60=Z*SIN60
+      DO 694 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 692 IJK=1,ILOT
+      A11=(A(IC+I)+A(IF+I))+(A(IB+I)+A(IE+I))
+      C(JA+J)=Z*((A(IA+I)+A(ID+I))+A11)
+      C(JC+J)=Z*((A(IA+I)+A(ID+I))-0.5*A11)
+      D(JC+J)=ZSIN60*((A(IC+I)+A(IF+I))-(A(IB+I)+A(IE+I)))
+      A11=(A(IC+I)-A(IF+I))+(A(IE+I)-A(IB+I))
+      C(JB+J)=Z*((A(IA+I)-A(ID+I))-0.5*A11)
+      D(JB+J)=ZSIN60*((A(IE+I)-A(IB+I))-(A(IC+I)-A(IF+I)))
+      C(JD+J)=Z*((A(IA+I)-A(ID+I))+A11)
+      I=I+INC3
+      J=J+INC4
+  692 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  694 CONTINUE
+      GO TO 900
+
+!     CODING FOR FACTOR 8
+!     -------------------
+  800 CONTINUE
+      IBAD=3
+      IF (ILA.NE.M) GO TO 910
+      IA=1
+      IB=IA+IINK
+      IC=IB+IINK
+      ID=IC+IINK
+      IE=ID+IINK
+      IF=IE+IINK
+      IG=IF+IINK
+      IH=IG+IINK
+      JA=1
+      JB=JA+ILA*INC2
+      JC=JB+2*M*INC2
+      JD=JC+2*M*INC2
+      JE=JD+2*M*INC2
+      Z=1.0/FLOAT(N)
+      ZSIN45=Z*SQRT(0.5)
+
+      DO 820 JL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 810 IJK=1,ILOT
+      C(JA+J)=Z*(((A(IA+I)+A(IE+I))+(A(IC+I)+A(IG+I)))+&
+         ((A(ID+I)+A(IH+I))+(A(IB+I)+A(IF+I))))
+      C(JE+J)=Z*(((A(IA+I)+A(IE+I))+(A(IC+I)+A(IG+I)))-&
+         ((A(ID+I)+A(IH+I))+(A(IB+I)+A(IF+I))))
+      C(JC+J)=Z*((A(IA+I)+A(IE+I))-(A(IC+I)+A(IG+I)))
+      D(JC+J)=Z*((A(ID+I)+A(IH+I))-(A(IB+I)+A(IF+I)))
+      C(JB+J)=Z*(A(IA+I)-A(IE+I))&
+         +ZSIN45*((A(IH+I)-A(ID+I))-(A(IF+I)-A(IB+I)))
+      C(JD+J)=Z*(A(IA+I)-A(IE+I))&
+         -ZSIN45*((A(IH+I)-A(ID+I))-(A(IF+I)-A(IB+I)))
+      D(JB+J)=ZSIN45*((A(IH+I)-A(ID+I))+(A(IF+I)-A(IB+I)))&
+         +Z*(A(IG+I)-A(IC+I))
+      D(JD+J)=ZSIN45*((A(IH+I)-A(ID+I))+(A(IF+I)-A(IB+I)))&
+         -Z*(A(IG+I)-A(IC+I))
+      I=I+INC3
+      J=J+INC4
+  810 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  820 CONTINUE
+
+!     RETURN
+!     ------
+  900 CONTINUE
+      IBAD=0
+  910 CONTINUE
+      IERR=IBAD
+      RETURN
+      END
Index: Source/Fortran/rpassm.f90
===================================================================
--- Source/Fortran/rpassm.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
+++ Source/Fortran/rpassm.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -0,0 +1,760 @@
+!c Copyright 1981-2016 ECMWF.
+!c
+!c This software is licensed under the terms of the Apache Licence 
+!c Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+!c
+!c In applying this licence, ECMWF does not waive the privileges and immunities 
+!c granted to it by virtue of its status as an intergovernmental organisation 
+!c nor does it submit to any jurisdiction.
+!c
+
+      SUBROUTINE RPASSM(A,B,C,D,TRIGS,INC1,INC2,INC3,INC4,ILOT,N,IFAC,ILA,IERR)
+      REAL A(*),B(*),C(*),D(*),TRIGS(*)
+!C
+!C     SUBROUTINE 'RPASSM' - PERFORMS ONE PASS THROUGH DATA AS PART
+!C     OF MULTIPLE REAL FFT (FOURIER SYNTHESIS) ROUTINE
+!C
+!C     A IS FIRST REAL INPUT VECTOR
+!C         EQUIVALENCE B(1) WITH A (ILA*INC1+1)
+!C     C IS FIRST REAL OUTPUT VECTOR
+!C         EQUIVALENCE D(1) WITH C(IFAC*ILA*INC2+1)
+!C     TRIGS IS A PRECALCULATED LIST OF SINES & COSINES
+!C     INC1 IS THE ADDRESSING INCREMENT FOR A
+!C     INC2 IS THE ADDRESSING INCREMENT FOR C
+!C     INC3 IS THE INCREMENT BETWEEN INPUT VECTORS A
+!C     INC4 IS THE INCREMENT BETWEEN OUTPUT VECTORS C
+!C     ILOT IS THE NUMBER OF VECTORS
+!C     N IS THE LENGTH OF THE VECTORS
+!C     IFAC IS THE CURRENT FACTOR OF N
+!C     ILA IS THE PRODUCT OF PREVIOUS FACTORS
+!C     IERR IS AN ERROR INDICATOR:
+!C              0 - PASS COMPLETED WITHOUT ERROR
+!C              1 - ILOT GREATER THAN 64
+!C              2 - IFAC NOT CATERED FOR
+!C              3 - IFAC ONLY CATERED FOR IF ILA=N/IFAC
+!C
+!C-----------------------------------------------------------------------
+!C
+      REAL A10(512),A11(512),&
+           A20(512),A21(512),&
+           B10(512),B11(512),&
+           B20(512),B21(512)
+!C
+      SAVE SIN36, SIN72, QRT5, SIN60
+!C
+      DATA SIN36/0.587785252292473/,SIN72/0.951056516295154/,&
+           QRT5/0.559016994374947/,SIN60/0.866025403784437/
+!C
+      M=N/IFAC
+      IINK=ILA*INC1
+      JINK=ILA*INC2
+      JUMP=(IFAC-1)*JINK
+      KSTOP=(N-IFAC)/(2*IFAC)
+!C
+      IBAD=1
+      IF (ILOT.GT.512) GO TO 910
+      IBASE=0
+      JBASE=0
+      IGO=IFAC-1
+      IF (IGO.EQ.7) IGO=6
+      IBAD=2
+      IF (IGO.LT.1.OR.IGO.GT.6) GO TO 910
+      GO TO (200,300,400,500,600,800),IGO
+!C
+!C     CODING FOR FACTOR 2
+!C     -------------------
+  200 CONTINUE
+      IA=1
+      IB=IA+(2*M-ILA)*INC1
+      JA=1
+      JB=JA+JINK
+!C
+      IF (ILA.EQ.M) GO TO 290
+!C
+      DO 220 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 210 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      C(JB+J)=A(IA+I)-A(IB+I)
+      I=I+INC3
+      J=J+INC4
+  210 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  220 CONTINUE
+      IA=IA+IINK
+      IINK=2*IINK
+      IB=IB-IINK
+      IBASE=0
+      JBASE=JBASE+JUMP
+      JUMP=2*JUMP+JINK
+      IF (IA.EQ.IB) GO TO 260
+      DO 250 K=ILA,KSTOP,ILA
+      KB=K+K
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      IBASE=0
+      DO 240 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 230 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      D(JA+J)=B(IA+I)-B(IB+I)
+      C(JB+J)=C1*(A(IA+I)-A(IB+I))-S1*(B(IA+I)+B(IB+I))
+      D(JB+J)=S1*(A(IA+I)-A(IB+I))+C1*(B(IA+I)+B(IB+I))
+      I=I+INC3
+      J=J+INC4
+  230 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  240 CONTINUE
+      IA=IA+IINK
+      IB=IB-IINK
+      JBASE=JBASE+JUMP
+  250 CONTINUE
+      IF (IA.GT.IB) GO TO 900
+  260 CONTINUE
+      IBASE=0
+      DO 280 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 270 IJK=1,ILOT
+      C(JA+J)=A(IA+I)
+      C(JB+J)=-B(IA+I)
+      I=I+INC3
+      J=J+INC4
+  270 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  280 CONTINUE
+      GO TO 900
+!C
+  290 CONTINUE
+      DO 294 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 292 IJK=1,ILOT
+      C(JA+J)=2.0*(A(IA+I)+A(IB+I))
+      C(JB+J)=2.0*(A(IA+I)-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  292 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  294 CONTINUE
+      GO TO 900
+!C
+!C     CODING FOR FACTOR 3
+!C     -------------------
+  300 CONTINUE
+      IA=1
+      IB=IA+(2*M-ILA)*INC1
+      IC=IB
+      JA=1
+      JB=JA+JINK
+      JC=JB+JINK
+!C
+      IF (ILA.EQ.M) GO TO 390
+!C
+      DO 320 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 310 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      C(JB+J)=(A(IA+I)-0.5*A(IB+I))-(SIN60*(B(IB+I)))
+      C(JC+J)=(A(IA+I)-0.5*A(IB+I))+(SIN60*(B(IB+I)))
+      I=I+INC3
+      J=J+INC4
+  310 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  320 CONTINUE
+      IA=IA+IINK
+      IINK=2*IINK
+      IB=IB+IINK
+      IC=IC-IINK
+      JBASE=JBASE+JUMP
+      JUMP=2*JUMP+JINK
+      IF (IA.EQ.IC) GO TO 360
+      DO 350 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      IBASE=0
+      DO 340 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 330 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+(A(IB+I)+A(IC+I))
+      D(JA+J)=B(IA+I)+(B(IB+I)-B(IC+I))
+      C(JB+J)=C1*((A(IA+I)-0.5*(A(IB+I)+A(IC+I)))-(SIN60*(B(IB+I)+B(IC+I))))&
+        -S1*((B(IA+I)-0.5*(B(IB+I)-B(IC+I)))+(SIN60*(A(IB+I)-A(IC+I))))
+      D(JB+J)=S1*((A(IA+I)-0.5*(A(IB+I)+A(IC+I)))-(SIN60*(B(IB+I)+B(IC+I))))&
+        +C1*((B(IA+I)-0.5*(B(IB+I)-B(IC+I)))+(SIN60*(A(IB+I)-A(IC+I))))
+      C(JC+J)=C2*((A(IA+I)-0.5*(A(IB+I)+A(IC+I)))+(SIN60*(B(IB+I)+B(IC+I))))&
+        -S2*((B(IA+I)-0.5*(B(IB+I)-B(IC+I)))-(SIN60*(A(IB+I)-A(IC+I))))
+      D(JC+J)=S2*((A(IA+I)-0.5*(A(IB+I)+A(IC+I)))+(SIN60*(B(IB+I)+B(IC+I))))&
+        +C2*((B(IA+I)-0.5*(B(IB+I)-B(IC+I)))-(SIN60*(A(IB+I)-A(IC+I))))
+      I=I+INC3
+      J=J+INC4
+  330 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  340 CONTINUE
+      IA=IA+IINK
+      IB=IB+IINK
+      IC=IC-IINK
+      JBASE=JBASE+JUMP
+  350 CONTINUE
+      IF (IA.GT.IC) GO TO 900
+  360 CONTINUE
+      IBASE=0
+      DO 380 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 370 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      C(JB+J)=(0.5*A(IA+I)-A(IB+I))-(SIN60*B(IA+I))
+      C(JC+J)=-(0.5*A(IA+I)-A(IB+I))-(SIN60*B(IA+I))
+      I=I+INC3
+      J=J+INC4
+  370 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  380 CONTINUE
+      GO TO 900
+!C
+  390 CONTINUE
+      SSIN60=2.0*SIN60
+      DO 394 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 392 IJK=1,ILOT
+      C(JA+J)=2.0*(A(IA+I)+A(IB+I))
+      C(JB+J)=(2.0*A(IA+I)-A(IB+I))-(SSIN60*B(IB+I))
+      C(JC+J)=(2.0*A(IA+I)-A(IB+I))+(SSIN60*B(IB+I))
+      I=I+INC3
+      J=J+INC4
+  392 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  394 CONTINUE
+      GO TO 900
+!C
+!C     CODING FOR FACTOR 4
+!C     -------------------
+  400 CONTINUE
+      IA=1
+      IB=IA+(2*M-ILA)*INC1
+      IC=IB+2*M*INC1
+      ID=IB
+      JA=1
+      JB=JA+JINK
+      JC=JB+JINK
+      JD=JC+JINK
+!C
+      IF (ILA.EQ.M) GO TO 490
+!C
+      DO 420 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 410 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+A(IC+I))+A(IB+I)
+      C(JB+J)=(A(IA+I)-A(IC+I))-B(IB+I)
+      C(JC+J)=(A(IA+I)+A(IC+I))-A(IB+I)
+      C(JD+J)=(A(IA+I)-A(IC+I))+B(IB+I)
+      I=I+INC3
+      J=J+INC4
+  410 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  420 CONTINUE
+      IA=IA+IINK
+      IINK=2*IINK
+      IB=IB+IINK
+      IC=IC-IINK
+      ID=ID-IINK
+      JBASE=JBASE+JUMP
+      JUMP=2*JUMP+JINK
+      IF (IB.EQ.IC) GO TO 460
+      DO 450 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      IBASE=0
+      DO 440 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 430 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+A(IC+I))+(A(IB+I)+A(ID+I))
+      D(JA+J)=(B(IA+I)-B(IC+I))+(B(IB+I)-B(ID+I))
+      C(JC+J)=C2*((A(IA+I)+A(IC+I))-(A(IB+I)+A(ID+I)))&
+        -S2*((B(IA+I)-B(IC+I))-(B(IB+I)-B(ID+I)))
+      D(JC+J)=S2*((A(IA+I)+A(IC+I))-(A(IB+I)+A(ID+I)))&
+        +C2*((B(IA+I)-B(IC+I))-(B(IB+I)-B(ID+I)))
+      C(JB+J)=C1*((A(IA+I)-A(IC+I))-(B(IB+I)+B(ID+I)))&
+        -S1*((B(IA+I)+B(IC+I))+(A(IB+I)-A(ID+I)))
+      D(JB+J)=S1*((A(IA+I)-A(IC+I))-(B(IB+I)+B(ID+I)))&
+        +C1*((B(IA+I)+B(IC+I))+(A(IB+I)-A(ID+I)))
+      C(JD+J)=C3*((A(IA+I)-A(IC+I))+(B(IB+I)+B(ID+I)))&
+        -S3*((B(IA+I)+B(IC+I))-(A(IB+I)-A(ID+I)))
+      D(JD+J)=S3*((A(IA+I)-A(IC+I))+(B(IB+I)+B(ID+I)))&
+        +C3*((B(IA+I)+B(IC+I))-(A(IB+I)-A(ID+I)))
+      I=I+INC3
+      J=J+INC4
+  430 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  440 CONTINUE
+      IA=IA+IINK
+      IB=IB+IINK
+      IC=IC-IINK
+      ID=ID-IINK
+      JBASE=JBASE+JUMP
+  450 CONTINUE
+      IF (IB.GT.IC) GO TO 900
+  460 CONTINUE
+      IBASE=0
+      SIN45=SQRT(0.5)
+      DO 480 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 470 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+A(IB+I)
+      C(JB+J)=SIN45*((A(IA+I)-A(IB+I))-(B(IA+I)+B(IB+I)))
+      C(JC+J)=B(IB+I)-B(IA+I)
+      C(JD+J)=-SIN45*((A(IA+I)-A(IB+I))+(B(IA+I)+B(IB+I)))
+      I=I+INC3
+      J=J+INC4
+  470 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  480 CONTINUE
+      GO TO 900
+!C
+  490 CONTINUE
+      DO 494 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 492 IJK=1,ILOT
+      C(JA+J)=2.0*((A(IA+I)+A(IC+I))+A(IB+I))
+      C(JB+J)=2.0*((A(IA+I)-A(IC+I))-B(IB+I))
+      C(JC+J)=2.0*((A(IA+I)+A(IC+I))-A(IB+I))
+      C(JD+J)=2.0*((A(IA+I)-A(IC+I))+B(IB+I))
+      I=I+INC3
+      J=J+INC4
+  492 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  494 CONTINUE
+      GO TO 900
+!C
+!C     CODING FOR FACTOR 5
+!C     -------------------
+  500 CONTINUE
+      IA=1
+      IB=IA+(2*M-ILA)*INC1
+      IC=IB+2*M*INC1
+      ID=IC
+      IE=IB
+      JA=1
+      JB=JA+JINK
+      JC=JB+JINK
+      JD=JC+JINK
+      JE=JD+JINK
+!C
+      IF (ILA.EQ.M) GO TO 590
+!C
+      DO 520 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 510 IJK=1,ILOT
+      C(JA+J)=A(IA+I)+(A(IB+I)+A(IC+I))
+      C(JB+J)=((A(IA+I)-0.25*(A(IB+I)+A(IC+I)))+QRT5*(A(IB+I)-A(IC+I)))&
+         -(SIN72*B(IB+I)+SIN36*B(IC+I))
+      C(JC+J)=((A(IA+I)-0.25*(A(IB+I)+A(IC+I)))-QRT5*(A(IB+I)-A(IC+I)))&
+         -(SIN36*B(IB+I)-SIN72*B(IC+I))
+      C(JD+J)=((A(IA+I)-0.25*(A(IB+I)+A(IC+I)))-QRT5*(A(IB+I)-A(IC+I)))&
+         +(SIN36*B(IB+I)-SIN72*B(IC+I))
+      C(JE+J)=((A(IA+I)-0.25*(A(IB+I)+A(IC+I)))+QRT5*(A(IB+I)-A(IC+I)))&
+         +(SIN72*B(IB+I)+SIN36*B(IC+I))
+      I=I+INC3
+      J=J+INC4
+  510 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  520 CONTINUE
+      IA=IA+IINK
+      IINK=2*IINK
+      IB=IB+IINK
+      IC=IC+IINK
+      ID=ID-IINK
+      IE=IE-IINK
+      JBASE=JBASE+JUMP
+      JUMP=2*JUMP+JINK
+      IF (IB.EQ.ID) GO TO 560
+      DO 550 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      KE=KD+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      C4=TRIGS(KE+1)
+      S4=TRIGS(KE+2)
+      IBASE=0
+      DO 540 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 530 IJK=1,ILOT
+!C
+      A10(IJK)=(A(IA+I)-0.25*((A(IB+I)+A(IE+I))+(A(IC+I)+A(ID+I))))&
+         +QRT5*((A(IB+I)+A(IE+I))-(A(IC+I)+A(ID+I)))
+      A20(IJK)=(A(IA+I)-0.25*((A(IB+I)+A(IE+I))+(A(IC+I)+A(ID+I))))&
+         -QRT5*((A(IB+I)+A(IE+I))-(A(IC+I)+A(ID+I)))
+      B10(IJK)=(B(IA+I)-0.25*((B(IB+I)-B(IE+I))+(B(IC+I)-B(ID+I))))&
+         +QRT5*((B(IB+I)-B(IE+I))-(B(IC+I)-B(ID+I)))
+      B20(IJK)=(B(IA+I)-0.25*((B(IB+I)-B(IE+I))+(B(IC+I)-B(ID+I))))&
+         -QRT5*((B(IB+I)-B(IE+I))-(B(IC+I)-B(ID+I)))
+      A11(IJK)=SIN72*(B(IB+I)+B(IE+I))+SIN36*(B(IC+I)+B(ID+I))
+      A21(IJK)=SIN36*(B(IB+I)+B(IE+I))-SIN72*(B(IC+I)+B(ID+I))
+      B11(IJK)=SIN72*(A(IB+I)-A(IE+I))+SIN36*(A(IC+I)-A(ID+I))
+      B21(IJK)=SIN36*(A(IB+I)-A(IE+I))-SIN72*(A(IC+I)-A(ID+I))
+!C
+      C(JA+J)=A(IA+I)+((A(IB+I)+A(IE+I))+(A(IC+I)+A(ID+I)))
+      D(JA+J)=B(IA+I)+((B(IB+I)-B(IE+I))+(B(IC+I)-B(ID+I)))
+      C(JB+J)=C1*(A10(IJK)-A11(IJK))-S1*(B10(IJK)+B11(IJK))
+      D(JB+J)=S1*(A10(IJK)-A11(IJK))+C1*(B10(IJK)+B11(IJK))
+      C(JE+J)=C4*(A10(IJK)+A11(IJK))-S4*(B10(IJK)-B11(IJK))
+      D(JE+J)=S4*(A10(IJK)+A11(IJK))+C4*(B10(IJK)-B11(IJK))
+      C(JC+J)=C2*(A20(IJK)-A21(IJK))-S2*(B20(IJK)+B21(IJK))
+      D(JC+J)=S2*(A20(IJK)-A21(IJK))+C2*(B20(IJK)+B21(IJK))
+      C(JD+J)=C3*(A20(IJK)+A21(IJK))-S3*(B20(IJK)-B21(IJK))
+      D(JD+J)=S3*(A20(IJK)+A21(IJK))+C3*(B20(IJK)-B21(IJK))
+!C
+      I=I+INC3
+      J=J+INC4
+  530 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  540 CONTINUE
+      IA=IA+IINK
+      IB=IB+IINK
+      IC=IC+IINK
+      ID=ID-IINK
+      IE=IE-IINK
+      JBASE=JBASE+JUMP
+  550 CONTINUE
+      IF (IB.GT.ID) GO TO 900
+  560 CONTINUE
+      IBASE=0
+      DO 580 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 570 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+A(IB+I))+A(IC+I)
+      C(JB+J)=(QRT5*(A(IA+I)-A(IB+I))+(0.25*(A(IA+I)+A(IB+I))-A(IC+I)))&
+         -(SIN36*B(IA+I)+SIN72*B(IB+I))
+      C(JE+J)=-(QRT5*(A(IA+I)-A(IB+I))+(0.25*(A(IA+I)+A(IB+I))-A(IC+I)))&
+         -(SIN36*B(IA+I)+SIN72*B(IB+I))
+      C(JC+J)=(QRT5*(A(IA+I)-A(IB+I))-(0.25*(A(IA+I)+A(IB+I))-A(IC+I)))&
+         -(SIN72*B(IA+I)-SIN36*B(IB+I))
+      C(JD+J)=-(QRT5*(A(IA+I)-A(IB+I))-(0.25*(A(IA+I)+A(IB+I))-A(IC+I)))&
+         -(SIN72*B(IA+I)-SIN36*B(IB+I))
+      I=I+INC3
+      J=J+INC4
+  570 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  580 CONTINUE
+      GO TO 900
+!C
+  590 CONTINUE
+      QQRT5=2.0*QRT5
+      SSIN36=2.0*SIN36
+      SSIN72=2.0*SIN72
+      DO 594 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 592 IJK=1,ILOT
+      C(JA+J)=2.0*(A(IA+I)+(A(IB+I)+A(IC+I)))
+      C(JB+J)=(2.0*(A(IA+I)-0.25*(A(IB+I)+A(IC+I)))&
+         +QQRT5*(A(IB+I)-A(IC+I)))-(SSIN72*B(IB+I)+SSIN36*B(IC+I))
+      C(JC+J)=(2.0*(A(IA+I)-0.25*(A(IB+I)+A(IC+I)))&
+         -QQRT5*(A(IB+I)-A(IC+I)))-(SSIN36*B(IB+I)-SSIN72*B(IC+I))
+      C(JD+J)=(2.0*(A(IA+I)-0.25*(A(IB+I)+A(IC+I)))&
+         -QQRT5*(A(IB+I)-A(IC+I)))+(SSIN36*B(IB+I)-SSIN72*B(IC+I))
+      C(JE+J)=(2.0*(A(IA+I)-0.25*(A(IB+I)+A(IC+I)))&
+         +QQRT5*(A(IB+I)-A(IC+I)))+(SSIN72*B(IB+I)+SSIN36*B(IC+I))
+      I=I+INC3
+      J=J+INC4
+  592 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  594 CONTINUE
+      GO TO 900
+!C
+!C     CODING FOR FACTOR 6
+!C     -------------------
+  600 CONTINUE
+      IA=1
+      IB=IA+(2*M-ILA)*INC1
+      IC=IB+2*M*INC1
+      ID=IC+2*M*INC1
+      IE=IC
+      IF=IB
+      JA=1
+      JB=JA+JINK
+      JC=JB+JINK
+      JD=JC+JINK
+      JE=JD+JINK
+      JF=JE+JINK
+!C
+      IF (ILA.EQ.M) GO TO 690
+!C
+      DO 620 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 610 IJK=1,ILOT
+      C(JA+J)=(A(IA+I)+A(ID+I))+(A(IB+I)+A(IC+I))
+      C(JD+J)=(A(IA+I)-A(ID+I))-(A(IB+I)-A(IC+I))
+      C(JB+J)=((A(IA+I)-A(ID+I))+0.5*(A(IB+I)-A(IC+I)))&
+         -(SIN60*(B(IB+I)+B(IC+I)))
+      C(JF+J)=((A(IA+I)-A(ID+I))+0.5*(A(IB+I)-A(IC+I)))&
+         +(SIN60*(B(IB+I)+B(IC+I)))
+      C(JC+J)=((A(IA+I)+A(ID+I))-0.5*(A(IB+I)+A(IC+I)))&
+         -(SIN60*(B(IB+I)-B(IC+I)))
+      C(JE+J)=((A(IA+I)+A(ID+I))-0.5*(A(IB+I)+A(IC+I)))&
+         +(SIN60*(B(IB+I)-B(IC+I)))
+      I=I+INC3
+      J=J+INC4
+  610 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  620 CONTINUE
+      IA=IA+IINK
+      IINK=2*IINK
+      IB=IB+IINK
+      IC=IC+IINK
+      ID=ID-IINK
+      IE=IE-IINK
+      IF=IF-IINK
+      JBASE=JBASE+JUMP
+      JUMP=2*JUMP+JINK
+      IF (IC.EQ.ID) GO TO 660
+      DO 650 K=ILA,KSTOP,ILA
+      KB=K+K
+      KC=KB+KB
+      KD=KC+KB
+      KE=KD+KB
+      KF=KE+KB
+      C1=TRIGS(KB+1)
+      S1=TRIGS(KB+2)
+      C2=TRIGS(KC+1)
+      S2=TRIGS(KC+2)
+      C3=TRIGS(KD+1)
+      S3=TRIGS(KD+2)
+      C4=TRIGS(KE+1)
+      S4=TRIGS(KE+2)
+      C5=TRIGS(KF+1)
+      S5=TRIGS(KF+2)
+      IBASE=0
+      DO 640 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 630 IJK=1,ILOT
+!C
+      A11(IJK)= (A(IE+I)+A(IB+I))+(A(IC+I)+A(IF+I))
+      A20(IJK)=(A(IA+I)+A(ID+I))-0.5*A11(IJK)
+      A21(IJK)=SIN60*((A(IE+I)+A(IB+I))-(A(IC+I)+A(IF+I)))
+      B11(IJK)= (B(IB+I)-B(IE+I))+(B(IC+I)-B(IF+I))
+      B20(IJK)=(B(IA+I)-B(ID+I))-0.5*B11(IJK)
+      B21(IJK)=SIN60*((B(IB+I)-B(IE+I))-(B(IC+I)-B(IF+I)))
+!C
+      C(JA+J)=(A(IA+I)+A(ID+I))+A11(IJK)
+      D(JA+J)=(B(IA+I)-B(ID+I))+B11(IJK)
+      C(JC+J)=C2*(A20(IJK)-B21(IJK))-S2*(B20(IJK)+A21(IJK))
+      D(JC+J)=S2*(A20(IJK)-B21(IJK))+C2*(B20(IJK)+A21(IJK))
+      C(JE+J)=C4*(A20(IJK)+B21(IJK))-S4*(B20(IJK)-A21(IJK))
+      D(JE+J)=S4*(A20(IJK)+B21(IJK))+C4*(B20(IJK)-A21(IJK))
+!C
+      A11(IJK)=(A(IE+I)-A(IB+I))+(A(IC+I)-A(IF+I))
+      B11(IJK)=(B(IE+I)+B(IB+I))-(B(IC+I)+B(IF+I))
+      A20(IJK)=(A(IA+I)-A(ID+I))-0.5*A11(IJK)
+      A21(IJK)=SIN60*((A(IE+I)-A(IB+I))-(A(IC+I)-A(IF+I)))
+      B20(IJK)=(B(IA+I)+B(ID+I))+0.5*B11(IJK)
+      B21(IJK)=SIN60*((B(IE+I)+B(IB+I))+(B(IC+I)+B(IF+I)))
+!C
+      C(JD+J)=C3*((A(IA+I)-A(ID+I))+A11(IJK))-S3*((B(IA+I)+B(ID+I))-B11(IJK))
+      D(JD+J)=S3*((A(IA+I)-A(ID+I))+A11(IJK))+C3*((B(IA+I)+B(ID+I))-B11(IJK))
+      C(JB+J)=C1*(A20(IJK)-B21(IJK))-S1*(B20(IJK)-A21(IJK))
+      D(JB+J)=S1*(A20(IJK)-B21(IJK))+C1*(B20(IJK)-A21(IJK))
+      C(JF+J)=C5*(A20(IJK)+B21(IJK))-S5*(B20(IJK)+A21(IJK))
+      D(JF+J)=S5*(A20(IJK)+B21(IJK))+C5*(B20(IJK)+A21(IJK))
+!C
+      I=I+INC3
+      J=J+INC4
+  630 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  640 CONTINUE
+      IA=IA+IINK
+      IB=IB+IINK
+      IC=IC+IINK
+      ID=ID-IINK
+      IE=IE-IINK
+      IF=IF-IINK
+      JBASE=JBASE+JUMP
+  650 CONTINUE
+      IF (IC.GT.ID) GO TO 900
+  660 CONTINUE
+      IBASE=0
+      DO 680 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 670 IJK=1,ILOT
+      C(JA+J)=A(IB+I)+(A(IA+I)+A(IC+I))
+      C(JD+J)=B(IB+I)-(B(IA+I)+B(IC+I))
+      C(JB+J)=(SIN60*(A(IA+I)-A(IC+I)))-(0.5*(B(IA+I)+B(IC+I))+B(IB+I))
+      C(JF+J)=-(SIN60*(A(IA+I)-A(IC+I)))-(0.5*(B(IA+I)+B(IC+I))+B(IB+I))
+      C(JC+J)=SIN60*(B(IC+I)-B(IA+I))+(0.5*(A(IA+I)+A(IC+I))-A(IB+I))
+      C(JE+J)=SIN60*(B(IC+I)-B(IA+I))-(0.5*(A(IA+I)+A(IC+I))-A(IB+I))
+      I=I+INC3
+      J=J+INC4
+  670 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  680 CONTINUE
+      GO TO 900
+!C
+  690 CONTINUE
+      SSIN60=2.0*SIN60
+      DO 694 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 692 IJK=1,ILOT
+      C(JA+J)=(2.0*(A(IA+I)+A(ID+I)))+(2.0*(A(IB+I)+A(IC+I)))
+      C(JD+J)=(2.0*(A(IA+I)-A(ID+I)))-(2.0*(A(IB+I)-A(IC+I)))
+      C(JB+J)=(2.0*(A(IA+I)-A(ID+I))+(A(IB+I)-A(IC+I)))&
+         -(SSIN60*(B(IB+I)+B(IC+I)))
+      C(JF+J)=(2.0*(A(IA+I)-A(ID+I))+(A(IB+I)-A(IC+I)))&
+         +(SSIN60*(B(IB+I)+B(IC+I)))
+      C(JC+J)=(2.0*(A(IA+I)+A(ID+I))-(A(IB+I)+A(IC+I)))&
+         -(SSIN60*(B(IB+I)-B(IC+I)))
+      C(JE+J)=(2.0*(A(IA+I)+A(ID+I))-(A(IB+I)+A(IC+I)))&
+         +(SSIN60*(B(IB+I)-B(IC+I)))
+      I=I+INC3
+      J=J+INC4
+  692 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  694 CONTINUE
+      GO TO 900
+!C
+!C     CODING FOR FACTOR 8
+!C     -------------------
+  800 CONTINUE
+      IBAD=3
+      IF (ILA.NE.M) GO TO 910
+      IA=1
+      IB=IA+ILA*INC1
+      IC=IB+2*ILA*INC1
+      ID=IC+2*ILA*INC1
+      IE=ID+2*ILA*INC1
+      JA=1
+      JB=JA+JINK
+      JC=JB+JINK
+      JD=JC+JINK
+      JE=JD+JINK
+      JF=JE+JINK
+      JG=JF+JINK
+      JH=JG+JINK
+      SSIN45=SQRT(2.0)
+!C
+      DO 820 IL=1,ILA
+      I=IBASE
+      J=JBASE
+!DIR$ IVDEP
+!VOCL LOOP,NOVREC
+      DO 810 IJK=1,ILOT
+      C(JA+J)=2.0*(((A(IA+I)+A(IE+I))+A(IC+I))+(A(IB+I)+A(ID+I)))
+      C(JE+J)=2.0*(((A(IA+I)+A(IE+I))+A(IC+I))-(A(IB+I)+A(ID+I)))
+      C(JC+J)=2.0*(((A(IA+I)+A(IE+I))-A(IC+I))-(B(IB+I)-B(ID+I)))
+      C(JG+J)=2.0*(((A(IA+I)+A(IE+I))-A(IC+I))+(B(IB+I)-B(ID+I)))
+      C(JB+J)=2.0*((A(IA+I)-A(IE+I))-B(IC+I))+SSIN45&
+          *((A(IB+I)-A(ID+I))-(B(IB+I)+B(ID+I)))
+      C(JF+J)=2.0*((A(IA+I)-A(IE+I))-B(IC+I))&
+         -SSIN45*((A(IB+I)-A(ID+I))-(B(IB+I)+B(ID+I)))
+      C(JD+J)=2.0*((A(IA+I)-A(IE+I))+B(IC+I))&
+         -SSIN45*((A(IB+I)-A(ID+I))+(B(IB+I)+B(ID+I)))
+      C(JH+J)=2.0*((A(IA+I)-A(IE+I))+B(IC+I))&
+         +SSIN45*((A(IB+I)-A(ID+I))+(B(IB+I)+B(ID+I)))
+      I=I+INC3
+      J=J+INC4
+  810 CONTINUE
+      IBASE=IBASE+INC1
+      JBASE=JBASE+INC2
+  820 CONTINUE
+!C
+!C     RETURN
+!C     ------
+  900 CONTINUE
+      IBAD=0
+  910 CONTINUE
+      IERR=IBAD
+      RETURN
+      END
Index: Source/Fortran/set99.f90
===================================================================
--- Source/Fortran/set99.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
+++ Source/Fortran/set99.f90	(revision 4f247989d67e9b1876b8e79ae2b81dd9450d0d47)
@@ -0,0 +1,65 @@
+!c Copyright 1981-2016 ECMWF.
+!c
+!c This software is licensed under the terms of the Apache Licence 
+!c Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+!c
+!c In applying this licence, ECMWF does not waive the privileges and immunities 
+!c granted to it by virtue of its status as an intergovernmental organisation 
+!c nor does it submit to any jurisdiction.
+!c
+
+      SUBROUTINE SET99(TRIGS,IFAX,N)
+      DIMENSION TRIGS(N),IFAX(*),JFAX(10),LFAX(7)
+!C
+!C     SUBROUTINE 'SET99' - COMPUTES FACTORS OF N & TRIGONOMETRIC
+!C     FUNCTIONS REQUIRED BY FFT99 & FFT991
+!C
+      DATA LFAX/6,8,5,4,3,2,1/
+      IXXX=1
+!C
+      DEL=4.0*ASIN(1.0)/FLOAT(N)
+      NIL=0
+      NHL=(N/2)-1
+      DO 10 K=NIL,NHL
+      ANGLE=FLOAT(K)*DEL
+      TRIGS(2*K+1)=COS(ANGLE)
+      TRIGS(2*K+2)=SIN(ANGLE)
+   10 CONTINUE
+!C
+!C     FIND FACTORS OF N (8,6,5,4,3,2; ONLY ONE 8 ALLOWED)
+!C     LOOK FOR SIXES FIRST, STORE FACTORS IN DESCENDING ORDER
+      NU=N
+      IFAC=6
+      K=0
+      L=1
+   20 CONTINUE
+      IF (MOD(NU,IFAC).NE.0) GO TO 30
+      K=K+1
+      JFAX(K)=IFAC
+      IF (IFAC.NE.8) GO TO 25
+      IF (K.EQ.1) GO TO 25
+      JFAX(1)=8
+      JFAX(K)=6
+   25 CONTINUE
+      NU=NU/IFAC
+      IF (NU.EQ.1) GO TO 50
+      IF (IFAC.NE.8) GO TO 20
+   30 CONTINUE
+      L=L+1
+      IFAC=LFAX(L)
+      IF (IFAC.GT.1) GO TO 20
+!C
+      WRITE(6,40) N
+   40 FORMAT('1N =',I4,' - CONTAINS ILLEGAL FACTORS')
+      RETURN
+!C
+!C     NOW REVERSE ORDER OF FACTORS
+   50 CONTINUE
+      NFAX=K
+      IFAX(1)=NFAX
+      DO 60 I=1,NFAX
+      IFAX(NFAX+2-I)=JFAX(I)
+   60 CONTINUE
+      IFAX(10)=N
+      RETURN
+      END
