dankelley / oce

R package for oceanographic processing
http://dankelley.github.io/oce/
GNU General Public License v3.0
142 stars 42 forks source link

freshwater data #1037

Closed mlcoleman closed 9 months ago

mlcoleman commented 8 years ago

I'd like to use oce for CTD data from freshwater. Are there any issues I should be aware of, e.g. functions that apply only to seawater? In particular I'm wondering about swRho and about calculations of pressure and salinity since my data files lack these values and they must be calculated. Thanks for developing such a useful tool!

dankelley commented 8 years ago

Hi. I am not exactly sure of the range over which the swRho formulae are meant to apply. If you do

library(oce)
?swRho

and look near the end, there are a lot of citations, and you might want to look at them. If you have any citations to freshwater formulae, please let me know, so I can do some tests and also put information on this topic into the documentation provided by ?swRho.

As for depth, I am quite sure that swDepth is going to give incorrect answers for fresh water. Actually, this function is not really very sophisticated, e.g. it does not take into account the depth variation of density, which might amount to 1 part in 1000 parts, I suppose. Then again, for freshwater applications it seems unlikely that depth will exceed 200m or so, and the rough estimate of error I've just stated would amount to a depth error of only 0.2m, and I doubt that the calibration of pressure for atmospheric pressure is much better than that. (I assume you know how a CTD is set up.)

As I'm sure you know, Seabird software uses constant density in the depth calculation for salt water. I bet they have a similar formula (with a different constant) for fresh water. So you could use that.

I cannot make a blanket statement about using oce for freshwater. I think whether it will work depends on what you want. For example, if you plot a temperature profile versus pressure (which I recommend, since that is what is measured ... and since this is the convention in oceanography) then it won't matter what the depth really is. Similarly, if you are calculating density to get a buoyancy frequency, I imagine the swRho formula for density at S=0 will be reasonably good, since calculating N^2 is tricky business, from the start, since it is based on derivatives of a noisy signal.

I think you should just go ahead and start trying oce, to see what works. And then post issues for anything that seems broken. Your first step, though, will likely be to send me some citations on documents for formulae for freshwater density, according to some international scientific body.

PS. of course, if you have some dissolved material in your lake, or sediments, etc., then you'll need to get more sophisticated. You might want to look up some papers by Rich Pawlowicz at the University of British Columbia, since he has worked in both lakes and oceans (and since he is software genius, in addition to being a science genius)

richardsc commented 8 years ago

Lots of good advice from @dankelley, particularly about looking at the work done by Rich Pawlowicz.

A few points I would add to the discussion. First, for salinity, my understanding is that the Practical Salinity Scale (PSS-78) is only considered to be valid for salinities in the range of 2 < SP < 42. There are some formulae in the literature which extend the range down to salinities close to zero, but it should be kept in mind that there is probably an assumption that this applies to diluted seawater, i.e. that the chemical composition is as it would be in the "deep ocean". Lakes with measurable conductivities rarely have chemical compositions that match the ocean. One hint is to check out the functions in the gsw package (on which oce depends), particularly the gsw_SP_from_C() function:

http://www.teos-10.org/pubs/gsw/html/gsw_SP_from_C.html

which states:

Note that the PSS-78 algorithm for Practical Salinity is only valid in the range 2 < SP < 42. If the PSS-78 algorithm produces a Practical Salinity that is less than 2 then the Practical Salinity is recalculated with a modified form of the Hill et al. (1986) formula. The modification of the Hill et al. (1986) expression is to ensure that it is exactly consistent with PSS-78 at SP = 2.

Hill, K.D., T.M. Dauphinee and D.J. Woods, 1986: The extension of the Practical Salinity Scale 1978 to low salinities. IEEE J. Oceanic Eng., OE-11, 1, 109 - 112.

Second, I'm a bit confused by your statement:

I'm wondering about swRho and about calculations of pressure and salinity since my data files lack these values and they must be calculated.

Salinity is always a derived quantity, but pressure should be present as an original measured quantity (there are no sensors that measure "depth" directly, as Dan discussed it needs to be derived from pressure, density, etc). I recommend if you have data files that only contain depth that you enquire as to the derivation employed, and perhaps try and obtain files that contain pressure instead. That way you can start from whatever formula for depth/density you find is best for freshwater, without having to re-calculate a derived quantity.

dankelley commented 8 years ago

I am pasting a screen snapshot, showing a page from the SBE processing manual, and a view of the oce::swDepth function. As you can see, the SBE formula for computing depth for freshwater is very simple. If your data file contains only depth and not pressure (which seems very weird to me, since a CTD does not measure depth but does measure pressure) then I suppose you could invert the formula. However, note that the SBE freshwater formula ignores latitude, which seems very strange to me. Perhaps SBE didn't bother trying to concoct a more reasonable formula because they didn't think their purchasers would be paying much attention to depth, as opposed to pressure.

screen shot 2016-08-03 at 6 53 56 am
dankelley commented 8 years ago

As a check, I ran the code

pressure <- 100
latitude <- seq(0, 90, 1)
depthFresh <- pressure * 1.019716
depthSalty <- swDepth(pressure, latitude)
plot(latitude, depthSalty/depthFresh,type='l')
grid()

and it gives as below. As you can see, the latitude correction amounts to about 0.5 percent variation in inferred depth. The salty/fresh alteration is about a 3 percent change.

screen shot 2016-08-03 at 7 05 01 am
richardsc commented 8 years ago

While none of this is currently coded into oce, I did a bit of looking at Rich's lim package for Matlab, which you can find here:

https://www.eoas.ubc.ca/~rich/#LIM

Much of what this package does (e.g. see limcond.m) is compute lake conductivity and salinity based on the chemical composition of the water. This in turn can be used with limstate.m to calculate density, according to: Chen and Millero, Limnol. Oceanogr. 31(3), 1986, 657-662.

The literature citation for Rich's package is:

http://www.aslo.org/lomethods/free/2008/0489.pdf

It's unclear to me how best to proceed in the case where the conductivity has been measured, but the chemical composition is not known. It looks like reading (and understanding) Rich's paper on page 493 is the right place:

Alternatively, it may be useful to use measurements of conductivity to estimate the salinity in a particular system. Rather than develop a second set of equations specifically for this purpose, we can perform this calculation by iteration. Note that salinity determination is not possible without making additional assumptions about the chemical composition of the system. In the simplest case, we assume that the system is composed of dilutions or concentrations of a known chemical composition

My guess is that in most cases, where the actual chemical composition of the lake is not known, that it is assumed to be similar enough to oceanic composition and the Hill et al extension to PSS-78 is used to estimate salinity and density. The screenshot below shows the difference between swSCTp() and gsw_SP_from_C() for a range of conductivities (with T=15, p=0), and between the densities inferred from those salinities (using swRho()): screen shot 2016-08-03 at 14 09 54

richardsc commented 8 years ago

Also, for reference, below is what the RBR software does for calculating density. Note that it doesn't make any distinction between freshwater or saltwater, but presumably that would be accounted for using a different density value (the units on that density must be g/ml, which is odd ... ). screen shot 2016-08-03 at 14 17 17

Note that there is the "Simplified" depth (no citation for this that I know of), and the "UNESCO" depth. I suspect that the UNESCO depth is the same as coded in swDepth(), however it's interesting to note that even that formulation doesn't specifically include the correction for depth variation of density, referred to in Unesco 1983, No. 44 as the "geopotential anomaly". They do state however that the correction based on density variations for p=10000 dbar is about 2m or less, significantly smaller than most pressure sensor accuracies.

dankelley commented 8 years ago

To @mlcoleman --

If you think the issue (as defined by its title) has been addressed, please close it. This is a generic message the developers send to reporters when it seems that an issue has been addressed and that discussion has ceased. Please bear in mind that we work on issues as defined by their titles. Narrow-scoped titles help us to focus our development efforts, and benefit other oce users.

If you think the issue has not been addressed, please look through all the comments and write a new one after this, stating what still needs to be done.

If you see that the issue has meandered from the topic defined by its title, please consider closing this and opening a new issue.

Thanks!

dankelley commented 8 years ago

Generic Message: The discussion of this issue seems to have ended, and the developers will likely close it soon, since it seems to have been addressed. This goes against the usual oce policy of hoping that reporters will close issues, so this is really a request to the original reporter to either (a) close the issue or (b) indicate that further action is desirable.

dankelley commented 8 years ago

Standardized comment. This issue has been closed, since it seems to have been addressed, and since discussion has stalled. The original reporter should feel free to reopen it, though, or to open new issues that might be related. Thanks.

slaeubli commented 10 months ago

Hi there, First of all: Thanks a lot for your great work! Oce is a fantastic package :)

I have a related question regarding freshwater data and depth calculations: I'm running into issues as I have freshwater CTD data from alpine lakes. I've tried with adjusting pressureAtmospheric in as.ctd to correct for the altitude but I'm trying out random numbers and see whether the pressure matches the known lake depths and the data. Because of the CTD calibration at sea level, i also get negative pressure values in the raw data. How do i correctly correct the pressure for the altitude in e.g., as.ctd?

dankelley commented 10 months ago

I am reopening this so that developers and other users can see it.

dankelley commented 10 months ago

As to your new comment, please post your code and your data. If the data are private, please email (or ftp, or dropbox, or google) the data to me at kelley dot dan at gmail.com and also to clark dot richards at gmail.com

Basically, you may need to apply an offset to your data. If you know that the instrument is in the water at the surface at a certain time, then the pressure at that time (plus perhaps 1m if you can see it 1m below the surface) will be your reference pressure. As for how to get z from p, you need to use the lake density and the hydrostatic equation $dp/dz=-\rho g$; see any textbook.

slaeubli commented 9 months ago

Ok, I see. But what does pressureAtmospheric do then in as.ctd?

dankelley commented 9 months ago

@slaeubli I think this is explained clearly in the help provided by typing ?as.ctd. If there is something about the description of pressureAtmospheric that is not clear, please let us know.

Below is an example.

library(oce)
#> Loading required package: gsw
as.ctd(35, 10, 3)@data
#> $scan
#> [1] 1
#> 
#> $salinity
#> [1] 35
#> 
#> $temperature
#> [1] 10
#> 
#> $pressure
#> [1] 3
as.ctd(35, 10, 3, pressureAtmospheric = 3)@data
#> $scan
#> [1] 1
#> 
#> $salinity
#> [1] 35
#> 
#> $temperature
#> [1] 10
#> 
#> $pressure
#> [1] 0

Created on 2023-11-15 with reprex v2.0.2

richardsc commented 9 months ago

As @dankelley has pointed out, pressureAtmospheric provides an offset to define the "zero pressure" -- i.e. the absolute (total) pressure that the CTD would experience just as it goes in the water. That should work fine for alpine lakes, but you haven't provided any data or examples to show the issue so I'm just guessing.

If your CTD records "sea pressure" (i.e. pressure adjusted for normal sealevel atmospheric pressure) then you will likely see a negative starting pressure at altitude. You should be able to adjust this with pressureAtmospheric, but again as was pointed out you should not rely on swDepth() to calculate depths because that algorithm is designed for seawater and normal sealevel starting pressures.