RAPID RADIOMETRIC ENHANCEMENT OF COLORED 3D POINT CLOUDS USING COLOR BALANCING

by

 

Modules for radiometric enhancement of colored 3D point clouds using color balancing

cbalance.m
function [enhanced] = cbalance(clouds,varargin)
%CBALANCE Balance colors of overlapping 3D point clouds.
%   [ENHANCED] = cbalance(CLOUDS,BASE,METHOD,...) balances the colors
%   of -preferably but not necessarily- overlapping 3D point clouds in such
%   a way that the color variances in overlapping parts of the point clouds
%   are reduced.
%
%   Point clouds must be stored in CLOUDS in a cell array of structs with
%   fields <points> and <colors>. <colors> will be reduced to an intensity
%   value by using the maximum color component.
%
%   BASE is the index of the point cloud which remains unchanged. BASE is
%   by default set to 1, which means the corresponding linear coefficients
%   are normalized using first point cloud in the collection.
%
%   METHOD can be one of the following:
%      1. 'none'    No overlap detection performed. All colors are of the
%                   point clouds are used for color balancing. If no METHOD
%                   is specified, 'none' is used.
%      2. 'box'     All points which fall into the intersection of bounding
%                   volumes of point clouds are used for color balancing.
%                   Bounding box boundary tolerance TOL can also given as
%                   input to the function. By default, TOL is 0.
%      3. 'random'  Randomly selected points which fall into the
%                   intersection of bounding volumes of point clouds are
%                   used for color balancing. For each candidate
%                   neighbouring points in the other point cloud are
%                   searched in a brute-force manner and matching pairs are
%                   used for color balancing. Matching tolerance TOL and
%                   number of samples NSAMPLES can also given as input to
%                   the function. By default, TOL is ? and NSAMPLES is 100.
%      4. 'precise' All points which have a neighbouring point in the other
%                    point cloud is used for color balancing. For each
%                    point the other point cloud is searched in a
%                    brute-force manner and matching pairs are extracted.
%                    Matching tolerance TOL can also given as input to the
%                    function. By default, TOL is ?.
%
%   See also KNNSEARCH, FIND
%
%   This implementation is based on "Rapid Radiometric Enhancement of
%   Colored 3D Point Clouds using Color Balacing", U. Yilmaz and O.
%   Hellwich, Proceedings of 3DTV-CON, pp. 1-4, 2010.
%   See www.cv.tu-berlin.de for further info.
%
% Copyright 2010  This file and its content belong to Ulas Yilmaz.
% You are welcome to use it for non-commercial purposes, such as
% student projects, research and personal interest. However,
% you are not allowed to use it for commercial purposes, without
% an explicit written and signed license agreement with the owner.
% Contact info at http://www.cv.tu-berlin.de
% Berlin University of Technology  Germany
% Last Modification: 06.08.2010
%

% Input
if nargin==0
    error('Too few input arguments.');
elseif nargin==1
    base = 1;
    method = 'none';
elseif nargin==2
    if ischar(varargin{1})
        base = 1;
        method = varargin{1};
        if isequal(method,'box') tol = 0; else tol = Inf; end;
        nsamples = 100;
    else
        base = varargin{1};
        method = 'none';
    end;
elseif nargin==3
    if ischar(varargin{1})
        base = 1;
        method = varargin{1};
        tol = varargin{2};
        nsamples = 100;
    else
        base = varargin{1};
        method = varargin{2};
        if isequal(method,'box') tol = 0; else tol = Inf; end;
        nsamples = 100;
    end;
elseif nargin==4
    if ischar(varargin{1})
        base = 1;
        method = varargin{1};
        tol = varargin{2};
        nsamples = varargin{3};
    else
        base = varargin{1};
        method = varargin{2};
        tol = varargin{3};
        nsamples = 100;
    end;
elseif nargin==5
    base = varargin{1};
    method = varargin{2};
    tol = varargin{3};
    nsamples = varargin{4};
elseif nargin>5
    error('Too many input arguments.');
end;

% Bounding boxes
nclouds = numel(clouds);
lcorners = zeros(nclouds,3);
ucorners = zeros(nclouds,3);
for n = 1:nclouds
    lcorners(n,:) = min(clouds{n}.points);
    ucorners(n,:) = max(clouds{n}.points);
end;

% Variances and means
S = zeros(nclouds);
M = zeros(nclouds);

% Overlaps
message = 'Please wait...';
h = waitbar(0,message,'Name','Finding overlaps...');
for n = 1:nclouds
    message = sprintf('Finding overlaps... %.0f complete!',100*n/nclouds);
    waitbar(n/nclouds,h,message);
    for m = (n+1):nclouds
        if isequal(lower(method),'none')
            % Without overlap detection
            b = 1:size(clouds{m}.points,1);
            a = 1:size(clouds{n}.points,1);
        else
            % With overlap detection
            lcorner = max([lcorners(n,:);lcorners(m,:)]);
            ucorner = min([ucorners(n,:);ucorners(m,:)]);
            if all(ucorner-lcorner>0)
                switch lower(method)
                    case {'box'}
                        b = find(inbox(clouds{m}.points,lcorner,ucorner,tol));
                        a = find(inbox(clouds{n}.points,lcorner,ucorner,tol));
                    case {'random','precise'}
                        nquery = size(clouds{m}.points,1);
                        if isequal(lower(method),'precise') nsamples = nquery; end;
                        s = randsample(nquery,nsamples);
                        [indices,distances] = knnsearch(clouds{n}.points,clouds{m}.points(s,:),'k',1,'distance','euclidean');
                        indices(distances>tol) = 0;
                        b = find(indices>0);
                        a = indices(b);
                        b = s(b);
                    otherwise
                        error('Unknown method.');
                end;
            else
                [a,b] = deal([]);
            end;
        end;
        % Variance and mean
        if ~isempty(a) && ~isempty(b)
            acolors = double(max(clouds{n}.colors(a,:),[],2));
            bcolors = double(max(clouds{m}.colors(b,:),[],2));
            S(n,m) = var(acolors);
            S(m,n) = var(bcolors);
            M(n,m) = mean(acolors);
            M(m,n) = mean(bcolors);
        end;
    end;
end;
delete(h);

% Variance coefficients
warning('off');
D = abs(null(S'-diag(sum(S,2))));
a = sqrt(D/D(base));
warning('on');

% Mean coefficients
warning('off');
A = ones(nclouds) - nclouds * eye(nclouds);
Y = diag(a) * M;
B = sum(Y,2)-sum(Y,1)';
b = A\B;
b(isnan(b)) = 0;
b = b-b(base);
warning('on');

% Linear transformations
enhanced = cell(nclouds,1);
for m = 1:nclouds
    enhanced{m}.points = clouds{m}.points;
    enhanced{m}.colors = uint8(a(m)*single(clouds{m}.colors) + b(m));
end;

end

%--------------------------------------------------------------------------
function [mask] = inbox(cloud,lcorner,ucorner,tol)
%--------------------------------------------------------------------------

if nargin==4
    lcorner = lcorner - tol;
    ucorner = ucorner + tol;
end;

mask = all([...
    cloud(:,1)>=lcorner(1),cloud(:,2)>=lcorner(2),cloud(:,3)>=lcorner(3),...
    cloud(:,1)<=ucorner(1),cloud(:,2)<=ucorner(2),cloud(:,3)<=ucorner(3)...
    ],2);

end

Contact us