!
!************************************************************************
!************************************************************************
!*                                                                     **
!*  PSS.F95    --- THIS PROGRAM TESTS THE POOLE-SOWELL-SPEAR           **
!*                 PROJECTION METHOD OF DETERMINING DIMENSIONALITY     **
!*                 USING DW-NOMINATE SCORES                            **
!*                                                                     **
!*                                                                     **
!*                           13 DECEMBER 2003                          **
!*                                                                     **
!*               FORTRAN 95  29 SEPTEMBER 2015                         **
!*                                                                     **
!************************************************************************
!************************************************************************
!
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      DIMENSION XLOG(200,20),PNTMAS(500,10),XLIKE(125000,7),&
                KROW(125000),KCOL(125000),KK(200)
!     
      OPEN(UNIT=20,FILE='tab2to8.dat',status='old')
      OPEN(UNIT=21,FILE='output21.dat')
      OPEN(UNIT=23,FILE='pntmass.dat',status='old')
      OPEN(UNIT=24,FILE='C:\DWNOMINATE_ISSUE_CODES\H01_113_CODES.TXT',&
                         STATUS='OLD')
  100 format(2i4,8f12.6)
  105 format(5x,9f8.6)
  200 FORMAT(I3,I5,30X,4I4)
  300 FORMAT(' NUMBER OF SCALABLE ROLL CALLS',I8)
  302 format(' Number of Rows in the Likelihood Matrix=',i8)
  303 format(' Number of Rows in the Point Mass Matrix=',i8)
  355 FORMAT(' NUMBER IN H01_113_CODES.TXT',2I7,7F12.4)
 1000 FORMAT(' START TIME OF JOB    ',3F12.7)
 1001 FORMAT(' END TIME OF JOB      ',3F12.7)
 1002 FORMAT(' ELAPSED TIME OF JOB  ',3F12.7)
!
      call cpu_time(result)
      write(*,1000)result
!
!
      NDIM=7
!
      DO 121 JJ=1,200
      DO 121 J=1,20
      KK(JJ)=0
      XLOG(JJ,J)=0.0
 121  CONTINUE
!
!  READ IN LIKELIHOOD TABLE
!
      I=0
  500 READ(20,100,END=400)KROW(I+1),KCOL(I+1),(XLIKE(I+1,J),J=1,NDIM)
      I=I+1
      GO TO 500
  400 WRITE(21,302)I
      write(*,302)i
!
!  READ IN POINT-MASS TABLE
!
      I=0
  525 READ(23,105,END=425)zulu,(PNTMAS(I,J),J=1,NDIM)
      I=I+1
      GO TO 525
  425 WRITE(21,303)I
      write(*,303)I
!
!  READ H01_113_CODES.TXT -- DW-NOMINATE CLASSIFICATIONS AND ISSUE CODES FILE
!
      I=0
  555 READ(24,200,END=455)ICONG,IRC,KYY,KYN,KNY,KNN
!      IF(ICONG.NE.113)GO TO 555
      I=I+1
      KYES=KYY+KNY
      KNO=KYN+KNN
      IF((KYES+KNO).EQ.0)GO TO 555
      KERR=KNY+KYN
      XERR=FLOAT(KNY+KYN)/FLOAT(KYES+KNO)
      IF(KNO.LT.KYES)THEN
         XMARG=FLOAT(KNO)/FLOAT(KYES+KNO)
      ENDIF
      IF(KNO.GE.KYES)THEN
         XMARG=FLOAT(KYES)/FLOAT(KYES+KNO)
      ENDIF
      IF(XMARG.GE.0.025)THEN
         KK(ICONG)=KK(ICONG)+1
      ENDIF
      IF(XMARG.LT.0.025)GO TO 555
!
!   KMARG IS THE MINORITY SIDE OF THE ROLL CALL -- M IN POOLE, SOWELL, AND SPEAR
!   KERR IS THE NUMBER OF ERRORS                -- E IN POOLE, SOWELL, AND SPEAR
! 
!       KMARG IS THE FIRST COLUMN OF TAB2TO8.DAT
!       KERR < KMARG IS THE SECOND COLUMN OF TAB2TO8.DAT
!
!       IF KMARG = KERR THEN USE PNTMASS.DAT WHERE THE ROW IS FLOAT(M)
!
      KMARG=INT(1000.0*XMARG)
      KERR=INT(1000.0*XERR)
      KEITH=(((KMARG-1)*KMARG)/2 - 300) + KERR + 1
      IF(ICONG.EQ.1)WRITE(21,1111)ICONG,I,KK(ICONG),KMARG,KEITH,KERR,&
                       XMARG,XERR
 1111 FORMAT(6I8,2F10.4)
 2222 FORMAT(6I8,F12.4)
      DO 90 JJJ=1,NDIM
      if(KMARG.ne.kerr)AA=XLIKE(KEITH,JJJ)
      IF(ABS(KMARG-KERR).LE..1)AA=PNTMAS(KMARG-24,JJJ)
      IF(AA.LE..0001)AA=0.0001
      XLOG(ICONG,JJJ)=XLOG(ICONG,JJJ)+DLOG(AA)
      XLOG(ICONG,JJJ+7)=XLOG(ICONG,JJJ+7)+AA
      IF(ICONG.EQ.1)WRITE(21,2222)ICONG,JJJ,KEITH,KMARG-24,KMARG,KERR,AA
  90  CONTINUE
      GO TO 555
  455 WRITE(*,300)I
!
      DO 999 ICONG=1,113
      WRITE(*,355)ICONG,KK(ICONG),&
             ((XLOG(ICONG,JJJ+7)/FLOAT(KK(ICONG))),JJJ=1,7)
      WRITE(21,355)ICONG,KK(ICONG),&
             ((XLOG(ICONG,JJJ+7)/FLOAT(KK(ICONG))),JJJ=1,7)
  999 CONTINUE
!
      call cpu_time(result1)
      write(*,1001)result1
      write(21,1001)result1
      TOTALTIME=RESULT1-RESULT
      WRITE(*,1002)TOTALTIME
      WRITE(21,1002)TOTALTIME
!
      STOP
      END
