%matplotlib inline
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np
import sys,os
ia898path = os.path.abspath('/etc/jupyterhub/ia898_1s2017/')
if ia898path not in sys.path:
    sys.path.append(ia898path)
import ia898.src as ia

# First case: unidimensional
# x: single value (single sample) or a row of values (many samples)
# mu and sigma are scalar
def fun1(x, mu, sigma):
    return (1./(np.sqrt(2 * np.pi) * sigma)) * np.exp(-1./2 * ((x-mu)/ sigma)**2)

print('Computing the Gaussian function at a single point')  
ex1 = "fun1( 10, 10, 5)"
print(ex1,"=>", eval(ex1))

ex2 = "fun1( 15, 10, 5)"
print(ex2,"=>", eval(ex2))

"""
Explanation: Demo Gaussian Generation
Illustrate the generation of a d-dimensional Gaussian image
Description
The sequence below shows a technique to a d-dimensional Gaussian image,
understanding the difficulties in computing an equation with vector and
matrix notation.
One dimensional case
The Gaussian function is a symmetric bell shaped function that is characterized by
two parameters: mean and variance. The one-dimensional Gaussian function at point
$x$ is given by the following equation, with mean $\mu$ and variance 
$\sigma^2$. The function is maximum at point $x=\mu$ and it falls by the
factor $e^{-\frac{1}{2}}$ (approx. 0.6) at point $x=\sigma$ away from the mean.
Equation
$$    f(x) = \frac{1}{\sqrt{2 \pi} \sigma} exp\left[ -\frac{1}{2}  \frac{\left(x - \mu\right)^2}{\sigma^2}  \right] $$
As this function is scalar, it is possible to compute this function on N samples represented
as a N x 1 vector ${\mathbf x} = [x_0, x_1, x_2, \ldots x_{N-1}]^\mathrm{T}$:
$$    f({\mathbf x}) = \frac{1}{\sqrt{2 \pi} \sigma} exp\left[ -\frac{1}{2}  \frac{\left({\mathbf x} - \mu\right)^2}{\sigma^2}  \right]$$
End of explanation
"""


ex3 = "fun1( np.array([[10,15,20]]).T, 10, 5)"
print(ex3,"=>\n", eval(ex3))

x = np.arange(-5,26).reshape(-1,1)
y = fun1(x, 10, 5)
plt.plot(x,y)

"""
Explanation: Computing the Gaussian function at many points, using the same code
End of explanation
"""


# Second case: d-dimensional, single sample
# x: single column vector (single sample with d characteristics)
# mu: column vector, 1 x d
# sigma: covariance matrix, square and symmetric, d x d
def funn(X, MU, COV):
    d = len(X)
    Xc = X - MU
    aux = np.linalg.inv(COV).dot(Xc)
    k = 1. * (Xc.T).dot(aux)
    return (1./((2 * np.pi)**(d/2.) * np.sqrt(np.linalg.det(COV)))) * np.exp(-1./2 * k)

print('\ncomputing the Gaussian function at a single 3-D sample')
X1 = np.array([[10],
               [5],
               [3]])
MU = X1
COV = np.array([[10*10, 0,   0],
                [0,     5*5, 0],
                [0,     0,   3*3]])
print('X1=',X1)
print('MU=',MU)
print('COV=',COV)
ex4 = "funn( X1, MU, COV)"
print(ex4,"=>", eval(ex4))

print('\nComputing the Gaussian function at two 3-D samples')
print('\nNote that it does not work')
X2 = 1. * X1/2
X = np.hstack([X1,X2])
print('X=',X)
ex5 = "funn( X, MU, COV)"
print(ex5,"=>", eval(ex5))

"""
Explanation: d-dimensional Case
If a sample point is a vector of dimension d: ${\mathbf x} = [x_0, x_1, \ldots x_{d-1}]^T$, 
the d-dimensional Gaussian function is characterized by the mean 
vector: ${\mathbf \mu} = [\mu_0, \mu_1, \ldots \mu_{d-1}]^T$ and the symmetric square 
covariance matrix:
$$     \Sigma_d = \left( 
                \begin{array}{cccc} 
                  \sigma_0^2           & \sigma_0\sigma_1     & \ldots & \sigma_0\sigma_{d-1} \
                  \sigma_1\sigma_0     & \sigma_1^2           & \ldots & \sigma_1\sigma_{d-1} \
                  \vdots               & \vdots               & \vdots & \vdots \
                  \sigma_{d-1}\sigma_0 & \sigma_{d-1}\sigma_1 & \ldots & \sigma_{d-1}^2
                \end{array}
                \right) $$
Equation
$$    f({\mathbf x}) = \frac{1}{(2 \pi)^{d/2}|\Sigma|^{1/2}} exp\left[ -\frac{1}{2}\left({\mathbf x} - {\mathbf \mu} \right)^\mathrm{T}\Sigma^{-1}\left({\mathbf x} - {\mathbf \mu} \right)\right] $$
End of explanation
"""


