'**************************************************************** '* Calculate the median value of an array without sorting * '* ------------------------------------------------------------ * '* REFERENCE: * '* "NUMERICAL RECIPES By W.H. Press, B.P. Flannery, * '* S.A. Teukolsky and W.T. Vetterling, Cambridge * '* University Press, 1986" [BIBLI 08]. * '* * '* Basic Version By J-P Moreau * '* (www.jpmoreau.fr) * '* ------------------------------------------------------------ * '* SAMPLE RUN: (see input file tmdian1.dat). * '* * '* N= 80 * '* * '* Input Table: * '* * '* 407.8 192.8 851.1 604.4 932.3 799.4 914.5 965.8 453.7 295.1 * '* 154.5 977.4 410.2 916.2 934.7 504.8 823.2 225.2 456.6 49.0 * '* 933.5 663.0 335.3 346.6 568.7 956.1 654.7 300.7 379.6 591.9 * '* 992.9 689.6 644.7 305.4 148.2 257.2 664.6 612.1 713.0 99.7 * '* 46.5 167.6 984.6 847.2 55.4 82.7 999.0 10.7 877.7 929.4 * '* 398.1 972.8 874.1 755.1 472.1 122.8 671.4 35.5 128.8 76.8 * '* 454.2 959.2 510.1 791.3 122.8 176.6 237.9 995.8 548.3 309.8 * '* 162.6 996.5 750.0 250.6 577.7 761.1 101.9 797.1 539.0 723.5 * '* * '* Median Value: 558.5 * '* * '* * '**************************************************************** DEFINT I-N CLS OPEN "c:\temp1\tmdian1.dat" FOR INPUT AS #1 'read size of table INPUT #1, N DIM A(N) 'read table of numbers FOR i = 1 TO N INPUT #1, A(i) NEXT i CLOSE #1 PRINT PRINT " N="; N PRINT PRINT " Input Table:" GOSUB 2000 'call TWRIT(N,A) 'call MDIAN1 subroutine GOSUB 1000 'call MDIAN1(A,N,AMED) PRINT PRINT " Median Value: "; AMED PRINT INPUT "", R\$ END '******************************************************* '* Given an array X of N numbers, returns their median * '* value XMED. The array X is not modified, and is * '* accessed sequentially in each consecutive pass. * '******************************************************* 1000 'Subroutine MDIAN1(A, N, AMED) BIG = 1.2E+16: AFAC = 1.5: AMP = 1.5 ' Here, AMP is an overconvergence factor: on each iteration, ' we move the guess by this factor. AFAC is a factor used to ' optimize the size of the "smoothing constant" EPS at each ' iteration. ' A1,AA,AM,AP,DUM,EPS,SUM,SUMX,XM,XP,XX:Real ' J,NM,NP:Integer A1 = .5 * (A(1) + A(N)) EPS = ABS(A(N) - A(1)) AP = BIG AM = -BIG 1 SUM = 0# SUMX = 0# NP = 0 NM = 0 XP = BIG XM = -BIG FOR J = 1 TO N XX = A(J) IF XX <> A1 THEN IF XX > A1 THEN NP = NP + 1 IF XX < XP THEN XP = XX ELSEIF XX < A1 THEN NM = NM + 1 IF XX > XM THEN XM = XX END IF DUM = 1 / (EPS + ABS(XX - A1)) SUM = SUM + DUM SUMX = SUMX + XX * DUM END IF NEXT J IF NP - NM >= 2 THEN 'guess is too low, make another pass AM = A1 AA = XP + MAX(0!, SUMX / SUM - A1) * AMP IF AA > AP THEN AA = .5 * (A1 + AP) EPS = AFAC * ABS(AA - A1) A1 = AA GOTO 1 ELSEIF NM - NP >= 2 THEN 'guess is too high AM = A1 AA = XM + MIN(0!, SUMX / SUM - A1) * AMP IF AA < AM THEN AA = .5 * (A1 + AM) EPS = AFAC * ABS(AA - A1) A1 = AA GOTO 1 ELSE 'guess is ok IF (N MOD 2) = 0 THEN 'for even N median is an average IF NP = NM THEN AMED = .5 * (XP + XM) ELSEIF NP > NM THEN AMED = .5 * (A1 + XP) ELSE AMED = .5 * (XM + A1) END IF ELSE IF NP = NM THEN AMED = A1 ELSEIF NP > NM THEN AMED = XP ELSE AMED = XM END IF END IF END IF RETURN 'write table of size N to standard output 2000 'SUBROUTINE TWRIT(N,A) F\$ = " ####.#" PRINT FOR i = 1 TO N PRINT USING F\$; A(i); IF INT(i / 10) = i / 10 THEN PRINT NEXT i RETURN 'end of file tmdian1.bas