-
Notifications
You must be signed in to change notification settings - Fork 42
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
tidem() and t_tide do not produce identical results #1653
Comments
I'll look into this. Please note that I created https://github.com/dankelley/oce-issues/tree/master/16xx/1653 and put the dataset there, along with some code I wrote to test, and your .pdf. If you could drop your Rmd into that directory, it might be helpful to me. The fact that only one constituent is affected tells me that it might be some kind of satellite / inference thing relating to grouped constituents. This is explained in the Foreman documents, referred to in Since I'll need to study the Foreman documents, and the tricky code that does some of this satellite/inference work, my guess is that this issue will take a while. However, your very detailed report will make that work a lot easier, so thanks for that!! Also, Clark, if you could put the matlab output file in the directory, that would be great, too. |
I'll post the Rmd and the matlab output -- likely tomorrow when I'm back in the office. |
Here is the Rmd: And a zip with all the files required (Matlab scripts and saved files): |
Further to our f2f regarding possible problems with phase inference for small-amplitude constituents, the code ## Test whether phase drifts for small-amplitude constituents
library(oce)
data(tidedata)
dt <- 15 * 60 # 15 minutes
t <- seq(as.POSIXct(Sys.time()), length.out=30*86400/dt, by=dt)
tt <- as.numeric(t)
## frequencies are in cycles/hour
M2 <- tidedata$const$freq[tidedata$const$name == "M2"] / 3600 * (2 * pi)
S2 <- tidedata$const$freq[tidedata$const$name == "S2"] / 3600 * (2 * pi)
M4 <- tidedata$const$freq[tidedata$const$name == "M4"] / 3600 * (2 * pi)
aM2 <- 1
aS2 <- 0.5
for (aM4 in 10^seq(1, -5)) {
eta <- aM2*sin(M2*tt) + aM2*cos(M2*tt) + aS2*sin(S2*tt) + aS2*cos(S2*tt) + aM4*sin(M4*tt) + aM4*cos(M4*tt)
sl <- as.sealevel(eta, t)
capture.output(m <- tidem(sl))
## summary(m)
phase <- m[["phase"]]
name <- m[["name"]]
cat(sprintf("aM4=%.5f yields phase=%.3f\n", aM4, phase[name=="S4"]))
} produces
I'm not sure I understand this, but thought it worth posting the results, to have them on the record. |
I dumped the CR zipfile into https://github.com/dankelley/oce-issues/tree/master/16xx/1653/cr and excerpted the table comparing results into https://github.com/dankelley/oce-issues/blob/master/16xx/1653/dk/1653dk_01.R so I could look in more detail. That R script has enough comments to explain, so I won't replicate what I did in words here. The output of https://github.com/dankelley/oce-issues/blob/master/16xx/1653/dk/1653dk_01.R is at Appendix 1 below. There's a lot I could look at, but my eye was drawn to the first big phase difference at MO3. I see big phase differences nearby (in freq) at MK3 and SK3. Now, take a look at Appendix 2, which is a table (figure, really) from Foreman (1981). This deals with frequencies that are near each other. This time-series is 768 hours long, so the bottom row of the table is irrelevant (i.e. the Rayleigh criterion means we do not fit for SO3). But look at the others. Note that M3 is related (in this sense) to MO3, MK3 and SK3. And (strokes 2-day old beard) this triplet, MO3, MK3 and SK3 have phase differences in the hundreds of degrees. I need to do some more reading on this topic. (I realize that I am not explaining this well here, but this is something @richardsc and I were discussing at our last f2f, so I think there's no real need to inarticulately state something here that is probably very clear in Foreman's document.) References
Appendix 1 Name Freq msnr phaseMismatch timeShiftHours
1 MM 0.001512 1.70 0.03 0.055
2 MSF 0.002822 0.41 0.02 0.020
3 *ALP1 0.034397 2.80 0.27 0.022
4 2Q1 0.035706 1.30 -0.30 -0.023
5 *Q1 0.037219 5.00 0.10 0.007
6 *O1 0.038731 150.00 -0.03 -0.002
7 NO1 0.040269 1.50 0.48 0.033
8 *K1 0.041781 220.00 0.06 0.004
9 J1 0.043293 0.87 -0.66 -0.042
10 OO1 0.044831 0.68 0.48 0.030
11 UPS1 0.046343 0.46 1.00 0.060
12 *EPS2 0.076177 2.40 0.09 0.003
13 *MU2 0.077689 4.60 -0.07 -0.003
14 *N2 0.078999 110.00 0.00 0.000
15 *M2 0.080511 1300.00 -0.03 -0.001
16 *L2 0.082024 15.00 -0.02 -0.001
17 *S2 0.083333 38.00 -0.04 -0.001
18 ETA2 0.085074 1.10 -0.22 -0.007
19 *MO3 0.119242 5.40 271.58 6.327
20 M3 0.120767 0.85 0.66 0.015
21 *MK3 0.122292 3.50 122.52 2.783
22 SK3 0.125114 0.12 109.09 2.422
23 *MN4 0.159511 2.60 23.88 0.416
24 *M4 0.161023 14.00 22.90 0.395
25 *SN4 0.162333 2.60 12.82 0.219
26 MS4 0.163845 1.60 11.58 0.196
27 S4 0.166667 0.52 0.56 0.009
28 2MK5 0.202804 1.00 135.69 1.859
29 *2SK5 0.208447 2.90 -249.15 -3.320
30 2MN6 0.240022 1.10 35.18 0.407
31 *M6 0.241534 4.00 34.47 0.396
32 *2MS6 0.244356 2.40 23.15 0.263
33 2SM6 0.247178 0.49 12.43 0.140
34 3MK7 0.283315 0.26 -216.01 -2.118
35 *M8 0.322046 2.60 45.74 0.395 Appendix 2 |
There's a lot of info in #1621 , which lead to the following exploration of
tidem()
andt_tide
on a particular dataset (but there's not reason to think that the results are specific to that dataset):compare_t_tide_with_tidem.pdf
Below is the result of a 10 constituent fit using
tidem()
:and below the same but for
t_tide
:Note that the fit amplitudes are all quite similar, as are the phases -- except for the
M4
phase (200.1 vs 177.15).Why is this?
The text was updated successfully, but these errors were encountered: