Function pH_total2nbs(S, T_C, pH_total) ' CONVERT pH (ON TOTAL SCALE) TO pH on NBS 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 T_C to degree K T = T_C + 273.15 ' SECOND: convert pH_total to pH_sws ' (CONVERT pH on TOTAL SCALE TO pH on SEAWATER SCALE) H_T = 10^(- pH_total) Hi = 10^(- (pH_total - 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_SWS = H_T + HF ' THIRD: convert pH_SWS to pH_NBS ' (CONVERT pH on SEAWATER SCALE TO pH on NBS SCALE) ' ' This is a two-step process: first adjust concentration units, ' then convert conc to activity... ' ' 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) H_L = H_SWS * rho_SW ' Now, we can convert H+ concentration to activity, aH ' aH = 10^(-pH_NBS) = fH*H_L ' 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 H+ concentration to pH_NBS fH = 1.2948 - .002036 * T fH = fH + (.0004607 - 1.475E-06 * T) * S^2 aH = H_L * fH pH_total2nbs = - log(aH) / log(10) End Function