A recent publication by Cogné et al. used a modified iolite 2/3 data reduction scheme to calculate ‘weighted’ U/Ca as part of an apatite fission track dating protocol. In this note, we will go through recreating that functionality as an iolite 4 user interface (UI) plugin.
Creating a new UI plugin
Creating a new UI plugin can be done by adding a python file to iolite’s UI plugins path with appropriate meta data and required function definitions. To make this process a bit easier, you can create a new plugin from the Tools → Plugins → Create menu item. In the dialog that pops up you can set the type to UI and enter reasonable things for the other parameters.
The first thing we need to do is define some meta data at the top of the python file. The only part that really matters is that the type is specified as UI. In doing so, we tell iolite to handle this script as a UI plugin and look for the associated functions (i.e. createUIElements below). The meta data that I wrote was as follows.
1
2
3
4
5
6
7
|
#/ Type: UI
#/ Name: WeightedUCa
#/ Description: Calculates weighted U/Ca
#/ Authors: J. Petrus
#/ References:
#/ Version: 1.0
#/ Contact: joe@iolite-software.com
|
Integration in iolite’s interface
As mentioned above, iolite looks for a function called createUIElements in UI plugins to establish hooks in the user interface that allow a user to interact with the plugin. In this example, we’ll just create a menu item that triggers the calculation by calling (a yet to be defined) run function. The way we do this is as follows.
1
2
3
4
5
6
7
|
from iolite.QtGui import QAction
def createUIElements():
action = QAction('Calculate Weighted U/Ca')
action.triggered.connect(run)
ui.setAction(action)
ui.setMenuName(['Tools'])
|
Calculating unweighted U/Ca
The first step in working out the weighted U/Ca is to calculate the unweighted U/Ca. The easiest way to do that is to use the U and Ca concentration channels previously calculated by a trace elements data reduction scheme. Recalling that you can get the data associated with a channel as data.timeSeries(channel_name).data() this is easily accomplished as follows.
1
2
3
|
# Create U/Ca from ppm channels
U_Ca = data.timeSeries('U238_ppm').data()/data.timeSeries('Ca43_ppm').data()
U_Ca_channel = data.createTimeSeries('U/Ca', data.Output, None, U_Ca)
|
Note that in the call to createTimeSeries we specify data.Output as the channel type and None as the time array. By specifying None as the time array it will default to whatever the current index time is.
Ensuring we have beam seconds
In the calculations that follow, we would like to use beam seconds and information about pit depth to work out the actual depth. One way we can do that is as follows.
1
2
3
4
|
try:
bs = data.timeSeries('BeamSeconds')
except:
bs = data.createBeamSeconds('log')
|
Here we have used a try-except pair to first try to get an existing beam seconds channel, and if that fails, create a beam seconds channel using the laser log. Note that the createBeamSeconds method is only available in iolite > 4.3.4.
Dialog and settings
One thing we would like to be able to do is specify a few input parameters to the calculation, such as pit depth, fission track length, and analysis duration. We can accomplish this by using QSettings and various components of QtGui as follows.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
|
# Start by importing the necessary modules
from iolite.QtGui import QDialog, QFormLayout, QLineEdit, QDialogButtonBox
from iolite.QtCore import QSettings
# Now use QSettings to retrieve previously saved values or use defaults
settings = QSettings()
pitDepth = settings.value("WeightedUCa/pitDepth", 12.5)
trackLength = settings.value("WeightedUCa/trackLength", 8.0)
duration = settings.value("WeightedUCa/duration", 30.0)
# Construct the dialog with a QFormLayout
d = QDialog()
d.setLayout(QFormLayout())
# Add an input for depth
depthLineEdit = QLineEdit(d)
depthLineEdit.text = str(pitDepth)
d.layout().addRow("Pit depth (μm)", depthLineEdit)
# Add an input for track length
lengthLineEdit = QLineEdit(d)
lengthLineEdit.text = str(trackLength)
d.layout().addRow("Track length (μm)", lengthLineEdit)
# Add an input for duration
durLineEdit = QLineEdit(d)
durLineEdit.text = str(duration)
d.layout().addRow("Duration (s)", durLineEdit)
# Add some buttons to bottom of the dialog
bb = QDialogButtonBox(QDialogButtonBox.Ok | QDialogButtonBox.Cancel, d)
d.layout().addRow(bb)
bb.accepted.connect(d.accept)
bb.rejected.connect(d.reject)
# Run the dialog requiring a user to click ok or cancel to dismiss it
# If cancelled, return (i.e. do nothing)
if d.exec_() == QDialog.Rejected:
return
# Get the values in the various inputs as floating point numbers
# (rather than strings) and store them in our settings
pitDepth = float(depthLineEdit.text)
settings.setValue("WeightedUCa/pitDepth", pitDepth)
trackLength = float(lengthLineEdit.text)
settings.setValue("WeightedUCa/trackLength", trackLength)
duration = float(durLineEdit.text)
settings.setValue("WeightedUCa/duration", duration)
|
Doing the calculation
Using the appendix of the Cogné et al. publication as a reference, the weighted U/Ca can be calculated as follows.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
|
import numpy as np
from math import pi
# Volume of a sphere
def VS(r):
return (4/3)*pi*r**3
# Partial volume of a sphere filled to h
def PVS(r, h):
return (1/3)*pi*h**2*(3*r - h)
def calcWtd(s, d, dur, l):
'''
s: selection
d: depth
dur: duration
l: track length
'''
# Get the beam seconds and unweighted U/Ca for the specified selection
# as numpy arrays
sbs = data.timeSeries('BeamSeconds').dataForSelection(s)
UCa = data.timeSeries('U/Ca').dataForSelection(s)
x = (d/dur)*sbs # Actual depth
dx = (d/dur)*(sbs[-1] - sbs[-2]) # Depth increment per time-slice
weight = 2*(PVS(l, l-x) - PVS(l, l-x-dx))/VS(l)
weight[weight<0] = 0 # Replace weights < 0 with 0 (i.e. when depth > track length)
# Return the weighted U/Ca
return np.sum(weight*UCa)/np.sum(weight)
|
Doing the calculation for all selections
Now we have the calculation in place taking 1 selection as a parameter, so we just need to iterate through the selections to calculate it for everything. One way to do that is as follows.
1
2
3
4
5
6
7
8
9
10
11
12
13
|
# Create an array of zeros matching the length of the unweighted U/Ca array
wtdUCa = np.zeros(len(U_Ca))
# Loop through groups that are RM or Sample type:
for group in data.selectionGroupList(data.ReferenceMaterial | data.Sample):
# Loop through selections in each group:
for s in group.selections():
# Calculate the weighted U/Ca for this selection and replace the part
# of the wtdUCa array corresponding to this selection with the value
wtdUCa[U_Ca_channel.selectionIndices(s)] = calcWtd(s, pitDepth, duration, trackLength)
# Lastly, create a new channel
data.createTimeSeries('Weighted U/Ca', data.Output, None, wtdUCa)
|
And that’s it! The newly created channel will have the weighted U/Ca and results for it can be exported along with your other data and/or viewed and manipulated in all the usual iolite ways.
Putting it all together
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
|
#/ Type: UI
#/ Name: WeightedUCa
#/ Description: Calculates weighted U/Ca
#/ Authors: J. Petrus
#/ References:
#/ Version: 1.0
#/ Contact: joe@iolite-software.com
from iolite.QtGui import QAction, QDialog, QFormLayout, QLineEdit, QDialogButtonBox
from iolite.QtCore import QSettings
import numpy as np
from math import pi
def createUIElements():
action = QAction('Calculate Weighted U/Ca')
action.triggered.connect(run)
ui.setAction(action)
ui.setMenuName(['Tools'])
def VS(r):
return (4/3)*pi*r**3
def PVS(r, h):
return (1/3)*pi*h**2*(3*r - h)
def calcWtd(s, d, dur, l):
sbs = data.timeSeries('BeamSeconds').dataForSelection(s)
UCa = data.timeSeries('U/Ca').dataForSelection(s)
x = (d/dur)*sbs
dx = (d/dur)*(sbs[-1] - sbs[-2])
weight = 2*(PVS(l, l-x) - PVS(l, l-x-dx))/VS(l)
weight[weight<0] = 0
return np.sum(weight*UCa)/np.sum(weight)
def run():
settings = QSettings()
pitDepth = settings.value("WeightedUCa/pitDepth", 12.5)
trackLength = settings.value("WeightedUCa/trackLength", 8.0)
duration = settings.value("WeightedUCa/duration", 30.0)
d = QDialog()
d.setLayout(QFormLayout())
depthLineEdit = QLineEdit(d)
depthLineEdit.text = str(pitDepth)
d.layout().addRow("Pit depth (μm)", depthLineEdit)
lengthLineEdit = QLineEdit(d)
lengthLineEdit.text = str(trackLength)
d.layout().addRow("Track length (μm)", lengthLineEdit)
durLineEdit = QLineEdit(d)
durLineEdit.text = str(duration)
d.layout().addRow("Duration (s)", durLineEdit)
bb = QDialogButtonBox(QDialogButtonBox.Ok | QDialogButtonBox.Cancel, d)
d.layout().addRow(bb)
bb.accepted.connect(d.accept)
bb.rejected.connect(d.reject)
if d.exec_() == QDialog.Rejected:
return
pitDepth = float(depthLineEdit.text)
settings.setValue("WeightedUCa/pitDepth", pitDepth)
trackLength = float(lengthLineEdit.text)
settings.setValue("WeightedUCa/trackLength", trackLength)
duration = float(durLineEdit.text)
settings.setValue("WeightedUCa/duration", duration)
# Create U/Ca from ppm channels
U_Ca = data.timeSeries('U238_ppm').data()/data.timeSeries('Ca43_ppm').data()
U_Ca_channel = data.createTimeSeries('U/Ca', data.Output, None, U_Ca)
# Make sure we have beam seconds
try:
bs = data.timeSeries('BeamSeconds')
except:
bs = data.createBeamSeconds('log')
# Loop through selections and calculate weighted U/Ca
wtdUCa = np.zeros(len(U_Ca))
for group in data.selectionGroupList(data.ReferenceMaterial | data.Sample):
for s in group.selections():
wtdUCa[U_Ca_channel.selectionIndices(s)] = calcWtd(s, pitDepth, duration, trackLength)
data.createTimeSeries('Weighted U/Ca', data.Output, None, wtdUCa)
|
Click here to discuss.
