Check Points within Rectangles

[vincentzack]

I have more than 40000 points coordinates and I need to determine whether the points are within about 100 named rectangles or not and get the name displaying on the sheet. The points may not lies in any rectangles.

Could anyone help to me to have the code to check it?

[Paul_Hossler]

Can you post a small workbook with sample data and rectangles?

[Kenneth Hobs]

It is hard to know just what you mean. A defined rectangle will have 4 sets of points. If you mean a rectangle as defined by an Excel shape, that will not relate to the same coordinate system as the point coordinates.

Duplicate post at: http://www.mrexcel.com/forum/excel-questions/952067-check-points-within-rectangles.html

[vincentzack]

Attached sample.

There are two sheets: 1. Point and 2. Area. I need to check whether the points in sheet (Point) fall in which area showing in sheet (Area) and display the name in the column C of Sheet (Point). The points may not fall in any rectangle.

Thank you very much!

[mikerickson]

Did my suggestion at the Cross post help?
http://www.mrexcel.com/forum/excel-questions/952067-check-points-within-rectangles.html

[mdmackillop]

Combined on one sheet for simplicity. BTW, your point positions don't relate to your table data.

[Paul_Hossler]

I think I really mis-understood the question

I thought all the rectangles were parallel to the axis

[vincentzack]

I think I really mis-understood the question

I thought all the rectangles were parallel to the axis

I'm sorry that I forget to mention the axis of the rectangles.

The axis of rectangles may not parallel to X-axis or Y-axis.

Could it possible to do it using VBA?

Thank you very much!

[Paul_Hossler]

Is the rectangle data correct?

Just looking a X1 it seems that P1 and P2 have the same X coordinates, some other issues also

16583

[vincentzack]

Is the rectangle data correct?

Just looking a X1 it seems that P1 and P2 have the same X coordinates, some other issues also

16583
Paul,

Sorry again! Please find the attached updated file for your reference.

I also added some data for your reference

[Kenneth Hobs]

This is what we call point in a polygon. There are several routines around that should work.


How do you want it shown? As a UDF, it would look like this for your first file. All the other areas are False.



Point
X
Y
Area X9


1
42.101
14.491
FALSE


2
38.664
14.491
FALSE


3
39.132
14.491
FALSE


4
39.132
14.982
FALSE


5
38.664
14.982
FALSE


6
38.664
14.491
FALSE


7
39.132
14.491
FALSE


8
39.132
14.982
FALSE


9
38.664
14.982
FALSE


10
39.132
14.491
FALSE


11
39.444
14.491
FALSE


12
39.444
14.982
FALSE


13
39.132
14.982
FALSE


14
39.132
14.491
FALSE


15
39.444
14.491
FALSE


16
39.444
14.982
FALSE


17
39.132
14.982
FALSE


18
39.444
14.491
TRUE


19
39.755
14.491
TRUE


20
12.131
21.536
TRUE


21
12.131
21.911
TRUE


22
11.667
21.911
TRUE


23
11.667
21.536
TRUE


24
12.131
21.536
TRUE


25
12.131
21.911
TRUE


26
11.667
21.911
TRUE


27
9.881
21.911
TRUE


28
10.310
21.911
TRUE


29
10.310
22.286
TRUE


30
9.881
22.286
TRUE


31
9.881
21.911
TRUE


32
10.310
21.911
TRUE


33
10.310
22.286
TRUE


34
9.881
22.286
TRUE


35
10.310
21.911
FALSE


36
10.739
21.911
FALSE


37
10.739
22.286
TRUE


38
10.310
22.286
TRUE


39
10.310
21.911
FALSE


40
10.739
21.911
FALSE


41
10.739
22.286
TRUE


42
10.310
22.286
FALSE


43
10.739
21.911
FALSE


44
11.203
21.911
FALSE


45
11.203
22.286
FALSE


46
10.739
22.286
FALSE


47
10.739
21.911
FALSE


48
11.203
21.911
FALSE


49
11.203
22.286
FALSE


50
10.739
22.286
FALSE


51
11.203
21.911
FALSE


52
11.667
21.911
FALSE


53
11.667
22.286
TRUE


54
11.203
22.286
TRUE


55
11.203
21.911
TRUE


56
11.667
21.911
TRUE


57
11.667
22.286
TRUE


58
11.203
22.286
TRUE


59
11.667
21.911
TRUE


60
12.131
21.911
TRUE


61
12.131
22.286
TRUE



Here is the UDF.

'Rick Rothstein, http://www.mrexcel.com/forum/excel-questions/713005-does-point-fall-within-polygon-visual-basic-applications-function-2.html
Public Function PtInPoly(Polygon As Range, Xcoord As Double, Ycoord As Double) As Boolean
Dim x As Long, m As Double, b As Double, Poly As Variant, NumSidesCrossed As Long
Poly = Polygon
For x = 1 To UBound(Poly) - 1
If Poly(x, 1) > Xcoord Xor Poly(x + 1, 1) > Xcoord Then
m = (Poly(x + 1, 2) - Poly(x, 2)) / (Poly(x + 1, 1) - Poly(x, 1))
b = (Poly(x, 2) * Poly(x + 1, 1) - Poly(x, 1) * Poly(x + 1, 2)) / (Poly(x + 1, 1) - Poly(x, 1))
If m * Xcoord + b > Ycoord Then NumSidesCrossed = NumSidesCrossed + 1
End If
Next
PtInPoly = NumSidesCrossed Mod 2
End Function

[Paul_Hossler]

More to think on

1. I did rearrange the P1-P4 to fit the link I found

2. Because of the floating point operations, there was some round off so I used epsilon = 0.00001 as "Close 'Nuff"

3. The code is very linear to follow the algorithm in the link

4. Only tested a little



Option Explicit
Const epsilon As Double = 0.00001

Sub test()
Dim rPoints As Range, rAreas As Range
Dim rPoint As Range, rArea As Range


Set rPoints = Worksheets("Point").Cells(1, 1).CurrentRegion
Set rPoints = rPoints.Cells(2, 1).Resize(rPoints.Rows.Count - 1, rPoints.Columns.Count)
Set rAreas = Worksheets("Area").Cells(1, 1).CurrentRegion
Set rAreas = rAreas.Cells(3, 1).Resize(rAreas.Rows.Count - 2, rAreas.Columns.Count)
For Each rPoint In rPoints.Rows

For Each rArea In rAreas.Rows
If Inside(rPoint, rArea) Then
If Len(rPoint.Cells(1, 4).Value) = 0 Then
rPoint.Cells(1, 4).Value = rArea.Cells(1, 1).Value
Else
rPoint.Cells(1, 4).Value = rPoint.Cells(1, 4).Value & ", " & rArea.Cells(1, 1).Value
End If
End If
Next
Next
End Sub

'http://math.stackexchange.com/questions/190111/how-to-check-if-a-point-is-inside-a-rectangle
' 1 2 3 4
'Point X Y Area
' 1 2 3 4 5 6 7 8 9
'name x1 y1 x2 y2 x3 y3 x4 y4
Private Function Inside(Pt As Range, Ar As Range) As Boolean
Dim a(1 To 4) As Double
Dim b(1 To 4) As Double
Dim u(1 To 4) As Double
Dim BigA(1 To 4) As Double
Dim X(1 To 4) As Double, Y(1 To 4) As Double
Dim i As Long
Dim pX As Double, pY As Double, RectArea As Double
'store area points
For i = LBound(a) To UBound(a)
X(i) = Ar.Cells(1, 2 * i).Value
Y(i) = Ar.Cells(1, 2 * i + 1).Value
Next i

'store point
pX = Pt.Cells(1, 2).Value
pY = Pt.Cells(1, 3).Value

'First we calculate the edge lengths:
a(1) = Sqr((X(1) - X(2)) ^ 2 + (Y(1) - Y(2)) ^ 2)
a(2) = Sqr((X(2) - X(3)) ^ 2 + (Y(2) - Y(3)) ^ 2)
a(3) = Sqr((X(3) - X(4)) ^ 2 + (Y(3) - Y(4)) ^ 2)
a(4) = Sqr((X(4) - X(1)) ^ 2 + (Y(4) - Y(1)) ^ 2)
'Next we calculate the lengths of the line segments:
b(1) = Sqr((X(1) - pX) ^ 2 + (Y(1) - pY) ^ 2)
b(2) = Sqr((X(2) - pX) ^ 2 + (Y(2) - pY) ^ 2)
b(3) = Sqr((X(3) - pX) ^ 2 + (Y(3) - pY) ^ 2)
b(4) = Sqr((X(4) - pX) ^ 2 + (Y(4) - pY) ^ 2)
'Then we calculate the areas using Heron's Formula:
u(1) = (a(1) + b(1) + b(2)) / 2
u(2) = (a(2) + b(2) + b(3)) / 2
u(3) = (a(3) + b(3) + b(4)) / 2
u(4) = (a(4) + b(4) + b(1)) / 2
RectArea = a(1) * a(2)

BigA(1) = Sqr(u(1) * (u(1) - a(1)) * (u(1) - b(1)) * (u(1) - b(2)))
BigA(2) = Sqr(u(2) * (u(2) - a(2)) * (u(2) - b(2)) * (u(2) - b(3)))
BigA(3) = Sqr(u(3) * (u(3) - a(3)) * (u(3) - b(3)) * (u(3) - b(4)))
BigA(4) = Sqr(u(4) * (u(4) - a(4)) * (u(4) - b(4)) * (u(4) - b(1)))
If (RectArea - (BigA(1) + BigA(2) + BigA(3) + BigA(4))) < epsilon Then
Inside = True
Else
Inside = False
End If
End Function

[p45cal]

I did some work on points in polygons with Aussiebear a couple of years ago (http://www.vbaexpress.com/forum/showthread.php?48478-Sorting-GPS-marks) and had a chance to assess the various udfs and found the PtInPoly function of Rick Rothstein to be the simplest and most reliable.
In the attached, I've used it to create a table on sheet1 which shows whether each point is in each area/rectangle. TRUE highlighted in green shows it does, FALSE otherwise.
I've also created a chart on sheet Area which always shows all the 106 points but shows each rectangle one at a time according to the position of the scrollbar in the vicinity of cell G24. Because many points are very close or duplicates I allowed individual points to be highlighted on the chart by adjusting the scrollbar in the vicinity of cell G29, and while doing that, cell G32 shows whether that point is in the rectangle that is showing. There are 2 more charts on that sheet which are just the same as the first chart but zoomed in to areas of interest. The graphic results agree with the sheet1 table

However… the results differ from both Kenneth Hobs' and Paul_Hossler's, so maybe I've done it wrong!


whether the points in sheet (Point) fall in which area showing in sheet (Area) and display the name in the column C of Sheet (Point). I'll try to do something in this regard soon.





ps. There have been some developments on this udf at http://www.excelfox.com/forum/showthread.php/1579-Test-Whether-A-Point-Is-In-A-Polygon-Or-Not

[Paul_Hossler]

@p45cal

Don't know if it matters, but I re-arranged the point number to fit the algorithm

I like how you made the graphs, really helps to see it

[p45cal]

@p45cal

Don't know if it matters, but I re-arranged the point number to fit the algorithmSo did I but I don't know if I did it in the same way; I swapped the positions of 3 and 4.

[Paul_Hossler]

@ken -- thanks for the link -- interesting stuff

In the link's comments, someone pointed out in #2



Just a small correction. When I used it in my project I got some false positives and after reviewing the original code where you took the algorithm I figured out why.

In their code the variable "j" is assigned the value "polyLines-1" because their arrays start at 0, in your case the arrays start at 1 so the correct assignment is "j=polyLines". After this fix it worked beautifully.


so I'm wondering if it shouldn't be



For x = LBound(Poly) To UBound(Poly) - 1

[p45cal]

Button in sheet Point of the attached.

[Kenneth Hobs]

Good work p45cal! Visuals are certainly worth viewing when validating a function such as this.

I did not view the excelfox link until now since chrome thinks it is a security risk. Firefox thinks that the excelfox link's file is a security risk. Now that I can see the concept that Rick used, double angle area theorem, I did not use Rick's short function properly.

If I get time, I may re-code Rick's routine to not require a two column two dimensional array/range with the last coordinate repeated. Code should be able to address coordinates laid out as the OP did, one row with no coordinate repeated.

Talk about repeating, looking at this thread and the concepts gives me deja vu.

[p45cal]

I did not view the excelfox link until now since chrome thinks it is a security risk. Firefox thinks that the excelfox link's file is a security risk.see msgs #8, 9 and 13 here: http://www.mrexcel.com/forum/excel-questions/951922-breakdown-option-code-automaticaly.html#post4573459 It's Google that the reporting link in FireFox goes to.





Code should be able to address coordinates laid out as the OP did, one row with no coordinate repeated.In attachment to msg#17 I use the OP's data layout in the code behind the Button (except for swapping points 3 and 4).

[snb]

I'd stick to arrays alltogether:


Sub M_snb()
With Sheets("Point").Range("A1").CurrentRegion
sn = .Resize(, .Columns.Count + 1)
End With
With Sheets("area").Cells(1).CurrentRegion
sp = .Resize(, .Columns.Count + 2)
End With

For j = 2 To UBound(sn)
For jj = 3 To UBound(sp)
sp(jj, 10) = sp(jj, 2)
sp(jj, 11) = sp(jj, 3)
If F_PtInPoly(Application.Index(sp, jj), sn(j, 2), sn(j, 3)) Then c00 = c00 & ";" & sp(jj, 1)
Next
sn(j, 4) = Mid(c00, 2)
c00 = ""
Next

Sheets("Point").Cells(1, 20).Resize(UBound(sn), UBound(sn, 2)) = sn
End Sub

Function F_PtInPoly(Poly, X, Y) As Boolean
For j = 2 To UBound(Poly) - 2 Step 2
If Poly(j) > X Xor Poly(j + 2) > X Then
m = (Poly(j + 3) - Poly(j + 1)) / (Poly(j + 2) - Poly(j))
b = (Poly(j + 1) * Poly(j + 2) - Poly(j) * Poly(j + 3)) / (Poly(j + 2) - Poly(j))
If m * X + b > Y Then p = p + 1
End If
Next

F_PtInPoly = p Mod 2
End Function

[vincentzack]

p45cal, Thank you very much! Your excel is excellent.

Could you add some notes (such as how to change the start row number and column number) on the code? So that I can amend the code by myself if necessary.

Thank you again! Your code make my calculation much faster.

[p45cal]

Instead of
Set CellsToFill = Sheets("Point").Range("A1").CurrentRegion
Set CellsToFill = CellsToFill.Offset(1).Resize(CellsToFill.Rows.Count - 1).Columns(1).Offset(, 3)
which in your case sets CellsToFill as cell D2:D107, you can manually work out your range and just have:
Set cellstofill = Sheets("Point").Range("D2:D107")Just remember that the code expects the x and y coordinates of a point to be 2 cells and 1 cell to the left respectively:
x = PtInPoly(Rect, cll.Offset(, -2), cll.Offset(, -1))

The only other range you need to know about is
For Each cel In Sheets("Area").Range("A3:A20").Cells '<< adjustAt the moment Sheets("Area").Range("A3:A20") contains the names of your rectangles. All you need to be aware of is that the code expects to see the 4 coordinates of the rectangle in the 8 cells to the right of each name, in order around the rectangle (doesn't matter if it's clockwise or anticlockwise - diagonally opposed corners of the rectangle next to each other is a no).

[snb]

I reduced the code so it's easier to understand.


Sub M_snb()
sn = Sheets("Point").Cells(1).CurrentRegion
sp = Sheets("area").Cells(1).CurrentRegion

For j = 2 To UBound(sn)
For jj = 3 To UBound(sp)
If F_PtInPoly(Application.Index(sp, jj), sn(j, 2), sn(j, 3)) Then c00 = c00 & ";" & sp(jj, 1)
Next
sn(j, UBound(sn, 2)) = Mid(c00, 2)
c00 = ""
Next

Sheets("Point").Cells(1, 10).Resize(UBound(sn), UBound(sn, 2)) = sn
End Sub

Function F_PtInPoly(Poly, X, Y) As Boolean
' based on Rick Rothstein, http://www.mrexcel.com/forum/excel-questions/713005-does-point-fall-within-polygon-visual-basic-applications-function-2.html

For j = 2 To UBound(Poly) Step 2
jj = j Mod (UBound(Poly) - 1) + 2
If Poly(j) > X Xor Poly(jj) > X Then
m = (Poly(jj+1) - Poly(j + 1)) / (Poly(jj) - Poly(j))
b = (Poly(j + 1) * Poly(jj) - Poly(j) * Poly(jj+1)) / (Poly(jj) - Poly(j))
If m * X + b > Y Then p = p + 1
End If
Next

F_PtInPoly = p Mod 2
End Function