Function pH_nbs2total(S, T_C, pH_NBS)
' CONVERT pH (ON NBS SCALE) TO pH on TOTAL SCALE
'
' This code was adapted from the work of others:
' CO2SYS.exe was originally written by Ernie Lewis and Doug Wallace.
' It was converted to Matlab by Richard E. Zeebe & Dieter A. Wolf-Gladrow.
' Parts were extracted and modified by Rachel M. Dunk.
'
' I am indebted to all of them for their initial development of this code.
'
' FIRST: convert pH_NBS to pH_SWS
' aH = 10^(-pH_NBS) = fH*H_SWS
' fH = activity coefficient of H+ ion & includes liquid junction effects
' The fit used in CO2SYS is valid for salinity 20-40 (from Takahashi et al,
' Chapter 3 in GEOSECS Pacific Expedition,v.3, 1982, p. 80).
' For future ref, CO2SYS assumes fH is independent of P
' Convert T to degree K
T = T_C + 273.15
' Convert pH to H+ concentration
aH = 10 ^ -pH_NBS
fH = 1.2948 - .002036 * T
fH = fH + (.0004607 - 1.475E-06 * T) * S^2
H = aH / fH
' On NBS scale, aH is in mol/dm3_H2O, whereas H_SWS is in mol/kg_SW
' to convert units we need the density of SW at given S,T
' Call DENSITY.bas (author=ETP), calculates rho(S,T,P=1atm) in g/cm3
' NOTE: in this function T is in deg C
' Equation of State is from Millero & Poisson (1981) DSR V28: 625-629.
rho_SW = DENSITY(S, T_C)
pH_SWS = -log(H / rho_SW) / log(10)
'SECOND: convert pH_SWS to pH_total
' CONVERT pH (on SEAWATER SCALE) TO pH on TOTAL SCALE
H_SWS = 10^(- pH_SWS)
Hi = 10^(- (pH_SWS + 0.01))
F_T = 0.00007 * S/35
lnK_F = 874/T - 9.68 + 0.111*S^0.5
K_F = exp(lnK_F)
HF = F_T / (1 + K_F/Hi)
H_T = H_SWS - HF
pH_nbs2total = - log(H_T) / log(10)
End Function