Discontinuities

A discontinuity is a plane that is typically used to mark a change in the physical or chemical properties of a rock mass. They are represented in the Maptek Python SDK as a set of coplanar points that are turned into a surface using facets similar to a surface. Unlike a surface, discontinuities do not support any kind of point or facet property. The example below demonstrates how to use a set of coplanar points to define a new discontinuity.

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.

Here is a an example of how to specify both the points and the dip, dip direction and location using code.

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]

If you run the above code, you should produce a discontinuity that looks like the image below.

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. Make sure to run all previous examples before running the example 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: