ArgoCanada / argoFloats

Tools for analyzing collections of oceanographic Argo floats
https://argocanada.github.io/argoFloats/index.html
17 stars 7 forks source link

does applyFlags work correctly? #594

Closed dankelley closed 1 year ago

dankelley commented 1 year ago

I won't bother explaining, since it's the weekend. Also, this is for my branch 591_plot_TS_control in a version that is not pushed to GH. In other words, this is a private note, just to remind me to take a look sometime.

library(oce)
library(argoFloats)
if (!exists("a")) {
    a <- getIndex() |> subset(ID="6990512") |> getProfiles() |> readProfiles()
    b <- a |> applyQC()
}
par(mfrow=c(1, 2))
plotArgoTS(a, cex=1, pch=1, debug=0)
usr <- par("usr")
plotArgoTS(b, cex=1, pch=1, xlim=usr[1:2], ylim=usr[3:4], debug=1, xaxsr="i", yaxs="i")

argo

dankelley commented 1 year ago

The discrepancy is colourization was done based on only S and T flags. But, if users prefer to use the new "gsw" / TEOS10 formulation, then pressure also has to have been flagged as OK., otherwise neither SA nor CT can be computed. (With the old unesco equation of state, salinity was not affected but temperature was done as potential temperature, so that was okay.)

I am finding that quite a lot of pressure data are flagged as bad. The pattern I'm seeing is that if p[i+1] equals p[i], then pflag[i+1] is set to a "bad" value. Frankly, I am not too sure I think this flag condition makes sense, for CT/SA computations anyway. Who cares if two numbers happen to agree to 10cm or whatever?

I'd like to discuss this with @richardsc at some point, because it affects the oce package also. It turns out that today is a holiday, so maybe Wednesday? (I noticed that it was very quiet walking in to work today. For 2 kilometres, I didn't see a single person or a single moving vehicle. Very creepy. I guess the fireworks last night might have been a clue, but wouldn't they be tonight?)

dankelley commented 1 year ago

I've resolved this in commit 8ba743b4bb49ca9b95b0868effa2394eb8e4bdcb of branch 591_plot_TS_control. As noted in the just-previous comment, the problem was that the new code (and possibly the old code) was only looking at the S and T flags, but ignoring the p flags. The following code and diagram illustrate. (I am also including as an attachment an archive that contains the code, the diagram, and an rda file that makes the code reproducible.)

In the diagram, I'm zooming in on a region of a particular cycle, and a particular zone of that cycle where I was seeing problems in my initial comment. The small symbols are what plot,argo-method() produces, as of the commit/branch stated above. red means bad data, black means good data. (There are no gray points, meaning unassessed data, in this zone.). The larger symbols are what I get by a test I did starting from scratch. In the larger symbols, green means good and blue means not-good (or bad, if you're a glass-half-empty person).

Note that the large and small symbols now agree. That's enough for me to state that the new code works as expected.

Note that I've rewritten a lot of the documentation, to be clearer. Also, note that I have modified all the documentation to refer to the 2022 version of the argo manual. They seem to not name authors anymore, which is kind of weird because I like the (author, date) format, but I decided to go with (reference 1) format.

I am not closing the issue. I want to wait until @richardsc gets the chance to weigh in (partly because this is a bit adjacent to the oce package). Also, I have not merged this into "develop", and don't intend to do so with a face-to-face with @richardsc, just to be sure. (I know @j-harbin works on other projects now, so I'm mainly pinging @richardsc now.)

argo594b

```R # investigate a profile in detail focusOnAPoint <- FALSE library(oce) library(argoFloats) options(digits=7) # helps in close examination of values # cache for speed if (file.exists("argo594.rda")) { load("argo594.rda") } else { a <- getIndex() |> subset(ID="6990512") |> getProfiles() |> readProfiles() b <- a |> applyQC() aOrig <- a bOrig <- b # number 5 has some bad data, so that's our test case a@data$argos <- list(aOrig@data$argos[[5]]) b@data$argos <- list(bOrig@data$argos[[5]]) save(a, b, file="argo594.rda") } if (!interactive()) png("argo594b.png", unit="in", width=7, height=5, res=200) par(mfrow=c(1, 2)) pch <- 1#0 cex <- 0.8 lwd <- 1.9 plotArgoTS(a, xlim=c(35.8, 36.2), ylim=c(19,21), cex=cex, pch=pch, debug=0, lwd=lwd) usr <- par("usr") #examine <- locator(1) # produces next (why is this removed?) if (focusOnAPoint) { SA0 <- 36.0664423618411 CT0 <- 19.9310114103092 points(SA0, CT0, cex=2, pch=3, col=4) SA <- unlist(a[["SA"]]) CT <- unlist(a[["CT"]]) i <- which.min((SA-SA0)^2/30^2 + (CT-CT0)^2/10^2) # 1142 } # Isolate data near SA0,CT0 for visual inspection aa <- a@data$argos[[1]] F <- as.data.frame(aa@metadata$flags) D <- as.data.frame(aa@data) span <- 10 if (focusOnAPoint) { FF <- F[i+seq(-span, span), c(1, 3, 5)] DD <- D[i+seq(-span, span), c(4, 10, 7, 1, 2)] } else { FF <- F[, c(1, 3, 5)] DD <- D[, c(4, 10, 7, 1, 2)] } names(FF) <- c("pFlag", "SFlag", "TFlag") names(DD) <- c("p", "S", "T", "lon", "lat") FFDD <- cbind(FF, DD) FFDD$SA <- gsw_SA_from_SP(FFDD$S, FFDD$p, FFDD$lon, FFDD$lat) FFDD$CT <- gsw_CT_from_t(FFDD$SA, FFDD$T, FFDD$p) FFDD good <- FFDD$pFlag==1 & FFDD$SFlag==1 & FFDD$TFlag==1 points(FFDD$SA[good], FFDD$CT[good], cex=2, col=3) # green for good points(FFDD$SA[!good], FFDD$CT[!good], cex=2, col=4) # blue for bad plotArgoTS(b, cex=cex, pch=pch, xlim=usr[1:2], ylim=usr[3:4], debug=1, xaxs="i", yaxs="i", lwd=lwd) points(FFDD$SA[good], FFDD$CT[good], cex=2, col=3) # green for good points(FFDD$SA[!good], FFDD$CT[!good], cex=2, col=4) # blue for bad #for (i in seq_len(a[["length"]])) { # cat("For i=",i,", table of salinity flags is:\n") # print(table(a[[i]]@metadata$flags$salinity)) #} if (!interactive()) dev.off() ``` [Archive.zip](https://github.com/ArgoCanada/argoFloats/files/12280169/Archive.zip)
github-actions[bot] commented 1 year ago

The Stale-bot has marked this issue as Stale, because no new comments have been added in the past 30 days. Unless a comment is added within the next 7 days, the Stale-bot will close the issue. The purpose of these automated actions is to prevent the developers from forgetting about unattended tasks. Note that adding a "pinned" label will turn this action off for a given issue.