#include "FlowSumCalc.h"
#include <math.h>
#include <fstream>

//sumRunoffToOutlet function sums the surface and subsurface runoff leaving the catchment and transfers water from GI devices to BulkArea
void FlowSumCalc::SumRunoffToOutlet(Inputs* input, DataFolder* folder, int timeStep)
{
	//Fraction of pervious area in folder
	double PerviousAreaOfFolder_frac;
	//Runoff (m) from saturation excess process
	double SatExQ;
	//Runoff (m) from infiltration excess process
	double infilExQ;
	//Runoff (m) from subsurface flow, TOPMODEL product
	double BQ;
	//Runoff (m) from pervious process
	double PAQ;
	//Runoff (m) from directly connected impervious areas
	double DCIAQ;
	//Runoff (m) sum of BQ, PAQ, DCIAQ
	double totalQ;

	//Set initial runoff sum values to zero
	double sumBQ = 0.0;
	double sumPAQ = 0.0;
	double sumDCIAQ = 0.0;
	double sumTotalQ = 0.0;

	//Fraction of pervious area in folder
	PerviousAreaOfFolder_frac = 0.0;
	//Runoff (m) from saturation excess process
	SatExQ = 0.0;
	//Runoff (m) from infiltration excess process
	infilExQ = 0.0;
	BQ = 0.0;
	//Runoff (m) from pervious process
	PAQ = 0.0;
	//Runoff (m) from directly connected impervious areas
	DCIAQ = 0.0;
	//Runoff sum of BQ, PAQ, DCIAQ
	totalQ = 0.0;

	//Note, multiplication of variables by fractional area below is to adjust depths from one analysis area to entire catchment area
	//	Rationale:	1. variables arrive as depths for an analysis area (e.g., pervious area), which are not conserved when transfered to the catchment area
	//				2. multiplying these depths by their analysis area gives a volume, which are conserved when transfered between different areas
	//				3. multiplying these depths by their fractional area is equivalent to multiplying by their analysis area and dividing by the catchment area, getting depth for catchment area

	//Fractional pervious area is from pervious under tree, short veg, and soil areas
	PerviousAreaOfFolder_frac = folder->ParamDict["TreeCanopyCover_overPervious_frac"] + folder->ParamDict["ShortVegCover_noTreeCanopy_frac"] + folder->ParamDict["SoilCover_noTreeCanopy_frac"];

	//compute the DataFolder area (m2) based on Type, noting GI folder needs to use giNumber
	//BulkArea DataFolder area
	if (folder->ParamStringDict["Type"] == "BulkArea" && timeStep == 0)
	{
		//Fraction of folder relative to total catchment area
		folder->VarDict["AreaOfCatchment_frac"] = folder->ParamDict["Area"] / input->SimulationNumericalParams["CatchmentArea_m2"];
	}
	//GI DataFolder area
	else if (folder->ParamStringDict["Type"] != "BulkArea" && timeStep == 0)
	{
		//Fraction of folder relative to total catchment area; this approach separates pervious area for GIstorageDrain_m3 and total area for GIsurfaceOutflow_m3 flows
		folder->VarDict["AreaOfCatchment_frac"] = (folder->ParamDict["Area"] * folder->ParamDict["giNumber"]) / input->SimulationNumericalParams["CatchmentArea_m2"];
		folder->VarDict["AreaOfAllGIDevices_m2"] = (folder->ParamDict["Area"] * folder->ParamDict["giNumber"]);
		folder->VarDict["AreaPerviousOfAllGIDevices_m2"] = (folder->ParamDict["TotalPerviousCover_frac"] * folder->ParamDict["giNumber"]);
	}

	//SatExQ (m) as saturation excess runoff depth (m) over pervious area
	SatExQ = folder->VarDict["Runoff_SatExcess_m"];
	//infilExQ (m) as infiltration excess runoff depth (m) over pervious area, and adjusted for fraction of catchment with infiltration excess
	infilExQ = folder->VarDict["Runoff_PervInfilEx_m"];

	//Note, multiplication of variables by fractional area below is to adjust depths from one analysis area to entire folder or catchment area
	//	Rationale:	1. variables arrive as depths for an analysis area (e.g., pervious area), which are not conserved when transfered to the catchment area
	//				2. multiplying these depths by their analysis area gives a volume, which are conserved when transfered between different areas
	//				3. multiplying these depths by their fractional area is equivalent to multiplying by their analysis area and dividing by the catchment area, getting depth for catchment area

	//DCIAQ (m) as folder area depth (e.g., bulk area) of impervious area (IA) runoff going to outlet, transformed in ImperviousDepressionStorageCalc.cpp from IA depth with * folder->ParamDict["ImperviousCover_DrainsToOutlet_frac"]
	//	Adjusting to folder area was done earlier in order to send IA water into pervious area for runon
	DCIAQ = folder->VarDict["Runoff_ImperviousToOutlet_m"];

	//PAQ (m) as folder area depth (e.g., bulk area) of pervious area runoff from saturation excess and infiltration excess, in FlowSumCalc.cpp with * PerviousAreaOfFolder_frac
	PAQ = (SatExQ + infilExQ) * PerviousAreaOfFolder_frac;

	//BQ (m) as folder area depth (e.g., bulk area) of subsurface runoff, draining catchment area underlying pervious area where Drainage_UnsatZone_m (m) supplies aveSMD (m) soil moisture defict
	BQ = folder->VarDict["Runoff_Subsurface_m"] * PerviousAreaOfFolder_frac;

	//cte 20210520 All GI transfers of water to Bulk Area likely need to be represented as depths for catchment area (GI+BulkArea), not just Bulk Area
	//cte 20210520 Finish work to ensure this is actually adjusting DCIAQ to catchment depth and contributing to its discharge
	//Transfer of GI pervious area drainage water and total area surface water to Bulk Area DCIA runoff; using a volume to uphold conservation of mass in transfer
	//DCIAQ (m) entering this statement is in folder area depth (e.g., BulkArea) and it is combined with GI runoff volumes, divided by folder area to get equivalent folder area depth
	if (input->RepoDict["GI_FolderOperational"] = 1 && folder->ParamStringDict["Type"] == "BulkArea") {
		//GIstorageDrain_m3 is from GI pervious area and GIsurfaceOutflow_m3 is from GI total area, for all GI units, and are divided by Bulk Area total area to get folder average depth
		//	GIsurfaceOutflow_m3 and GIstorageDrain_m3 are from prior timestep, distributed across Bulk Area pervious area but should likely be distributed across Catchment pervious area
		DCIAQ = DCIAQ + (input->RepoDict["GIsurfaceOutflow_m3"] + input->RepoDict["GIstorageDrain_m3"]) / folder->ParamDict["Area"];
	}

	// totalQ (m) as folder area depth (e.g., bulk area), pervious area flow (m), DCIA fraction of impervious flow (m), subsurface baseflow (m) as catchment depths
	totalQ = PAQ + DCIAQ + BQ;

	//Note, multiplication of variables by fractional area below is to adjust depths from one analysis area to entire catchment area
	//Runoff as catchment depth (m) transformed from DataFolder fractional area depth, to maintain balance in BulkArea DataDrawer 
	//Impervious runof (m)
	DCIAQ *= folder->VarDict["AreaOfCatchment_frac"];
	//Pervious runoff (m)
	PAQ *= folder->VarDict["AreaOfCatchment_frac"];
	//Subsurface runoff (m)
	BQ *= folder->VarDict["AreaOfCatchment_frac"];
	//Total runoff (m) = BQ + PAQ + DCIAQ
	totalQ *= folder->VarDict["AreaOfCatchment_frac"];

	// create catchment runoff depths (m) by summing all DataFolder runoff within DataDrawer
	sumDCIAQ += DCIAQ;
	sumPAQ += PAQ;
	sumBQ += BQ;
	sumTotalQ += totalQ;

	// store catchment depth of non-routed (NR) runoff depth (m) terms for water balance output
	folder->VarDict["ImperviousQ_NR_m"] = sumDCIAQ;
	folder->VarDict["PerviousQ_NR_m"] = sumPAQ;
	folder->VarDict["BaseQ_NR_m"] = sumBQ;
	folder->VarDict["TotalQ_NR_m"] = sumTotalQ;

	//Store catchment runoff depth (m) as non-time distributed vector values (later distributed in FlowRouteCalc.cpp)
	//Subsurface runoff
	input->BQi_NR.push_back(sumBQ);
	//Pervious runoff 
	input->PAQi_NR.push_back(sumPAQ);
	//Directly connected impervious area runoff
	input->DCIAQi_NR.push_back(sumDCIAQ);

	//simulation totals are stored in subsequent StoreTimeStepTotalsCalc.cpp 
}
