SUBROUTINE FLOWV3(IFLAG,DLAT,DLON,RATIO)
      implicit none
C     ****
C     ****     TRANSFORM TO AURORAL CIRCLE COORDINATES
C     ****
      include "params.h"
      include "cflowv3.h"
      include "cons.h"
      integer :: istar,ihem
      real :: theta0,offa,offc,dskofa,dskofc,phid,phin,phidp0,
     |  phidm0,phinp0,phinm0,psim,psie,pcen,r1,dduumm,thetac
      COMMON /IONCR/ ISTAR,IHEM,THETA0(2),OFFA(2),OFFC(2),DSKOFA(2),
     1  DSKOFC(2),PHID(2),PHIN(2),PHIDP0(2),PHIDM0(2),PHINP0(2),
     2  PHINM0(2),PSIM(2),PSIE(2),PCEN(2),R1(2)
     3  ,DDUUMM(IMAXMP*JMX0+IMX0*JMX0+2)
C    3  ,THETAC(2),R2(2)       !   FOR OLD HEELIS POLAR CAP POTENTIAL
      real :: cosofa,sinofa,aslona,pi,pi2
      COMMON /OVALPOS/ COSOFA(2),SINOFA(2),ASLONA(2),PI,PI2
!
! Args:
      integer,intent(in) :: IFLAG(ZIMX)
      real,intent(in) :: DLAT(ZIMX),DLON(ZIMX),RATIO(ZIMX)
!
! Local:
      integer :: n,imaxd2,i
      real :: PSI(8),e,pih,ofda,ofdc,sinth0,sinthr1
      real :: PHDPMX(2),PHDMMX(2),PHNPMX(2),PHNMMX(2),
     1  COSOFC(2),SINOFC(2),ASLONC(2)
!
      DATA E/1.E-6/,PI/3.1415927/,PI2/6.2831853/,PIH/1.5707963/
C     ****
C     ****     CONVERT PARAMETERS TO RADIANS ETC.
      IF(ISTAR.EQ.0)THEN
        ISTAR=1
        PI=4.*ATAN(1.)
        PI2=2.*PI
        PIH=.5*PI
        DO 50 N=1,2
          OFDA = SQRT(OFFA(N)**2 + DSKOFA(N)**2)
          COSOFA(N) = COS(OFDA)
          SINOFA(N) = SIN(OFDA)
          ASLONA(N) = ASIN(DSKOFA(N)/OFDA)
          OFDC = SQRT(OFFC(N)**2 + DSKOFC(N)**2)
          COSOFC(N) = COS(OFDC)
          SINOFC(N) = SIN(OFDC)
          ASLONC(N) = ASIN(DSKOFC(N)/OFDC)
          IF(PHIN(N).LT.PHID(N))PHIN(N)=PHIN(N)+PI2
          PHDPMX(N)=.5*AMIN1(PI,(PHIN(N)-PHID(N)))
          PHNPMX(N)=.5*AMIN1(PI,(PHID(N)-PHIN(N)+PI2))
          PHNMMX(N)=PHDPMX(N)

!       write(6,"(/'flowv3: n=',i2,' offa(n)=',e12.4,' dskofa(n)=',
!    |    e12.4,' ofda=',e12.4)") n,offa(n),dskofa(n),ofda
!       write(6,"('  cosofa(n)=',e12.4)") cosofa(n)

   50     PHDMMX(N)=PHNPMX(N)
      ENDIF
C     ****               SET IHEM = 1,2  FOR  S,N HEMISPHERE
      IMAXD2 = MAX0(1,IMAX/2)
      IHEM = IFIX(DLAT(IMAXD2)*2./3.1416 + 2.)
      SINTH0=SIN(THETA0(IHEM))
C     ****
C     ****     AVERAGE AMIE RESULTS SHOW R1=-2.6 FOR 11.3 DEGREES
C     ****       (0.1972 RAD) BEYOND THETA0.
C     ****
      SINTHR1 = SIN(THETA0(IHEM)+0.1972)
C     G0=(SIN(THETAC(IHEM))/SINTH0)**R2(IHEM)    !   FOR OLD HEELIS
C     A2=1./((SIN(THETA0(IHEM)+THETAC(IHEM))/SINTH0)**R2(IHEM)-G0)
      PSI(1)=PSIE(IHEM)
      PSI(3)=PSIM(IHEM)
      DO 60 N=2,4,2
        PSI(N)=PSI(N-1)
   60 CONTINUE
      DO 70 N=1,4
        PSI(N+4)=PSI(N)
   70 CONTINUE
C     ****
C     ****     TRANSFORM TO AURORAL CIRCLE CORDINATES
C     ****       I.E. CALCULATE COLAT AND ALON FOR POLHT AND AURHT
C     ****
      DO 8 I=1,IMAX
        SINLAT(I)=SIN(ABS(DLAT(I)))
        COSLAT(I)=COS(DLAT(I))
        SINLON(I) = SIN(DLON(I)+ASLONA(IHEM))
        COSLON(I) = COS(DLON(I)+ASLONA(IHEM))
        COLAT(I) = COSOFA(IHEM)*SINLAT(I) - SINOFA(IHEM)*COSLAT(I)*
     1    COSLON(I)
        COLAT(I) = ACOS(COLAT(I)) * SQRT(RATIO(I))
        ALON(I) = AMOD(ATAN2(+SINLON(I)*COSLAT(I),SINLAT(I)*
     1    SINOFA(IHEM)+COSOFA(IHEM)*COSLAT(I)*COSLON(I)) -
     2    ASLONA(IHEM)+3.*PI,PI2) - PI

!       write(6,"(/'flowv3: i=',i2,' ihem=',i2,' dlat(i)=',e12.4,
!    |    ' dlon(i)=',e12.4)") i,ihem,dlat(i),dlon(i)
!       write(6,"('  aslona(ihem)=',e12.4,' cosofa(ihem)=',e12.4)")
!    |    aslona(ihem),cosofa(ihem)

    8 CONTINUE
      RETURN
      END
C