IsInside points on boundaries

ROOT
1

When creating a TCutG are points on the boundaries considered inside or outside?

I created some tcutg and then use the points defining the tcutg to test whether they are inside or out and it seems that the condition is marked as true only in some cases. Looking at the source for IsInside I believe that the if statement checking the x position should have a less than or equal sign not just less than for the case that y[i] == yp.

  if (x[i]+(yp-y[i])/(y[j]-y[i])*(x[j]-x[i])<=xp) {
template <typename T> Bool_t TMath::IsInside(T xp, T yp, Int_t np, T *x, T *y)
{
   // Function which returns kTRUE if point xp,yp lies inside the
   // polygon defined by the np points in arrays x and y, kFALSE otherwise.
   // Note that the polygon may be open or closed.

   Int_t i, j = np-1 ;
   Bool_t oddNodes = kFALSE;

   for (i=0; i<np; i++) {
      if ((y[i]<yp && y[j]>=yp) || (y[j]<yp && y[i]>=yp)) {
         if (x[i]+(yp-y[i])/(y[j]-y[i])*(x[j]-x[i])<=xp) {
            oddNodes = !oddNodes;
         }
      }
      j=i;
   }

   return oddNodes;
}
2

Any comment?

3

Most of the time floating point equality isn’t useful since computers don’t work with exact representations of floating point numbers. For most realistic use cases (e.g. when the points are from measurements or results of some calculations) what happens exactly on the boundary doesn’t make a difference.

4

The boundaries can be very important in a measurement if the statistics are low. I agree there is some ambiguity about whether the boundary should be considered or not, but a consistent choice should be made.

5

The boundary status does not seems to be very clear as shown with the following example

#include <TMath.h>

void testinside()
{
   Double_t x[5] = {0,0,1,1,0};
   Double_t y[5] = {0,1,1,0,0};

   printf("%d\n",TMath::IsInside(0.5, 0.5, 5, x, y));
   printf("%d\n",TMath::IsInside(0.0, 0.5, 5, x, y));
   printf("%d\n",TMath::IsInside(1.0, 0.5, 5, x, y));
   printf("%d\n",TMath::IsInside(0.5, 0.0, 5, x, y));
   printf("%d\n",TMath::IsInside(0.5, 1.0, 5, x, y));
}

it gives:

root [0] .x testinside.C++
Info in <TUnixSystem::ACLiC>: creating shared library /Users/couet/roottest/./testinside_C.so
1
0
1
0
1
root [1]
6

Does this mean its a bug or a “feature”?

7

It is implemented that way. I guess you can call it a feature. What would you expect ?

8

I would have expected it to return either all 1s or all 0s. Why should have the points on the boundary be considered inside why the others are not? I’m just looking for a consistent treatment.

9

I will check the algorithm. I suspect the vertical line might be special.
I will let you know

10

Even in your example it treats the left and right vertical lines differently. First result was zero and second was one.

[quote=“couet”]The boundary status does not seems to be very clear as shown with the following example


   printf("%d\n",TMath::IsInside(0.0, 0.5, 5, x, y));
   printf("%d\n",TMath::IsInside(1.0, 0.5, 5, x, y));

[/quote]

11

Yes I noticed that

12

The algorithm used in isinside is fast but If the test point is on the border of the polygon, this algorithm will deliver unpredictable results;

alienryderflex.com/polygon/

I will see if there is a better one fast and with predictable edges status.

13

I’ve also noticed that points close to the boundary, but outside are sometimes considered inside. This is even worse than those points on the boundaries.

For example:

[code]#include <float.h>
#include <TMath.h>

void testinside()
{
Double_t x[5] = {0,0,1,1,0};
Double_t y[5] = {0,1,1,0,0};

printf("Inside Points\n");

printf("%d\n",TMath::IsInside(0.5, 0.5, 5, x, y));
printf("%d\n",TMath::IsInside(+DBL_MIN, 0.5, 5, x, y));
printf("%d\n",TMath::IsInside(1.0-DBL_MIN, 0.5, 5, x, y));
printf("%d\n",TMath::IsInside(0.5, +DBL_MIN, 5, x, y));
printf("%d\n",TMath::IsInside(0.5, 1.0-DBL_MIN, 5, x, y));

printf("Outside Points\n");

printf("%d\n",TMath::IsInside(-DBL_MIN, 0.5, 5, x, y));
printf("%d\n",TMath::IsInside(1.0+DBL_MIN, 0.5, 5, x, y));
printf("%d\n",TMath::IsInside(0.5, -DBL_MIN, 5, x, y));
printf("%d\n",TMath::IsInside(0.5, 1.0+DBL_MIN, 5, x, y));
}
[/code]

Results:

root [0] .x testinside.C+ Inside Points 1 1 1 1 1 Outside Points 0 1 0 1 root [1]

14

This test fails on your side. Since a double can only hold a limited number of significant digits you have a predictable rounding error.

assert(1.0 - DBL_MIN == 1.0 + DBL_MIN);
assert(1.0 - DBL_MIN == 1.0);

http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html

15

You are correct, I mistakenly used DBL_MIN when i meant to use DBL_EPSILON. With DBL_EPSILON the point is still incorrectly considered inside.

I discovered that TMath::IsInside expects floats and not doubles and therefore one should not use it if this kind of precision is necessary.

16

In fact it is a templated function. It accept anything:

template <typename T> Bool_t TMath::IsInside(T xp, T yp, Int_t np, T *x, T *y)
17

I should have made one more click through the documentation. I had stopped right at the list of members when it said it expected floats. Is there anyway to update the methods list to reflect that it is templated? (I should have already known it was templated from my first post.)

18

[quote=“ksmith”]With DBL_EPSILON the point is still incorrectly considered inside.
[/quote]

This is not true, I must have mixed my types and CINT let me get away with it. IsInside works as expected for very small numbers near the boundaries.

Sorry for the confusion.