#include "dims.h"
      SUBROUTINE POTM
      use cons_module,only: pi
      implicit none
C     ****
C     ****     CALCULATE HEELIS POTENTIAL IN GEOMAGNETIC COORDINATES
C     ****
#include "params.h"
#include "consdyn.h"
#include "dynphi.h"
C     ****
C     ****     WORK SPACE
C     ****
      integer :: iflag,ifn
      real :: dlat,dlon,ratio,psi,phi,phifun,colat,alon,sinlat,coslat,
     |  sinlon,coslon,wk1,wk2,wk3
      COMMON/flwv32_com/IFLAG(IMAXM), DLAT(IMAXM), DLON(IMAXM), 
     +  RATIO(IMAXM), PSI(8), PHI(IMAXM,8), IFN(IMAXM), PHIFUN(IMAXM), 
     +  COLAT(IMAXM), ALON(IMAXM), SINLAT(IMAXM),
     +  COSLAT(IMAXM), SINLON(IMAXM), COSLON(IMAXM), WK1(IMAXM),
     +  WK2(IMAXM),WK3(IMAXM)
!
! Local:
      integer :: i,j
C     ****
C     ****     /TRGM/ CONTAINS MAGNETIC INFO
! (dlons is updated at each iter by sun)
C     ****
      real :: rlatm,rlonm,dipmag,decmag,sndec,csdec,sn2dec,sncsdc,
     |  rlatmp,rlonmp,dlons
      COMMON/TRGM/RLATM(ZIMXP,ZJMX),RLONM(ZIMXP,ZJMX),
     1  DIPMAG(ZIMXP,ZJMX),DECMAG(ZIMXP,ZJMX),SNDEC(ZIMXP,ZJMX),
     2  CSDEC(ZIMXP,ZJMX),SN2DEC(ZIMXP,ZJMX),SNCSDC(ZIMXP,ZJMX),
     3  RLATMP(ZJMX),RLONMP(ZJMX),DLONS(ZJMX)
C     ****
C     ****     LOOP OVER MAGNETIC LATITUDE
C     ****
      DO 1 J = 1,JMAXM
C     ****
C     ****     FILL ARRAYS IFLAG, DLAT, DLON, RATIO
C     ****
        DO 2 I = 1,IMAXM
          IFLAG(I) = 1
          DLAT(I) = YLATM(J)
          DLON(I) = YLONM(I)-DLONS(1)
          RATIO(I) = 1.
    2   CONTINUE
C       ****
C       ****     CALCULATE HEELIS POTENTIAL IN PHIHM(IMAXMP,JMAXM)
C       ****
        IF(ABS(YLATM(J)).GT.pi/6.)THEN
          CALL FLWV32(PHIHM(1,J))
        ELSE
          DO 3 I = 1,IMAXM
            PHIHM(I,J) = 0.
    3     CONTINUE
        ENDIF
    1 CONTINUE
C     ****
C     ****    PERIODIC POINTS
C     ****
      DO 4 J = 1,JMAXM
        PHIHM(IMAXMP,J) = PHIHM(1,J)
    4 CONTINUE
      RETURN
      END