
# plotGamma.py
# info: https://en.wikipedia.org/wiki/Gamma_function
#       https://scipython.com/book/chapter-8-scipy/examples/the-gamma-function/

import math
import matplotlib.pyplot as plt
from frange import *

def ngamma(x):
  try:
    return math.gamma(x)
  except ValueError:
    return math.inf

def ngamma2(x):
  if x == 0:
    return math.inf
  elif x > 0:
    return math.gamma(x)
  else:  # use Euler reflection
    return math.pi / (math.sin(math.pi*x) * math.gamma(1-x))


# some tests
print("gamma(-5):", ngamma(-4))
print("gamma(5):", ngamma(5))  # 4!
print("gamma(0.5):", ngamma(0.5))
print("gamma(-0.5):", ngamma(-0.5))
print("gamma(-1.5):", ngamma(-1.5))


xs = linspace(-4, 4, 1000)
ys = [ngamma(x) for x in xs]

plt.plot(xs, ys, color="red")
plt.ylim(-50, 50)   #  to cut off infinities at x=0 and neg ints

plt.axhline(y=0, color='k')
plt.axvline(x=0, color='k')
plt.title('Gamma(x) for reals')

plt.show()