code 34502; procedure BOUNDS(N,A,RE,IM,RELE,ABSE,RECENTRE,IMCENTRE,BOUND); value N, RELE, ABSE; integer N; real RELE, ABSE; array RE, IM, A, RECENTRE, IMCENTRE, BOUND; begin integer I, J, K, L, INDEX1, INDEX2; boolean GOON; real H, MIN, RECENT, IMCENT, GIA, XK, YK, ZK, CORR; array RC, C, RCE[0:N], CLUST[1:N]; real procedure G(RAD, RECENT, IMCENT, K, M); value RAD, RECENT, IMCENT, K, M; real RAD, RECENT, IMCENT; integer K, M; begin real S, H1, H2; integer I; S:= SQRT(RECENT * RECENT + IMCENT * IMCENT) + RAD; H1:= RC[1]; H2:= RC[0]; for I:= 2 step 1 until N do H1:= H1*S + RC[I]; for I:= 1 step 1 until M-1, M+K step 1 until N do H2:= H2 * ABS(SQRT((RE[I]-RECENT)**2 + (IM[I]-IMCENT)**2) - RAD); G:= if H1=0 then 0 else if H2=0 then -10 else H1 / H2 end; procedure KCLUSTER(K, M); value K, M; integer K, M; begin integer I, J, STOP, L; boolean NONZERO; real RECENT, IMCENT, D, PROD, RAD, GR, R; array DIST[M: M+K-1]; RECENT:= RE[M]; IMCENT:= IM[M]; STOP:= M+K-1; L:= SIGN(IMCENT); NONZERO:= L ^= 0; for I:= M+1 step 1 until STOP do begin RECENT:= RECENT+RE[I]; if NONZERO then begin NONZERO:= L = SIGN(IM[I]); IMCENT:= IMCENT+IM[I] end end; RECENT:= RECENT/K; IMCENT:= if NONZERO then IMCENT/K else 0; D:= 0; RAD:= 0; for I:= M step 1 until STOP do begin RECENTRE[I]:= RECENT; IMCENTRE[I]:= IMCENT; DIST[I]:= SQRT((RE[I] -RECENT)**2 + (IM[I]-IMCENT)**2); if D < DIST[I] then D:= DIST[I] end; GR:= ABS(G(0, RECENT, IMCENT, K, M)); if GR > 0 then begin for J:= 1, 1 while PROD <= GR do begin R:= RAD; RAD:= D + EXP(LN(1.1*GR)/K); if RAD = R then RAD:= EXP(LN(1.1)/K) * RAD; GR:= G(RAD, RECENT, IMCENT, K, M); PROD:= 1; for I:= M step 1 until STOP do PROD:= PROD*(RAD-DIST[I]) end end; for I:= M step 1 until STOP do begin BOUND[I]:= RAD; CLUST[I]:= K end; procedure SHIFT(INDEX, NEW); value INDEX, NEW; integer INDEX, NEW; begin integer J, PLACE, CLUSTIN; real BOUNDIN, IMCENT, RECENT; real array WA1, WA2[1:CLUST[INDEX]]; CLUSTIN:= CLUST[INDEX]; BOUNDIN:= BOUND[INDEX]; IMCENT:= IMCENTRE[INDEX]; RECENT:= RECENTRE[INDEX]; for J:= 1 step 1 until CLUSTIN do begin PLACE:=INDEX+J-1; WA1[J]:= RE[PLACE]; WA2[J]:= IM[PLACE]; end; for J:= INDEX-1 step -1 until NEW do begin PLACE:= J+CLUSTIN; RE[PLACE]:= RE[J]; IM[PLACE]:= IM[J]; CLUST[PLACE]:= CLUST[J]; BOUND[PLACE]:= BOUND[J]; RECENTRE[PLACE]:= RECENTRE[J]; IMCENTRE[PLACE]:= IMCENTRE[J] end; for J:= NEW+CLUSTIN-1 step -1 until NEW do begin PLACE:= J+1-NEW; RE[J]:= WA1[PLACE]; IM[J]:= WA2[PLACE]; BOUND[J]:= BOUNDIN; CLUST[J]:= CLUSTIN; RECENTRE[J]:= RECENT; IMCENTRE[J]:= IMCENT end end; GIA:= GIANT; RC[0]:= C[0]:= A[N]; RCE[0]:= ABS(C[0]); K:= 0; for I:= 1 step 1 until N do begin RC[I]:= RCE[I]:= 0 ; C[I]:= A[N-I] end; for I:= 0 while K < N do begin K:= K+1; XK:= RE[K]; YK:= IM[K]; ZK:= XK*XK+YK*YK; for J:= K step -1 until 1 do RCE[J]:= RCE[J]+RCE[J-1]*SQRT(ZK); if YK = 0 then begin for J:= K step -1 until 1 do RC[J]:= RC[J]-XK*RC[J-1] end else begin K:= K+1; if K <= N & XK = RE[K] & YK = -IM[K] then begin XK:= -2*XK; for J:= K step -1 until 1 do RCE[J]:= RCE[J]+RCE[J-1]*SQRT(ZK); for J:= K step -1 until 2 do RC[J]:= RC[J]+XK*RC[J-1]+ZK*RC[J-2]; RC[1]:= RC[1]+XK*RC[0] end end end; RC[0]:= RCE[0]; CORR:= 1.06*ARREB; for I:= 1 step 1 until N-1 do RC[I]:= ABS(RC[I]-C[I])+RCE[I]*CORR*(N+I-2)+RELE*ABS(C[I])+ABSE; RC[N]:= ABS(RC[N]-C[N])+RCE[N]*CORR*(N-1)+RELE*ABS(C[N])+ABSE; for I:= 1 step 1 until N do KCLUSTER(1, I); GOON:= true; for L:= 1 while GOON do begin INDEX1:= INDEX2:= 0; MIN:= GIANT; I:= N-CLUST[N]+1; for I:= I while I >= 2 do begin J:= I; RECENT:= RECENTRE[I]; IMCENT:= IMCENTRE[I]; for J:= J while J >= 2 do begin J:= J-CLUST[J-1]; H:= SQRT((RECENT-RECENTRE[J])**2 + (IMCENT-IMCENTRE[J])**2); if H < BOUND[I] + BOUND[J] & H <= MIN then begin INDEX1:= J; INDEX2:= I; MIN:= H end end; I:= I-CLUST[I-1] end; if INDEX1 = 0 then GOON:= false else begin if IMCENTRE[INDEX1] = 0 then begin if IMCENTRE[INDEX2] ^= 0 then CLUST[INDEX2]:= 2*CLUST[INDEX2] end else if IMCENTRE[INDEX2] = 0 then CLUST[INDEX1]:= 2*CLUST[INDEX1]; K:= INDEX1+CLUST[INDEX1]; if K ^= INDEX2 then SHIFT(INDEX2, K); K:= CLUST[INDEX1]+CLUST[K]; KCLUSTER(K, INDEX1) end end end; eop