program State_By_State_Sim

!     Program to produce state-by-state simulation results
!
!
!     First two data cards contain title, a72/a72.
!     Third data card contains KAPPA1, KAPPA2, INPUT1, INPUT2, NJOB, 5i10.
!       KAPPA1 and KAPPA2 control subroutine VOTCOR.
!       If INPUT1 is less than, equal, or greater than 1, the TPRIMES
!          are read in on cards, input tape 9, or tape 9,
!          respectively, according to variable format FORMT1.
!       If Input2 is less than, equal, or greater than 1, the Kennedy
!          and Nixon votes for each cell are read in on cards, input
!          tape 10, or tape 10, respectively, according to FORMT2.
!       The program reads in data, computes and writes output for
!          NJOB sets of Kennedy, Nixon votes.  If NJOB is less than 1,
!          it is set equal to 1.
!     Fourth data card contains the variable format FORMT1, a72.
!     Fifth data card contains the variable format FORMT2, a72.
!     Next four data cards contain VKACT, actual votes for Kennedy,
!       for the 48 states, 12f6.5.
!     Next four data cards contain weights for 48 states, 12f6.5.
!     Next four data cards contain names of the states, 12a6
!     If INPUT1 is less than 1, deck of TPRIMES follows, FORMT1.
!     If INPUT2 is less than 1, NJOB decks of Kennedy, Nixon votes for
!       each cell follow, FORMT2.

  use VotCor_m, only: VotCor

  implicit NONE

  real :: ADJ
  logical :: BadInput     ! LCOUNT is not zero and not -111
  real :: CSHIFT
  character(72) :: FORMT1 ! Input format for TPRIMES
  character(72) :: FORMT2 ! Input format Kennedy, Nixon votes
  integer :: I, IALL, II(1628), III(480), INPUT1, INPUT2, IS, IV
  integer :: IKN             ! 1 = Kennedy, 2 = Nixon
  integer, parameter :: ItemTest(5) = [ 5, 9, 11, 12, 11 ]
  integer :: J, JJ(1628), JJJ(480)
  integer :: K
  !  KAPPA1, KAPPA2 control correlation calculation and printing:
  integer :: KAPPA1 ! < 2 -> All 48 states including south
                    ! = 2 -> 32 Non-southern or border states
                    ! > 2 -> Both
  integer :: KAPPA2 ! < 1 -> Each weighted equally
                    ! = 1 -> Each state weighted by size
                    ! > 1 -> Both
  integer :: KK(1628), KKK(480)
  character, parameter :: KN(2) = [ "K", "N" ] ! Kennedy or Nixon
  integer :: L, LLL(480)
  integer :: LCOUNT(5,10,10) ! For testing correct reading of VOTEK, VOTEN,
                             ! LCOUNT is either -111 (no input), or zero
                             ! (inputs correctly indexed).
                             ! Later, LCOUNT(:,i,j) should be ItemTest(i).
  integer :: M, MCOUNT(48), MM(1628)
  character(6) :: NAME(48)
  integer :: NJOB
  real :: SUMS(48,4)
  character(72) :: Title(2)
  real :: TPRIME(1628,3)
  real :: VK(48,4), VKACT(48), VKNORM(48,4) ! Kennedy
  real :: VKN(48,4,2)                       ! Both
  real :: VN(48,4), VNNORM(48,4)            ! Nixon
  equivalence ( VKN(1,1,1), VK ), ( VKN(1,1,2), VN )
  real :: VOTE(48)  ! Predicted vote for Kennedy in each state
  real :: VOTEK(480)       ! Kennedy
  real :: VOTEN(480)       ! Nixon
  real :: VOTEKN(480,2)    ! Both
  equivalence ( VOTEKN(1,1), VOTEK ), (VOTEKN(1,2), VOTEN )
  real :: WEIGHT(48)
  real :: XK(5,10,10,3)    ! Kennedy
  real :: XN(5,10,10,3)    ! Nixon
  real :: XKN(5,10,10,3,2) ! Both
  equivalence ( XKN(1,1,1,1,1), XK ), ( XKN(1,1,1,1,2), XN )

! Read parameter cards and print title

  read "(a)", title
  print 300, title
300 format ("1SIMULMATICS Project, program to produce state-by-state ",&
     & "simulation results, now in operation."//a/a/1h1)
  read 201, kappa1, kappa2, input1, input2, njob,    &
          & formt1, formt2, vkact(1:48), weight(1:48)
201 format ( 5i10/a/a/(12f6.5) )
  read "(12a6)", name(1:48)

! Read in T-PRIMES for all cells

  select case ( input1 )
  case ( :0 )
    read ( *, formt1 ) ( ii(i), jj(i), kk(i), mm(i), tprime(i,:), i=1,1628 )
  case ( 1 )
    read ( 9, formt1 ) ( ii(i), jj(i), kk(i), mm(i), tprime(i,:), i=1,1628 )
  case ( 2: )
    read ( 9 )         ( ii(i), jj(i), kk(i), mm(i), tprime(i,:), i=1,1628 )
  end select

! Initialize sums and counters

  do is = 1, max(njob,1)
    lcount = -111 ! "No Input for this cell" sentinel
    mcount = 0
    vk = 0.0
    vn = 0.0
    sums = 0.0

!   Read in and index Kennedy and Nixon votes

    select case ( input2 )
    case ( :0 )
      read (*, formt2 )  cshift, ( iii(i), jjj(i), kkk(i), lll(i), &
                       & votek(i), voten(i), i = 1, 480 )
    case ( 1 )
      read (10, formt2 ) cshift, ( iii(i), jjj(i), kkk(i), lll(i), &
                       & votek(i), voten(i),i = 1, 480 )
    case ( 2: )
      read ( 10 )        cshift, ( iii(i), jjj(i), kkk(i), lll(i), &
                       & votek(i), voten(i),i = 1, 480 )
    end select

    do iv = 1, 480
      i = iii(iv) + 1
      j = jjj(iv) + 1
      k = kkk(iv) + 1
      l = lll(iv)
      xkn(i,j,k,l,:) = votekn(iv,:)
      ! Bump one digit of LCOUNT.  Each digit should be bumped exactly once,
      ! or no digit should be bumped, so the result is either zero of -111.
      lcount(i,j,k) = lcount(i,j,k) + 10 ** (3 - l) 
    end do ! IV

!   Check for errors in reading in VOTEK, VOTEN.  Elements of LCOUNT should
!   be -111 (no input), or zero (inputs correctly indexed).

    badInput = .false.
    do i = 1, 5
      do j = 1, 10
        do k = 1, 10
          if ( lcount(i,j,k) /= 0 .and. lcount(i,j,k) /= -111 ) then
            print 301, i - 1, j - 1, k - 1, lcount(i,j,k)
        301 format ( " Error in reading VOTEK, VOTEN." &
                   & " LCOUNT(", i1, ",", i1, ",", i1, ") =",i4, &
                   & ". Run terminated.")
            badInput = .true.
          end if
        end do ! K
      end do ! J
    end do ! I
    if ( badInput ) stop

!   Compute state-by-state results

    do iall = 1, 1628
      i = ii(iall) + 1
      j = jj(iall) + 1
      k = kk(iall) + 1
      m = mm(iall)
      if ( lcount(i,j,k) < 0 ) then ! LCOUNT is -111 here
        print 302, ii(iall), jj(iall), kk(iall), m, name(m), lcount(i,j,k)
    302 format ( " For cell ", 3i1," in state ", i0, " (",a,"), lcount =", i5, &
               & ". VK, VN, SUMS computed without this cell."//)
      else ! LCOUNT is zero here
        vk(m,1:3) = xk(i,j,k,:) * tprime(iall,:) + vk(m,1:3)
        vn(m,1:3) = xn(i,j,k,:) * tprime(iall,:) + vn(m,1:3)
        vk(m,4) = sum(xk(i,j,k,:) * tprime(iall,:)) + vk(m,4)
        vn(m,4) = sum(xn(i,j,k,:) * tprime(iall,:)) + vn(m,4)
        sums(m,1:3) = tprime(iall,:) + sums(m,1:3)
        sums(m,4) = sum(tprime(iall,:)) + sums(m,4)
        lcount(i,j,k) = lcount(i,j,k) + 1
      end if
      mcount(m) = mcount(m) + 1
    end do ! IALL

    vk = vk / sums
    vn = vn / sums
    vknorm = vk / ( vk + vn )
    vnnorm = 1 - vknorm
    vote = vknorm(:,4)

!   Check for errors in input or computation

    do m = 1, 48
      i = merge ( 16, 36, m <= 5 ) ! 16 if m <= 5, else 36
      if ( mcount(m) /= i ) print 303, m, name(m), mcount(m), i
    303 format (" MCOUNT(",i0,"=",a,") =", i3,". It should be ", i2,"."// )
      do l = 1, 4
        do ikn = 1, 2 ! 1 = Kennedy, 2 = Nixon
          if ( vkn(m,l,ikn) < 0 ) then
            print 305, kn(ikn), m, name(m), l, vkn(m,l,ikn)
        305 format ( " V", a, "(", i0, "=",a,",", i1, ") = ", f8.5, &
                   & ". It will be set equal to 0.0" // )
            vkn(m,l,ikn) = 0.0
          end if
        end do
        ! Might print bogus messages with sum very close to 1.0 due to
        ! floating-point round-off.
        ! if ( vk(m,l) + vn(m,l) > 1.0 ) then
        ! Allow one round-off digit relative error (about 1.0e-6) per TPRIME:
        if ( abs ( vk(m,l) + vn(m,l) - 1.0 ) > &
           & 10 * size(tprime,1) * epsilon(1.0) ) then
          print 307, m, l, vk(m,l), m, l, vn(m,l)
      307 format ( " VK(", i2, ",", i1, ") = ",f8.5, " and VN(",I2,",",i1,") =", F8.5, &
                 & ". They will be adjusted to sum to 1.0"// )
          adj = ( vk(m,l) + vn(m,l) - 1.0 ) / 2.0
          vkn(m,l,:) = vkn(m,l,:) - adj
        end if
      end do ! l
      ! There was originally no message printed using format 308.
      ! Floating-point arithmetic is not associative.
      ! Floating-point rounding might have caused bogus messages.
  !   if (sums(m,4)-sums(m,1)-sums(m,2)-sums(m,3) /= 0 ) then
      ! Allow one round-off digit relative error (about 1.0e-6) per TPRIME:
      if ( abs(sums(m,4)-sum(sums(m,1:3))) > &
         & 10 * size(tprime,1) * epsilon(sums(m,4)) * sums(m,4) ) then
      print 308, sums(m,:)
  308 format ( " SUMS(",i2,",1:3) are ", 3f7.4, &
             & ", Which do not add up to SUMS(",i2,",4) = ",f7.0//)
      end if
    end do ! M

    do i = 1, 5
      do j = 1, 10
        do k = 1, 10
          if ( lcount(i,j,k) /= 0 .and. lcount(i,j,k) /= itemTest(i) ) then
            print 309, i - 1, j - 1, k - 1, lcount(i,j,k), itemTest(i)
        309 format (" LCOUNT(",I1,",",I1,",",I1,") =",I3,". It should equal", i3//)
          end if
        end do
      end do
    end do

!   Compute correlation and write output

    !     Predicted vote for Kennedy in each state
    !     PUNCH 9999,VOTE ! Debugging output?  Or input to another program?
    write ( 7, "(12f6.5)" ) vote
    print 310, is, cshift
310 format ( "1Correlations follow for cycle", i3, &
           & ", with the shift parameter equal to", f8.5/ )
    call votcor ( kappa1, kappa2, weight, vote, vkact )
    print 315, ( name(m), vkact(m), vote(m), m = 1, 48 ) 
315 format ( " Actual and predicted Kennedy vote for each state" / &
           &     " STATE   ACTUAL   PREDICTED" / ( 1x, a6, f8.5, f9.5 ) )
    write ( 2, 312 ) title, is, cshift, "eligible electorate", &
                   & ( name(m), ( vk(m,l), vn(m,l), l = 1, 4 ), m = 1, 48 )
312 format ( "1",a / a / &
           & 21x, "Output cycle", i3,".  Shift parameter =", f8.5 // &
           & 29x, "Kennedy, Nixon votes as percentages of ", a, "." // &
           & 14x, "Democrats", 14x, "Republicans", 13x, "Independents", &
           & 27x, "Total" / " State" / &
           & "  NO.", 3x, 3 ( 6x,  "Kennedy", 6x, "Nixon" ), &
           &                  18x, "Kennedy", 6x, "Nixon" /// &
           & ( 1x, a6, 6f12.5, 12x, 2f12.5 ) )
    write ( 2, 312 ) title, is, cshift, "two-party vote", &
                   & ( name(m), ( vknorm(m,l), vnnorm(m,l), l = 1, 4 ), m = 1, 48 )
  end do

  print "( '1',40X,'END OF RUN.')"

end program State_By_State_Sim

! File Name: box1batch2.bin, Fingerprint: 70c5bb-358f00
