二項分布

############################
# 二項分布の計算
#      coded by Y.Suganuma
###########################

############
# 組合せ nCr
############

def comb(n, r)

	c = 1.0
	x = 1.0
	y = 1.0

	if r > 0
		if r > n-r
			r = n - r
		end
		i1 = n
		while i1 > n-r
			x  *= i1
			i1 -= 1
		end
		for i1 in 2 ... r+1
			y *= i1
		end
		c = x / y
	end

	return c
end

######################################
# 二項分布の計算（P(X = x), P(X < x)）
#      x : データ
#      n,p : パラメータ
#      pr : P(X = x)
#      return : P(X < x)
#      coded by Y.Suganuma
######################################

def binomial(x, n, p, pr)

	q     = 1.0 - p
	pr[0] = comb(n, x) * (p ** x) * (q ** (n-x))
	f     = pr[0]
	for i1 in 0 ... x
		f += comb(n, i1) * (p ** i1) * (q ** (n-i1))
	end

	return f
end

print("n は？ ")
n = Integer(gets())
print("p は？ ")
p = Float(gets())
print("グラフ出力？（=1: yes,  =0: no） ")
sw = Integer(gets())
pr = Array.new(1)
			# 密度関数と分布関数の値
if sw == 0
	print("   データは？ ")
	x = Integer(gets())
	f = binomial(x, n, p, pr)
	print("P(X = " + String(x) + ") = " + String(pr[0]) + ",  P( X < " + String(x) + ") = " + String(f) + " (n = " + String(n) + ", p = " + String(p) + ")\n")
			# グラフ出力
else
	print("   密度関数のファイル名は？ ")
	file1 = gets().strip()
	print("   分布関数のファイル名は？ ")
	file2 = gets().strip()
	out1  = open(file1, "w")
	out2  = open(file2, "w")
	for x in 0 ... n+1
		f = binomial(x, n, p, pr)
		out1.print(String(x) + " " + String(pr[0]) + "\n")
		out2.print(String(x) + " " + String(f) + "\n")
	end
	out1.close()
	out2.close()
end