PROGRAM freq2d ! Generate pairs of random numbers and put into 256 x 256 bins. IMPLICIT NONE INTEGER :: i, j, k, l, lower, ndf, nrows, ncols, upper REAL :: x(2), aver, stdev, chisq INTEGER, ALLOCATABLE :: freq(:,:), seed(:) CHARACTER (LEN=200) :: text WRITE(*, *) 'Enter number of rows & columns (e.g. 24 64): ' READ(*, *) nrows, ncols ALLOCATE( freq(0:nrows-1,0:ncols-1) ) ! Set the random number seed. CALL RANDOM_SEED(size=k) ALLOCATE (seed(k)) CALL RANDOM_SEED(get=seed) WRITE(*, *)'Old random number seeds: ', seed WRITE(*, '(a, i4, a)') ' Enter ', k, ' integers as random number seeds: ' READ(*, *) seed CALL RANDOM_SEED(put=seed) ! Generate 4096 x 4096 pairs of random numbers freq = 0 DO i = 1, 4096 DO j = 1, 4096 CALL RANDOM_NUMBER(x) k = nrows * x(1) l = ncols * x(2) freq(k,l) = freq(k,l) + 1 END DO END DO ! Average & variance of number in each cell = 4096 x 4096 / (nrows x ncols) aver = 4096. * 4096. / (nrows * ncols) stdev = SQRT(aver) upper = aver + 3.*stdev + 0.5 lower = aver - 3.*stdev + 0.5 chisq = 0.0 DO k = 0, nrows-1 text = ' ' DO l = 0, ncols-1 IF (freq(k,l) >= lower .AND. freq(k,l) <= upper) THEN text(l+3:l+3) = '+' ELSE IF (freq(k,l) > upper) THEN text(l+3:l+3) = CHAR(219) END IF chisq = chisq + (freq(k,l) - aver)**2 END DO WRITE(*, '(a)') TRIM(text) END DO WRITE(*, *) ' If most of the cells are NOT pluses, its a bad generator' chisq = chisq / aver ndf = nrows * ncols - 1 WRITE(*, '(a, f10.1, a, i6, a)') & ' Chi-squared =', chisq, ' with', ndf, ' deg. of freedom' chisq = ndf + 2.0*SQRT(2.0*ndf) WRITE(*, '(a, f10.1)') ' Chi-squared should be less than', chisq STOP END PROGRAM freq2d