{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Reflectance FTIR Spectra\n", "\n", "- This Jupyter notebook applies the Newtonian inversion to the molar absorptivity data. \n", "\n", "- The Jupyter notebook and data can be accessed here: https://github.com/SarahShi/PyIRoGlass/blob/main/docs/examples/inversion/. \n", "\n", "- You need to have the PyIRoGlass PyPi package on your machine once. If you have not done this, please uncomment (remove the #) symbol and run the cell below. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "#!pip install PyIRoGlass" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Load Python Packages" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "tags": [] }, "outputs": [], "source": [ "# Import packages\n", "\n", "import numpy as np\n", "import pandas as pd\n", "\n", "import PyIRoGlass as pig\n", "\n", "from matplotlib import pyplot as plt\n", "\n", "%matplotlib inline\n", "%config InlineBackend.figure_format = 'retina'\n", "\n", "pig.__version__" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Load Data \n", "\n", "Load all compiled molar absorptivity data and assign uncertainties to the peaks. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_5200 = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='NIRRegress')\n", "tau_5200 = df_5200['Tau']\n", "sigma_tau_5200 = tau_5200 * 0.025\n", "epsilon_5200 = df_5200['Epsilon_5200']\n", "sigma_epsilon_5200 = epsilon_5200 * 0.10\n", "\n", "df_4500 = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='NIRRegress')\n", "tau_4500 = df_4500['Tau']\n", "sigma_tau_4500 = tau_4500 * 0.025\n", "epsilon_4500 = df_4500['Epsilon_4500']\n", "sigma_epsilon_4500 = epsilon_4500 * 0.20\n", "\n", "df_3550 = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='3550Regress')\n", "tau_3550 = df_3550['Tau']\n", "sigma_tau_3550 = tau_3550 * 0.025\n", "epsilon_3550 = df_3550['Epsilon_3550']\n", "sigma_epsilon_3550 = epsilon_3550 * 0.10\n", "\n", "df_1635 = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='1635Regress')\n", "tau_1635 = df_1635['Tau']\n", "sigma_tau_1635 = tau_1635 * 0.025\n", "epsilon_1635 = df_1635['Epsilon_1635']\n", "sigma_epsilon_1635 = epsilon_1635 * 0.05\n", "\n", "df_carbonate = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='CarbonateRegress')\n", "eta = df_carbonate['Eta']\n", "sigma_eta = eta * 0.025\n", "epsilon_carbonate = df_carbonate['Epsilon_Carbonate']\n", "sigma_epsilon_carbonate = epsilon_carbonate * 0.10" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# ε5200 Dataset\n", "\n", "Display `df_5200`, the DataFrame of calibration data for the $\\mathrm{H_2O_{m,5200}}$ peak." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_5200" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# ε4500 Dataset\n", "\n", "Display `df_4500`, the DataFrame of calibration data for the $\\mathrm{OH^-_{4500}}$ peak." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_4500" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# ε3550 Dataset\n", "\n", "Display `df_3550`, the DataFrame of calibration data for the $\\mathrm{H_2O_{t, 3550}}$ peak." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_3550" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# ε1635 Dataset\n", "\n", "Display `df_1635`, the DataFrame of calibration data for the $\\mathrm{H_2O_{m,1635}}$ peak." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_1635" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# εCarbonate Dataset\n", "\n", "Display `df_carbonate`, the DataFrame of calibration data for the $\\mathrm{CO_{3}^{2-}}$ peak." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df_carbonate" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "We're ready to use the `pig.inversion`, `pig.least_squares`, `pig.inversion_fit_errors`, and `pig.inversion_fit_errors_plotting` functions now. We input the arguments: \n", "\n", "- `tau` or `eta`: Compositional parameter\n", "- `epsilon`: Molar absorptivity \n", "- `sigma_tau` or `sigma_eta`: Uncertainty on compositional parameter \n", "- `sigma_epsilon`: Uncertainty on absorption coefficient \n", "\n", "and output: \n", "\n", "- `mest`: Best-fit inversion parameters \n", "- `covm`: Covariance matrix for inversion parameters\n", "- `covepsilon`: Covariance on fit molar absorptivity\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "mest_5200, covm_est_5200, covepsilon_5200 = pig.inversion(tau_5200, epsilon_5200, sigma_tau_5200, sigma_epsilon_5200, intercept_zero=False)\n", "mls_5200, covls_5200 = pig.least_squares(tau_5200, epsilon_5200, sigma_epsilon_5200)\n", "E_calib_5200, see_5200, r2_5200, rmse_5200, rrmse_5200, ccc_5200 = pig.inversion_fit_errors(tau_5200, epsilon_5200, mest_5200, covepsilon_5200)\n", "\n", "mest_4500, covm_est_4500, covepsilon_4500 = pig.inversion(tau_4500, epsilon_4500, sigma_tau_4500, sigma_epsilon_4500, intercept_zero=False)\n", "mls_4500, covls_4500 = pig.least_squares(tau_4500, sigma_tau_4500, sigma_epsilon_4500)\n", "E_calib_4500, see_4500, r2_4500, rmse_4500, rrmse_4500, ccc_4500 = pig.inversion_fit_errors(tau_4500, epsilon_4500, mest_4500, covepsilon_4500)\n", "\n", "mest_3550, covm_est_3550, covepsilon_3550 = pig.inversion(tau_3550, epsilon_3550, sigma_tau_3550, sigma_epsilon_3550, intercept_zero=False)\n", "mls_3550, covls_3550 = pig.least_squares(tau_3550, epsilon_3550, sigma_epsilon_3550)\n", "E_calib_3550, see_3550, r2_3550, rmse_3550, rrmse_3550, ccc_3550 = pig.inversion_fit_errors(tau_3550, epsilon_3550, mest_3550, covepsilon_3550)\n", "\n", "mest_1635, covm_est_1635, covepsilon_1635 = pig.inversion(tau_1635, epsilon_1635, sigma_tau_1635, sigma_epsilon_1635, intercept_zero=False)\n", "mls_1635, covls_1635 = pig.least_squares(tau_1635, epsilon_1635, sigma_epsilon_1635)\n", "E_calib_1635, see_1635, r2_1635, rmse_1635, rrmse_1635, ccc_1635 = pig.inversion_fit_errors(tau_1635, epsilon_1635, mest_1635, covepsilon_1635)\n", "\n", "mest_carbonate, covm_est_carbonate, covepsilon_carbonate = pig.inversion(eta, epsilon_carbonate, sigma_eta, sigma_epsilon_carbonate, intercept_zero=False)\n", "mls_carbonate, covls_carbonate = pig.least_squares(eta, epsilon_carbonate, sigma_epsilon_carbonate)\n", "E_calib_carbonate, see_carbonate, r2_carbonate, rmse_carbonate, rrmse_carbonate, ccc_carbonate = pig.inversion_fit_errors(eta, epsilon_carbonate, mls_carbonate, covepsilon_carbonate)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Assessing uncertainty of prediction:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "tau_arr_5200, epsilon_5200_arr, conf_lower_5200, conf_upper_5200, pred_lower_5200, pred_upper_5200 = pig.inversion_fit_errors_plotting(tau_5200, epsilon_5200, mest_5200)\n", "tau_arr_4500, epsilon_4500_arr, conf_lower_4500, conf_upper_4500, pred_lower_4500, pred_upper_4500 = pig.inversion_fit_errors_plotting(tau_4500, epsilon_4500, mest_4500)\n", "tau_arr_3550, epsilon_3550_arr, conf_lower_3550, conf_upper_3550, pred_lower_3550, pred_upper_3550 = pig.inversion_fit_errors_plotting(tau_3550, epsilon_3550, mest_3550)\n", "tau_arr_1635, epsilon_1635_arr, conf_lower_1635, conf_upper_1635, pred_lower_1635, pred_upper_1635 = pig.inversion_fit_errors_plotting(tau_1635, epsilon_1635, mest_1635)\n", "eta_arr, epsilon_carbonate_arr, conf_lower_carbonate, conf_upper_carbonate, pred_lower_carbonate, pred_upper_carbonate = pig.inversion_fit_errors_plotting(eta, epsilon_carbonate, mest_carbonate)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "We can now make Figure 5 from the paper. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "sz = 150\n", "\n", "fig, ax = plt.subplots(3, 2, figsize=(14, 19))\n", "ax = ax.flatten()\n", "\n", "epsilon_5200_mandeville = -2.463 + 4.899*tau_arr_5200\n", "fuego_idx = np.where((tau_arr_5200 > 0.653) & (tau_arr_5200 < 0.715))\n", "legend_5200 = r'$\\mathregular{ƐH_2O_{m, 5200}}$ = ' + f'{round(mest_5200[0],3)}(±{round(np.sqrt(np.diag(covm_est_5200))[0],3)}) + {round(mest_5200[1],3)}(±{round(np.sqrt(np.diag(covm_est_5200))[1],3)})'+ '·' + r'$\\tau$'+ f', N={len(tau_5200)}'\n", "ax[0].plot(tau_arr_5200, epsilon_5200_arr, 'k', lw=2, zorder=0, label=legend_5200)\n", "mand, = ax[0].plot(tau_arr_5200, epsilon_5200_mandeville, 'k-.', lw=2, zorder=0, label='Mandeville et al., 2002')\n", "mand.set_dashes([1.5, 1, 3, 1])\n", "ax[0].fill_between(tau_arr_5200, conf_lower_5200, conf_upper_5200, color='k', alpha=0.20, edgecolor=None,\n", " zorder=-5, label='68% Confidence Interval')\n", "ax[0].plot(tau_arr_5200, pred_upper_5200, 'k--', lw=0.5, zorder=0, dashes=(16, 10))\n", "ax[0].plot(tau_arr_5200, pred_lower_5200, 'k--', lw=0.5, zorder=0, dashes=(16, 10), label='68% Prediction Interval')\n", "ax[0].fill_between(tau_arr_5200[fuego_idx], conf_lower_5200[fuego_idx], conf_upper_5200[fuego_idx], color='r', alpha=0.30,\n", " edgecolor=None, zorder=-5, label='Fuego Interval')\n", "ax[0].errorbar(tau_5200, epsilon_5200, yerr = sigma_epsilon_5200, xerr = sigma_tau_5200, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[0].scatter(tau_5200, epsilon_5200, s=sz, c='#0C7BDC', edgecolors='black', linewidth=0.5, zorder=15)\n", "ax[0].set_xlim([ 0.5, 1.0])\n", "ax[0].set_ylim([-0.5, 3.5])\n", "xlabel_5200 = r'$\\tau$=(Si+Al)/Total Cations'\n", "ax[0].set_xlabel(xlabel_5200)\n", "ax[0].set_ylabel(r'$\\mathregular{ƐH_2O_{m, 5200}}$')\n", "ax[0].legend(loc='upper left', labelspacing=0.4, handletextpad=0.5, handlelength=1.50, prop={'size': 12}, frameon=False)\n", "ax[0].tick_params(axis=\"x\", direction='in', length=5, pad=6.5)\n", "ax[0].tick_params(axis=\"y\", direction='in', length=5, pad=6.5)\n", "\n", "\n", "epsilon_4500_mandeville = -2.026+4.054*tau_arr_4500\n", "fuego_idx = np.where((tau_arr_4500 > 0.653) & (tau_arr_4500 < 0.715))\n", "legend_4500 = r'$\\mathregular{ƐOH^{-}_{4500}}$ = ' + f'{round(mest_4500[0],3)}(±{round(np.sqrt(np.diag(covm_est_4500))[0],3)}) + {round(mest_4500[1],3)}(±{round(np.sqrt(np.diag(covm_est_4500))[1],3)})'+ '·' + r'$\\tau$'+ f', N={len(tau_4500)}'\n", "ax[1].plot(tau_arr_4500, epsilon_4500_arr, 'k', lw=2, zorder=0, label=legend_4500)\n", "mand, = ax[1].plot(tau_arr_4500, epsilon_4500_mandeville, 'k-.', lw=2, zorder=0, label='Mandeville et al., 2002')\n", "mand.set_dashes([1.5, 1, 3, 1])\n", "ax[1].fill_between(tau_arr_4500, conf_lower_4500, conf_upper_4500, color='k', alpha=0.20, edgecolor=None,\n", " zorder=-5, label='68% Confidence Interval')\n", "ax[1].plot(tau_arr_4500, pred_upper_4500, 'k--', lw=0.5, zorder=0, dashes=(16, 10))\n", "ax[1].plot(tau_arr_4500, pred_lower_4500, 'k--', lw=0.5, zorder=0, dashes=(16, 10), label='68% Prediction Interval')\n", "ax[1].fill_between(tau_arr_4500[fuego_idx], conf_lower_4500[fuego_idx], conf_upper_4500[fuego_idx], color='r', alpha=0.30,\n", " edgecolor=None, zorder=-5, label='Fuego Interval')\n", "ax[1].errorbar(tau_4500, epsilon_4500, yerr = sigma_epsilon_4500, xerr = sigma_tau_4500, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[1].scatter(tau_4500, epsilon_4500, s=sz, c='#0C7BDC', edgecolors='black', linewidth=0.5, zorder=15)\n", "ax[1].set_xlim([ 0.5, 1.0])\n", "ax[1].set_ylim([-0.5, 3.5])\n", "xlabel_4500 = r'$\\tau$=(Si+Al)/Total Cations'\n", "ax[1].set_xlabel(xlabel_4500) \n", "ax[1].set_ylabel(r'$\\mathregular{ƐOH^{-}_{4500}}$')\n", "ax[1].legend(loc='upper left', labelspacing=0.4, handletextpad=0.5, handlelength=1.50, prop={'size': 12}, frameon=False)\n", "ax[1].tick_params(axis=\"x\", direction='in', length=5, pad=6.5)\n", "ax[1].tick_params(axis=\"y\", direction='in', length=5, pad=6.5)\n", "\n", "\n", "fuego_idx = np.where((tau_arr_3550 > 0.653) & (tau_arr_3550 < 0.715))\n", "legend_3550 = r'$\\mathregular{ƐH_2O_{t, 3550}}$ = ' + f'{round(mest_3550[0],3)}(±{round(np.sqrt(np.diag(covm_est_3550))[0],3)}) + {round(mest_3550[1],3)}(±{round(np.sqrt(np.diag(covm_est_3550))[1],3)})'+ '·' + r'$\\tau$'+ f', N={len(tau_3550)}'\n", "ax[2].plot(tau_arr_3550, epsilon_3550_arr, 'k', lw=2, zorder=0, label=legend_3550)\n", "mand.set_dashes([1.5, 1, 3, 1])\n", "ax[2].fill_between(tau_arr_3550, conf_lower_3550, conf_upper_3550, color='k', alpha=0.20, edgecolor=None,\n", " zorder=-5, label='68% Confidence Interval')\n", "ax[2].plot(tau_arr_3550, pred_upper_3550, 'k--', lw=0.5, zorder=0, dashes=(16, 10))\n", "ax[2].plot(tau_arr_3550, pred_lower_3550, 'k--', lw=0.5, zorder=0, dashes=(16, 10), label='68% Prediction Interval')\n", "ax[2].fill_between(tau_arr_3550[fuego_idx], conf_lower_3550[fuego_idx], conf_upper_3550[fuego_idx], color='r', alpha=0.30,\n", " edgecolor=None, zorder=-5, label='Fuego Interval')\n", "ax[2].errorbar(tau_3550, epsilon_3550, yerr = sigma_epsilon_3550, xerr = sigma_tau_3550, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[2].scatter(tau_3550, epsilon_3550, s=sz, c='#0C7BDC', edgecolors='black', linewidth=0.5, zorder=15)\n", "ax[2].set_xlim([0.4, 1.0])\n", "ax[2].set_ylim([20, 110])\n", "xlabel_3550 = r'$\\tau$=(Si+Al)/Total Cations'\n", "ax[2].set_xlabel(xlabel_3550) \n", "ax[2].set_ylabel(r'$\\mathregular{ƐH_2O_{t, 3550}}$')\n", "ax[2].legend(loc='upper left', labelspacing=0.4, handletextpad=0.5, handlelength=1.50, prop={'size': 12}, frameon=False)\n", "ax[2].tick_params(axis=\"x\", direction='in', length=5, pad=6.5)\n", "ax[2].tick_params(axis=\"y\", direction='in', length=5, pad=6.5)\n", "ax[2].scatter(tau_3550[epsilon_3550==63.027], epsilon_3550[epsilon_3550==63.027], s=sz, c='#0C7BDC', edgecolors='black', linewidth=2, zorder=15)\n", "\n", "epsilon_1635_mandeville = -57.813+131.94*tau_arr_1635\n", "fuego_idx = np.where((tau_arr_1635 > 0.653) & (tau_arr_1635 < 0.715))\n", "legend_1635 = r'$\\mathregular{ƐH_2O_{m, 1635}}$ = ' + f'{round(mest_1635[0],3)}(±{round(np.sqrt(np.diag(covm_est_1635))[0],3)}) + {round(mest_1635[1],3)}(±{round(np.sqrt(np.diag(covm_est_1635))[1],3)})'+ '·' + r'$\\tau$'+ f', N={len(tau_1635)}'\n", "ax[3].plot(tau_arr_1635, epsilon_1635_arr, 'k', lw=2, zorder=0, label=legend_1635)\n", "mand, = ax[3].plot(tau_arr_1635, epsilon_1635_mandeville, 'k-.', lw=2, zorder=0, label='Mandeville et al., 2002')\n", "mand.set_dashes([1.5, 1, 3, 1])\n", "ax[3].fill_between(tau_arr_1635, conf_lower_1635, conf_upper_1635, color='k', alpha=0.20, edgecolor=None,\n", " zorder=-5, label='68% Confidence Interval')\n", "ax[3].plot(tau_arr_1635, pred_upper_1635, 'k--', lw=0.5, zorder=0, dashes=(16, 10))\n", "ax[3].plot(tau_arr_1635, pred_lower_1635, 'k--', lw=0.5, zorder=0, dashes=(16, 10), label='68% Prediction Interval')\n", "ax[3].fill_between(tau_arr_1635[fuego_idx], conf_lower_1635[fuego_idx], conf_upper_1635[fuego_idx], color='r', alpha=0.30,\n", " edgecolor=None, zorder=-5, label='Fuego Interval')\n", "ax[3].errorbar(tau_1635, epsilon_1635, yerr = sigma_epsilon_1635, xerr = sigma_tau_1635, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[3].scatter(tau_1635, epsilon_1635, s=sz, c='#0C7BDC', edgecolors='black', linewidth=0.5, zorder=15)\n", "ax[3].set_xlim([0.5, 1.0])\n", "ax[3].set_ylim([0, 90])\n", "xlabel_1635 = r'$\\tau$=(Si+Al)/Total Cations'\n", "ax[3].set_xlabel(xlabel_1635) \n", "ax[3].set_ylabel(r'$\\mathregular{ƐH_2O_{m, 1635}}$')\n", "ax[3].legend(loc='upper left', labelspacing=0.4, handletextpad=0.5, handlelength=1.50, prop={'size': 12}, frameon=False)\n", "ax[3].tick_params(axis=\"x\", direction='in', length=5, pad=6.5)\n", "ax[3].tick_params(axis=\"y\", direction='in', length=5, pad=6.5)\n", "\n", "epsilon_carbonate_dixonpan = 451-342*eta_arr\n", "fuego_idx = np.where((eta_arr > 0.389) & (eta_arr < 0.554))\n", "df_carbonate = pd.read_excel('./EpsilonRegression.xlsx', sheet_name='CarbonateRegress')\n", "low_df = df_carbonate[df_carbonate.Epsilon_Location == 'Low']\n", "high_df = df_carbonate[df_carbonate.Epsilon_Location == 'High']\n", "shi_low = low_df[low_df.Compilation=='Shi']\n", "shi_high = low_df[low_df.Compilation=='Shi']\n", "ax[4].errorbar(low_df['Eta'], low_df['Epsilon_Carbonate'], yerr = low_df['Epsilon_Carbonate']*0.1, xerr = low_df['Eta']*0.025, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[4].scatter(low_df['Eta'], low_df['Epsilon_Carbonate'], s=sz, c='#0C7BDC', edgecolors='black', linewidth=0.5, zorder=15, label=r'$\\mathregular{CO_{3, 1430}^{2-}}$, N='+str(len(low_df)))\n", "ax[4].errorbar(high_df['Eta'], high_df['Epsilon_Carbonate'], yerr = high_df['Epsilon_Carbonate']*0.10, xerr = high_df['Eta']*0.025, ls='none', elinewidth=0.5, ecolor='k')\n", "ax[4].scatter(high_df['Eta'], high_df['Epsilon_Carbonate'], s=sz, c='#E42211', marker='s', edgecolors='black', linewidth=0.5, zorder=15, label=r'$\\mathregular{CO_{3, 1515}^{2-}}$, N='+str(len(high_df)))\n", "dixonpan, = ax[4].plot(eta_arr, epsilon_carbonate_dixonpan, 'k-.', lw=1.5, zorder=0, label='Dixon and Pan, 1995')\n", "dixonpan.set_dashes([1.5, 1, 3, 1])\n", "legend_carbonate = r'$\\mathregular{ƐCO_3^{2-}}$= ' + f'{round(mest_carbonate[0],3)}(±{round(np.sqrt(np.diag(covm_est_carbonate))[0],3)})-{round(mest_carbonate[1],3)*-1}(±{round(np.sqrt(np.diag(covm_est_carbonate))[1],3)})' + '·' + r'Na/(Na+Ca)'\n", "ax[4].plot(eta_arr, epsilon_carbonate_arr, 'k', lw=2, zorder=0, label=legend_carbonate)\n", "ax[4].fill_between(eta_arr, conf_lower_carbonate, conf_upper_carbonate, color='k', alpha=0.20, edgecolor=None,\n", " zorder=-5, label='68% Confidence Interval')\n", "ax[4].plot(eta_arr, pred_upper_carbonate, 'k--', lw=0.5, zorder=0, dashes=(16, 10))\n", "ax[4].plot(eta_arr, pred_lower_carbonate, 'k--', lw=0.5, zorder=0, dashes=(16, 10), label='68% Prediction Interval')\n", "ax[4].fill_between(eta_arr[fuego_idx], conf_lower_carbonate[fuego_idx], conf_upper_carbonate[fuego_idx], color='r', alpha=0.30, edgecolor=None, zorder=-5, label='Fuego Interval')\n", "ax[4].set_xlim([0.1, 0.9])\n", "ax[4].set_ylim([0, 500])\n", "ax[4].set_xlabel(r'$\\eta$ = Na/(Na+Ca)') \n", "ax[4].set_ylabel(r'$\\mathregular{ƐCO_3^{2-}}$')\n", "ax[4].legend(loc='lower left', labelspacing=0.4, handletextpad=0.5, handlelength=1.50, prop={'size': 12}, frameon=False)\n", "ax[4].tick_params(axis=\"x\", direction='in', length=5, pad=6.5)\n", "ax[4].tick_params(axis=\"y\", direction='in', length=5, pad=6.5)\n", "ax[4].scatter(shi_high['Eta'], shi_high['Epsilon_Carbonate'], s=sz, c='#E42211', marker='s', edgecolors='black', linewidth=2, zorder=15)\n", "ax[4].scatter(shi_low['Eta'], shi_low['Epsilon_Carbonate'], s=sz, c='#0C7BDC', edgecolors='black', linewidth=2, zorder=15)\n", "\n", "fig.delaxes(ax[5])\n", "plt.tight_layout()\n", "plt.savefig('AllEpsilonRegress.pdf')" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Apply Inversion Parameters to determine Molar Absorptivities with Uncertainties " ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Load in our glass composition data with the `pig.SampleDataLoader` class. We input: \n", "\n", "- `chemistry_thickness_path`: String pointing to CSV file with glass chemistry and thickness data\n", "\n", "and use the method `load_chemistry_thickness` to return: \n", "\n", "- `chemistry`: DataFrame of chemical data\n", "- `thickness`: DataFrame of thickness data" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "chemistry_thickness_path = 'ChemThick.csv'\n", "loader = pig.SampleDataLoader(chemistry_thickness_path=chemistry_thickness_path)\n", "chemistry, thickness = loader.load_chemistry_thickness()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We'll use the `pig.calculate_epsilon` function, which takes in the parameters: \n", "\n", "- `chemistry`: DataFrame of compositions \n", "- `T`: Room temperature at time of FTIR analysis, given the sensitivity of density to temperature\n", "- `P`: Room pressure at time of FTIR analysis, given the sensitivity of density to pressure\n", "\n", "and returns: \n", "\n", "- `epsilon`: DataFrame of the appropriate molar absorptivity and the corresponding uncertainty, following Equation 8 from the paper. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "T = 25 # C\n", "P = 1 # Bar\n", "\n", "epsilon = pig.calculate_epsilon(chemistry, T, P)\n", "epsilon" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.18" } }, "nbformat": 4, "nbformat_minor": 4 }