determine whether two points with known normals are facing each other or not (matlab)

Question

I am trying to find a solution to the following problem but it is not clear how to solve it. Imagine that I have the following points in the space as in the image below:

If I consider that my only known information is the point positions and their normals I would like to determine whether two points (considering as a reference the position of the first point) are facing each other or not. For example from the above image for points a, b, c, d and e I have:

Point a faces points c and e but not points b and d.

Point b faces points d and e but not points a and c.

Point c faces point a but not points b, d and e.

Point d faces points b and e but not points a and c.

and finaly

Point e faces points a, b and d but not point c.

My first though was to play with the signed angles between the two normal vectors of each pair by using the solutions proposed here but this works for some pairs while not for others. The idea regarding what two point are facing each other is that if we consider a point as the origin then it faces another point if the other point is within the origin point's 180 degrees field of view and its normal vector is going inwards (kind of "towards") the origin point.

Any ideas what could help.

Thanks.

---

Update:

to try to be a bit more clear and answer some of the comments below. In principle its point in the space corresponds to the centroid of a face. However, I do not have this information beforehand (i.e. that each point corresponds to the center of a face, or the list of faces and their vertices). So in a higher level, if we were dealing with faces the problem would be how to determine if two faces are visible to each other or not but as I said the only information that I have now are the actual points in the space and their normals.

Sample points:

a = [26415.3720833199 11986.0504166605 739];
na = [0 0 1];

b = [27263.8100000023 11103.1983333336 1512.50000000021];
nb = [0.102791963903622 -0.994702876318771 0];

c = [28059.5700000001 11185.4316666667 962.499999999998];
nc = [-0.102791963903623 0.994702876318771 -9.06557542353252e-16];

d = [26606.7112499615 10390.7487916521 739];
nd = [0 0 1];

e = [27792.4499999996 9225.36499999984 2782];
ne = [0 0 -1];

Answer 1

You can solve your problem with a few simple dot products...

Based on your description, a point b is within the field of view (FOV) of another point a if the angle between the normal of a (i.e. na) and a vector going from a to b is less than or equal to 90 degrees. As described here, the angle can be found by taking the dot product of b-a and na, dividing by the length of b-a (and assuming the length of na is already 1), and taking the inverse cosine of the result. Putting it into an anonymous function, you have:

isInFOV = @(b, a, na) (acosd(dot(b-a, na)./norm(b-a)) <= 90);

You can then define a point b as "pointing toward" another point a if the component of nb (the normal of b) running along the vector going from b to a is positive. As described here, the component can be found by taking the dot product of a-b and nb and dividing by the length of a-b (and assuming the length of nb is already 1). Putting it into an anonymous function, you have:

isPointingToward = @(b, nb, a) (dot(a-b, nb)./norm(a-b) > 0);

We can then define whether a point a is "facing" another point b as:

isFacing = @(a, na, b, nb) (isInFOV(b, a, na) && isPointingToward(b, nb, a));

Note that I used the logical short circuit AND operator && since isPointingToward doesn't need to be evaluated if isInFOV already evaluates to false.

Vectorizing

You can reformulate the above equations to vectorize the operations, using functions like bsxfun or replacing a call to dot with standard matrix operations. This will allow you to check which points in a set a given point is facing. A vectorized version of the function isFacing is given below:

function index = isFacing(a, na, b, nb)

  V = bsxfun(@minus, b, a);                     % Compute b-a for all b
  V = bsxfun(@rdivide, V, sqrt(sum(V.^2, 2)));  % Normalize each row
  index = (acosd(V*na.') <= 90);                % Find points in FOV of a
  index(index) = (sum(V(index, :).*nb(index, :), 2) < 0);  % Of those points in FOV,
                                                           %   find those pointing
                                                           %   towards a

end

Example

Using the sample data in the question:

pointMat = [26415.3720833199 11986.0504166605 739; ...               % Point a
            27263.8100000023 11103.1983333336 1512.50000000021; ...  % Point b
            28059.5700000001 11185.4316666667 962.499999999998; ...  % Point c
            26606.7112499615 10390.7487916521 739];                  % Point d
normalMat = [0 0 1; ...
             0.102791963903622 -0.994702876318771 0; ...
             -0.102791963903623 0.994702876318771 -9.06557542353252e-16; ...
             0 0 1];
p = [27792.4499999996 9225.36499999984 2782];  % Point e
np = [0 0 -1];

>> isFacing(p, np, pointMat, normalMat)

ans =

  4×1 logical array

   1    % Facing a
   1    % Facing b
   0    % Not facing c
   1    % Facing d