'Two dimensional Gaussian function'

@staticmethod
def evaluate(x, y, amplitude, x_mean, y_mean, x_stddev, y_stddev, theta):

    cost2 = (np.cos(theta) ** 2)
    sint2 = (np.sin(theta) ** 2)
    sin2t = np.sin((2.0 * theta))
    xstd2 = (x_stddev ** 2)
    ystd2 = (y_stddev ** 2)
    xdiff = (x - x_mean)
    ydiff = (y - y_mean)
    a = (0.5 * ((cost2 / xstd2) + (sint2 / ystd2)))
    b = (0.5 * ((sin2t / xstd2) - (sin2t / ystd2)))
    c = (0.5 * ((sint2 / xstd2) + (cost2 / ystd2)))
    return (amplitude * np.exp((- (((a * (xdiff ** 2)) + ((b * xdiff) * ydiff)) + (c * (ydiff ** 2))))))
