sunAngle function not recognizing the length of useRefraction

[rroday]

Hi! I'm new to the oce package and I am using the sunAngle function for a dataframe with ~40,000 different time x latitude/longitude combinations. When I use the function on my dataframe, it doesn't like that the lengths of the different parameters don't match (which they do!).

Here is a reproducible bit of code:

library(oce)
times <- as.POSIXct(rep("2011-03-03 06:49:00", 4))    
sunAngle(times, latitude = c(39, 38, 34, 35), longitude =  c(-75, -74, -73, -72), useRefraction = length(times))

Here is my error message: Error in sunAngle(times, latitude = c(39, 38, 34, 35), longitude = c(-75, : lengths of longitude and useRefraction must match

Someone on stack exchange recommended using Vectorize in the sunAngle function to fix this issue, but I thought I would repost my issue on your git since I am not the original poster on this link:

Session Info:

R version 4.3.0 (2023-04-21 ucrt) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 19045)  Matrix products: default   locale: [1] LC_COLLATE=English_United States.utf8  LC_CTYPE=English_United States.utf8    LC_MONETARY=English_United States.utf8 LC_NUMERIC=C                           [5] LC_TIME=English_United States.utf8      time zone: America/New_York tzcode source: internal  attached base packages: [1] grid      stats     graphics  grDevices utils     datasets  methods   base       other attached packages:  [1] oce_1.8-2               gsw_1.1-1               chillR_0.75             lunar_0.2-1             multcompView_0.1-9      rstatix_0.7.2            [7] cowplot_1.1.1           tm_0.7-11               NLP_0.2-1               rcompanion_2.4.34       scales_1.2.1            ggtext_0.1.2            [13] dataRetrieval_2.7.12    ggpubr_0.6.0            ggsignif_0.6.4          FSA_0.9.4               vcd_1.4-11              epitools_0.5-10.1       [19] lawstat_3.6             ggforce_0.4.1           plotrix_3.8-2           gifski_1.12.0           gganimate_1.0.8         mapdata_2.3.1           [25] sf_1.0-13               maps_3.4.1              ggrepel_0.9.3           ggspatial_1.1.8         rnaturalearthdata_0.1.0 rnaturalearth_0.3.3     [31] viridis_0.6.3           viridisLite_0.4.2       tidygraph_1.2.3         ggraph_2.1.0            igraph_1.4.3            network_1.18.1          [37] lubridate_1.9.2         forcats_1.0.0           stringr_1.5.0           purrr_1.0.1             readr_2.1.4             tidyr_1.3.0             [43] tibble_3.2.1            tidyverse_2.0.0         dplyr_1.1.2             data.table_1.14.8       glatos_0.4.0            ggplot2_3.4.2            loaded via a namespace (and not attached):   [1] libcoin_1.0-10       rstudioapi_0.14      jsonlite_1.8.5       magrittr_2.0.3       TH.data_1.1-2        modeltools_0.2-23    farver_2.1.1           [8] fields_15.2          vctrs_0.6.2          memoise_2.0.1        RCurl_1.98-1.14      progress_1.2.2       broom_1.0.4          cellranger_1.1.0      [15] GenSA_1.1.14         KernSmooth_2.23-20   plyr_1.8.8           sandwich_3.0-2       cachem_1.0.8         rootSolve_1.8.2.3    zoo_1.8-12            [22] ecmwfr_1.5.0         lifecycle_1.0.3      pkgconfig_2.0.3      fastmap_1.1.1        Matrix_1.5-4         R6_2.5.1             rbibutils_2.2.13      [29] digest_0.6.31        Exact_3.2            colorspace_2.1-0     patchwork_1.1.3      metR_0.15.0          fansi_1.0.4          timechange_0.2.0      [36] httr_1.4.6           polyclip_1.10-4      abind_1.4-5          compiler_4.3.0       proxy_0.4-27         remotes_2.5.0        withr_2.5.0           [43] backports_1.4.1      carData_3.0-5        DBI_1.1.3            R.utils_2.12.2       MASS_7.3-58.4        classInt_0.4-9       gld_2.6.6             [50] tools_4.3.0          units_0.8-2          lmtest_0.9-40        R.oo_1.25.0          glue_1.6.2           gridtext_0.1.5       checkmate_2.3.1       [57] generics_0.1.3       gtable_0.3.3         nortest_1.0-4        tzdb_0.4.0           R.methodsS3_1.8.2    class_7.3-21         lmom_3.0              [64] hms_1.1.3            sp_1.6-1             xml2_1.3.4           car_3.1-2            coin_1.4-3           utf8_1.2.3           pillar_1.9.0          [71] spam_2.10-0          splines_4.3.0        Kendall_2.2.1        tweenr_2.0.2         lattice_0.21-8       survival_3.5-5       tidyselect_1.2.0      [78] gridExtra_2.3        stats4_4.3.0         pls_2.8-3            expm_0.999-7         graphlayouts_1.0.0   matrixStats_1.0.0    stringi_1.7.12        [85] statnet.common_4.9.0 boot_1.3-28.1        codetools_0.2-19     cli_3.6.1            DescTools_0.99.50    Rdpack_2.4           munsell_0.5.0         [92] Rcpp_1.0.10          readxl_1.4.2         coda_0.19-4          XML_3.99-0.16.1      parallel_4.3.0       assertthat_0.2.1     prettyunits_1.1.1     [99] dotCall64_1.1-1      bitops_1.0-7         slam_0.1-50          mvtnorm_1.2-1        e1071_1.7-13         crayon_1.5.2         rlang_1.1.1          [106] multcomp_1.4-25

[richardsc]

Hi @rroday -- thanks for using oce!

@dankelley is more up-to-speed on the various sun functions, as they were written by him, but I think that to make your example work you just need to follow the hint in the error. Note that length(times) returns a vector of length 1. Also, the help page for sunAngle() indicates that useRefraction= is a boolean, i.e. TRUE/FALSE. So, the following should work:

library(oce)
#> Loading required package: gsw

times <- as.POSIXct(rep("2011-03-03 06:49:00", 4))
sunAngle(times, latitude = c(39, 38, 34, 35), longitude = c(-75, -74, -73, -72), 
                useRefraction = rep(TRUE, length(times)))
#> $time
#> [1] "2011-03-03 10:49:00 UTC" "2011-03-03 10:49:00 UTC"
#> [3] "2011-03-03 10:49:00 UTC" "2011-03-03 10:49:00 UTC"
#> 
#> $azimuth
#> [1] 91.75708 92.52231 93.62700 94.07240
#> 
#> $altitude
#> [1] -8.776681 -7.957257 -6.934286 -6.179128
#> 
#> $diameter
#> [1] 0.5378894 0.5378894 0.5378894 0.5378894
#> 
#> $distance
#> [1] 0.9912818 0.9912818 0.9912818 0.9912818
#> 
#> $declination
#> [1] -6.865559 -6.865559 -6.865559 -6.865559
#> 
#> $rightAscension
#> [1] -16.12483 -16.12483 -16.12483 -16.12483

^{Created on 2024-08-02 with reprex v2.1.0}

[richardsc]

This does look related to #2178, which arose from that SO post, and it looks like we maybe haven't done a CRAN release since then. My example above was using the develop branch of oce. @rroday, are you able to install oce from source following the instructions at:

https://dankelley.github.io/oce/#installing-oce

and see if that fixes your issue?

(@dankelley; we should think about a CRAN release)

[dankelley]

I just altered the "develop" branch (commit 9e408adcd3ab7aa84fbe99f11fd06fc043847e79) so that it now detects if useRefraction is not a logical value. See reprex below.

library(oce)
#> Loading required package: gsw
times <- as.POSIXct(rep("2011-03-03 06:49:00", 4))
sunAngle(times,
    latitude = c(39, 38, 34, 35),
    longitude = c(-75, -74, -73, -72),
    useRefraction = length(times)
)
#> Error in sunAngle(times, latitude = c(39, 38, 34, 35), longitude = c(-75, : useRefraction must be a logical value

^{Created on 2024-08-02 with reprex v2.1.1}

[dankelley]

I see that the docs say that longitude, latitude and useRefraction are made to have the same length as time. This is not what the code was doing. I will alter the code to do that, because we try to make code follow docs.

I will also add a note that this is useful only in for one use, namely that of constructing time-series of the angles. For other uses, the user ought to set up the arguments to be of equal length (if that makes sense) or to call the function repeatedly. I don't see a clever way of handling the range of possible things a person might want. Say one person wants to determine angles at a certain and certain latitude, but a bunch of longitudes ... and another person wants to use a certain time, some number N of latitudes and some other number M of latitudes ... and on and on. Handling all possibilities would be tricky. This is what expand.grid() is for ...

[rroday]

Hi Dr. Kelley and Dr. Richards,

Thank you immensely for your timely and thorough responses. I really appreciate you addressing my comments so quickly! I realized that my initial comment had an error. Rather than using "useRefraction = length(times)", I was using "useRefraction = rep(TRUE, length(times))", and was receiving the following error message: Error in sunAngle(as.POSIXct(dtc_final3$detection_timestamp_utc, tz = "UTC"), :

lengths of longitude, latitude and useRefraction must be 1 or the length of time

I'm not sure if this changes things that you have worked on this morning, but I will reinstall the "oce" package and rerun the code. My apologies for the error.

Best, Rachel Roday

On Fri, Aug 2, 2024 at 11:25 AM Dan Kelley @.***> wrote:

I see that the docs say that longitude, latitude and useRefraction are made to have the same length as time. This is not what the code was doing. I will alter the code to do that, because we try to make code follow docs.

I will also add a note that this is useful only in for one use, namely that of constructing time-series of the angles. For other uses, the user ought to set up the arguments to be of equal length (if that makes sense) or to call the function repeatedly. I don't see a clever way of handling the range of possible things a person might want. Say one person wants to determine angles at a certain and certain latitude, but a bunch of longitudes ... and another person wants to use a certain time, some number N of latitudes and some other number M of latitudes ... and on and on. Handling all possibilities would be tricky. This is what expand.grid() is for ...

— Reply to this email directly, view it on GitHub https://github.com/dankelley/oce/issues/2225#issuecomment-2265642224, or unsubscribe https://github.com/notifications/unsubscribe-auth/ASOQWUC3FPHAUVC2DBLTQYDZPOQHFAVCNFSM6AAAAABL4PU7ZCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDENRVGY2DEMRSGQ . You are receiving this because you were mentioned.Message ID: @.***>

[dankelley]

@rroday -- please note that you'll want to install oce from github (as @richardsc noted). I just made some modifications a moment ago, improving the documentation for sunAngle() to suggest that only 3 possibilities are available:

1. all parameters are of length 1
2. time has length exceeding 1, in which case rep() is used on the other arguments to make them of the same length.
3. for other cases, either call sunAngle() repeatedly, or use expand.grid() to make same-length parameters

Hopefully the docs are clear on that now. I think expand.grid() might be the way you want to go, depending on your purpose.

[dankelley]

@rroday if you find it helpful to format the text in comments, check https://docs.github.com/en/get-started/writing-on-github/getting-started-with-writing-and-formatting-on-github/basic-writing-and-formatting-syntax which gives hints. I took the liberty of reformatting your initial comment so I could read it a bit more easily.

[rroday]

@dankelley Fantastic, thanks so much!! And thanks for this link, this is the first issue I've opened on github.

[dankelley]

@rroday, in case it's of interest, this function was used a lot in the following paper. I can't remember, but I might have written it just for that occasion.

Hamelin, Kayla M., Daniel E. Kelley, Christopher T. Taggart, and Michael C. James. “Water Mass Characteristics and Solar Illumination Influence Leatherback Turtle Dive Patterns at High Latitudes.” Ecosphere 5, no. 2 (February 2014): 1–20. https://doi.org/10.1890/ES13-00158.1.

[dankelley]

@rroday, for fun, I made the following. It shows how the sunrise line (I'm sure there's a name for that) sweeps along Nova Scotia at this time of year. The squiggly contours go away if you set useRefraction to FALSE. I made this for the "official" sunrise time for Halifax, but the source I have does not give the seconds part, nor does it say what they mean by Halifax (city centre? centroid of environs boundaries?). Also, I'm not 100% sure what sunrise means on that site. Time of seeing edge of sun? Time when centre of sun is at horizon? Etc. So I'm not too worried about why the zero line does not go through my apartment.

```R # Demonstrate how to get sun angle at a collection of locations # and times. png("sunrise.png") library(oce) data(coastlineWorldFine, package = "ocedata") plot(coastlineWorldFine, clongitude = -63, clatitude = 45, span = 800) usr <- par("usr") lon <- seq(usr[1], usr[2], 0.1) lat <- seq(usr[3], usr[4], 0.1) # Halifax sunrise on Aug 6 is 06:06 local (09:06 UTC) t <- as.POSIXct("2024-08-06 09:06:00", tz = "UTC") g <- expand.grid(t = t, lon = lon, lat = lat) sa <- sunAngle( t = g$t, longitude = g$lon, latitude = g$lat, useRefraction = TRUE ) contour(lon, lat, matrix(sa$altitude, nrow = length(lon)), add = TRUE, labcex = 1 ) mapScalebar() mtext(t) ```

[dankelley]

@rroday do you think anything more needs to be done on this issue? In the oce project, we ask reporters to close issues, with developers only closing them if discussion has gone silent after having asked the reporters if things are ok. This is a way to avoid developers saying "not an issue" and silencing reporters.

Since I plan to release oce to CRAN on August 16th (when CRAN will again be accepting submissions), I'm keen to fix up this issue, if possible.

Thanks!

Dan