function varargout = hmode(varargin) %#codegen

%% DEFAULT PARAMETERS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Input processing

if nargin == 0,
  varargout{1} = hmode('gradient'); % load default params
  return %empty call, probably to get default structures
elseif nargin == 1
  % get default parameters for each case
  def_params.modeltype = varargin{1};
  switch def_params.modeltype
    case 'imposed'
      def_params.params = Hmode_imposed;
    case 'gradient'
      def_params.params = Hmode_gradient;
    case 'feedback'
      def_params.params = Hmode_feedback;
    case 'none'
      def_params.params.active = false;
    otherwise
      error('invalid hmode model type ''%s''',def_params.modeltype)
  end
  varargout{1} = def_params; % return appropriate parameter structure
  return
else
  stap = varargin{1};
  geop = varargin{2};
  trap = varargin{3};
  model = varargin{4};
  pp = varargin{5};
  channel = varargin{6};
end

%% Parse channel
[species,eq] = parsechannel(channel);

%% Compute output
switch model.hmode.modeltype
  case 'imposed'
    switch eq
      case 'temperature'
        % no hmode modification of chie
        varargout{1} = varargin{7};
        % no hmode modification of dchie_dx
        varargout{2} = varargin{8};
        % no hmode modification of dchie_dxdot
        varargout{3} = varargin{9};
        % no hmode modification of dchie_du
        varargout{4} = varargin{10};
      case 'density'
        % no hmode modification of dne
        varargout{1} = varargin{7};
        % no hmode modification of vne
        varargout{2} = varargin{9};
        % no hmode modification of ddne_dx
        varargout{3} = varargin{8};
        % no hmode modification of dvne_dx
        varargout{4} = varargin{10};
    end
  case 'gradient'
    switch eq
      case 'temperature'
        chi0 = varargin{7};
        dchi0_dx = varargin{8};
        dchi0_dxdot = varargin{9};
        dchi0_du = varargin{10};
        [varargout{1},varargout{2},varargout{3},varargout{4}] = Hmode_gradient(stap,geop,trap,model,pp,species,eq,chi0,dchi0_dx,dchi0_dxdot,dchi0_du);
      case 'density'
        dne0 = varargin{7};
        ddne0_dx = varargin{8};
        vne0 = varargin{9};
        dvne0_dx = varargin{10};
        [varargout{1},varargout{2},varargout{3},varargout{4}] = ...
          Hmode_gradient(stap,geop,trap,model,pp,species,eq,dne0,ddne0_dx,vne0,dvne0_dx);
    end
  case 'feedback'
    switch eq
      case 'temperature'
        chi0 = varargin{6};
        dchi0_dx = varargin{7};
        [varargout{1},varargout{2}] = Hmode_feedback(stap,geop,model,pp,species,chi0,dchi0_dx);
      case 'density'
        error('Hmode_feedback for density equation not supported');
    end
  otherwise
    error('invalid setting for echcd model');
end


end

function [species,equation] = parsechannel(channel)

switch channel
  case 'te';
    species = 'e'; equation='temperature';
  case 'ti';
    species = 'i'; equation='temperature';
  case 'ne';
    species = 'e'; equation='density';
  case 'ni';
    error('not tested yet');
    %species = 'i'; equation='n';
  otherwise
    error('unknown channel');
end
end