! N. R. BADNELL  XTRCT            UoS v1.8                      08/01/19
!
      PROGRAM MAIN
C
      IMPLICIT REAL*8 (A-H,O-Z)
C
      LOGICAL BFORM
!
! EXTRACT OMEGA IS,JS OR KS
!
! DEFAULT ASSUMES NEUTRAL OMEGAS INCLUDE ELASTIC, IONS NOT.
! SET IS INDEX = -IS IF NEUTRALS DO NOT CONTAIN ELASTIC OR IONS DO.
!
! DEFAULT WRITES OMEGA FOR ALL CHANNELS OPEN ONLY.
! SET JS INDEX = -JS TO GET OMEGA FROM THRESHOLD.
!
      PARAMETER (ZERO=0.0D0)
!
      ALLOCATABLE :: ENAT(:),LAT(:),ISAT(:),OMEGA(:),NRD(:)
!
      call cpu_time(timei)
!
      INQUIRE (FILE='OMEGA',EXIST=BFORM)
      IF(BFORM)THEN                           !USE FORMATTED OMEGA
        WRITE(*,*)'USING FORMATTED OMEGA FILE'
        OPEN(7,FILE='OMEGA',STATUS='OLD',FORM='FORMATTED')   !INPUT FILE
        OPEN(9,STATUS='SCRATCH',FORM='FORMATTED')
      ELSE
        WRITE(*,*)'USING UNFORMATTED OMEGA FILE'
        OPEN(7,FILE='OMEGAU',STATUS='OLD',FORM='UNFORMATTED')!INPUT FILE
      ENDIF
!
      OPEN(8,FILE='xout',STATUS='UNKNOWN')          !OUTPUT ENERGY(RYD),
!                                                    COLLISION STRENGTH
      WRITE(*,*)'INPUT TRANSITION TO BE EXTRACTED: I,J OR K,0'
      READ(*,*)IS,JS
      KS=0
      IF(JS.EQ.0)KS=ABS(IS)
!
!  READ-IN BASIC DATA
!
      IF(BFORM)READ(7,*)NZED,NELC
      IF(.NOT.BFORM)READ(7)NZED,NELC
!
      TAZ=NZED-NELC
      IF(NZED.EQ.NELC)TAZ=1.0
      TAZ=TAZ*TAZ
!
      IF(BFORM)READ(7,*)NAST,MXE0,NOMWRT
      IF(.NOT.BFORM)READ(7)NAST,MXE0,NOMWRT
!
      IF(IS.LT.0)THEN
        IF(NZED.EQ.NELC)THEN
          NELC=-NELC
        ELSE
          NELC=NZED
        ENDIF
      ENDIF
!
      IF(NZED.EQ.NELC)THEN
        NOMT=(NAST*(NAST+1))/2
        IONE=0
        NOT=(NAST*(NAST-1))/2
      ELSE
        NOMT=(NAST*(NAST-1))/2
        IONE=1
        NOT=(NAST*(NAST+1))/2
      ENDIF
      IF(NOMWRT.NE.0)THEN
        NOMT=IABS(NOMWRT)
        IF(NOMT.EQ.NOT)                                                 &
     &  WRITE(*,*)'*** POSSIBLE USER FLAG ERROR OF ELASTIC TRANSITIONS?'
      ENDIF
!
      MXE=IABS(MXE0)                                    !CASE AS-DW FLAG
!
      ALLOCATE (ENAT(NAST),LAT(NAST),ISAT(NAST),OMEGA(NOMT),NRD(MXE)    &
     &         ,STAT=IERR)
      IF(IERR.NE.0)THEN                              
        WRITE(0,*)'ALLOCATION FAILS FOR NAST,NOMT,MXE=',NAST,NOMT,MXE
        STOP 'ABORT'                                         
      ENDIF                                          
!
      IF(BFORM)READ(7,*)(ISAT(I),LAT(I),I=1,NAST)
      IF(BFORM)READ(7,*)(ENAT(I),I=1,NAST)
      IF(.NOT.BFORM)READ(7)(ISAT(I),LAT(I),I=1,NAST)
      IF(.NOT.BFORM)READ(7)(ENAT(I),I=1,NAST)
!
      IF(NOMWRT.GE.0)THEN                            !AND CASE MXE0.LT.0
        DO I=1,MXE
          NRD(I)=NOMT
        ENDDO
        IF(MXE0.GT.0)THEN
          EINCR=ZER0
        ELSE
          EINCR=2.D-7
        ENDIF
      ELSE
        DO I=1,NOMT
          OMEGA(I)=ZERO
        ENDDO
        IF(BFORM)READ(7,*)E0                    !CHECK FIRST E
        IF(.NOT.BFORM)READ(7)E0                 !CHECK FIRST E
        BACKSPACE(7)
        I=1
        IF(E0.GE.ENAT(NAST))GO TO 1    !ALL OPEN - EARLY EXIT
        IF(BFORM)THEN
          I0=1
          NO=NOPEN(E0,ENAT,NAST,IONE,I0)
          NO=MIN(NO,NOMT)
          READ(7,*)E0,(DUM,N=1,NO)
          READ(7,*)E                   !=E0+EINCR
          EINCR=E-E0
          E0=E0-EINCR
          EINCH=EINCR*0.5
          NREC=2+(NO-1)/6
          DO N=1,NREC
            BACKSPACE(7)
          ENDDO
        ENDIF
        I0=1
        DO I=1,MXE
          IF(BFORM)THEN            !CONSTANT MESH
            E=E0+I*EINCR
            WRITE(9,700)E
            BACKSPACE(9)
            READ(9,700)E
            BACKSPACE(9)
          ELSE
            READ(7)E               !VARIABLE MESH
          ENDIF
          IF(E.GE.ENAT(NAST))GO TO 1
          NRD(I)=NOPEN(E,ENAT,NAST,IONE,I0)
          NRD(I)=MIN(NRD(I),NOMT)
        ENDDO
        I=MXE+1
   1    I0=I
        DO I=I0,MXE
          NRD(I)=NOMT
        ENDDO
        IF(.NOT.BFORM)THEN
          REWIND(7)
          NREC=4
          DO N=1,NREC
            READ(7)
          ENDDO
        ENDIF
      ENDIF
!
      IF(JS.LT.0)THEN
!
! ONLY WRITE OMEGA'S ABOVE THRESHOLD
!
        ETEST=(ENAT(-JS)-EINCR/2)*TAZ
        JJS=-1
!
      ELSE
!
! ONLY WRITE OMEGAS ABOVE ALL THRESHOLDS (NO RESONANCES)
!
        ETEST=ENAT(NAST)*TAZ
        JJS=1
!
      ENDIF
!
! DETERMINE TRANSITIONS CONTAINED IN OMEGA FILE
!
      IS=ABS(IS)
      JS=ABS(JS)
      DO I=1,NAST
        ENAT(I)=ENAT(I)*TAZ
      ENDDO
!
      IF(NOMWRT.GE.0)THEN
        N=0
        DO I=1,NAST
          DO J=I+IONE,NAST
            N=N+1
            IF(I.EQ.IS.AND.J.EQ.JS)KS=N
            IF(JS.EQ.0.AND.KS.EQ.N)THEN
              IS=I
              JS=J
              GO TO 5
            ENDIF
          ENDDO
        ENDDO
      ELSE
        N=0
        DO J=1+IONE,NAST
          DO I=1,J-IONE
            N=N+1
            IF(I.EQ.IS.AND.J.EQ.JS)KS=N
            IF(JS.EQ.0.AND.KS.EQ.N)THEN
              IS=I
              JS=J
              GO TO 5
            ENDIF
          ENDDO
        ENDDO
      ENDIF
!
   5  IF(KS.GT.NOMT)THEN
        WRITE(*,*)'REQUESTED TRANSITION NOT PRESENT IN OMEGA FILE'
        STOP
      ENDIF
!
! READ/WRITE OMEGA FILE
!
      WRITE(8,100)IS,JS,ENAT(JS)-ENAT(IS)
      IF(JJS.LT.0)WRITE(8,200)ENAT(JS),ENAT(JS)-ENAT(IS),ZERO,ZERO
      DO I=1,MXE
        IF(BFORM)THEN
          READ(7,*,END=999)E,(OMEGA(N),N=1,NRD(I)) !E RELATIVE TO GROUND
        ELSE
          IF(KS.LE.NRD(I))THEN
            READ(7,END=999)E,(OMEGA(N),N=1,KS)     !E RELATIVE TO GROUND
          ELSE
            READ(7,END=999)E
          ENDIF
        ENDIF
        IF(NOMWRT.LT.0.AND.BFORM)THEN
          ET=E0+I*EINCR
          IF(ABS(E-ET).GT.EINCH)THEN
            IF(NRD(I).LT.NOMT)THEN
              WRITE(*,*)'*** NON-CONSTANT ENERGY MESH?'
              WRITE(*,*)I,NRD(I),ET,E
              STOP'*** MIS-MATCH DURING READ OF FORMATTED OMEGAS'
            ENDIF
          ENDIF
        ENDIF
        IF(MXE0.GT.0)THEN                    !E RELATIVE TO GROUND
          IF(E.EQ.ZERO)E=ENAT(JS)/TAZ
          E1=E*TAZ                           !E RELATIVE TO GROUND
          E2=E1-ENAT(IS)                     !E RELATIVE TO INITIAL
          E3=E1-ENAT(JS)                     !E RELATIVE TO FINAL
        ELSE                                 !E RELATIVE TO FINAL
          E3=E*TAZ                           !E RELATIVE TO FINAL
          E1=E3+ENAT(JS)                     !E RELATIVE TO GROUND
          E2=E1-ENAT(IS)                     !E RELATIVE TO INITIAL
        ENDIF
        IF(E1.GE.ETEST)WRITE(8,200)E1,E2,E3,OMEGA(KS)
      ENDDO
!
   10 WRITE(*,*)'OUTPUT IN FILE xout'
!
      DEALLOCATE (ENAT,LAT,ISAT,OMEGA,NRD,STAT=IERR)
      IF(IERR.NE.0)THEN                              
        WRITE(0,*)'DE-ALLOCATION FAILS FOR NAST,NOMT,MXE...'
        STOP 'ABORT'                                         
      ENDIF                                          
!
      call cpu_time(timef)
!
      time=(timef-timei)/60.
      write(*,777)time
 777  format(' CPU TIME=',f9.3,' MIN')
!
      STOP 'NORMAL FINISH TO EXTRACTION'
!
 999  BACKSPACE (7)
      BACKSPACE (7)
      WRITE(*,300)I
      IF(BFORM)THEN
        READ(7,*,END=10)E
      ELSE
        READ(7,END=10)E
      ENDIF
      WRITE(*,310)E
      GO TO 10
!
 100  FORMAT('#',I3,' -',I4,' TRANSITION, ENERGY=',F12.6)
 200  FORMAT(1X,4(1PE14.6))
 300  FORMAT(/' UNEXPECTED END-OF-LINE/FILE AT ENERGY INDEX=',I6)
 310  FORMAT(37X,' ENERGY=',1PE16.8)
 700  FORMAT(1PE14.8,6(1PE11.3)/(14X,6(E11.3)))
!
      END
!
!************************************************
!
      INTEGER FUNCTION NOPEN(E,ENAT,NAST,IONE,I0)
      IMPLICIT REAL*8 (A-H,O-Z)
!
      DIMENSION ENAT(NAST)
!
      DO I=I0,NAST
        IF(E.LT.ENAT(I))GO TO 1
      ENDDO
      I=NAST+1
   1  I0=I-1
      NOPEN=(I0*(I0-2*IONE+1))/2
!
      RETURN
      END
