From b49752b15663b2a3a2674e6de863463cbc0e79b0 Mon Sep 17 00:00:00 2001 From: Elian Dib <elian.dib@inria.fr> Date: Thu, 27 Jun 2019 16:48:40 +0200 Subject: [PATCH] Fixed missing triangle intersections when using pointLocation due to floating point errors. The intersections are evaluated after shifting the reference and query values by the same offset in different directions which reduces the number of errors. --- gridCoords.m | 29 ++++++++++++++++++++++++++--- triangulation.m | 35 +++++++---------------------------- 2 files changed, 33 insertions(+), 31 deletions(-) diff --git a/gridCoords.m b/gridCoords.m index 69ec322..f1dac78 100644 --- a/gridCoords.m +++ b/gridCoords.m @@ -2,10 +2,34 @@ function [Xq,Yq,X,Y,varargout] = gridCoords(X,Y,Xq,Yq,varargin) %GRIDCOORDS Summary of this function goes here % Detailed explanation goes here +I = (1:numel(X))'; +offset = 1e-2; + +p = inputParser(); +p.addOptional('I' , I , @isnumeric); +p.addOptional('offset', offset, @isnumeric); + +p.parse(varargin{:}); +I = p.Results.I; +offset = p.Results.offset; + sz = size(Xq); -[t,I] = utils.triangulation(X,Y,varargin{:}); -[ind,b] = t.pointLocation(Xq(:),Yq(:)); +ind = nan(numel(Xq),1); +b = nan(numel(Xq),3); +nanmask = true(numel(Xq),1); + +[dX,dY] = ndgrid([-offset,0,offset]); +sigma = [5,4,6,2,8,1,3,7,9]; + +for r = 1:numel(sigma) + dx = dX(sigma(r)); + dy = dY(sigma(r)); + + t = utils.triangulation(X+dx,Y+dy,I); + [ind(nanmask),b(nanmask,:)] = t.pointLocation(Xq(nanmask)+dx,Yq(nanmask)+dy); + nanmask = isnan(ind); +end ind = reshape(ind,sz); u = reshape(b(:,2),sz); @@ -13,7 +37,6 @@ v = reshape(b(:,3),sz); nanind = isnan(ind); ind(nanind) = 1; -ind = I(ind); ut = ~mod(ind,2); ind(nanind) = nan; diff --git a/triangulation.m b/triangulation.m index ab50c80..7ce095f 100644 --- a/triangulation.m +++ b/triangulation.m @@ -1,16 +1,13 @@ -function varargout = triangulation(X,Y,varargin) +function t = triangulation(X,Y,varargin) %TRIANGULATION Summary of this function goes here % Detailed explanation goes here +I = (1:numel(X))'; p = inputParser(); -p.addOptional('Z' , [] , @isnumeric); -p.addOptional('order', 'none', @ischar); -p.addOptional('fun' , @(x) min(x,[],2)); +p.addOptional('I', I, @isnumeric); p.parse(varargin{:}); -Z = p.Results.Z; -order = p.Results.order; -fun = p.Results.fun; +I = p.Results.I; sz = size(X); @@ -20,26 +17,8 @@ TY = reshape( (Y_(:)+[0,0,1,1,1,0])' ,3,[])'; T = sub2ind(sz,TX,TY); -switch order - case 'none' - I = (1:size(T,1))'; - otherwise - Z_ = Z(T); - Z_ = fun(Z_); - [~,I] = sort(Z_,order); - T = T(I,:); -end - -if nargout>0 - varargout{1} = triangulation(T,X(:),Y(:)); -end - -if nargout>1 - varargout{2} = I; -end - -if nargout>2 - varargout{3} = triangulation(T,X(:),Y(:),Z(:)); -end +Ii(I) = 1:numel(I); + +t = triangulation(Ii(T),X(I(:)),Y(I(:))); end \ No newline at end of file -- GitLab