function [q,dq_dx] = eval_q(x,g,v,model,gauss)
% function [q,dq_dx] = eval_q(x,g,~,model,whichgrid)
% #codegen

if ~gauss;
    ngrid = model.rgrid.nrho;
else
    ngrid = model.rgrid.nrhogauss;
end

if nargout>1
        der = true; % compute derivative
        dq_dx = zeros(ngrid,numel(x)); % init
else 
        der = false;
end

% q=1/iota
% dq/dx = -1/iota.^2 * diota_dx = -q.^2 * diota_dx

if der
    [iota,diota_dx] = eval_iota(x,g,v,model,gauss);
    q = 1./iota;
    dq_dx = bsxfun(@times,-q.^2,diota_dx);
else
    [iota] = eval_iota(x,g,v,model,gauss);
    q = 1./iota;
end

return
