function [params_vec,fval,exitflag,output] = gate_fit(alpha,beta,g_ss,tau,V,T)
%% Fit bond graph parameters to model
% params: [kf, zf, kr, zr];
error_func_alpha = @(params) square_error(alpha - calc_alpha(params,V/1000,T));
error_func_beta = @(params) square_error(beta - calc_beta(params,V/1000,T));
error_func_gss = @(params) square_error(g_ss - p2gss(params,V,T));
error_func_tau = @(params) square_error(tau- p2tau(params,V,T));
error_func = @(params) error_func_alpha(params) + error_func_beta(params) + error_func_gss(params) + error_func_tau(params);
A = [];
b = [];
Aeq = [];
beq = [];
lb = [0; -10; 0; -10];
ub = [Inf; 10; Inf; 10];
options_ps = optimoptions('particleswarm','UseParallel',true,'HybridFcn',@fmincon);
[params_vec,fval,exitflag,output] = particleswarm(error_func,4,lb,ub,options_ps);