#pragma once
#include "../MFMTwoCurrency.h"
#include "UsesHullWhiteJ0.h"

//////////////////////////////////////////////////////////////////////////////////////
// A 3-process Markov Functional model, with two currencies and their FX rate,
// employing hull white method, both mean reversions set to 0
//////////////////////////////////////////////////////////////////////////////////////
template <class CalibratorType>
class HWMFMTwoCurrency : public MFMTwoCurrency<CalibratorType>
{

private:

	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	// Create domestic and foreign currency processes - delegate final work to next constructor
	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	HWMFMTwoCurrency(
		size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times,
		const std::array<double, 2>& discount_factors_t1,
		const std::array<gelib::d_vector, 2>& accrual_factors, const std::array<gelib::d_vector, 2>& forward_libors,
		const std::array<gelib::d_vector, 3>& process_variances, const std::array<gelib::d_vector, 3>& process_covariances,
		const std::array<gelib::d_vector, 3>& process_volatilities_squared,
		gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance, double spot_fx,
		gelib::d_vector const& mu_X, gelib::d_vector const& mu_Y, gelib::d_vector const& mu_Zhat, gelib::d_vector const& mu_Z,
		gelib::d_vector const& mu_R, gelib::d_vector const& mu_Rhat, gelib::d_vector const& mu_y, gelib::d_vector const& mu_z,
		gelib::d_vector const& C_XX, gelib::d_vector const& C_ZZ, gelib::d_vector const& C_YY, gelib::d_vector const& C_zz,
		gelib::d_vector const& C_yy, gelib::d_vector const& C_RR, gelib::d_vector const& C_RhatRhat, gelib::d_vector const& C_XR,
		gelib::d_vector const& C_YR, gelib::d_vector const& C_ZR, gelib::d_vector const& C_XRhat, gelib::d_vector const& C_YRhat,
		gelib::d_vector const& C_ZRhat, gelib::d_vector const& C_ZhatRhat, gelib::d_vector const& C_Xy, gelib::d_vector const& C_Ry,
		gelib::d_vector const& C_Rhaty, gelib::d_vector const& C_RhatR, gelib::d_vector const& C_Yz, gelib::d_vector const& C_Zz,
		gelib::d_vector const& C_XY, gelib::d_vector const& C_XZ, gelib::d_vector const& C_YZ, gelib::d_vector const& C_ZhatZhat,
		std::unique_ptr<CalibratorType>(*ad_constructor) (void))
		:
		// delegate to next constructor
		HWMFMTwoCurrency(gridpoints_per_slice, sd_spread, poly_order, times, discount_factors_t1,
			process_variances,
			process_covariances,
			process_volatilities_squared,
			y_larrow_mean, y_larrow_variance, spot_fx,

			// Construct an array interest rate processes (one per time slice) to represent domestic currency,
			// using the customised Hull White constructor which sets up the appropriate
			// libor in arrears calculation
			std::move(InterestRateTimeSliceProcess::construct(
				sd_spread, gridpoints_per_slice, this, 0.0, times, forward_libors[0],
				process_volatilities_squared[0], process_variances[0], accrual_factors[0], discount_factors_t1[0],
				mu_R, mu_X, C_RR, C_XX, C_XR)),

			// Construct an array interest rate processes (one per time slice) to represent foreign currency,
			// using the customised Hull White constructor which sets up the appropriate
			// libor in arrears calculation
			std::move(InterestRateTimeSliceProcess::construct(
				sd_spread, gridpoints_per_slice, this, 0.0, times, forward_libors[1],
				process_volatilities_squared[2], process_variances[2], accrual_factors[1], discount_factors_t1[1],
				mu_Rhat + C_YRhat, mu_Z + C_YZ, C_RhatRhat, C_ZhatZhat, C_ZhatRhat)),

			mu_X, mu_y,
			C_XX, C_yy, C_Xy,
			mu_z, mu_Y, mu_Z,
			C_zz, C_YY, C_ZZ,
			C_YZ, C_Yz, C_Zz,
			ad_constructor)
	{

		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		// uncomment two lines below to set up calibrator to replace domestic libors with analytic values
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		/*
		for (size_t i = 0; i < this->m_number_of_slices - 1; ++i)
			(this->m_time_slices.get(i)->get_calibrator()).prepare_analytic_xlibor(mu_R[i], C_RR[i]);
		*/

		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		// uncomment two lines below to set up calibrator to replace foreign libors with analytic values
		// TODO WARNING 
		// Note that currently, in correlated case with alpha non-zero, the foreign libors
		// aren't matching up to these values in the Hull White model. 
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		/*
		for (size_t i = 0; i < this->m_number_of_slices - 1; ++i)
			(this->m_time_slices.get(i)->get_calibrator()).prepare_analytic_zlibor(mu_Rhat[i], C_RhatRhat[i]);
		*/
	}

	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	// Use the more general MFMTwoCurrency constructor here to finalise the model which is then ready to be calibrated.
	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	HWMFMTwoCurrency(
		size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times, const std::array<double, 2>& discount_factors_t1,
		const std::array<gelib::d_vector, 3>& process_variances, const std::array<gelib::d_vector, 3>& process_covariances,
		const std::array<gelib::d_vector, 3>& process_volatilities_squared,
		gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance, double spot_fx,
		std::vector<std::unique_ptr<InterestRateTimeSliceProcess>> dom_processes,
		std::vector<std::unique_ptr<InterestRateTimeSliceProcess>> fgn_processes,
		gelib::d_vector const& mu_X, gelib::d_vector const& mu_y,
		gelib::d_vector const& C_XX, gelib::d_vector const& C_yy, gelib::d_vector const& C_Xy,
		gelib::d_vector const& mu_z, gelib::d_vector const& mu_Y, gelib::d_vector const& mu_Z,
		gelib::d_vector const& C_zz, gelib::d_vector const& C_YY, gelib::d_vector const& C_ZZ,
		gelib::d_vector const& C_YZ, gelib::d_vector const& C_Yz, gelib::d_vector const& C_Zz,
		std::unique_ptr<CalibratorType>(*ad_constructor) (void))
		:
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		// initialise the two currency model
		// provide the processes and their covariances at each time slice, and the calbrator class they will need to calculate
		// j0 on each grid
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		MFMTwoCurrency<CalibratorType>(
			times.size() + 1, poly_order, std::array<double, 2>({ 0.0, 0.0 }), discount_factors_t1,
			y_larrow_mean, y_larrow_variance, gridpoints_per_slice, sd_spread, times,
			ad_constructor, &SliceProcessCorrelations3d::construct,
			process_variances,
			process_covariances,

			std::move(dom_processes),
			std::move(fgn_processes),

			// Construct FX processes, one at each time slice, note that we need the domestic and foreign currencies from the previous constructor
			// to do so.
			std::move(FXTimeSliceProcess::construct(sd_spread, mu_y, gridpoints_per_slice, this, times, process_volatilities_squared[1],
				process_variances[1], this->m_mean_reversion_param_domestic, this->m_mean_reversion_param_foreign, spot_fx, std::move(dom_processes), std::move(fgn_processes),
				mu_y, mu_X, C_yy, C_XX, C_Xy))
			)

	{
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		// Set up the calibrator to be able to be able to calculate J0y and J0z
		//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
		for (size_t i = 0; i < this->m_time_slices.size() - 1; ++i)
			this->m_time_slices[i]->get_calibrator().pre_calibrate(
				mu_y[i], mu_X[i], C_yy[i], C_XX[i], C_Xy[i],
				mu_z[i], mu_Y[i] - mu_Z[i], C_zz[i], C_YY[i] + C_ZZ[i] - 2 * C_YZ[i], C_Yz[i] - C_Zz[i]
			);
	}

public:

	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	// Note this is the only public constructor - it can be used from outside the class to create a HWMFMTwoCurrency
	// object which is ready for calibration.
	// The input data is extracted here and all the work is delegated to next constructor
	//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
	HWMFMTwoCurrency(
		size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times,
		const std::array<double, 2>& discount_factors_t1, const std::array<gelib::d_vector, 2>& accrual_factors,
		const std::array<gelib::d_vector, 2>& forward_libors, const std::array<gelib::d_vector, 3>& process_variances,
		const std::array<gelib::d_vector, 3>& process_covariances, gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance,
		double spot_fx, gelib::d_matrix const& mu_data, gelib::d_matrix const& C_data, std::unique_ptr<CalibratorType>(*ad_constructor) (void))
		:
		HWMFMTwoCurrency(
			gridpoints_per_slice, sd_spread, poly_order, times, discount_factors_t1, accrual_factors, forward_libors,
			process_variances,
			process_covariances,
			// Get piecewise process volatilities from process variances
			MFMTwoCurrency<CalibratorType>::variances_to_squared_volatilites(0.0, times, process_variances),
			y_larrow_mean, y_larrow_variance, spot_fx,
			// here we just extract the appropriate data from the input matrices, assuming the ordering is correct
			mu_data.extract_vector(-1, 0), mu_data.extract_vector(-1, 1), mu_data.extract_vector(-1, 2), mu_data.extract_vector(-1, 3),
			mu_data.extract_vector(-1, 4), mu_data.extract_vector(-1, 5), mu_data.extract_vector(-1, 6), mu_data.extract_vector(-1, 7),
			C_data.extract_vector(-1, 0), C_data.extract_vector(-1, 1), C_data.extract_vector(-1, 2), C_data.extract_vector(-1, 3),
			C_data.extract_vector(-1, 4), C_data.extract_vector(-1, 5), C_data.extract_vector(-1, 6), C_data.extract_vector(-1, 7),
			C_data.extract_vector(-1, 8), C_data.extract_vector(-1, 9), C_data.extract_vector(-1, 10), C_data.extract_vector(-1, 11),
			C_data.extract_vector(-1, 12), C_data.extract_vector(-1, 13), C_data.extract_vector(-1, 14), C_data.extract_vector(-1, 15),
			C_data.extract_vector(-1, 16), C_data.extract_vector(-1, 17), C_data.extract_vector(-1, 18), C_data.extract_vector(-1, 19),
			C_data.extract_vector(-1, 20), C_data.extract_vector(-1, 21), C_data.extract_vector(-1, 22), C_data.extract_vector(-1, 1),
			ad_constructor
		)
	{}

	static void constructor_check(size_t gridpoints_per_slice, double sd_spread, size_t poly_order,
		gelib::d_vector const& times, const std::array<double, 2>& discount_factors_t1,
		const std::array<gelib::d_vector, 2>& accrual_factors, const std::array<gelib::d_vector, 2>& forward_libors,
		const std::array<gelib::d_vector, 3>& process_variances, const std::array<gelib::d_vector, 3>& process_covariances,
		gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance, double spot_fx,
		gelib::d_vector const& mu_X, gelib::d_vector const& mu_Y, gelib::d_vector const& mu_Zhat, gelib::d_vector const& mu_Z,
		gelib::d_vector const& mu_R, gelib::d_vector const& mu_Rhat, gelib::d_vector const& mu_y, gelib::d_vector const& mu_z,
		gelib::d_vector const& C_XX, gelib::d_vector const& C_ZZ, gelib::d_vector const& C_YY, gelib::d_vector const& C_zz,
		gelib::d_vector const& C_yy, gelib::d_vector const& C_RR, gelib::d_vector const& C_RhatRhat, gelib::d_vector const& C_XR,
		gelib::d_vector const& C_YR, gelib::d_vector const& C_ZR, gelib::d_vector const& C_XRhat, gelib::d_vector const& C_YRhat,
		gelib::d_vector const& C_ZRhat, gelib::d_vector const& C_ZhatRhat, gelib::d_vector const& C_Xy, gelib::d_vector const& C_Ry,
		gelib::d_vector const& C_Rhaty, gelib::d_vector const& C_RhatR, gelib::d_vector const& C_Yz, gelib::d_vector const& C_Zz,
		gelib::d_vector const& C_XY, gelib::d_vector const& C_XZ, gelib::d_vector const& C_YZ, gelib::d_vector const& C_ZhatZhat);

	static void constructor_check(size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times, const std::array<double, 2>& discount_factors_t1, const std::array<gelib::d_vector, 2>& accrual_factors, const std::array<gelib::d_vector, 2>& forward_libors, const std::array<gelib::d_vector, 3>& process_variances, const std::array<gelib::d_vector, 3>& process_covariances, gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance, double spot_fx, gelib::d_matrix const& mu_data, gelib::d_matrix const& C_data);
};

template<class CalibratorType>
void HWMFMTwoCurrency<CalibratorType>::constructor_check(
	size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times,
	const std::array<double, 2>& discount_factors_t1, const std::array<gelib::d_vector, 2>& accrual_factors,
	const std::array<gelib::d_vector, 2>& forward_libors, const std::array<gelib::d_vector, 3>& process_variances,
	const std::array<gelib::d_vector, 3>& process_covariances, gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance,
	double spot_fx,
	gelib::d_vector const& mu_X, gelib::d_vector const& mu_Y, gelib::d_vector const& mu_Zhat,
	gelib::d_vector const& mu_Z, gelib::d_vector const& mu_R, gelib::d_vector const& mu_Rhat,
	gelib::d_vector const& mu_y, gelib::d_vector const& mu_z, gelib::d_vector const& C_XX,
	gelib::d_vector const& C_ZZ, gelib::d_vector const& C_YY, gelib::d_vector const& C_zz,
	gelib::d_vector const& C_yy, gelib::d_vector const& C_RR, gelib::d_vector const& C_RhatRhat,
	gelib::d_vector const& C_XR, gelib::d_vector const& C_YR, gelib::d_vector const& C_ZR,
	gelib::d_vector const& C_XRhat, gelib::d_vector const& C_YRhat, gelib::d_vector const& C_ZRhat,
	gelib::d_vector const& C_ZhatRhat, gelib::d_vector const& C_Xy, gelib::d_vector const& C_Ry,
	gelib::d_vector const& C_Rhaty, gelib::d_vector const& C_RhatR, gelib::d_vector const& C_Yz,
	gelib::d_vector const& C_Zz, gelib::d_vector const& C_XY, gelib::d_vector const& C_XZ,
	gelib::d_vector const& C_YZ, gelib::d_vector const& C_ZhatZhat)
{
	MFMGeneral<CalibratorType>::constructor_check(gridpoints_per_slice, sd_spread, poly_order, times);

	std::array<size_t, 45> sizes = {
		times.size(), accrual_factors[0].size(), accrual_factors[1].size(),
		forward_libors[0].size(), forward_libors[1].size(),
		process_variances[0].size(), process_variances[1].size(), process_variances[2].size(),
		process_covariances[0].size(), process_covariances[1].size(), process_covariances[2].size(),
		y_larrow_mean.size(), y_larrow_variance.size(),
		mu_X.size(),mu_Y.size(),  mu_Zhat.size(),  mu_Z.size(),
		mu_R.size(),  mu_Rhat.size(), mu_y.size(),  mu_z.size(),
		C_XX.size(),  C_ZZ.size(), C_YY.size(), C_zz.size(),
		C_yy.size(),  C_RR.size(),  C_RhatRhat.size(),  C_XR.size(),
		C_YR.size(), C_ZR.size(), C_XRhat.size(),  C_YRhat.size(),
		C_ZRhat.size(), C_ZhatRhat.size(), C_Xy.size(),  C_Ry.size(),
		C_Rhaty.size(), C_RhatR.size(), C_Yz.size(), C_Zz.size(),
		C_XY.size(),  C_XZ.size(), C_YZ.size(), C_ZhatZhat.size() };

	if (!std::equal(sizes.begin() + 1, sizes.end(), sizes.begin()))
		throw std::runtime_error("Every column vector input must be the same length, (= number of slices T_1 -> T_n))");
}

template<class CalibratorType>
void HWMFMTwoCurrency<CalibratorType>::constructor_check(size_t gridpoints_per_slice, double sd_spread, size_t poly_order, gelib::d_vector const& times, const std::array<double, 2>& discount_factors_t1, const std::array<gelib::d_vector, 2>& accrual_factors, const std::array<gelib::d_vector, 2>& forward_libors, const std::array<gelib::d_vector, 3>& process_variances, const std::array<gelib::d_vector, 3>& process_covariances, gelib::d_vector const& y_larrow_mean, gelib::d_vector const& y_larrow_variance, double spot_fx, gelib::d_matrix const& mu_data, gelib::d_matrix const& C_data)
{
	if (mu_data.cols() != 8)
		throw std::runtime_error("mu_data should contain 8 columns");
	if (C_data.cols() != 23)
		throw std::runtime_error("C_data should contain 23 columns");

	HWMFMTwoCurrency::constructor_check(
		gridpoints_per_slice, sd_spread, poly_order, times, discount_factors_t1, accrual_factors, forward_libors,
		process_variances, process_covariances, y_larrow_mean, y_larrow_variance, spot_fx,
		mu_data.extract_vector(-1, 0), mu_data.extract_vector(-1, 1), mu_data.extract_vector(-1, 2), mu_data.extract_vector(-1, 3),
		mu_data.extract_vector(-1, 4), mu_data.extract_vector(-1, 5), mu_data.extract_vector(-1, 6), mu_data.extract_vector(-1, 7),
		C_data.extract_vector(-1, 0), C_data.extract_vector(-1, 1), C_data.extract_vector(-1, 2), C_data.extract_vector(-1, 3),
		C_data.extract_vector(-1, 4), C_data.extract_vector(-1, 5), C_data.extract_vector(-1, 6), C_data.extract_vector(-1, 7),
		C_data.extract_vector(-1, 8), C_data.extract_vector(-1, 9), C_data.extract_vector(-1, 10), C_data.extract_vector(-1, 11),
		C_data.extract_vector(-1, 12), C_data.extract_vector(-1, 13), C_data.extract_vector(-1, 14), C_data.extract_vector(-1, 15),
		C_data.extract_vector(-1, 16), C_data.extract_vector(-1, 17), C_data.extract_vector(-1, 18), C_data.extract_vector(-1, 19),
		C_data.extract_vector(-1, 20), C_data.extract_vector(-1, 21), C_data.extract_vector(-1, 22), C_data.extract_vector(-1, 1)
	);
}
