!
! Test for allelic association following Aston and Wilson.
! 2 locus version
!
      integer i, j, k, l, nall1, nall2 , cells(20,20,20,20), gam, mod
      integer a1(20),a2(20),p1(20,20),p2(20,20)
      integer e(20,20),d(20,20)
      open(3,file='assoc.in')
      write(*,'(1x,a,$)') ' Enter number of alleles for Locus A (<21): '
      read(*,*) nall1  
      write(*,'(1x,a,$)') ' Enter number of alleles for Locus B (<21): '
      read(*,*) nall2  
      write(*,*) nall1,nall2
      gam=nall1*nall1*nall2*nall2
      mod=nall1+nall2-1+
     2	    (nall1-1)*(nall1-1)+(nall2-1)*(nall2-1)+
     3	    2*(nall1-1)*(nall2-1)
      k=0
!
! produce all A**2 * B**2 gametes and map onto genotypes
!
      do 10 i=1,nall1
      do 10 i2=1,i
      do 10 j=1,nall2
      do 10 j2=1,j
	k=k+1
	cells(i,i2,j,j2)=k
	cells(i,i2,j2,j)=k
	cells(i2,i,j,j2)=k
	cells(i2,i,j2,j)=k
   10 continue
      write(3,'(3(1x,a/),1x,a,i2,a,i2,a/1x,a/1x,a,i4//1x,a,i4)') 
     2	 '!', '! Test for allelic association', '!',
     3	 '! Locus A: ',nall1,' alleles; Locus B: ',nall2,' alleles',
     4   '!','data ',k,'cells ',gam
      do 20 i=1,nall1
      do 20 i2=1,nall1
   20 write(3,'(20i3,:)') ((cells(i,i2,j,j2),j=1,nall2),j2=1,nall2)
      write(3,*) 'design ',mod
      do 30 i=1,nall1
      do 30 i2=1,nall1
      do 30 j=1,nall2
      do 30 j2=1,nall2
	do 40 k=2,nall1
	  a1(k)=0
	  if (i.eq.k) a1(k)=a1(k)+1
	  if (i2.eq.k) a1(k)=a1(k)+1
	  do 45 l=2,nall1
	    p1(k,l)=0
	    if (i.eq.k.and.i2.eq.l) p1(k,l)=1
   45	  continue
	  do 47 l=2,nall2
	    e(k,l)=0
	    if (i.eq.k.and.j.eq.l) e(k,l)=e(k,l)+1
	    if (i2.eq.k.and.j2.eq.l) e(k,l)=e(k,l)+1
   47	  continue
	  do 49 l=2,nall2
	    d(k,l)=0
	    if (i.eq.k.and.j2.eq.l) d(k,l)=d(k,l)+1
	    if (i2.eq.k.and.j.eq.l) d(k,l)=d(k,l)+1
   49	  continue
   40	continue
	do 50 k=2,nall2
	  a2(k)=0
	  if (j.eq.k) a2(k)=a2(k)+1
	  if (j2.eq.k) a2(k)=a2(k)+1
	  do 55 l=2,nall2
	    p2(k,l)=0
	    if (j.eq.k.and.j2.eq.l) p2(k,l)=1
   55	  continue
   50	continue
	write(3,*) ' ',1,(a1(k),k=2,nall1),(a2(k),k=2,nall2),
     2             ((p1(k,l),k=2,nall1),l=2,nall1),
     3		   ((p2(k,l),k=2,nall2),l=2,nall2),
     4             ((e(k,l),k=2,nall1),l=2,nall2),
     5             ((d(k,l),k=2,nall1),l=2,nall2) 
   30 continue	
      write(3,*) '!',('--',i=1,mod)
      write(3,*) '!',(i,i=1,mod)
      write(3,*) '! i',(' a',k=2,nall1),(' a',k=2,nall2),
     2           ((' p',k=2,nall1),l=2,nall1),
     3           ((' p',k=2,nall2),l=2,nall2),
     4           ((' e',k=2,nall1),l=2,nall2),
     5           ((' d',k=2,nall1),l=2,nall2) 
      close(3)
      end