Re: [OctDev] signed distance of point from polygon

[JuanPi <aj...@gm...>]

On Tue, Sep 11, 2012 at 8:46 AM, c. <car...@gm...> wrote:
>
> On 10 Sep 2012, at 13:20, JuanPi wrote:
>
>> On Mon, Sep 10, 2012 at 1:01 PM, c. <car...@gm...> wrote:
>>> Hi,
>>>
>>> I need to compute the "signed distance" of a point P from a polygon Q,
>>> i.e., the distance of P from the polyline including Q with a "+" sign if
>>> the point is inside the polygon and a "-" sign if it is outside.
>>>
>>> I see "distancePointPolygon" in the geometry package computes the distance
>>> but always returns a positive number. Is there already a function to do
>>> what I need?
>>>
>>> If not, any idea where I can look-up the algorithm to modify "distancePointPolygon"
>>> according to my needs?
>>>
>>> Thanks,
>>> c.
>>
>> Hi,
>>
>> I would not modify the function, but use it together with inpolygon
>> (in the core geomtry package).
>>
>> I guess this is what you want
>> http://agora.octave.org/snippet/FRXH/
>>
>> Let me know
>
> I found that the function "distancePointPolyline" contains a for loop that can be easily vectorized away.
> As I work with a large number of points this has a strong impact on running time for me.
> Do you mind if I commit this patch?
> c.
>
> Index: inst/polygons2d/distancePointPolyline.m
> ===================================================================
> --- inst/polygons2d/distancePointPolyline.m     (revision 10796)
> +++ inst/polygons2d/distancePointPolyline.m     (working copy)
> @@ -1,5 +1,6 @@
> ## Copyright (C) 2003-2011 David Legland <dav...@gr...>
> ## Copyright (C) 2012 Adapted to Octave by Juan Pablo Carbajal <car...@if...>
> +## Copyright (C) 2012 Carlo de Falco (Speed up by vectorization)
> ## All rights reserved.
> ##
> ## Redistribution and use in source and binary forms, with or without
> @@ -51,17 +52,14 @@
>   # allocate memory for result
>   minDist = inf * ones(Np, 1);
>
> -  # process each point
> -  for p = 1:Np
> -      # construct the set of edges
> -      edges = [poly(1:end-1, :) poly(2:end, :)];
> +  ## construct the set of edges
> +  edges = [poly(1:end-1, :) poly(2:end, :)];
>
> -      # compute distance between current each point and all edges
> -      dist = distancePointEdge(point(p, :), edges);
> +  ## compute distance between current each point and all edges
> +  dist = distancePointEdge(point, edges);
>
> -      # update distance if necessary
> -      minDist(p) = min(dist);
> -  end
> +  ## update distance if necessary
> +  minDist = min(dist, [], 2);
>
>   # process output arguments
>   if nargout<=1
>
>

Please go ahead!. Can you leave the original code commented at the end
of the function? If you can't do not worry.
Thank you!

-- 
JuanPi Carbajal
-----
"The bad economist pursues a small present good, which will be
followed by a great evil to come, while the true economist pursues a
great good to come, at the risk of a small present evil." - Frédéric
Bastiat
-----
http://ailab.ifi.uzh.ch/carbajal/