This page was generated from
docs/Examples/Liquid_Ol_Liq_Themometry/Liquid_only_Thermometry.ipynb.
Interactive online version:
.
Liquid-only thermometry
This script shows how to use the various options for liquid-only thermometry
This tutorial is also a great starting point when using other thermobarometers, as the choice of inputs is the same as for other phases
You can download the Excel spreadsheet you need here: https://github.com/PennyWieser/Thermobar/blob/main/docs/Examples/Liquid_Ol_Liq_Themometry/Liquid_only_Thermometry.xlsx
You need to install Thermobar once on your machine, if you haven’t done this yet, uncomment the line below (remove the #)
[1]:
#!pip install Thermobar
[2]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import Thermobar as pt
Load in data
The function “import_excel” allows you to specify the name of the excel file and the sheet name.
Data should have the headings SiO2_Liq, TiO2_Liq etc. etc. The order of headings doesn’t matter.
As in the 5 min intro example, you could also have no _Liq in the headings, and load them in the functoin using suffix=”_Liq”
You can also have any other columns, e.g., estimate of pressure from any other proxy (melt inclusions, geophysics), and anything else you might want to plot (e.g., latitude, longitude)
[3]:
out=pt.import_excel('Liquid_only_Thermometry.xlsx', sheet_name="Liquid_only")
my_input=out['my_input']
myLiquids1=out['Liqs']
Inspect inputs to check you didn’t have any funny column headings.
It is always useful to inspect the outputs from import_excel, sometimes your column headings may have funny characters due to use of spaces, subscripts etc. in journal pdf tables. Check that all the columns you entered have numbers. If, say your SiO2_Liq heading had funny characters, this column will be filled with zeros when you inspect it.
[4]:
myLiquids1.head()
[4]:
| SiO2_Liq | TiO2_Liq | Al2O3_Liq | FeOt_Liq | MnO_Liq | MgO_Liq | CaO_Liq | Na2O_Liq | K2O_Liq | Cr2O3_Liq | P2O5_Liq | H2O_Liq | Fe3Fet_Liq | NiO_Liq | CoO_Liq | CO2_Liq | Sample_ID_Liq | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 57.023602 | 0.623106 | 16.332899 | 4.36174 | 0.103851 | 4.19180 | 6.94858 | 3.59702 | 0.896895 | 0.000000 | 0.226584 | 5.59 | 0.2 | 0.0 | 0.0 | 0.0 | 0 |
| 1 | 57.658600 | 0.654150 | 17.194799 | 3.90621 | 0.084105 | 2.86892 | 5.91538 | 3.85948 | 1.018600 | 0.000000 | 0.214935 | 6.55 | 0.2 | 0.0 | 0.0 | 0.0 | 1 |
| 2 | 60.731201 | 0.862054 | 17.144199 | 4.07781 | 0.077488 | 2.50867 | 5.22075 | 4.45556 | 1.414160 | 0.000000 | 0.319638 | 3.14 | 0.2 | 0.0 | 0.0 | 0.0 | 2 |
| 3 | 61.532799 | 0.440860 | 16.508801 | 3.32990 | 0.037520 | 1.64150 | 4.34294 | 4.40860 | 1.407000 | 0.000000 | 0.215740 | 6.20 | 0.2 | 0.0 | 0.0 | 0.0 | 3 |
| 4 | 52.969101 | 0.803412 | 17.563000 | 5.93217 | 0.149472 | 3.78351 | 7.65110 | 3.80219 | 0.551178 | 0.037368 | 0.196182 | 6.58 | 0.2 | 0.0 | 0.0 | 0.0 | 4 |
If at any point you need help with what functions inputs to use, you can type help(pt.function_name)
[5]:
# Options for different thermometers laid out in the "equationT" =
help(pt.calculate_liq_only_temp)
Help on function calculate_liq_only_temp in module Thermobar.liquid_thermometers:
calculate_liq_only_temp(*, liq_comps, equationT, P=None, H2O_Liq=None, print=False)
Liquid-only thermometery. Returns a temperature in Kelvin.
Parameters
-------
liq_comps: pandas.DataFrame
liquid compositions with column headings SiO2_Liq, MgO_Liq etc.
equationT: str
If has _sat at the end, represents the saturation surface of that mineral.
Equations from Putirka et al. (2016).
| T_Put2016_eq3_amp_sat (saturation surface of amphibole)
Equations from Putirka (2008) and older studies:
| T_Put2008_eq13
| T_Put2008_eq14
| T_Put2008_eq15
| T_Put2008_eq16
| T_Put2008_eq34_cpx_sat
| T_Put2008_eq28b_opx_sat
| T_Put1999_cpx_sat
* Following 3 thermometers are adaptations of olivine-liquid thermometers with DMg calculated using Beattie 1993,
This means you can use them without knowing an olivine composition. ocan be applied when you haven't measured an olivine composiiton.
| T_Put2008_eq19_BeattDMg
| T_Put2008_eq21_BeattDMg
| T_Put2008_eq22_BeattDMg
Equations from Sugawara (2000):
| T_Sug2000_eq1
| T_Sug2000_eq3_ol
| T_Sug2000_eq3_opx
| T_Sug2000_eq3_cpx
| T_Sug2000_eq3_pig
| T_Sug2000_eq6a
| T_Sug2000_eq6b
Equations from Helz and Thornber (1987):
| T_Helz1987_MgO
| T_Helz1987_CaO
Equation from Molina et al. (2015)
| T_Molina2015_amp_sat
Equation from Montrieth 1995
| T_Montierth1995_MgO
Equation from Beattie (1993)
| T_Beatt1993_opx
P: float, int, pandas.Series, str ("Solve")
Pressure in kbar
Only needed for P-sensitive thermometers.
If enter P="Solve", returns a partial function
Else, enter an integer, float, or panda series
H2O_Liq: optional.
If None, uses H2O_Liq column from input.
If int, float, pandas.Series, uses this instead of H2O_Liq Column
Returns
-------
pandas series
Temperature in K
Example 1 - P and H2O-independent thermometers
The thermometer of Helz and Thornber 1987 is not dependent on pressure or H2O.
Thus, the function only requires 2 inputs, the name of pandas dataframe with liquid compositions (here, myLiquids1) and the name of the equationT you have choosen
All functions output in Kelvin and Kbar, so to get celcius, we subtract 273.15
[6]:
T_Helz1987=pt.calculate_liq_only_temp(liq_comps=myLiquids1, equationT="T_Helz1987_MgO")-273.15
# In jupyter lab/notebooks, outputs don't automatically display,
# but it prints the last line in each cell
T_Helz1987
[6]:
0 1098.255182
1 1071.665294
2 1064.424269
3 1046.994150
4 1090.048550
5 1079.746499
6 1067.772325
7 1099.844089
8 1077.533086
9 1121.514497
Name: MgO_Liq, dtype: float64
Example 2 - Pressure-dependent thermometer
If you select an equation which is Pressure-dependent, and don’t specify a pressure, the code returns an error
You can either select P=“Solve” which returns a partial function. This means you can evaluate it at any pressure you want easily
Or you can specify a fixed pressure, e.g., P=5 (which runs all calculations at 5 kbar)
Or you can specify pressure based on a column in your original spreadsheet containing pressure
Example 2a - Don’t specify a pressure
returns an error, telling you to input a pressure
[7]:
Teq15_2kbar_6wt=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15")-273.15 # Convert to Celcius
Teq15_2kbar_6wt
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
Input In [7], in <cell line: 1>()
----> 1 Teq15_2kbar_6wt=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
2 equationT="T_Put2008_eq15")-273.15 # Convert to Celcius
3 Teq15_2kbar_6wt
File g:\my drive\postdoc\pymme\mybarometers\thermobar_outer\src\Thermobar\liquid_thermometers.py:526, in calculate_liq_only_temp(liq_comps, equationT, P, H2O_Liq, print)
524 if sig.parameters['P'].default is not None:
525 if P is None:
--> 526 raise ValueError(f'{equationT} requires you to enter P, or specify P="Solve"')
528 #else:
529 # print('gothere2')
530 # if P is not None:
(...)
533 # P=None
534 # print('got here')
538 kwargs = {name: anhyd_cat_frac[name] for name, p in sig.parameters.items() if p.kind == inspect.Parameter.KEYWORD_ONLY}
ValueError: T_Put2008_eq15 requires you to enter P, or specify P="Solve"
Example 2b - specify P=“Solve” to return a partial function
[8]:
Teq15_partial=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15", P="Solve")# Output is in Kelvin
[9]:
# You can then evaluate this partial at any pressure you want (here, 10 kbar).
# This can save computation time
Teq15_10kbar=Teq15_partial(10)-273.15
Teq15_10kbar
[9]:
0 1082.230315
1 1031.095340
2 1065.719538
3 991.644728
4 1049.453769
5 1033.581813
6 1021.013571
7 946.501846
8 919.539200
9 1108.248454
dtype: float64
Example 2c - evaluate at constant pressure
[10]:
# Here P=5 kbar
Teq15_5kbar=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15", P=5)-273.15 # Output is in Kelvin, convert to C
Teq15_5kbar
[10]:
0 1062.650315
1 1011.515340
2 1046.139538
3 972.064728
4 1029.873769
5 1014.001813
6 1001.433571
7 926.921846
8 899.959200
9 1088.668454
dtype: float64
Example 2d - Evaluate at values given by P_kbar column in inputted excel spreadsheet
your column can be called anything, you just need state P=my_input[‘column name’]. If you stored pressure in GPa say, you could do P=10*my_input[‘P_GPa’]
[11]:
Teq15_input=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15", P=my_input['P_kbar'])-273.15
# Output is in Kelvin, convert to c.
Teq15_input
[11]:
0 1050.902315
1 1003.683340
2 1046.139538
3 972.064728
4 1029.873769
5 1014.001813
6 1001.433571
7 926.921846
8 899.959200
9 1088.668454
dtype: float64
Example 3 - Specifying water content
The dataframe for liquids in this example has a H2O column it read from the excel spreadsheet. By default, H2O-dependent thermometers will use the values entered in this column
If you don’t have this column, it’ll be filled with zeros, so calculations will be done at H2O=0
However, you can also override this by specifying H2O_Liq=value in the function.
[12]:
Teq15_2H2O=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15", P=5, H2O_Liq=2)-273.15 # Output is in Kelvin
Teq15_2H2O
[12]:
0 1108.710017
1 1069.891842
2 1060.765739
3 1025.950726
4 1088.635168
5 1075.714113
6 1058.655371
7 1083.961043
8 1058.281402
9 1129.211252
dtype: float64
We can also investigate how sensitive our calculations are to water, say we think maybe water ranges from 2-4 wt%, we can do the same calculation at 4 wt%
[13]:
Teq15_4H2O=pt.calculate_liq_only_temp(liq_comps=myLiquids1,
equationT="T_Put2008_eq15", P=5, H2O_Liq=4)-273.15 # Output is in Kelvin
Teq15_4H2O
[13]:
0 1083.050017
1 1044.231842
2 1035.105739
3 1000.290726
4 1062.975168
5 1050.054113
6 1032.995371
7 1058.301043
8 1032.621402
9 1103.551252
dtype: float64
We can now calculate the difference between these by subtracting one from the other. Here, as H2O_Liq is just a constant on the function, the difference is constant for all liquids
[14]:
Teq15_4H2O-Teq15_2H2O
[14]:
0 -25.66
1 -25.66
2 -25.66
3 -25.66
4 -25.66
5 -25.66
6 -25.66
7 -25.66
8 -25.66
9 -25.66
dtype: float64
Exporting data
There are a number of ways you can export data to excel
The first example appends these new calculations onto a copy of the “my_input” dataframe. This means you can export a spreadsheet that looks exactly like your input spreadsheet, but with the values for the new calculations
[15]:
# First make a copy of the dataframe so you don't overwrite the original
my_input_c=my_input.copy()
# To make a new column, you specify my_input_c['new column name']= value
my_input_c['T_Helz1987']=T_Helz1987
# You can add as many new columns as you wish
my_input_c['Teq15_4H2O']=Teq15_4H2O
# You can give your columns in the new dataframe any name you wish...
my_input_c['Temp eq 15 4 wt% H2O']=Teq15_4H2O
# Now, you ue the "to_excel" panda function to save this new appended dataframe to excel. The name of the dataframe you've been appended to goes on the let, and the name of the
# file you want to make goes bewteen the ' ' signs
my_input_c.to_excel('Thermometry_out1.xlsx')
[ ]: