logNsamplin.m

% function [xp yp] = logNsamp(x,f,N)
% x is the x-data with even linear spacing
% f is the function value for the x values
%     NOTE: f should be on linear scale, because logNsample takes log(f)
% N is the number of divisions
%
% xp is evenly spaced on log
% yp is function value on xp
% NOTE: Size of xp and yp is 1:N+1 which includes the endpoints
% xp(1) = log(x(1))  and xp(N+1) = log(x(szx))
%
% See also logevendist.m


function [xp yp] = logNsamplin(x,f,N)

[dum sz] = size(f);
n = 1:N+1;      %

logx = log(x);
%logf = log(f);

delx = x(2) - x(1);
delogx = diff(logx);

delxp = (log(x(sz)) - log(x(1)))/N;

xp = log(x(1)) + (n-1).*delxp;      %
yp = zeros(size(xp));

%keyboard

[X, I] = changesign(delxp - 2*delogx);

I;      %cross-over location between needing to interpolate versus needing to average
%ind = floor(logx(I)/delxp);
ind = floor((logx(I)-logx(1))/delxp);

%keyboard

yplo = intersamp(logx(1:I+1),f(1:I+1),xp(1:ind+1));
yphi = avgsamp(logx(I+1:sz),f(I+1:sz),xp(ind+1:N));

%keyboard

yp(1:ind+1) = (yplo(1:ind+1));
yp(ind+2) = (0.5*(yplo(ind+1)+ yphi(2)));
yp(ind+3:N+1) = (yphi(2:(N-ind)));

%keyboard



function [X, I] = changesign(a)

[dum sz1] = size(a);
b(1:(sz1-1)) = a(2:sz1);
tst = a(1:(sz1-1)).*b;

[X, Itmp] = min(tst);
I = Itmp + 1;



function y = intersamp(x,f,xp)

y = interp1(x,f,xp,'linear','extrap');



function y = avgsamp(x,f,xp)

[dum szx] = size(x);
[dum szp] = size(xp);
delxp = xp(2) - xp(1);

%keyboard

fcum = zeros(1,szp);
cnt = zeros(1,szp);
for loop = 1:szx
    bin = ceil((x(loop)-(xp(1) + 0.5*delxp))/delxp);  %
    %keyboard
    if bin>0
        fcum(bin) = fcum(bin) + f(loop);
        cnt(bin) = cnt(bin) + 1;
    end
end
for loop = 1:szp
    if cnt(loop) ~=0
        fav(loop) = fcum(loop)/cnt(loop);
    else
        fav(loop) = 0;
    end
end

y = fav;