APSG tutorial - Part 1
APSG defines several new python classes to easily manage, analyze and visualize orientation structural geology data. There are several classes to work with orientation data, namely Vector3 for vectorial data and Lineation, Foliation for axial data
Basic usage
APSG module could be imported either into its own namespace or into the active one for easier interactive work.
[1]:
from apsg import *
Basic operations with vectors in 2D and 3D
Instance of vector object Vector3 could be created by passing three arguments corresponding to 3 components to function vec:
[2]:
u = vec(1, -2, 3)
v = vec(-2, 1, 1)
Alternative ways to create a vector are to pass a single iterable object as list, tuple or array, or to provide geological orientation according to APSG notation (default is direction of dip and dip angle for planar and trend and plunge for linear features)
[3]:
coords = (-2, 2, 3)
a = vec(coords)
b = vec(120, 60)
print(a, b)
Vector3(-2, 2, 3) Vector3(-0.25, 0.433, 0.866)
The instance of a 2D vector Vector2 could be created passing 2 arguments corresponding to 2 components or a single iterable of components to the function vec2:
[4]:
k = vec2(1, 1)
coords = (-1, 2)
l = vec2(coords)
print(k, l)
Vector2(1, 1) Vector2(-1, 2)
Alternatively, single argument is interpreted as direction in degrees (0-top, 90-right, 180-down, 270-left)
[5]:
m = 5 * vec2(120)
m
[5]:
Vector2(-2.5, 4.33)
For common vector operations we can use standard mathematical operators or special methods using dot notation
[6]:
u + v
[6]:
Vector3(-1, -1, 4)
[7]:
u - v
[7]:
Vector3(3, -3, 2)
[8]:
3*u - 2*v
[8]:
Vector3(7, -8, 7)
Its magnitude or length is most commonly defined as its Euclidean norm and could be calculated using abs
[9]:
abs(v)
[9]:
2.449489742783178
[10]:
abs(u + v)
[10]:
4.242640687119285
For dot product we can use the dot method or the @ operator
[11]:
u.dot(v)
[11]:
-1
[12]:
u @ v
[12]:
-1.0
For cross product we can use the cross method
[13]:
u.cross(v)
[13]:
Vector3(-5, -7, -3)
To project vector u onto vector v we can use method proj
[14]:
u.proj(v)
[14]:
Vector3(0.333, -0.167, -0.167)
To find angle (in degrees) between to vectors we use method angle
[15]:
u.angle(v)
[15]:
96.26395271992722
The rotate method rotates a vector around another vector. For example, to rotate vector u around vector v for 45°
[16]:
u.rotate(v, 45)
[16]:
Vector3(2.248, 0.558, 2.939)
Classes Lineation and Foliation
To work with orientation data in structural geology, APSG provides two classes: Foliation to represent planar features and Lineation to represent linear features. Both classes support all Vector3 methods and operators, but it should be noted, that dot and angle respect their axial nature, i.e. angle between two lineations can’t be bigger than 90 degrees.
To create instance of Lineation or Foliation, we can use functions lin and fol. Arguments have similar syntax to vec.
[17]:
lin(120, 60), fol(216, 62)
[17]:
(L:120/60, S:216/62)
We can also cast Vector3 instance to Foliation or Lineation
[18]:
lin(u), fol(u)
[18]:
(L:297/53, S:117/37)
Vector methods for Lineation and Foliation
To find angle between two linear or planar features we can use method angle
[19]:
l1 = lin(110, 40)
l2 = lin(160, 30)
l1.angle(l2)
[19]:
41.59741268003547
[20]:
p1 = fol(330, 50)
p2 = fol(250, 40)
p1.angle(p2)
[20]:
54.696399321975335
We can use cross product to construct a planar feature defined by two linear features
[21]:
l1.cross(l2)
[21]:
S:113/40
or to construct linear feature defined by intersection of two planar features
[22]:
p1.cross(p2)
[22]:
L:278/36
Cross product of planar and linear features could be used to construct a plane defined by a linear feature and the normal of a planar feature
[23]:
l2.cross(p2)
[23]:
S:96/53
or to find the perpendicular linear feature on a given plane
[24]:
p2.cross(l2)
[24]:
L:276/37
To rotate structural features we can use method rotate
[25]:
p2.rotate(l2, 45)
[25]:
S:269/78
Classes Pair and Fault
To work with paired orientation data like foliations and lineations or fault data in structural geology, APSG provides the Pair base class and the derived Fault class. Both classes are instantiated providing dip direction and dip of planar and linear measurements, which are automatically orthogonalized. If misfit is too high, warning is raised. The Fault class expects one more argument providing sense of movement information, either 1 or -1 for normal/reverse movement.
To create instance of Pair, we have to pass two arguments for planar and two arguments for linear features following geological notation to function pair:
[26]:
p = pair(120, 40, 162, 28)
p
[26]:
P:118/39-163/30
[27]:
p.misfit
[27]:
3.5623168411508175
Planar and linear features are accessible using fol and lin properties
[28]:
p.fol, p.lin
[28]:
(S:118/39, L:163/30)
To rotate a Pair instance we can use the rotate method
[29]:
p.rotate(lin(45, 10), 60)
[29]:
P:314/83-237/61
Instantiation of Fault class is similar, we just have to provide an argument to define the sense of movement
[30]:
f = fault(120, 60, 110, 58, 1) # 1 for normal fault
f
[30]:
F:120/59-110/59 N
Note the change in sense of movement after Fault rotation
[31]:
f.rotate(lin(45, 10), 60)
[31]:
F:312/62-340/59 N
For simple fault analyses the Fault class also provides p, t, m, and d properties to get PT-axes, kinematic plane and dihedra separation plane
[32]:
f.p, f.t, f.m, f.d
[32]:
(L:315/75, L:116/14, S:27/85, S:290/31)
Class Cone
General feature type to store small or great circles as a cone. It could be defined by axis, secant and revolving angle or by axis and apical angle. The revolving angle is by default 360° defining full cone. For segments of small circles, the revolving angle could be different.
[33]:
c = cone(lin(140, 50), lin(140, 75))
c
[33]:
C:140/50 [25]
To create small circle segment, you can provide revolving angle
[34]:
c = cone(lin(90, 70), lin(45,30), 115)
The tips of small circle segments could be obtained from the cone properties secant and rotated_secant
[35]:
lin(c.secant), lin(c.rotated_secant)
[35]:
(L:45/30, L:137/30)
To define cone using apical angle, use axis and number arguments. Note, that revolving angle is 360 by default
[36]:
c = cone(lin(140, 50), 25)
c
[36]:
C:140/50 [25]
[37]:
c.revangle
[37]:
360.0
Feature sets
APSG provides several classes to process, analyze, and visualize sets of data. There are e.g. vecset, linset and folset classes to store groups of vec, lin, and fol objects. All these feature sets are created from homogeneous list of data with optional name attribute.
[38]:
v = vecset([vec(120,60), vec(116,50), vec(132,45), vec(90,60), vec(84,52)], name='Vectors')
v
[38]:
V3(5) Vectors
[39]:
l = linset([lin(120,60), lin(116,50), lin(132,45), lin(90,60), lin(84,52)], name='Lineations')
l
[39]:
L(5) Lineations
[40]:
f = folset([fol(120,60), fol(116,50), fol(132,45), fol(90,60), fol(84,52)], name='Foliations')
f
[40]:
S(5) Foliations
To simplify interactive group creation, you can use function G
[41]:
g = G([lin(120,60), lin(116,50), lin(132,45), lin(90,60), lin(84,52)], name='L1')
g
[41]:
L(5) L1
Method len returns number of features in group
[42]:
len(v)
[42]:
5
Most of the vec, lin and fol methods could be used for feature sets as well. For example, to measure angles between all features in group and another feature, we can use method angle:
[43]:
l.angle(lin(110,50))
[43]:
array([11.49989817, 3.85569115, 15.61367789, 15.11039885, 16.3947936 ])
To rotate all features in group around another feature, we can use method rotate
[44]:
lr = l.rotate(lin(150, 30), 45)
To show data in list you can use data property
[45]:
l.data
[45]:
(L:120/60, L:116/50, L:132/45, L:90/60, L:84/52)
[46]:
lr.data
[46]:
(L:107/35, L:113/26, L:126/30, L:93/26, L:94/18)
Method R returns the resultant of all features in set. Note that Lineation and Foliation are axial in nature, so resultant vector is not reliable. Check the orientation tensor analysis below.
[47]:
v.R()
[47]:
Vector3(-0.941, 2.649, 3.993)
[48]:
lin(v.R())
[48]:
L:110/55
There are several methods to infer spherical statistics as spherical variance, Fisher’s statistics, confidence cones on data etc.
[49]:
l.var()
[49]:
0.02337168447438498
[50]:
v.fisher_statistics()
[50]:
{'mu': Vector3(-0.193, 0.542, 0.818),
'k': np.float64(42.786817573888605),
'alpha': 13.26402990511733,
'csd': np.float64(12.383118421320239),
'uniform': False}
[51]:
v.fisher_cone()
[51]:
C:110/55 [13.264]
[52]:
v.fisher_cone_csd()
[52]:
C:110/55 [12.3831]
[53]:
v.delta()
[53]:
12.411724720740516
[54]:
v.rdegree()
[54]:
95.32566310512297
To calculate orientation tensor of all features in group, we can use ortensor method.
[55]:
v.ortensor()
[55]:
OrientationTensor3
[[ 0.074 -0.096 -0.143]
[-0.096 0.284 0.421]
[-0.143 0.421 0.642]]
(S1:0.977, S2:0.201, S3:0.0758)