view general/numerical/scalar/lambertw.m @ 61:eff6bddf82e3 tip

Finally implemented perceptual brightness thing.

function W = lambertw(x)
% lambertw - Lambert's W function
%
% lambertw :: real -> real.
%
% evaluates W(y) for -1 branch of Lambert W function where
% y = -exp(-x)
%
% The W function satisfies W(t)*exp(W(t)) = t

if x<1, error('lambertw: argument must be >= 1'); end
W=fixpoint(@(w)(-x-log(abs(w))),-x,'its',60);