logo Sign In

Post #265267

Author
Laserman
Parent topic
Colour matching for fan edits
Link to post in topic
https://originaltrilogy.com/post/id/265267/action/topic#265267
Date created
11-Jan-2007, 7:01 AM
OK, I found some sample matlab code(not my code, but the base I used - I'll post credit it for it when I find whre I got it from originally), I'm still digging through my CDs to find the whitepapers, but this should give you an overview, you can see that the computation load is minimal.
Matlab is of course the perfect way to show this as the process stays linear and is a matrix fest, you could say matlab was written to do this

Edit: Tried to make smileys go away but didn't work.

%
% simple implementation of N-Dimensional PDF Transfer
%
% [DR] = pdf_transferND(D0, D1, rotations);
%
% D0, D1 = NxM matrix containing N-dimensional features
% rotations = { {R_1}, ... , {R_n} } with R_i PxN
%
% note that we can use more than N projection axes. In this case P > N
% and the inverse transformation is done by least mean square.
% Using more than N axes leads to a more stable (but also slower)
% convergence.
%
function [DR] = pdf_transferND(D0, D1, Rotations)

nb_dims = size(D0,1);

relaxation = 1;
DR = D0;

for it=1:length(Rotations)

R = Rotations{it};
nb_projs = size(R,1);

%% apply rotation

D0R = R * D0;
D1R = R * D1;

%% find data range
for i=1:nb_projs
datamin(i) = min([D0R(i, D1R(i,]);
datamax(i) = max([D0R(i, D1R(i,]);
end

%% get the marginals
for i=1:nb_projs
step = (datamax(i) - datamin(i))/300;
p0R{i} = hist(D0R(i,, datamin(i):stepatamax(i));
p1R{i} = hist(D1R(i,, datamin(i):stepatamax(i));
end

%% match the marginals
for i=1:nb_projs
f{i} = pdf_transfer1D(p0R{i}, p1R{i});
scale = (length(f{i})-1)/(datamax(i)-datamin(i));
D0R_(i, = interp1(0:length(f{i})-1, f{i}', (D0R(i, - datamin(i))*scale)/scale + datamin(i);
end

D0 = relaxation * R \ (D0R_ - D0R) + D0;
end

end
%%
%% 1D - PDF Transfer
%%
function f = pdf_transfer1D(pX,pY)
nbins = max(size(pX));

PX = cumsum(pX);
PX = PX/PX(end);

PY = cumsum(pY);
PY = PY/PY(end);

% inversion
PX = [0 PX nbins] + (0:nbins+1)*1e-10;
PY = [0 PY nbins] + (0:nbins+1)*1e-10;

f = interp1(PY, [0 ((0:nbins-1)+1e-16) (nbins+1e-10)], PX,'linear');
f = f(2:end-1);
end