/*
 * SUMMARY:      MassBalance.c - calculate basin-wide mass balance
 * USAGE:        Part of DHSVM
 *
 * AUTHOR:       Mark Wigmosta
 * ORG:          Battelle - Pacific Northwest National Laboratory
 * E-MAIL:       ms_wigmosta@pnl.gov
 * ORIG-DATE:    Oct-96
 * Last Change: Thu Dec  5 10:34:17 2002 by Scott Waichler <waichler@tuff.pnl.gov>
 * DESCRIPTION:  Calculate water mass balance errors
 *               
 * DESCRIP-END.
 * FUNCTIONS:    FinalMassBalance()
 * COMMENTS:
 */

#ifndef lint
static char vcid[] = "$Id: FinalMassBalance.c,v 1.4 2001/05/15 17:30:25 waichler Exp $";
#endif /* lint */

#include <stdio.h>
#include <stdlib.h>
#include "settings.h"
#include "data.h"
#include "DHSVMerror.h"
#include "functions.h"
#include "constants.h"

/*****************************************************************************
  Aggregate()
  
  Calculate the average values for the different fluxes and state variables
  over the basin.  Only the runoff is calculated as a total volume instead
  of an average.  In the current implementation the local radiation
  elements are not stored for the entire area.  Therefore these components
  are aggregated in AggregateRadiation() inside MassEnergyBalance().

  The aggregated values are set to zero in the function RestAggregate,
  which is executed at the beginning of each time step.
*****************************************************************************/
void FinalMassBalance(FILES *Out, AGGREGATED *Total, WATERBALANCE *Mass)
{
  float NewWaterStorage;        /* water storage at the end of the time step */
  float Input;					/* total water flux entering the basin;  */
  float Output;                 /* total water flux leaving the basin;  */
  float MassError;              /* mass balance error m  */

  NewWaterStorage = Total->Runoff + Total->CanopyWater + Total->SoilWater + 
    Total->Snow.Swq + Total->Soil.SatFlow + Total->GwDepth;

  /* Note for Input that SnowVaporFlux is positive for condensation on snow
  Remember, the watershed control volume doesn't include the channels.
  That's why culvert return flow is considered an input.  SRW 11-17-98   */
  Input = Mass->CumPrecipIn + Mass->CumSnowVaporFlux + Mass->CumCulvertReturnFlow
          + Mass->CumChannelLoss;
  /*Output = Mass->CumChannelInt + Mass->CumRoadInt + Mass->CumET;*/
  Output = Mass->CumChannelInt + Mass->CumRoadInt + Mass->CumET + Mass->CumGwDeepLoss;
  
  /* in - out = change in storage.  */
  MassError = Input - Output - (NewWaterStorage - Mass->StartWaterStorage);

  printf( "\nFinal Mass Balance\n");
  printf( " StartWaterStorage (mm): %f\n", Mass->StartWaterStorage * 1000.); 
  printf( " NewWaterStorage (mm): %f\n", NewWaterStorage * 1000.); 
  printf( " Input (mm): %f\n", Input * 1000.); 
  printf( " Output (mm): %f\n", Output * 1000.);
  printf( " Mass Error (mm): %f\n", MassError * 1000.); 
  printf( " Mass Error (percent of Input): %.3f\n\n", (MassError/Input) * 100.); 

  printf( "Cumulative fluxes\n");
  printf( " Precip (input) (mm): %f\n", Mass->CumPrecipIn * 1000.); 
  printf( " SnowVaporFlux (input) (mm): %f\n", Mass->CumSnowVaporFlux * 1000.);   
  printf( " ChannelInt (output) (mm): %f\n", Mass->CumChannelInt * 1000.);
  printf( " ET (output) (mm): %f\n", Mass->CumET * 1000.); 
  printf( " Soil Evap (mm): %f\n", Mass->CumSoilEvap * 1000);
  printf( " Recharge to Deep Groundwater Zone (mm): %f\n", Mass->CumGwRecharge * 1000 );
  printf( " Loss from Deep Groundwater Zone (ultimate recharge) (mm): %f\n", 
            Mass->CumGwDeepLoss * 1000 );
  printf( " Runoff ratio (ChannelInt/(Precip+SnowVapor)): %.2f\n", 
          (Mass->CumChannelInt/(Mass->CumPrecipIn + Mass->CumSnowVaporFlux)) );
  printf( " Channel Infiltration Loss (mm): %f\n", Mass->CumChannelLoss * 1000 );

  printf( " \nCulvert and Runoff Information\n");
  printf( " RoadInt (output) (mm): %f\n", Mass->CumRoadInt * 1000.); 
  printf( " CulvertReturnFlow (input) (mm): %f\n", Mass->CumCulvertReturnFlow * 1000.);
  printf( " CulvertToChannel (mm): %f\n", Mass->CumCulvertToChannel * 1000.);
  printf( " Average Runoff over all cells and time steps (mm): %f\n", Mass->CumRunoff * 1000. / 365); 
  printf( " RunoffToChannel (mm): %f\n", Mass->CumRunoffToChannel * 1000.);
  printf( " Final Runoff Depth (mm): %f\n", Total->Runoff * 1000.); 
  printf( " Volume to QuickFlow (mm): %f\n", Mass->CumQuickVol * 1000.); 
  printf( " Final Deep Groundwater Depth (mm): %f\n", Total->GwDepth * 1000.); 
  printf( " Total Deep Groundwater to Channel (mm): %f\n", Mass->CumGwToChannel * 1000.); 

  printf( "\n Final SatFlow (mm): %f\n", Total->Soil.SatFlow * 1000.); 
  printf( " Final CanopyWater (mm): %f\n", Total->CanopyWater * 1000.); 
  printf( " Final SoilWater (mm): %f\n", Total->SoilWater * 1000.); 
  printf( " Final Snow.Swq (mm): %f\n", Total->Snow.Swq * 1000.); 
}

