code 41001; real procedure BININV(PROB, N, P, LEFT); value PROB, N, P, LEFT; real PROB, N, P; Boolean LEFT; begin integer X; real PX, PCUM; if PROB ≤ 0 ∨ PROB ≥ 1 then STATAL3 ERROR(“BININV”, 1, PROB) else if N > ENTIER(N) ∨ N < 0 then STATALS ERROR(“BININV”, 2, N) else if P < 0 ∨ P > 1 then STATAL3 ERROR(“BININV”, 3, P); if P = 0 ∨ N = 0 then BININV:= (if LEFT then -1 else 1) else if P = 1 then BININV:= (if LEFT then N - 1 else N + 1) else if LEFT then begin X:= PHINV(PROB) × SQRT(N × P × (1 - P)) - 0.5 + N × P; if X < 0 then X:= 0 else if X > N then X:= N; if PROB < (1 - P) ⭡ N then BININV:= -1 else begin PX:= BINPROB(X, N, P); PCUM:= BIN(X, N, P); if PCUM > PROB then begin for PCUM:= PCUM - PX while PCUM > PROB do begin PX:= PX × X × (1 - P) / (N - X + 1) / P; X:= X - 1 end; X:= X - 1 end else begin for PX:= PX × (N - X) / (X + 1) × P / (1 - P) while PCUM + PX < PROB do begin X:= X + 1; PCUM:= PCUM + PX end end; BININV:= X end end else begin X:= PHINV(1 - PROB) × SQRT(N × P × (1 - P)) + 0.5 + N × P; if X < 0 then X:= 0 else if X > N then X:= N; if PROB < P ⭡ N then BININV:= N + 1 else begin PCUM:= 1 - BIN(X - 1, N, P); PX:= BINPROB(X, N, P); if PCUM < PROB then begin for PX:= PX × X × (1 - P) / (N - X + 1) /P while PCUM + PX < PROB do begin X:= X - 1; PCUM:= PCUM + PX end end else begin for PCUM:= PCUM - PX while PCUM > PROB do begin PX:= PX × (N - X) × P / (X + 1) / (1 - P); X:= X + 1 end; X:= X + 1 end; BININV:= X end end end BININV; eop