!-------------------------------------------------------------------
      SUBROUTINE FOF2INT(XNE,HTS,NCOL,HMF2OUT,FOF2OUT,IPRINT,I,J)
      DIMENSION XNE(NCOL),HTS(NCOL),XXIN(3),YYIN(3),FOF2IN(ncol)
!
      IF (IPRINT.EQ.1) THEN
        write(6,"(' fof2int: ncol=',i3)") ncol
        WRITE(6,"(' FOF2INT: HTS=',/(5F9.2))") HTS
        WRITE(6,"(' FOF2INT: XNE=',/(5E9.2))") XNE
      ENDIF
      DO 50 K=1,NCOL
 50     FOF2IN(K) = SQRT(XNE(K)/1.24E4)
      IF (IPRINT.EQ.1)
     +  WRITE(6,"(' FOF2INT: FOF2 FROM XNE=',/(5E9.2))") FOF2IN
      FOF2OUT = -1.E30
      DO 100 K=1,NCOL
        FOF2 = SQRT(XNE(K)/1.24E+4)
!
! 4/29/98: below change in conditional requested by Ganglu:
!       IF (FOF2.GT.FOF2OUT) THEN
        IF (FOF2.GT.FOF2OUT .and. K.GT.10) THEN
          FOF2OUT = FOF2
          KMHT = K
        ENDIF
 100  CONTINUE
!
! Occaisionally, XNE will increase steadily to top of column, in which
! case, we do no interp, and set fof2max and hmf2max at top of column.
! (This happened with ETS18A1, UT=0 (1st UT) at 7.5deg lat (j=20),
!  and 155-180 deg lon (i=68-73))
!
      IF (KMHT.EQ.NCOL) THEN
        FOF2OUT = SQRT(XNE(KMHT)/1.24E+4)
        HMF2OUT = HTS(KMHT)
!       WRITE(6,"(' NOTE FOF2INT: MAX FOF2 AT TOP OF COLUMN: I J=',2I3
!    +    )") I,J
!       WRITE(6,"(' >>> FOF2INT: MAX FOF2 AT TOP OF COLUMN: I J=',2I3,
!    +    ' XNE=',/(5E12.4))") I,J,XNE
!       WRITE(6,"(' I J=',2I3,' HTS=',/(5F12.4))") I,J,HTS
!       WRITE(6,"(' WILL USE VALUES AT TOP OF COLUMN: FOF2OUT=',E12.4,
!    +    ' HMF2OUT=',E12.4)") FOF2OUT,HMF2OUT
        RETURN
      ENDIF
      XXIN(1)=HTS(KMHT-1)
      XXIN(2)=HTS(KMHT)
      XXIN(3)=HTS(KMHT+1)
      YYIN(1)=SQRT(XNE(KMHT-1)/1.24E+4)
      YYIN(2)=SQRT(XNE(KMHT)/1.24E+4)
      YYIN(3)=SQRT(XNE(KMHT+1)/1.24E+4)

      if ((yyin(1)<=yyin(2)+1.e-10.and.yyin(1)>=yyin(2)-1.e-10).and.
     +    (yyin(1)<=yyin(3)+1.e-10.and.yyin(1)>=yyin(3)-1.e-10)) then
        fof2out = yyin(1)
        hmf2out = xxin(2)
      else
        CALL HINTPO(XXIN,YYIN,HMF2OUT,FOF2OUT)
      endif

      IF (IPRINT.EQ.1)
     +  WRITE(6,"('KMHT=',I3,' HMF2OUT=',F9.4,
     +  ' FOF2OUT=',F9.4,/'  XXIN=',3F9.4,' YYIN=',3F9.4)")
     +  KMHT,HMF2OUT,FOF2OUT,XXIN,YYIN
      RETURN
      END
      SUBROUTINE HINTPO(X,Y,XOUX,YOUX)
!
! Interp routine written by Mihail Codrescu:
!
      DIMENSION X(3),Y(3)
      X12=X(1)*X(1)
      X22=X(2)*X(2)
      X32=X(3)*X(3)
      DELTX=X12*X(2)+X22*X(3)+X32*X(1)-X32*X(2)-X12*X(3)-X22*X(1)
      ATX=Y(1)*X(2)+Y(2)*X(3)+Y(3)*X(1)-X(2)*Y(3)-X(3)*Y(1)-X(1)*Y(2)
      AX=ATX/DELTX
      BTX=X12*Y(2)+X22*Y(3)+X32*Y(1)-X32*Y(2)-X12*Y(3)-X22*Y(1)
      BX=BTX/DELTX
      CTX=X12*X(2)*Y(3)+X22*X(3)*Y(1)+X32*X(1)*Y(2)-
     *    X32*X(2)*Y(1)-X12*X(3)*Y(2)-X22*X(1)*Y(3)
      CX=CTX/DELTX
      XOUX=-(BX/(2.*AX))
      YOUX=-((BX*BX-4.*AX*CX)/(4.*AX))
      RETURN
      END
!-------------------------------------------------------------------

