function M=expmap(p1,p2,p3)
% expmap - Create exponentially spaced discretisation map
%
% expmap :: Min:real, Max:real, N:natural -> dmap(N).
% expmap :: [Min,Max]:[[1,2]], N:natural -> dmap(N).
%
% Parameters work almost but not quite like logspace: the
% third parameter is the number of GAPS, not the number of
% points, defining the map, So for a map with edges at 0:10,
% use linmap(0,10,10). To get the same edges from linspace,
% you would need linspace(0,10,11).

% unpack parameters
if nargin==2, min=log(p1(1)); max=log(p1(2)); bins=p2;
else min=log(p1); max=log(p2); bins=p3; end

dx=(max-min)/bins;
M=dmap(bins,@(t)map(bins,min,1/dx,t),@(t)revmap(min,dx,t));

function I=map(N,min,k,X)
	U=log(X);
	I=1+floor(k*(U-min));
	I(I<1)=-inf;
	I(I>N)=inf;

function X=revmap(min,ik,I)
	I=shiftdim(shiftdim(I),-1);
	U1=min+(I-1)*ik;
	U2=min+(I)*ik;
	X=exp(cat(1,U1,U2));
	
