oscarbranson / cbsyst

Python module for calculating carbon and boron solution chemistry.
MIT License
29 stars 17 forks source link

Separate T/S/P input & output conditions. #3

Open oscarbranson opened 7 years ago

oscarbranson commented 7 years ago

Currently, input and output conditions are the same. It would be good to be able to separate input/output T, S, P conditions.

This should be relatively simple and inexpensive, and involve a second call to MyAMI_get_Ks and re-calculating the final speciations. Suggest wrapping this in a separate function which takes Ks, and outputs all conditions from a single parameter combination (e.g. DIC and H for Csys).

tompc35 commented 7 years ago

@oscarbranson I could work on implementing this, as a step toward comparing with CO2SYS. Is separate input/output needed for S? I would expect T and P to differ between lab and in situ conditions, but not S.

oscarbranson commented 7 years ago

Probably not in practice - but it's no extra work, so might as well add it for completeness?

Just pushed an updated cbsyst.py to the co2sys_compare branch (https://github.com/oscarbranson/cbsyst/commit/9110fd31994b8be43cd3c438722139a67c2c25ba), which contains my notes on this so far.

I'm also heading towards re-factoring the CBsys function, as it's rather inefficient as-is.

oscarbranson commented 7 years ago

I think I've completed this.

Just need to test. Best way of testing probably against GLODAPv2 data, which have both in-situ ('phtsinsitutp') and 25 C 0 Bar ('pht25p0') variables.

Preliminary test looks good, but there's a slight offset in the residuals: image

Method: calculate GLODAPv2 carbon system at in-situ conditions using 'standardised' pH ('pht25p0'), DIC, Temp and P. Compare output pH to GLODAPv2 in-situ pH ('phtsinsitutp'):

import pandas as pd
import cbsyst as cb  # co2sys_compare branch

gd = pd.read_csv('./GLODAPv2_pH_DIC_ALK_subset.csv')  # obtained by running get_GLODAP_data.py
calcpH = cb.CBsys(pHtot=gd.phts25p0, DIC=gd.tco2,
                  T_in=25., P_in=0., S_in=gd.salinity,
                  TP=gd.phosphate, TSi=gd.silicate,
                  T_out=gd.temperature, P_out=gd.pressure)

Note: Can't find any sensible patterns in the residuals. Could it be caused by the methods (esp. Constants) used when GLODAPv2 converted 'pht25p0' to 'phtsinsitutp'?

Still to do: