#ifndef CHISQUARED_H
#define CHISQUARED_H

#include <cmath>

namespace Distribution {

	/**
	 * https://en.wikipedia.org/wiki/Chi-squared_distribution
	 * @brief The ChiSquared class
	 */
	template <typename Scalar> class ChiSquared {

	private:

		// degrees of freedom
		int k;

	public:

		/** ctor */
		ChiSquared(const int k) : k(k) {
			;
		}

		Scalar get(const Scalar val) const {

			const Scalar k2 = k/((Scalar)2);
			const Scalar k2_1 = k2 - 1;
			const Scalar gamma = std::tgamma(k2);

			return 1.0 / (std::pow(2, k2)*gamma) * std::pow(val, k2_1) * std::exp(-val/2);

		}


		Scalar getInvCDF(const Scalar val) const {

			// brute-force get the inverse CDF...

			const Scalar ss = 0.002;
			Scalar sum = 0;

			for (float t = 0; t < 20; t += ss) {
				sum += get(t) * ss;
				if (sum >= val) {
					return t;
				}
			}

			throw "failed";

		}

	};

}

#endif // CHISQUARED_H