% lightens or darkens image based on range-range of my empirical law
%
% parameters based on mannfamily coon series
%
% steve mann, 1997
%
% g=lightdark(f,4,2.2*2,.3/2)
% g=lightdark(f,4,2.2*3,.3/3) % is sharper
%
% currently i am still trying to figure out why NaN and 0 in vector doesnt work
% e.g. why is it that:
%>>[NaN Inf].^(1/2)
%>>ans =
%      NaN +    NaNi      NaN +    NaNi

function g = f(f,k,a,c,method)

if nargin<5
  method=2;
  dis('method defaulting to %d',method)
end%if

if nargin<4
  if method==1 % exp
    c=0.1423;
  end%if
  if method==2 % linear
    c=.3/2;
  end%if
  dis('c defaulting to %d',c)
end%if

if nargin<3
  if method==1
    a=3.8783;
  end%if
  if method==2
    a=2.2*2;
  end%if
  dis('a defaulting to %d',a)
end%if

if nargin<2
  k=4;
  dis('k defaulting to %d',k)
end%if

%g = exp(1../( 1+1/(k^a)*( log(f).^(1/c) - 1 ) ));

%g=h(1../( 1+1/(k^a)*( (1../hinv(f)).^(1/c) - 1 ) ));

% let h(x)=log(x+1) and hinv(x)=exp(x)-1

f=f./255;
%whydoessqrtofnanandzerogivenanandnan=(f~=0);
whicharenan=isnan(f);
f=chnan(f,127);
fl = (1../( 1+1/(k^a)*( (1../(f)).^(1/c) - 1 ) )).^c; % f limited
%fl=whydoessqrtofnanandzerogivenanandnan.*fl;
fl=fl + chval(whicharenan,1,NaN);

g=fl.*255;  % scale back up again

if ~isreal(g)
  disp('lightdark.m: warning: output image is complex; taking real part')
  g=real(g); % in case there are NaN in f, then output would be complex NaN
end%if

%%keyboard