Skip to content

Commit 57e6229

Browse files
committed
alpha release
1 parent 335d49a commit 57e6229

11 files changed

+564
-93
lines changed

README.md

+89-3
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,92 @@
1-
## intersystems-objectscript-template
2-
This is a template for InterSystems ObjectScript Github repository.
3-
The template goes also with a few files which let you immediately compile your ObjectScript files in InterSystems IRIS Community Edition in a docker container
1+
# SpatialIRIS
2+
3+
SpatialIRIS is a library of basic geographic features -- points, lines and polygons -- that implement common spatial operations, such as those involving distances between points, and containment of a geometry within a polygon.
4+
5+
Many of these operations match those found in PostGIS, but the selection here is much more limited, and there are no spatial indices applied, so they will execute fairly slowly. However, this library brings useful functionality to InterSystems IRIS, such as identifying all patients with 10km of a hospital, or within a zip code or county.
6+
7+
## Classes
8+
9+
### geo.model.Point
10+
11+
A Point is a 0-dimensional geometry that represents a single location in coordinate space.
12+
13+
### geo.model.LineString
14+
15+
A LineString is a 1-dimensional line formed by a contiguous sequence of line segments. Each line segment is defined by two points, with the end point of one segment forming the start point of the next segment.
16+
17+
### geo.model.Polygon
18+
19+
A Polygon is a 2-dimensional planar region, delimited by an exterior boundary. Each boundary is a LineString for which the first and last points must be equal, and the line must not self-intersect.
20+
21+
## Functions
22+
23+
### Implemented by all geometry types
24+
25+
#### dimension() returns %Integer
26+
27+
Return the topological dimension of this Geometry object -- 0 for POINT, 1 for LINESTRING, 2 for POLYGON.
28+
29+
#### geometryType() returns %String
30+
31+
Returns the type of the geometry as a string. Eg: 'LINESTRING', 'POLYGON', 'POINT'.
32+
33+
#### numPoints() returns %Integer
34+
35+
Return the number of points in a geometry
36+
37+
#### asText() returns %String
38+
39+
Returns the OGC Well-Known Text (WKT) representation of the geometry/geography
40+
41+
#### asGeoJSON() returns %JSON
42+
43+
Returns a geometry as a GeoJSON "feature" (See the [GeoJSON specifications RFC 7946](https://tools.ietf.org/html/rfc7946)).
44+
45+
### Implemented by geo.model.Point
46+
47+
#### pointInsideCircle(centerX As %Float, centerY As %Float, radius As %Float) returns %Boolean
48+
49+
Returns true if the point is inside the circle with center `center_x,center_y` and radius (in kilometers) `radius`.
50+
51+
#### distance(geo.model.Point) returns %Integer
52+
53+
Returns the minimum 2D Cartesian (planar) distance between two Points, in kilometers
54+
55+
### Implemented by geo.model.Line
56+
57+
#### startPoint() returns geo.model.Point
58+
59+
Returns the first point of a LineString.
60+
61+
#### endPoint() returns geo.model.Point
62+
63+
Returns the last point of a LineString.
64+
65+
### Implemented by geo.model.Polygon
66+
67+
#### isClosed() returns %Boolean
68+
69+
Returns TRUE if the Polygon's start and end points are coinciden
70+
71+
#### contains(geo.model.AbstractGeometry) returns %Boolean
72+
73+
Returns TRUE if and only if no points of the passed geometry lie in the exterior of Polygon, and at least one point of the interior of the geometry lies in the interior of Polygon.
74+
75+
#### intersects(geo.model.AbstractGeometry) returns %Boolean
76+
77+
Returns TRUE if any point of the passed geometry falls within the polygon.
78+
79+
#### within(geo.model.Polygon) returns %Boolean
80+
81+
Returns TRUE if all points of the passed geometry fall within the polygon.
82+
83+
## Installation
84+
85+
TBD
86+
87+
## Example
88+
89+
TBD
490

591
## Prerequisites
692
Make sure you have [git](https://git-scm.com/book/en/v2/Getting-Started-Installing-Git) and [Docker desktop](https://www.docker.com/products/docker-desktop) installed.

data/kendall.json

+48
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
{
2+
"type": "FeatureCollection",
3+
"features": [
4+
{
5+
"type": "Feature",
6+
"properties": {},
7+
"geometry": {
8+
"coordinates": [
9+
[
10+
[
11+
-71.07296971120546,
12+
42.36751524646158
13+
],
14+
[
15+
-71.08872828096423,
16+
42.36897060744346
17+
],
18+
[
19+
-71.08982906341102,
20+
42.36605985176101
21+
],
22+
[
23+
-71.09226237197649,
24+
42.36717280368296
25+
],
26+
[
27+
-71.09434806503295,
28+
42.36511810781076
29+
],
30+
[
31+
-71.08814892178214,
32+
42.358268636220345
33+
],
34+
[
35+
-71.0776045846639,
36+
42.36019512557601
37+
],
38+
[
39+
-71.07296971120546,
40+
42.36751524646158
41+
]
42+
]
43+
],
44+
"type": "Polygon"
45+
}
46+
}
47+
]
48+
}

data/kendall_bbox.json

+31
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,31 @@
1+
{
2+
"type": "FeatureCollection",
3+
"features": [
4+
{
5+
"type": "Feature",
6+
"properties": {},
7+
"geometry": {
8+
"coordinates": [
9+
[
10+
[
11+
-71.09792960473557,42.36934185151014
12+
],
13+
[
14+
-71.09792960473557,42.35841621918186
15+
],
16+
[
17+
-71.07704919078053,42.35841621918186
18+
],
19+
[
20+
-71.07704919078053,42.36934185151014
21+
],
22+
[
23+
-71.09792960473557,42.36934185151014
24+
]
25+
]
26+
],
27+
"type": "Polygon"
28+
}
29+
}
30+
]
31+
}

src/geo/Testing.cls

+77
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,77 @@
1+
Class geo.Testing
2+
{
3+
4+
ClassMethod createSimpleTestPolygon() As geo.model.Polygon
5+
{
6+
set l = ##class(geo.model.LineString).%New()
7+
do l.addPoint(##class(geo.model.Point).%New(0.01,0.01))
8+
do l.addPoint(##class(geo.model.Point).%New(2.01,0.01))
9+
do l.addPoint(##class(geo.model.Point).%New(2.01,2.01))
10+
do l.addPoint(##class(geo.model.Point).%New(0.01,2.01))
11+
do l.addPoint(##class(geo.model.Point).%New(0.01,0.01))
12+
set polygon = ##class(geo.model.Polygon).%New()
13+
set polygon.line = l
14+
Return polygon
15+
}
16+
17+
ClassMethod createCambridgePolygon() As geo.model.Polygon
18+
{
19+
set l = ##class(geo.model.LineString).%New()
20+
do l.addPoint(##class(geo.model.Point).%New(-71.07296971120546, 42.36751524646158))
21+
do l.addPoint(##class(geo.model.Point).%New(-71.08872828096423, 42.36897060744346))
22+
do l.addPoint(##class(geo.model.Point).%New(-71.08982906341102, 42.36605985176101))
23+
do l.addPoint(##class(geo.model.Point).%New(-71.09226237197649, 42.36717280368296))
24+
do l.addPoint(##class(geo.model.Point).%New(-71.09434806503295, 42.36511810781076))
25+
do l.addPoint(##class(geo.model.Point).%New(-71.08814892178214, 42.358268636220345))
26+
do l.addPoint(##class(geo.model.Point).%New(-71.0776045846639, 42.36019512557601))
27+
do l.addPoint(##class(geo.model.Point).%New(-71.07296971120546, 42.36751524646158))
28+
set polygon = ##class(geo.model.Polygon).%New()
29+
set polygon.line = l
30+
Return polygon
31+
}
32+
33+
ClassMethod createLineIntersectingKendall() As geo.model.LineString
34+
{
35+
set l = ##class(geo.model.LineString).%New()
36+
do l.addPointFromCoords(-71.09009501700031,42.36703809)
37+
do l.addPointFromCoords(-71.08100073568535,42.3670380723421)
38+
do l.addPointFromCoords(-71.085000,42.3670380723)
39+
Return l
40+
}
41+
42+
ClassMethod createLineInsideKendall() As geo.model.LineString
43+
{
44+
set l = ##class(geo.model.LineString).%New()
45+
do l.addPointFromCoords(-71.08100073568535,42.3670380723421)
46+
do l.addPointFromCoords(-71.085000,42.3670380723)
47+
Return l
48+
}
49+
50+
ClassMethod createLineOutsideKendall() As geo.model.LineString
51+
{
52+
set l = ##class(geo.model.LineString).%New()
53+
do l.addPointFromCoords(-71.09009501700031,42.3670380723421)
54+
do l.addPointFromCoords(-71.09009501700031,42.60000)
55+
Return l
56+
}
57+
58+
ClassMethod createPointOutsideKendall() As geo.model.Point
59+
{
60+
Return ##class(geo.model.Point).%New(-71.09009501700031,42.3670380723421)
61+
}
62+
63+
ClassMethod createPointInsideKendall() As geo.model.Point
64+
{
65+
Return ##class(geo.model.Point).%New(-71.08100073568535,42.3670380723421)
66+
}
67+
68+
ClassMethod test()
69+
{
70+
set polygon = ..createSimpleTestPolygon()
71+
set isin = ##class(geo.Tools).PointInPolygon(polygon, ##class(geo.model.Point).%New(1.01,1.01))
72+
w isin,!
73+
set isout = ##class(geo.Tools).PointInPolygon(polygon, ##class(geo.model.Point).%New(3.01,3.01))
74+
w isout,!
75+
}
76+
77+
}

src/geo/Tools.cls

+9-50
Original file line numberDiff line numberDiff line change
@@ -4,67 +4,29 @@ Class geo.Tools Extends %RegisteredObject
44
/// Earth's radius in km
55
Parameter RADIUSEARTH = 6371;
66

7-
ClassMethod createPolygon() As geo.model.Polygon
8-
{
9-
set l = ##class(geo.model.Line).%New()
10-
do l.addPoint(##class(geo.model.Point).%New(0.01,0.01))
11-
do l.addPoint(##class(geo.model.Point).%New(2.01,0.01))
12-
do l.addPoint(##class(geo.model.Point).%New(2.01,2.01))
13-
do l.addPoint(##class(geo.model.Point).%New(0.01,2.01))
14-
do l.addPoint(##class(geo.model.Point).%New(0.01,0.01))
15-
set polygon = ##class(geo.model.Polygon).%New()
16-
set polygon.line = l
17-
Return polygon
18-
}
19-
20-
ClassMethod test()
21-
{
22-
set polygon = ..createPolygon()
23-
set isin = ..PointInPolygon(polygon, ##class(geo.model.Point).%New(1.01,1.01))
24-
w isin,!
25-
set isout = ..PointInPolygon(polygon, ##class(geo.model.Point).%New(3.01,3.01))
26-
w isout,!
27-
}
28-
29-
ClassMethod test2()
30-
{
31-
set polygon = ..createPolygon()
32-
set isin = ..PointInPolygon2(polygon, ##class(geo.model.Point).%New(1.01,1.01))
33-
w isin,!
34-
set isout = ..PointInPolygon2(polygon, ##class(geo.model.Point).%New(3.01,3.01))
35-
w isout,!
36-
}
37-
387
ClassMethod PointInPolygon(poly As geo.model.Polygon, pt As geo.model.Point) As %Boolean
398
{
409
#; https://stackoverflow.com/questions/36399381/whats-the-fastest-way-of-checking-if-a-point-is-inside-a-polygon-in-python
4110
#; https://matplotlib.org/stable/api/path_api.html
4211
#; import numpy as np
4312
set mp = ##class(%SYS.Python).Import("matplotlib")
4413
set mpltPath = mp.path
45-
set polylist = poly.getAsList()
14+
set polylist = poly.asList()
4615
set nothing = ##class(Py.Helper).toPyListOrString(polylist,.poly2)
4716
zw poly2
48-
set point = ##class(Py.Helper).toPyListOrString(pt.getAsList())
17+
set point = ##class(Py.Helper).toPyListOrString(pt.asList())
4918
set path = ..invoke(mpltPath,poly2)
5019
Return path."contains_point"(point)
5120
}
5221

53-
ClassMethod PointInPolygon2(poly As geo.model.Polygon, pt As geo.model.Point) As %Boolean
22+
ClassMethod invoke(ByRef mpltPath As %SYS.Python, ByRef poly2) As %SYS.Python [ Language = python ]
5423
{
55-
set mp = ##class(%SYS.Python).Import("matplotlib")
56-
set mpltPath = mp.path
57-
set polylist = poly.getAsList()
58-
set nothing = ##class(Py.Helper).toPyListOrString(polylist,.poly2)
59-
zw poly2
60-
set path = mpltPath.Path(poly2)
61-
set point = ##class(Py.Helper).toPyListOrString(pt.getAsList())
62-
Return path."contains_point"(point)
24+
return mpltPath.Path(poly2)
6325
}
6426

65-
/// Computes the great-circle distance between two points on a sphere
27+
/// Computes the great-circle distance between a point and the closest point to it on a line
6628
/// returns: distance in kilometers
67-
ClassMethod DistancePointandLine(pt As geo.model.Point, line As geo.model.Line) As %Float
29+
ClassMethod DistancePointandLine(pt As geo.model.Point, line As geo.model.LineString) As %Float
6830
{
6931
set dist = 0
7032
for i = 2:1:line.points.Count() {
@@ -74,12 +36,9 @@ ClassMethod DistancePointandLine(pt As geo.model.Point, line As geo.model.Line)
7436
Return dist
7537
}
7638

77-
ClassMethod invoke(ByRef mpltPath As %SYS.Python, ByRef poly2) As %SYS.Python [ Language = python ]
78-
{
79-
return mpltPath.Path(poly2)
80-
}
81-
82-
ClassMethod Distance(p1 As geo.model.Point, p2 As geo.model.Line) As %Float
39+
/// Computes the great-circle distance between two points on a sphere
40+
/// returns: distance in kilometers
41+
ClassMethod Distance(p1 As geo.model.Point, p2 As geo.model.Point) As %Float
8342
{
8443
Set lat1Radianos = (p1.latitude * $ZPI) / 180
8544
Set lng1Radianos = (p1.longitude * $ZPI) / 180

src/geo/model/AbstractGeometry.cls

+46
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,46 @@
1+
Class geo.model.AbstractGeometry Extends %SerialObject
2+
{
3+
4+
// Returns the geometry as a $LIST
5+
6+
Method asList() [ Abstract ]
7+
{
8+
}
9+
10+
// Returns a geometry as a GeoJSON "geometry"
11+
12+
Method asGeoJSON() As %String [ Abstract ]
13+
{
14+
}
15+
16+
// Returns the OGC Well-Known Text (WKT) representation of the geometry
17+
18+
Method asText() As %String [ Abstract ]
19+
{
20+
}
21+
22+
// Return the number of points in a geometry.
23+
24+
Method numPoints() As %Integer [ Abstract ]
25+
{
26+
}
27+
28+
// Returns the type of the geometry as a string. Eg: 'LINESTRING', 'POLYGON', 'POINT'
29+
30+
ClassMethod geometryType() As %String [ Abstract ]
31+
{
32+
}
33+
34+
// Return the topological dimension of this Geometry object, which must be less than or equal to the coordinate dimension. OGC SPEC s2.1.1.1 - returns 0 for POINT, 1 for LINESTRING, 2 for POLYGON
35+
36+
ClassMethod dimension() As %Integer [ Abstract ]
37+
{
38+
}
39+
40+
Storage Default
41+
{
42+
<StreamLocation>^geo.model.AbstractGeometryS</StreamLocation>
43+
<Type>%Storage.Serial</Type>
44+
}
45+
46+
}

0 commit comments

Comments
 (0)