-
Notifications
You must be signed in to change notification settings - Fork 16
SHRIMP: Sub ThUfromA1A2
The subroutine evaluates, on an analysis-by-analysis basis, the 232Th/238U, by ratioing the 206Pb/238U calibration constant and the 208Pb/232Th calibration constant. Note that this subroutine only applies when Switch.DirectAltPD = TRUE, and piNumDauPar = 2 (i.e. it is essential that calibration constants for both U-Pb and Th-Pb systems have been calculated: Perm2 and Perm4).
ThUfromA1A2(Std, SpotRow, Only1)
Std: Boolean input which dictates, for each specific invocation of ThUfromA1A2, whether the analysis is that of a 'primary-ratio' reference material, or not.
SpotRow: Index number of the output-row, to which the value calculated by ThUfromFormula is to be written.
Only1: Boolean input which dictates, for each specific invocation of ThUfromA1A2, whether the 232Th/238U value column is to be populated with the dummy value '1', or not. DEFAULT value = FALSE.
According to Ludwig, Only1 = TRUE happens when "the two WtdMeanA formulae and values have not yet been calculated, [it is necessary to] pass first time & put dummy values of 1. Put in real formulae later." With context from both invocations, it might be possible to rework this later.
Values of type Boolean
Std, Only1
Values of type Integer
p, SpotRow, Th2U8col, Th2U8ecol
Values of type String
Exp232, Exp238, t1, t2
The function of the subroutine is paraphrased as follows:
--Recalling that Std = -1 if TRUE, 0 if FALSE, use p to define
--sheet/column for output. Note also t2 is explicitly cleared:
p = -Std
Th2U8col = piaTh2U8col[p]
Th2U8ecol = piaTh2U8ecol[p]
t2 = ""
As per subroutine Tot68_82_fromA, there is a slight complication here, arising from the more restricted function of SQUID 2.50 relative to that envisaged/implemented in SQUID 3.0. You will recall that for SQUID 3.0, I requested all the permissible index-isotope ratio-combinations to be calculated for each permutation (e.g. the pairs of calibration constants calculated in Perm2 and Perm4 can have 204corr- and 207corr- variants: see the set of Cases in the Intro of the Procedural Framework Part 3).
The issue is that SQUID 2.50 calculates only one (pre-specified by 'index isotope') variant of each calibration constant, so there is never any ambiguity when performing calculations on those values, and using the outputs of those calculations in other expressions downstream. Obviously, such ambiguity does exist in our SQUID 3.0 implementation, so we need to resolve it (at least in a preliminary way) via a similar user-specified control (which would live on the same screen as the controls for SBM-normalisation and SpotAvg vs LinReg) called 'Preferred index isotope', with options 204Pb, 207Pb, 208Pb (strictly, the availability of 208Pb on this list should hinge on whether the selected Task is Perm1-type or not; something to look at later). Down the track, it would be good if there was scope for the user to switch 'preferred index isotope' after the StandardData and SampleData sheets are produced, but we can talk about that later too.
This is only a minor change - I still want you to "calculate everything" as you have been doing. The aim of the control is solely to specify which specific permutation should be propagated to downstream calculations, in a bid to minimise the number of pointless permutations we generate downstream.
The (paraphrased) code continues:
If Only1 = TRUE
t1 = "1" --note that this is a string (in VBA), not a number
ElseIf Th2U8col > 0 --i.e. a destination exists for 232Th/238U values
If (Std = TRUE) AND (piStdRad86col > 0) --the second part means that if
--'preferred index isotope' = 204Pb, then there exists a column named
--["4-corr208Pb*/206Pb*"], or alternatively, if 'preferred index isotope'
--= 207Pb, then there exists a column named ["7-corr208Pb*/206Pb*"]. In
--SQUID 3.0, I think piStdRad86col > 0 can simply be assumed TRUE.
--Recall that Lambda238Ma and Lambda232Ma are the decay constants of 238U
--and 232Th respectively, expressed in units of "Ma^-1" i.e. 1.55125e-04
--and 4.9475e-05 respectively.
If 'preferred index isotope' = 204Pb
Exp238 = "( EXP ( Lambda238Ma * ["4-corr206Pb/238UAge(Ma)"] ) - 1 )"
Exp232 = "( EXP ( Lambda232Ma * ["4-corr208Pb/232ThAge(Ma)"] ) - 1 )"
t1 = "=["4-corr208Pb*/206Pb*"] * Exp238 / Exp232"
Else --'preferred index isotope' = 207Pb
Exp238 = "( EXP ( Lambda238Ma * ["7-corr206Pb/238UAge(Ma)"] ) - 1 )"
Exp232 = "( EXP ( Lambda232Ma * ["7-corr208Pb/232ThAge(Ma)"] ) - 1 )"
t1 = "=["7-corr208Pb*/206Pb*"] * Exp238 / Exp232"
End If
ElseIf Std = FALSE --sample spots only; use uncorrected values
t1 = ["208/206"] * ["Total206Pb/238U"] / ["Total208Pb/232Th"]
End If --(Std = TRUE) AND (piStdRad86col > 0)
Now evaluate the associated uncertainty by quadratic addition. The unfinished ElseIf (Th2U8col > 0) continues:
If Th2U8ecol > 0 --i.e. a destination exists for 232Th/238U uncertainties
If Std = TRUE
If 'preferred index isotope' = 204Pb
t2 = '=SQRT( ["4-corr208Pb*/206Pb* %err"]^2 +
["4-corr206Pb/238Ucalibr.const %err"]^2 +
["4-corr208Pb/232Thcalibr.const %err"]^2 )
Else --'preferred index isotope' = 207Pb
t2 = '=SQRT( ["7-corr208Pb*/206Pb* %err"]^2 +
["7-corr206Pb/238Ucalibr.const %err"]^2 +
["7-corr208Pb/232Thcalibr.const %err"]^2 )
End If
Else --i.e. sample spots; Ludwig wrote
--"Elseif 'preferred index isotope' is not 208Pb", but this seems redundant:
--entry into ThUfromA1A2 precludes 208Pb as preferred index isotope
t2 = '=SQRT( ["208/206 %err"]^2 + ["Total206Pb/238U" %err]^2 +
["Total208Pb/232Th" %err]^2 )
End If --(Std = TRUE)
End If --(Th2U8ecol > 0)
End If --(ElseIf Th2U8col > 0)
Finally, place the analysis-by-analysis formulae according to the expression-strings developed above. First, the value:
PlaceFormulae t1, SpotRow, Th2U8col
And secondly, the uncertainty (if the expression has been developed):
If t2 <> ""
PlaceFormulae t2, SpotRow, Th2U8ecol
End If
End Sub