Discontinuities

A discontinuity is a plane that is typically used to mark a change in the physical or chemical properties of a rock mass. A discontinuity is represented in the Maptek Python SDK as a set of coplanar points that are visualised as a surface-like object. Unlike a surface however, discontinuities do not support any kind of point or facet properties.

The example below demonstrates how to use a set of coplanar points to define a new Discontinuity object.

from mapteksdk.project import Project
from mapteksdk.data import Discontinuity

points = [[0, 0, -0.5], [1, 0, 0], [0, 1, -2], [1, 1, -1.5]]

project = Project()

with project.new("geotechnical/-x+3y+2z+1=0", Discontinuity) as plane:
  plane.planar_points = points

The plane of the discontinuity is defined using the first three points. Any points beyond the first are projected to be on the plane defined by the first three points. This is not necessary for the discontinuity created in the above example because all four points used lie on the plane -x + 3y + 2z + 1 = 0. The discontinuity is shown below.

Alternatively, you can create a discontinuity by defining the dip, dip direction and location of the discontinuity, where:

• Dip is the angle the discontinuity is rotated by in the dip direction.

• Dip direction is the angle about the z-axis that the discontinuity is rotated in.

• Location is the coordinate of the centre of the discontinuity.

When you use this method to create a discontinuity, the planar points are automatically generated. However, because no points were provided when creating the discontinuity, querying the length or area of the discontinuity results in NaN (“Not a Number”). The example below shows creating a discontinuity in this way.

import math

from mapteksdk.project import Project
from mapteksdk.data import Discontinuity

project = Project()

with project.new("geotechnical/angled_discontinuity", Discontinuity) as plane:
  plane.dip_direction = math.radians(32.5)
  plane.dip = math.radians(16.35)
  plane.location = [-12, 42.8, 32.7]

The image below shows the discontinuity created by the above example. Note that the annotations match with the angles in degrees used to create the discontinuity.

The most robust method for creating a discontinuity requires specifying both the points and the dip, dip direction and location. If you use this method to create a discontinuity, the points are projected onto the plane defined by the dip and dip direction and translated to be centred around the specified location. This allows for the shape of the discontinuity and its properties to be specified independently.

The example below shows how to create a Discontinuity by specifying both the points and the dip, dip direction and location.

import math

from mapteksdk.project import Project
from mapteksdk.data import Discontinuity

project = Project()

# These are the points of a unit hexagon.
sqrt3_2 = math.sqrt(3) / 2
points = [[0, 0, 0], [1, 0, 0], [0.5, sqrt3_2, 0], [-0.5, sqrt3_2, 0],
          [-1, 0, 0], [-0.5, -sqrt3_2, 0], [0.5, -sqrt3_2, 0]]

with project.new("geotechnical/both_discontinuity", Discontinuity) as plane:
  plane.planar_points = points
  plane.dip_direction = math.radians(32.5)
  plane.dip = math.radians(16.35)
  plane.location = [-12, 42.8, 32.7]

Running this script produces the following discontinuity:

Discontinuities possess other properties aside from those shown in the previous examples. The following example queries many of these properties for the discontinuities created in the previous examples. Be sure to run all previous example scripts before running the one below.

import math

from mapteksdk.project import Project

project = Project()

discontinuities = ("-x+3y+2z+1=0", "angled_discontinuity", "both_discontinuity")

for discontinuity in discontinuities:
  path = f"geotechnical/{discontinuity}"
  print(f"Properties for: {discontinuity}")
  with project.read(path) as plane:
    print("Dip:", math.degrees(plane.dip))
    print("Dip direction:", math.degrees(plane.dip_direction))
    print("Strike:", math.degrees(plane.strike))
    print("Plunge", math.degrees(plane.plunge))
    print("Trend", math.degrees(plane.trend))
    print("Location:", plane.location)
    print("Area:", plane.area)
    print("Length:", plane.length)
    print("Colour", plane.planar_colour)
  print("-" * 25)

Polarity

Vulcan GeologyCore 2023 introduced the concept of discontinuity polarity. This allows for a discontinuity to define whether it is UPRIGHT, OVERTURNED, or UNKNOWN. The discontinuity polarity property allows for setting and reading this information in Python. The following script demonstrates creating a discontinuity with each possible polarity:

from mapteksdk.project import Project
from mapteksdk.data import Discontinuity, Polarity, Text2D, HorizontalAlignment
from mapteksdk.operations import open_new_view

if __name__ == "__main__":
    with Project() as project:
        objects_to_view = []
        for i, polarity in enumerate(Polarity):
            with project.new(f"geotechnical/examples/label_{polarity.name}", Text2D) as text:
                text.location = [1, i * 2 - 2, 0]
                text.text = polarity.name
                text.horizontal_alignment = HorizontalAlignment.CENTRED
                objects_to_view.append(text.id)
            with project.new(f"geotechnical/examples/{polarity.name}", Discontinuity) as plane:
                plane.location = [0, i * 2 - 2, 0]
                plane.polarity = polarity
                plane.planar_colour = [0, 255, 0, 255]
                plane.back_colour = [100, 100, 100]
                objects_to_view.append(plane.id)
        open_new_view(objects_to_view)

The following animation shows the three discontinuities created by the above script, animated by switching between view from above and view from below: