function out = mm_efficiency_measure(measured,simulated,efficiency,varargin)
%MM_EFFICIENCY_MEASURE Compute efficiency measures
% Use this function to compute efficiency measures used in hydrology
%
% Input:
%   measures    ... vector with measured/observed values
%   simulated   ... vector with simulated/modelled values
%   efficiency  ... efficiency abbreviation (string):
%                   'NSE': Nash_and_Sutcliffe (1970) River flow forecasting
%                   through conceptual models part A discussion of principles 
%                   'KGE2009': Gupta et al (2009) Decomposition of the mean
%                   squared error and NSE performance criteria:
%                   Implications for improving hydrological modelling 
%   varargin(1) ... scales used for KGE. Vector with 3 values. By default
%                   [1 1 1] = [scale r, scale alpha, scale beta].
%
% Output:
%   out         ... computed efficiency coefficient
%
% Example1:
%   out = mm_efficiency_measure(measured,simulated,'NSE')
%
% Example2:
%   out = mm_efficiency_measure(measured,simulated,'KGE2009')
%
% Example2:
%   out = mm_efficiency_measure(measured,simulated,'KGE2009',[0.4 0.3 0.3])
%
%
%                                                   M.Mikolaj, GFZ Potsdam

% Check for NaNs = remove NaNs
r = find(isnan(measured+simulated));
measured(r) = [];
simulated(r) = [];

% Computed efficiency based on user input
switch efficiency
    case 'NSE'
        out = 1 - sum((simulated-measured).^2)./sum((measured-mean(measured)).^2);
    case 'KGE2009'
        r = cov(simulated,measured);
        r = r(1,2)/(std(simulated)*std(measured));
        alpha = std(simulated)/std(measured);
        beta = mean(simulated)/mean(measured);
        % Use either default [1 1 1] scales or user input
        if nargin == 4
            s = varargin{1};
        else
            s = [1 1 1];
        end
        out = 1 - sqrt((s(1)*(r - 1))^2 + (s(2)*(alpha - 1))^2 + (s(3)*(beta - 1))^2);
end
    