Utility Modules

segLibrary

Operations over lists of line segments.

pyspatiotemporalgeom.utilities.segLibrary.calcIntersectionPoints(segs, laList=None, lbList=None, indexList=None)[source]

computes the intersection points between segments in the list segs. labels are updated for overlapping segs, but not intersecting segs.

Note

Labels are preserved for all segs that do not overlap with another segment. If two segments overlap, their labels are summed.

The constructSegsWithLabels() function will preserve these labels and thier sums.

If the desire is to keep track of which areas are covered by cycles with particular labels, use the mapLibrary versions of these functions to construct a map

segs is in the usual format: [((x1,y1),(x2,y2)),((x3,y3),(x4,y4)),...]

output: a list of indexed points. An indexed point has the format: (index,(x,y)) where index is the index of the seg in which the intersetion point occurred.

Assumes, that snapEndPoints has already been called

Use constructSegsWithLabels to build back the actual segs broken at intersection points

pyspatiotemporalgeom.utilities.segLibrary.calcIntersectionPointsRedBlue(segs1, segs2, indexList1=None, indexList2=None)[source]

computes intersection points between red and blue segs.

  • assumes seg snapping has taken place
  • assumes segs are ALL floats for end points
  • Index lists are used if strip decomposition has been performed.
  • Otherwise, pass None and everything will be taken care of

output:

(pois1,pois2): a tuple. pois1 is a list of intersection points occuring in segs1. pois2 is a list of intersection points occurrin in segs2. A point in pois1 or pois2 has the strucutre: intersectPoi = (0,(0.0,0.0)). The first integer is the index of the line segment in which the intersection point occurs.

pyspatiotemporalgeom.utilities.segLibrary.calcNonIntersectingSegs(segs, laList=None, lbList=None)[source]

break segs at intersection points. Preserve labels if they are present.

Uses a strip decomposition scheme to speed things up

input:

segs: a list of segments in the usual format [((x1,y2),(x2,y2)),((x3,y3),(x4,y4)),...]

laList, lbList: parallel arrays with segs, indicating label values.

output: a list of segs in the usual format.

if label lists are passed in, the output is a tuple: (segs, laList, lblist)

pyspatiotemporalgeom.utilities.segLibrary.calcNonIntersectingSegsIntEndPoints(segs)[source]

take a set of line segments.

  1. Snap all end points to intergers
  2. Break all segs at intersection points.

the result may have floating point values at intersection points. This is not such a useful function. Typically, you would use on of the other line segment intersection functions.

input: segs: a list of line segments of the format [((x1,y1),(x2,y2)),((x3,y3),(x4,y4)),...]

output: a list of line segments, same structure as input.

pyspatiotemporalgeom.utilities.segLibrary.checkNumRobustnes(segs)[source]

DO NOT USE!!!

This function was intended to remove errors introduced by floating point rounding, but it takes a long time to compute, and doesn’t actually fix the main problems.

INSTEAD.. use snapEndPoints()

takes labeld segs (same structure as left halfsegs)

pyspatiotemporalgeom.utilities.segLibrary.collinear(L, R, tolerance=1e-07)[source]

Test if two line segments are collinear within a tolerance.

input:

L: a line segment: ((x1,y1),(x2,y2))

R: a line segment like L

tolerance: a tolerance value. Set to 0 for exact collinearity test, higher for greater tolerance.

output:

boolean value

pyspatiotemporalgeom.utilities.segLibrary.collinearPoints(p1, p2, p3, tolerance=1e-07)[source]

uses left hand turn test to indicate if three points are collinear.

input

p1,p2,p3
3 points of the form (x,y)
tolerance
If the absolute value of the value computed by the left hand turn test is smaller that tolerance, it is considered colinear

output

boolean value

pyspatiotemporalgeom.utilities.segLibrary.collinearValue(p1, p2, p3)[source]

uses the left hand turn test based on the sign of the area of the trianngle defiened by points p1,p2, and p3, except that the actual value is returned. This is useful if you want to see how close the three points are to being collinear.

Left hand turn test assumes you are traveling from p1 to p2 and then on to p3, and indicates if you must make a left or right turn at p2 to reach p3.

input: 3 points, p1, p2, p3 of the format (x,y)

output: the signed area of the square enclosing the triangle defined by the given points.

  • will return a 0 if p1, p2,p3 are exactly collinear
  • will return positive value if it forms a left turn
  • will return neg value if it forms a right turn
pyspatiotemporalgeom.utilities.segLibrary.constructSegs(pois, segs)[source]

Red/blue line intersection returns a list of indexed points such that the points indicate where line segment intersections occurred, and the index indicates in which line segment the intersection occurred. This function jsut breaks the line segments at those points

*input:

pois: a list of indexed poins of the format (index,(x,y)) where index is the index of a seg in segs

segs: a list of segments in the usual format [((x1,y1),(x2,y2)),((x3,y3),(x4,y4)),...]

output:

a list of segs in the usual format.

pyspatiotemporalgeom.utilities.segLibrary.constructSegsWithLabels(pois, segs, laList, lbList)[source]

See documentation for constructSegs( pois, segs ).

This function works the same as that one, but preserves labelling.

laList and lbList are parallel lists with segs. When a seg is split at an intersection point, the labels are adjusted accordingly.

pyspatiotemporalgeom.utilities.segLibrary.countInteriors(segs, laList, lbList)[source]

marks each segment with ALL interiors it bounds.

assumes a list of segs describing multiple regions that intersect only at endpoints, and no duplicates (these are already taken care of).

Assumes left end point less than right endpoint.

output: a tuple (segs, laList, lbList)

pyspatiotemporalgeom.utilities.segLibrary.createRandomSegs(numSegs, minVal=0, maxVal=10000)[source]

Create a bunch of random line segments within a bounding box. Segments may or may not intersect, but collinear and overlapping segments are not allowed.. Duplicates will be removed.

input:

numSegs: The number of random line segments to create

minVal: the lowest value allowed as an X or Y value to construct a point

maxVal: the largest value allowed as an X or Y value to construct a point

pyspatiotemporalgeom.utilities.segLibrary.indexedPoiCmp(p1, p2)[source]

comparison of indexed points for sorting.

order is based on index, then x values the y values.

input: p1 and p2, two indexed points of the format (index,(x,y))

output:

-1 if p1 < p2

0 if p1 == p2

1 if p1 > p2

pyspatiotemporalgeom.utilities.segLibrary.isCollinearAndOverlapping(L, R, tolerance=1e-07)[source]

test if two segments are both collinear AND overlapping. Collinearity within a tolerance is allowed.

input:

L: a line segment with format ((x1,y1),(x2,y2))

R: line segment with structure as L

tolerance: set to 0 for exact collinearity test. Otherwise, allows almost collinear to count as collinear

returns:

True if L and R are collinear within tolerance and overlapping

False otherwise

pyspatiotemporalgeom.utilities.segLibrary.isLeftTurn(p1, p2, p3)[source]

tests if when traveling from point p1 to point p2, and then on to p3, whether the traveler must make a left or right turn at p2 in order to reach p3.

uses sign of the area of the triangle computation.

input:

p1,p2,p3: points of the format (x,y)

output:

-1 if a left turn is needed

0 if the points are collinear

1 if a right turn is needed

pyspatiotemporalgeom.utilities.segLibrary.keepIn(origSegs1, origSegs2)[source]

point in polygon test for each seg

test if segs from seg1 are in seg2

Assumes non interior intersecting segs.

output: a list of segs from origsegs1 that lie inside the region defined by origSegs2

pyspatiotemporalgeom.utilities.segLibrary.keepInnerBoundary(origSegs1, laList1, lbList1, origSegs2, laList2, lbList2, keepEqualWithSameLA=True)[source]

input segs are Labeled segs!

point in polygon test for each seg

test if segs from seg1 are out of seg2

  • Assumes non interior intersecting segs (except for equivalent segs).
  • Assumes all segs have floating point end points
  • Assumes end point snapping is done

Note that the use of segIntersection(origSegs1, origSegs2) satisfies all assumptions

Used to compute set operations between regions! keep the segments from origSegs1 that are inside or on origSegs2, but take into account the labeling.

returns a list of segments in the usual format.

pyspatiotemporalgeom.utilities.segLibrary.keepOnOrIn(origSegs1, origSegs2)[source]

point in polygon test for each seg

test if segs from seg1 are in or on seg2

Assumes non interior intersecting segs.

output: a list of segs from origsegs1 that lie inside or ON the region defined by origSegs2

pyspatiotemporalgeom.utilities.segLibrary.keepOuterBoundary(origSegs1, laList1, lbList1, origSegs2, laList2, lbList2, keepEqualWithSameLA=True)[source]

input segs are Labeled segs!

point in polygon test for each seg

test if segs from seg1 are out of seg2

  • Assumes non interior intersecting segs (except for equivalent segs).
  • Assumes all segs have floating point end points
  • Assumes end point snapping is done

Note that the use of segIntersection(origSegs1, origSegs2) satisfies all assumptions

Used for region set operations!

output a list of segments in the usual format.

pyspatiotemporalgeom.utilities.segLibrary.removeRandSegs(segs, percentToRemove=0.05)[source]

Remove some segments at random from a list of segs.

This is used by the generateRandomRegion functions to increase the amount of randomness.

input:

segs: a list of line segments, with the usual format: [((x1,y1),(x2,y2)),((x3,y3),(x4,y4)),...]

percentToRemove: A real number n such that 0 <= n <= 1 indicating the percentage of segments to remove.

output:

a list of segs.

pyspatiotemporalgeom.utilities.segLibrary.segIntersection(segs1, segs2)[source]

intersect segs1 and segs2 . returns a tuple containing 2 lists: 1 list has segs1 broken at intersection points with segments in segs2, and vice versa.

thus calls the reb/blue intersection function, but actually returns the resulting line segs (the arrangment), instead of just the intersection points.

There is code use a multi-threaded scheme based on strip decomposition of the embedding space, but it is commented out at the moment.

The call to snapEndPoints call will make sure all segs are left, and snap endpoints that are very close to be the slame point.

input:

segs1: a list of line segments.

segs2: a list of line segments.

output: (resgs1, rsegs2): a tuple. resgs1 is segs1 with segments broken at intersection points with segments in segs2. same for rsegs2.

pyspatiotemporalgeom.utilities.segLibrary.snapEndPoints(segs)[source]

Does 2 things:

  1. if any segs are too short (below a hard coded threshold, they are removed. Too short segs cause issues in sorting and intersection computations
  2. if two segs have end points that are super close together, but that are not identical, one of those end points is shifted to the other. A rather clever algorithm using hash tables does this very quickly

This function does NOT identify segment end points that are super close to another segment interior. In other words, this is not a full on snap rounding implementation, it jsut handles very close end points.

input: A list of segments in the usual format: [((x1,y1),(x2,y2)),((x3,y3),(x4,y4)),...]

output: a list of segments that have been snapped at end points.

pyspatiotemporalgeom.utilities.segLibrary.snapEndPointsWithLabels(segs)[source]

takes a list of LABELED segs, removes too short segs, and snaps too close end points together ASSUMES segs end points are ALL floats

pyspatiotemporalgeom.utilities.segLibrary.stripDecomposition(segs, laList=None, lbList=None, minx=None, maxx=None, forceNumStrips=None)[source]

See docstring for stripDecompositionWithBounds

pyspatiotemporalgeom.utilities.segLibrary.stripDecompositionWithBounds(segs, laList=None, lbList=None, minx=None, maxx=None, forceNumStrips=None)[source]

returns a list of lists. each list will contain copies of the segs that fall into that strip. Stips are equal length, and will be based on the min and max xvalues of the segs, unless specific min and max values are passed.

is the laList and lbList are == None, no label lists are returned

returns the lists of strip segs, the lists of their labels (if lists are not none) the list of indexes mapping to the original seg in newly sorted list segs the list of newly sorted above labels (if labels are passed in ) the list of newly sorted below labels (if labels are passed in ) the list of segment boundaries containin

hsegLibrary

Operations over lists of halfsegments (typically that form a region)

pyspatiotemporalgeom.utilities.hsegLibrary.addVerticalConnectorsPS(hsegs, highestCycleNum=None)[source]

Adding vertical connectors (for use in regionInterpolator module)

pyspatiotemporalgeom.utilities.hsegLibrary.assignCycleLabelsForConcavityMapping(hsegs)[source]

first we need to relabel everything for our new mapping scheme assume connectors have both negative labels, segs have a pos an nevg labelling described in getIndexOfNextOuterCycleWalkSkipProcessed function

pyspatiotemporalgeom.utilities.hsegLibrary.cmpForIndexedPoints(p1, p2)[source]

comparison function for indexed points. stucture is (index, (x,y))

pyspatiotemporalgeom.utilities.hsegLibrary.collinear(L, R)[source]

returns true if L and R are EXACTLY collinear. Check the collinearity functions in segLibrary if a tolerance is allowed.

pyspatiotemporalgeom.utilities.hsegLibrary.convertSegsToHsegs(segs)[source]

convert a list of segments to a list of halfsegments.

Labels are all set to 0.

input: a list of segs in the usual structure

output: a sorted list of halfsegments.

pyspatiotemporalgeom.utilities.hsegLibrary.doWeSwitchSides(currVal, currSeg, nextSeg)[source]

When walking a cycle boundary and assigning labels to halfsegments, this function determines when we should switch the labeling of the hsegs (swith the LA and LB values).

input:

currVal: a bool indicating if we are currently labeling the interior of the region as ABOVE the current hseg

currSeg: The hseg that was labeled according to the currVal value

nextSeg: The next hseg in cyclic order around a cycle

output:

a boolean value indicating if the nextSeg should have the interior of the region as ABOVE (true) or BELOW (False)

pyspatiotemporalgeom.utilities.hsegLibrary.extractAllLargeValidCycles(segs)[source]

Extracts cycles from an input seg set.

segs shoule be a list of line segments. a segments is a tuple ((x1,y1),(x2,y2))

This algorithm favors larger cycles, due to CCW rotations around hseg endpoints.

will not return any cycles that touch another cycle from the interior

Returns a valid region, correctly labeled.

Each cycle get is own unique label

pyspatiotemporalgeom.utilities.hsegLibrary.getConnectedCycleMapping(hsegs)[source]

returns a vector where the value V of the Ith Item indicates that cycle I either directly touches cycle V or touches a sequence of 1 or more cycles that touch V. Cycle V cill be the first cycle in hseg order that is touched by all other cycles whose value in the vector is V.

ASSUMES that each cycle is labeled with a unique label.

The basic idea is to find connected components (cycles that share points)

pyspatiotemporalgeom.utilities.hsegLibrary.getConnectedCycleMappingConcav(hsegs)[source]

Assumes that assignCycleLabelsForConcavityMapping is already called Create a list such that the position at a cycle’s index will have the cycle it is connected to Notes:

  • Exactly 1 cycle will be the outer cycle for the entire region
  • The first hseg will have the cycle number for the outer cycle as its above label (In this implementation the first cycle number is always 2 )
  • Each cycle will map to exactly 1 other cycle

Once we have the list, make lists of segs for each connected cycle chain remove the connected cycle chain from the hsegs list. return the map list, the smaller hsegs list, and the connected cycle segment lists, and a dict mapping a point to all cycles connected at that point

pyspatiotemporalgeom.utilities.hsegLibrary.getIndexOfNextOuterCycleWalkSkipProcessed(searchHseg, hsegs, indexLookupDict)[source]

this walk will assume the following labelling scheme after processed la*lb*-1 will be positive, with -1 being the exterior label processed connectors will have 1,-1 processed segs will have x, -1 where x>1

connectors: segs: processed: x,-x x, -1 : x>1 (-1 is exterior) unprocessed: -1,-1 1,0 (0 is exterior)

pyspatiotemporalgeom.utilities.hsegLibrary.getListOfCycleLabels(hsegs)[source]

get a list of all cycle labels used in a region.

input: an hseg list representing a region

output: a list if inegers

pyspatiotemporalgeom.utilities.hsegLibrary.getNextInnerCycleWalkIndex(searchHseg, hsegs, indexLookupDict)[source]

when traversing a cycle in cyclic order, find the next halfsegment by rotating the current segment through the INTERIOR of the cycle.

This will find smaller cycles, but will interpret a well formed region correctly.

input:

searchHseg. The current hseg from which the next hseg in cyclic order is found

hsegs: a hseg list forming a region

indexLookupDict: a hash table mapping segments to their corresponding halfsegment index in the hsegs list.

output: the next hsed in cyclic order.

pyspatiotemporalgeom.utilities.hsegLibrary.getNextOuterCycleWalk(searchHseg, hsegs, indexLookupDict)[source]

get the next hseg in cyclic order from searchHseg along the outer cycle of a region structure.

input

searchHseg: the hseg from which the next hseg in cyclic order will be computed

hsegs: an hseg list (typically a well formed region, but not required)

indexLookupDict: A dictionary that maps segs to the index of thier equivalent hseg in the hsegs list

output: a hseg

pyspatiotemporalgeom.utilities.hsegLibrary.getNextOuterCycleWalkIndex(searchHseg, hsegs, indexLookupDict)[source]

identical to getNextOuterCycleWalk, except the index of the hseg in the hsegs list is returned instead of the actual hseg.

pyspatiotemporalgeom.utilities.hsegLibrary.getOuterWalkPointSequence(hsegs)[source]

given a list oh hsegs, get the liset of points, in cyclic order, that define the boundary of outer cycle of the leftmost face of the region.

This function is meant for regions with a single face. If there are multiple faces, only the leftmost face is porcessed.

input: a list of hsegs

output: a list of points.

pyspatiotemporalgeom.utilities.hsegLibrary.getStopHsegForOuterCycleWalkSkipProcessed(startIndex, startHseg, hsegs)[source]

helper function for building well-formed regions from possibly non-well-formed input

pyspatiotemporalgeom.utilities.hsegLibrary.getYvalAtX(x, seg)[source]

given a segment, and an x value that lies on the segment, return the y value on the segment at that x value.

pyspatiotemporalgeom.utilities.hsegLibrary.hsegCmpForActiveList(h1, h2)[source]

comparison function for hsegs in an active list.

Used in PLANE SWEEP algorithms

pyspatiotemporalgeom.utilities.hsegLibrary.hsegComp(L, R)[source]

Halfsegment comparison function. Used to provide a total ordering on halfsegments

input: two halfsegments

output: True if L < R, False otherwise

pyspatiotemporalgeom.utilities.hsegLibrary.hsegCompForSorted(L, R)[source]

wrapper for halfsegment comparison function used for the sorted library.

pyspatiotemporalgeom.utilities.hsegLibrary.hsegTriangulateStripCompForSorted(L, R)[source]

comparison of halfsegments used in triangulating regions.

this a special case of halfsegment comparisons.

pyspatiotemporalgeom.utilities.hsegLibrary.isLeft(hseg)[source]

returns true if the input halfsegment is a left halfsegment, false otherwise

A halfsegment has the structure: (((x1,y1),(x2,y2)),la,lb) where la,lb are ints

pyspatiotemporalgeom.utilities.hsegLibrary.isLeftTurn(p1, p2, p3)[source]

Another left hand turn test.

When traveling from point p1 to pont p2, do we take a left or right turn to travel on th point p3?

input: 3 points, p1, p2, p3 in (x,y) format

output: -1 if its a left turn, 0 if the points ar collinear, 1 if its a right turn.

pyspatiotemporalgeom.utilities.hsegLibrary.labelUniqueCycles(segs, onlyReturnOuterCycle=False)[source]

will do interior walks around cycles.

Each cycle will get a unique label.

Label numbers will not necessarily start at 2, and may skip numbers

Will also remove sticks!!! a full service solution!

Cycles will have unique labels, but holes in particular will have labelling flipped.

Use: def switchLabelsForCorrectCycleLabelling( hsegs ): to finalize labels

This function is used to create well-formed regions from possibly non-wellformed input.

input:

segs: a list of segs.

onlyReturnOuterCycle: if you want a simple region, set this to true!

output: a hseg list representing a well formed region. Agian, call def switchLabelsForCorrectCycleLabelling( hsegs ): to finalize labels!!! Everything is labeled as an outercycle after this call.

pyspatiotemporalgeom.utilities.hsegLibrary.mergeListSets(currCycle, meetingSets, fullSet)[source]

a helper function to merge cycle lists.

pyspatiotemporalgeom.utilities.hsegLibrary.mergeLists(currCycle, meetingSets, cyleTofirstConnectedCycleMapList, firstCycleInChain)[source]

helper function identifying connected components in a region (cycles that meet at a point)

pyspatiotemporalgeom.utilities.hsegLibrary.mergeSortedHsegLists(L1, L2)[source]

merge two sorted hseg lists into a single hseg list.

pyspatiotemporalgeom.utilities.hsegLibrary.relabelTouchingCyclesToFirstCycleNum(hsegs, cycleMapList)[source]

relabel segments in touching cycles to have the same cycle label

pyspatiotemporalgeom.utilities.hsegLibrary.slope(seg)[source]

return the slope of a line segment.

input: A SEGMENT, not a halfsegment.

output: a float

pyspatiotemporalgeom.utilities.hsegLibrary.splitSegsOnIndexedPoints(hsegs, indexedPoints)[source]

break apart halfsegments at the points indicated in indexedPoints. Similar functions are in the segLibrary.

In this case, we must break the hseg, and its brother, and maintain labels appropriately

pyspatiotemporalgeom.utilities.hsegLibrary.switchLabelsForCorrectCycleLabelling(hsegs)[source]

Assumes all cycles have a unique label. Labels can start at any number and can skip numbers.

Use labelUniqueCycles() to guarantee that cycles have unique label numbers.

Assumes that cycles are labeled consistently, but possibly with the exterior on the incorrect side.

This function will flip labels for all cycles that need it.

Approach: find the first seg in each cycle. Do a point in polygon test for the midpoint of that seg.

this will tell us if it is a hole or outer cycle. then flip labels accordingly

input: hsegs: a list of hsegs indicating a well-formed region, but with hole labels labeled as outer cycles.

output: a well formed region with appropriate labels.

pyspatiotemporalgeom.utilities.hsegLibrary.triangulate(hsegs)[source]

polygon triangulator. handles arbitrary polygons – holes, nested holes, convex, etc

pass it an ordered list of hsegs, get a bunch of triangles

Not coded for maximum efficiency, but should be robust.

input: hsegs – an ordered list of hsegs. for example, call region.createRegionFromSegs( ) on a bunch of segs

output: A list of triangles. a triangle is a tuple containing 3 points. triangles are not in a particular CW or CCW order

pyspatiotemporalgeom.utilities.hsegLibrary.writeHsegListToFile(theRegion, theFileName)[source]

write halfsegments to a file.

1 hseg per line

format: x1 y1 x2 y2 la lb

triangleLibrary

Operations involving triangles in 3D

pyspatiotemporalgeom.utilities.triangleLibrary.crossProduct(v, u)[source]

returns cross product of u and v

u and v are 3-tuples representing a 3D point: (x,y,z)

This function is used to compute triangle/triangle intersections

Returns a 3-tuple of numbers (a 3D point): (x,y,z)

pyspatiotemporalgeom.utilities.triangleLibrary.dot(u, v)[source]

returns dot product of vectors u and v

u and v are 3-tuples representing a 3D point: (x,y,z)

This function is used to compute triangle/triangle intersections

Returns a number

pyspatiotemporalgeom.utilities.triangleLibrary.segment3DTriangle3DIntersection(seg, tri)[source]

Test if a triangle and a 3D segment intersect. Whereas some of the other functions assume a vertical segment in the Z dimension, this function assumes NO orientation of the segment. It can be pointing in any direction.

To convert this algorithm to use rays, change the line:

if r <= 0.0 or r >= 1.0:

to:

if r <= 0.0:

Input:

  • seg is a lien segment of the usual format ((x1,y1),(x2,y2))
  • tri is a triangle a tuple of 3, 3D points. Triangle vertices in arbitrary order

Returns: a tuple containing a value indicating the type of intersection, and the actual intersection point (if it exists, None otherwise):

(value, point) where value is:

  • -1: triangle is degenerate (seg or a point)
  • 0: disjoint
  • 1: intersect in unique point
  • 2: intersects in the same plane. (the intersection is a line segment)

and point is :

  • None if there is no intersection point (either no intersection or intersection is a line)
  • (x,y,z) the point of intersection

This function is modeled from triangle ray intersection algorithm at http://geomalgorithms.com/body_a06-_intersect-2.html#intersect3D_RayTriangle(), although this algorithm is different to handle specific intersection cases differently. Comments are included that identify portions specific to our intersection needs The code at http://geomalgorithms.com/body_a06-_intersect-2.html#intersect3D_RayTriangle() has the following copyright:

//Copyright 2001 softSurfer, 2012 Dan Sunday
// This code may be freely used and modified for any purpose
// providing that this copyright notice is included with it.
// SoftSurfer makes no warranty for this code, and cannot be held
// liable for any real or imagined damage resulting from its use.
// Users of this code must verify correctness for their application.
pyspatiotemporalgeom.utilities.triangleLibrary.triangle3DTriangle3DIntersection(tri1, tri2)[source]

Input: two triangles of the form: ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3))

Warning

triangles always have points [0] and [1] from the same slice. In other words, for the representation listed above, z1 = z2 ALWAYS.

Both triangles cover the same time interval. So, the the Z values will be the same among both intput triangles.

So ... always make segs with [0][2] and [1][2] for tri tri intersection

Triangles that share an edge are not considered intersecting. This function does return true if two triangles meet at a point and that point lies on the boundary of the triangles.

Input:

tri1, tri2

Two triangles of the form:

tri1 = ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3))

tri1 = ((x4,y4,z4),(x5,y5,z5),(x6,y6,z6))

Where z1 == z2 and z4==z5 and ((z1 == z4 and z3 == z6) or (z1 == z6 and z3 == z4))

Also, the triangles cover the same interval in the Z dimension. so a Z value that appears in tri1 also appears in tri2

Returns:

True
if the triangles intersect in their interior or at a point along their boundaries
False
if the triangle’s interiors do not intersect.

Warning

triangle with shared edges but disjoint interiors return False

This function uses the general intersection test that makes no assumptions about the orientations of triangles.

pyspatiotemporalgeom.utilities.triangleLibrary.triangle3DTriangle3DIntersectionPoint(tri1, tri2)[source]

Input: two triangles of the form: ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3))

Warning

triangles always have points [0] and [1] from the same slice. In other words, for the representation listed above, z1 = z2 ALWAYS.

Both triangles cover the same time interval. So, the the Z values will be the same among both intput triangles.

So ... always make segs with [0][2] and [1][2] for tri tri intersection

Triangles that share an edge are not considered intersecting. This function does return true if two triangles meet at a point and that point lies on the boundary of the triangles.

This function returns a single point of intersection. If a segment from one triangle is coplanar with the other triangle, this function will not return an intersection point. triangle3DTriangle3DIntersection( tri1, tri2 ) will consider such a pair of triangles to be overlapping and will return true. However, THIS function only returns an intersection point if there is a SINGLE intersection point between a line segment from one triangle and the other triangle.

Input:

tri1, tri2

Two triangles of the form:

tri1 = ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3))

tri1 = ((x4,y4,z4),(x5,y5,z5),(x6,y6,z6))

Where z1 == z2 and z4==z5 and ((z1 == z4 and z3 == z6) or (z1 == z6 and z3 == z4))

Also, the triangles cover the same interval in the Z dimension. so a Z value that appears in tri1 also appears in tri2

Returns:

A list of points where a point is defined as (x,y)

[(x1,y1,zn), ..., (xn,yn, zn)]
For every point at which a triangle edge intersects the opposing triangle at a single point, that point will appear in the list.
[]
An empty list, otherwise.

Warning

triangle with shared edges but disjoint interiors return False

This function uses the general intersection test that makes no assumptions about the orientations of triangles.

pyspatiotemporalgeom.utilities.triangleLibrary.verticalSegment3DAndSegment3DIntersection(verticalSegProjectedTo2D, seg3D)[source]

find the intersection point of a vertical segment (projected out of the Z dimension) and 3D segment (if it exists.). If there is no intersection, a None value is returned.

Input:

  • verticalSegProjectedTo2D: a point tuple (x,y) defining a vertical line in 3D
  • seg3D: a 3d line segment tuple: ((x1,y1,z1),(x2,y2,z2))

Returns:

  • A list() containing 1 or 2 intersection points. A point is a tuple (x,y,z). The x and y values of the intersection point will always be equal to the values in lineProjected if the intersection point exists. Or...
  • None if there is no intersection
pyspatiotemporalgeom.utilities.triangleLibrary.verticalSegment3DAndTriangle3DIntersection(ray, tri)[source]

This function returns the points of intersection between a triangle in 3D and vertical 3D segment.

If a ray and a traingle are co-planar and intersect, they will intersect at 2 points (unless the ray grazes the tip of the triangle, in which case they will intersect at a single point). If the ray is collinear and overlapping a segment, this will return the segment end points.

Warning

Again, as with all the functions in this file, we are assuming rays that are vertical in the z (temporal) dimension and triangles in the usual format: tri: ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3)) where z1 == z2. So... the triangle always has a segment that is in the Z plane.

The function will check explicitly for intersection between the segment and each triangle edge (using a leftHandTurn test). The approach for this is to project the triangle segments to 2d space; then, find which projected segs s1,s2 contain the point p induced by projecting the vertical segment out of the Z dimension. Finally, the z dimension of the intersection of s1 and p are computed, and likewise for s2.

Becuase 1 triangle seg is always planar in the Z dimension, we refer to the triangle edges as the base (the one planar in the Z dimension) and stem1 and stem2. The base may be at the top (highest Z value) or bottom of the triangle. The stems always share one end point.

Input:

  • tri: a triangle of the form: ((x1,y1,z1),(x2,y2,z2),(x3,y3,z3))
  • (triangles always have points [0] and [1] from the same slice. In other words, for the representation listed above, z1 = z2 ALWAYS.)
  • ray: a line segment of the form ((x1,y1,z1),(x2,y2,z2)) that is vertical in the Z dimension and that spans the Z dimension of the triangle

Returns:

  • a list containing 0, 1, or 2 3D points: [], [((x,y,z),(x2,y2,z2))] or [((x,y,z))]. The segment may graze the triangle at a point or along a segment, it may enter and then leave a co-planar triangle, or it may miss the triangle.

mapLibrary

pyspatiotemporalgeom.utilities.mapLibrary.breakSegsAndPreserveSegLabels(segs, la, lb)[source]

segs is a list of segments that defines cycles. It is of the usual format: [((x1,y1),(x2,y2)), ... , ((xn,yn),(xm,ym))]

la and lb are lists of integers, are parallel to segs, and carry label infromation for the area above (la) and below (lb) its associated seg.

  1. Each cycle is assumed to be a well-formed, simple cycle.
  2. Each cyle is assumed to carry a unique label.
  3. Segments can be in ANY order in segs, just as long as the list is parallel to la and lb
  4. la and lb are parallel to segs. If the interior of the cycle lies above the segment, then la will have the cycle’s labeland lb will have a value of 0. Vice versa if the interior of the cycle lies below the segment. Note that a cycle may correspond to a hole geometry within a region, in which case the cycle will be labeled such that the cycle label lies on the exterior of the cycle.
  5. Segments from one cycle can intersect segments from another cycle, but indeividual cycles are simple (non-self intersecting)

The basic idea is to label a bunch of cycles using, for instance hsegLibrary.labelUniqueCycles(), then just append them all to a segment list and label lists, and call this function.

Note that duplicate segments are maintained with their original labels. Therefore, this function can be used to break segments intersecting segments between cycles, then simply retrieve the individual cycles by label number.

Input:

  1. segs: a list of segments. It is assumed that the segments form simple cycles. Each simple cycle is well formed. Segments from different cycles can intersect.
  2. la: a list of labels indicating if the interior of a cycle lies above the corresponding segment. la is parallel to lb and segs
  3. lb: same as la, but indicates if the interior of a cycle lies below the corresponding segment.

Output:

  1. a list of segments. No two segments in the list will intersect in thier interiors. The only exception is that overlapping segments from the input will be duplicated. One instance of the overlapping seg will exist for each cycle that contains the seg.
  2. a list of above labels
  3. a list of below labels
pyspatiotemporalgeom.utilities.mapLibrary.createMapFromCycles(segs, la, lb)[source]

segs is a list of segments that defines cycles. It is of the usual format: [((x1,y1),(x2,y2)), ... , ((xn,yn),(xm,ym))]

la and lb are lists of integers, are parallel to segs, and carry label infromation for the area above (la) and below (lb) its associated seg.

The purpose of this function is to create a map from a collection of individual cycles, so two cycles may overlap.

Because we are creating a map, the idea is to treat the input cycles as imposing a spatial partition on the embedding space. This function will assign to each region in that partition the labels of all cycles that cover that region. Thus, an area covered by 2 cycles will carry the labels of both of those cycles.

It is assumed that input segments have correct labeling.

Input:

  1. segs a list of segs in the format [ ((x1,y1),(x2,y2)), ... ,((xn,yn),(xm,ym))]. The segs form cycles. Segments cannot intersect each other, excet for overlapping segs where two cycles share a boundary.
  2. la a list of labels that is a parallel list to the segments (item 0 in la corresponds to item 0 in lb corresponds to item 0 in segs, etc).
  3. lb a list of labels that is a parallel list to the segments (item 0 in la corresponds to item 0 in lb corresponds to item 0 in segs, etc).

la and lb are lists of integers. It is assumed each cycle simply has a single integer as a label for this function

Returns:

  1. A list of segs forming a map
  2. Parallel la and lb arrays with label IDs corresponding to the segs
  3. A dictionary to translate label indexes in the la and lb arrays to a set indicating all input labels for all input cycles that cover the corresponding region.

Warning

This function assumes that the exterior is labeled with 0

convexHull

Convex hull computations.

pyspatiotemporalgeom.utilities.convexHull.angleFromVertical(seg)[source]

Computes the angle (in degrees) formed by extending a ray vertically through the leftmost end point of the segment, and rotating it counter-clockwise until overlaps with the segment.

input

seg: A line segment in the format ((x1,y1),(x2,y2))

output

the angle in degrees

pyspatiotemporalgeom.utilities.convexHull.getHullFromSequence(pois)[source]

Returns the convex hull of a list of points. Uses the Graham’s Scan technique (with associated floating point rounding problems)

input:

pois: A list of points of the format: [(x1,y1),(x2,y2),...(xn,yn)]

output:

a list of points that, traversed in order, form the convex hull of the input points
pyspatiotemporalgeom.utilities.convexHull.getHullFromSequenceMelkman(pois)[source]

!!! This function does not quite work, but the general algorithm is there. Use as a reference only!!

Returns the convex hull of a list of points. Uses the Melkmans’ Algorithm technique (with associated floating point rounding problems)

input:

pois: A list of points of the format: [(x1,y1),(x2,y2),...(xn,yn)]

output:

a list of points that, traversed in order, form the convex hull of the input points

pyspatiotemporalgeom.utilities.convexHull.isLeftTurn(p1, p2, p3)[source]

Indicates if the person traveling from point p1 to point p2 and then on to point p3 must take a left or right turn at p2 in order to reach p3.

The approach is to use the sign of the area of the triangle formed by the points (actually the square, since the area is not divided by 2)

input:

p1,p2,p3: points of the form (x,y)

output:

-1 – if a left turn was made 0 – if the points are collinear 1 – if a right turn was made

regionInterpolator

Region interpolation functions.

pyspatiotemporalgeom.utilities.regionInterpolator.angleFromVertical(seg)[source]

A utility to compute the angle a segment forms with a vertical line emanating downwards from the least end point of the segment.

Input: A line segment ((x1,y1),(x2,y2))

Returns: an angle in degrees

pyspatiotemporalgeom.utilities.regionInterpolator.appendToAnimationFile(triList, fileObject, numFrames=100)[source]

For an triangle list defining a moving region over a single time interval (for instance, one element of the tuple returned by interpolateRegions()), write snapshots of that region as it moves over the time interval. In essence, create files that show an animation of the moving region.

Input –

triList: A list of triangles in the usual format: ( (x0,y0,z0), (x1,y1,z1), (x2,y2,z2) )

fileObject: An open file object that the snapshot will be written to.

numFrames: The number of snapshots to generate (the number of frames in thr resulting animation if one were to create a movie out of them).

pyspatiotemporalgeom.utilities.regionInterpolator.checkBoundaryExistence(tris)[source]

Becuase of the way the interpolation algorithm works, some segments do not actually exist at the temporal extrema of the interval. This function appends a boolean to each triangle indicating if the segment portion of the trinagle should be printed (exists) at the temporal extrema of the interval.

pyspatiotemporalgeom.utilities.regionInterpolator.createMotionPlan(r1, r1Hull, r1ConCycles, r1PoiToCycLabelDict, r2, r2Hull, r2ConCycles, r2PoiToCycLabelDict, startTime, endTime)[source]

Used in the interpolateRegions function to actuallt create trinagles.

create motion triangles from r1 to r2. start time must be earlier than end time. r1 must be the earlier (source) region, r2 the later (destination) region

walk the outer cycles (cycles with a 2 label from the prep functions above) and convex hulls concurrently any concavities map to a single point on the opposing region. Any cycles are treated as concavities. connectors other than connectors with a 2,-2 label will disappear...poof

pyspatiotemporalgeom.utilities.regionInterpolator.interpolate(r1, r2, startTime, endTime, noTriIntersectionChecks=False)[source]

This is where the magic happens. Create an interpolation between two well-formed regions over a time interval (defined by startTime and endTime) such that at every instant within that time interval, the region resulting from the interpolation at that instant conforms to the type definition of complex regions as defined in [1]. Note that the various region generators and region creation functions int he region.py file create well formed regions according to [1]. In otherwords, the moving region resulting from this function conforms to the type definition of moving regions in [2].

This function is an extension of the algorithm in [3] to handle both simple (1 simple cycle with no holes) regions and complex regions.

[1] Markus Schneider and Thomas Behr. 2006. Topological relationships between complex spatial objects. ACM Trans. Database Syst. 31, 1 (March 2006), 39-81. DOI=10.1145/1132863.1132865 http://doi.acm.org/10.1145/1132863.1132865

[2] Ralf Hartmut Guting, Michael H. Bohlen, Martin Erwig, Christian S. Jensen, Nikos A. Lorentzos, Markus Schneider, and Michalis Vazirgiannis. 2000. A foundation for representing and querying moving objects. ACM Trans. Database Syst. 25, 1 (March 2000), 1-42. DOI=10.1145/352958.352963 http://doi.acm.org/10.1145/352958.352963

[3] Mark McKenney and James Webb. 2010. Extracting moving regions from spatial data. In Proceedings of the 18th SIGSPATIAL International Conference on Advances in Geographic Information Systems (GIS ‘10). ACM, New York, NY, USA, 438-441. DOI=10.1145/1869790.1869856 http://doi.acm.org/10.1145/1869790.1869856

Input:

  1. r1, r2: two well formed regions represented as lists of hlafsegments. Any of the region creation functions in region.py will do.
  2. startTime, endTime: two numbers defining a time interval. These numbers are used as the 3D dimension when extrapolating into 3D space.
  3. noTriIntersectionChecks. See paper [3]. The algorithm first creates an interpolation between the input regions. It is possible that the interpolation will result in a self-intersecting region at some point. The triangle/triangle intersection test is then performed. This test is very computationally intensive (especially for python) and can take a LONG time to compute. If you pass a True for this argument, the tri/tri intersection test is skipped, and the interpolation returned AS-IS (possibly with self-intersections). This makes the algorithm O(n \lg n) instead of O(n^2).

Output:

A 3-tuple. See [3]. The algorithm will create at MOST, 3 interval regions to describe the interpolation of r1 to r2 over the defined time interval. Not all 3 interval regions are always required, so 1 or 2 of the values in the tuple may be None, but a 3 tuple is ALWAYS returned. If the noTriIntersectionChecks argument is set to True, or the original interpolation succeeds, then the return value will look like this: (None, triList, None).

Any non-None value in the return tuple will be a list of trinagles describing the movinement of line segments in r1 as they travel across the defined interval to r2 (between intermediate states of r1 and r2 if necessary).

pyspatiotemporalgeom.utilities.regionInterpolator.prepRegionForInterpolation(hsegs)[source]

This function is used to prepare regions for interpolation.

It is expected that this function is only useful for the interpolateRegions() function.

First, connected cycles in a region are identified (cycles that touch). Then connected cycles are relabeled so the interior label is identical for each group of connected cycles. Then vertical connector segments are added. The point sequence of an outercycle walk is computed. The convex hull of that sequence is computed. Finally, labels are assigned for concavity mappings.

ASSUMPTIONS:

  1. hsegs is a valid region in hseg order returned from region.createRegionFromSegs() or any of the region creation functions.
  2. Each cycle in the region has its own unique label, and -1 for exterior labels.

RETURNS teh following as a tuple:

  1. Returns the region with cycles labeled for interoplation
  2. Returns the convex hull points in CCW sequence
  3. Returns the cycle mapping.
  4. Returns list of connected cycles
  5. Returns a dict that maps a point to a list of cycles connected to that point
pyspatiotemporalgeom.utilities.regionInterpolator.writeHsegsToTextFile(fileName, hsegs)[source]

Used in debugging. Open a file named filename and clobber it if it exists. Then write all the halfsegments in the list hsegs to that file. 1 hseg per line, each point and label seperated by a space.

pyspatiotemporalgeom.utilities.regionInterpolator.writeTrisToFile(triTup, fileName)[source]

Write triangles to a file.

This was used during debugging.

expects an iterable containing iterables ( a list of lists).

Each iterable contained in the first contains triangles.

A trianlge is a 3 tuple of points (p0, p1, p2) == ( (x0,y0,z0), (x1,y1,z1), (x2,y2,z2) ) or 3d points

A point is a 3 tuple.

Input:

triTup = (listofTris1, listofTris2,...)

where listofTrisX = [t1, t2, t3,...]

where tX = (p1,p2,p3)

where pX=(x0, y0, z0,)

Prints the contents to a file named fileName. will clobber that file

regionOverlayGrid

Compute region overlays using a grid decomposition scheme.

pyspatiotemporalgeom.utilities.regionOverlayGrid.assignCells(r1, cellSize, r1Dict, usedCellSet)[source]

find cells using left hand turn test

pyspatiotemporalgeom.utilities.regionOverlayGrid.assignCellsOld(r1, cellSize, r1Dict, usedCellSet)[source]

This is DEPRECATED. do not use. left over as a reference.

Assign segments to cells. Uses a modified Bresenham’s algorithm

r1
list of segments defining regions
cellSize
size of a cell
r1Dict

Output value: Modified by the function

a dictionay of cell->list of segs. cell is just an x,y. Segs are ((x1,y1),(x2,y2))

usedCellsSet
A set of cell addresses used. Output value: Modified by the function.
Notes:
If a segment is placed in 2 adjacent cells in the Y direction, it will be placed in the uppermost cell once, and twice in the lower cells. This is so the ray shooting test will still work.
pyspatiotemporalgeom.utilities.regionOverlayGrid.breakSegments(r1, cellsPerSeg=3)[source]

given a bunch of line segments, break them at intersection points. Handles overlapping segments appropriately (duplicates are removed).

Input:

r1
A list of segments of the form ((x1,y1),(x2,y2))
cellsPerSeg
an integer indicating how many cells (roughly) a segment should fall into. This is computed based on the average length of line segments in r1.

Output:

A list of segments such that segments only overlap at end points.

pyspatiotemporalgeom.utilities.regionOverlayGrid.breakSegmentsRedBlue(r1, r2, cellsPerSeg=3)[source]

Red Blue line intersection. Self intersections in the red set or blue set are not discovered

Given a bunch of line segments, break them at intersection points. Handles overlapping segments appropriately (duplicates are removed).

Input:

r1
A list of segments of the form ((x1,y1),(x2,y2))
cellsPerSeg
an integer indicating how many cells (roughly) a segment should fall into. This is computed based on the average length of line segments in r1.

Output:

A list of segments such that segments only overlap at end points.

pyspatiotemporalgeom.utilities.regionOverlayGrid.computeTopology(r1, usedCells2, cellDict2, cellSize)[source]

This is meant to be called from the overlay function

Computes the topolgy of r1 with respect to a region r2. r2 is not passed directly; rather, a cell dictionary that maps cells to the segments that lie in that cell are passed. Also, a set of cells that appear in cellDict2 (the keys of the dict) are passed.

NOTE: cellSize should be the size used to create the cellDict2

INPUT:

r1
a list of segs where each seg is has an index, and a interiorBelow label: ((x1,y1),(x2,y2),index,interiorBelow)
usedCells2
They keys of the cellDict2
cellDict2
A dict mapping a cell to the segs in that cell
pyspatiotemporalgeom.utilities.regionOverlayGrid.overlay(r1, r2, cellsPerSeg=3)[source]

Perform an overlay of r1 and r2

Warning

We assume the first end point in the tuple is less than the second. Otherwise, the algorithm will not return correct results.

Input:

r1, r2

A list of labeled segments. Segments have the following structure:

((x1,y1),(x2,y2),interiorBelow) where interiorBelow > 0 if the interior of the region lies below the segment, < 0 otherwise

cellsPerSeg
an integer indicating how many cells (roughly) a segment should fall into. This is computed based on the average length of line segments in r1 and r2.

Output:

(resultR1, resultR2)
A tuple containing two lists. Each list is the respective input regions line segments broken where they intersect line segments from the other region. Furthermore, each line segment carries a label consisting of a 1 or a 0 and boolean flag. The boolean flag indicates if a line segment is equal to a line segment in the opposing result region. The label indicates if the seg lies in the interior of the opposing region, or if it is equivalent to a line seg in the opposing region, if the interior of the opposing region lies below the seg.
pyspatiotemporalgeom.utilities.regionOverlayGrid.overlayHalfsegmentInput(r1, r2, cellsPerSeg=3)[source]

A wrapper for def overlay(r1, r2, cellsPerSeg = 3 ):

This function accepts regions defined by a list of halfsegments.

A halfsegment has the structure ((x1, y1),(x2,y2),la, lb) where la and lb indicate if the interior if the region lies above or below the shalfsegment. An interior label is greater than 0. An exterior label is 0 or less.

temporalAggregate

Functions to compute temporal aggregate operations over moving regions.

pyspatiotemporalgeom.utilities.temporalAggregate.createIntervalRegionAndComputeTemporalAggregate(R1, R2)[source]

Interpolate regions R1,R2. Use 1 of the interval regions returned from the interpolate function to call the aggregate function.

Get the points with maximal coverage.Find if it contains cycles.

Extracts points with maximum coverage,segments that are not part of cycles,cycles with maximal temporal coverage for the given regions if present

input:

R1,R2: two regions which are an ordered list of half segments with valid labels.

A half segment is a tuple with the following:
((x,y)(x,y), labelAbove, labelBelow)

labels are integral values

an interior label will be positive. An exterior label will be -1. The side of the hseg on which the interior (exterior) of the region lies will have an interior (exterior) label.

hsegs will be ordered in half segment order

Each cycle in hsegs will have its own unique label number

returns:

tuple: tuple consists of the :
  1. points with maximum duration :

    [((x,y),[z]),((x,y),[z])..] if present else [‘None’]

  2. segments that are not part of cycles :

    [((x1,y1),(x2,y2)),((x1,y1),(x2,y2))..] if present else [‘None’]

  3. cycles with maximum duration.

    [((x,y)(x,y), labelAbove, labelBelow)..] if present else [‘None’]

pyspatiotemporalgeom.utilities.temporalAggregate.createProjectedSegsWithIDs(triangles)[source]

Extracts segments out of the given triangles and assigns a numeric value to each segment

Input:

  • triangles:

Output:

  • a list of segments that have a unique identifier