ポアソン分布

# -*- coding: UTF-8 -*-
from math import *
import numpy as np

#############
# 階乗の計算
#############

def fact(n) :
	x = 1.0
	for i1 in range(2, n+1) :
		x *= i1
	return x

##########################################
# ポアソン分布の計算（P(X = x), P(X < x)）
#      x : データ
#      ram : パラメータ
#      pr : P(X = x)
#      return : P(X < x)
#      coded by Y.Suganuma
##########################################

def Poisson(x, ram, pr) :

	f     = 0.0
	pr[0] = (ram ** x) * exp(-ram) / fact(x)
	f     = pr[0]
	for i1 in range(0, x) :
		f += (ram ** i1) * exp(-ram) / fact(i1)

	return f

----------------------------------

# -*- coding: UTF-8 -*-
import numpy as np
import sys
from math import *
from function import Poisson

###########################
# ポアソン分布の計算
#      coded by Y.Suganuma
###########################

s   = input("λ は？ ")
ram = float(s)
s   = input("グラフ出力？（=1: yes,  =0: no） ")
sw  = int(s)
pr  = np.empty(1, np.float)
			# 密度関数と分布関数の値
if sw == 0 :
	s = input("   データは？ ")
	x = int(s)
	f = Poisson(x, ram, pr)
	print("P(X = " + str(x) + ") = " + str(pr[0]) + ",  P( X < " + str(x) + ") = " + str(f) + " (λ = " + str(ram) + ")")
			# グラフ出力
else :
	file1 = input("   密度関数のファイル名は？ ")
	file2 = input("   分布関数のファイル名は？ ")
	s     = input("   データの上限は？ ")
	up    = int(s)
	out1  = open(file1, "w")
	out2  = open(file2, "w")
	for x in range(0, up+1) :
		f = Poisson(x, ram, pr)
		out1.write(str(x) + " " + str(pr[0]) + "\n")
		out2.write(str(x) + " " + str(f) + "\n")
	out1.close()
	out2.close()