---
title: "skylight use"
author: "Koen Hufkens"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{skylight use}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

```{r}
# check version number
process <- as.numeric(version$major) + as.numeric(version$minor)/10 > 4.2
```

### Single date/time and location

skylight values can be calculated for a single point and date using the below call. This will generate a data frame with model values.

```{r}
# load the library
library(skylight)

# calculate sky illuminance values for
# a single date/time and location
df <- skylight(
      longitude = -135.8,
      latitude = -23.4,
      date = as.POSIXct("1986-12-18 21:00:00", tz = "GMT"),
      sky_condition = 1
    )

print(df)
```

### Multiple dates/times and/or locations

The skylight function is vectorized, so you can provide vectors of input parameters instead of using a loop and the above function call.

```{r}
# Generate a dataset with 15 minute values
# for approximately two months
input <- data.frame(
  longitude = 0,
  latitude = 50,
  date =  as.POSIXct("2020-06-18 00:00:00", tz = "GMT") + seq(0, 29*24*3600, 1800),
  sky_condition = 1
)

# calculate sky illuminance values for
# a single date/time and location
df <- skylight(
      longitude = input$longitude,
      latitude = input$latitude,
      date = input$date,
      sky_condition = 1
    )

print(head(df))

# previous results are of the same dimension (rows)
# as the input data and can be bound together
# for easy plotting
input <- cbind(input, df)
```

```{r}
library(ggplot2)

ggplot(input) +
  geom_tile(
    aes(
      as.Date(date),
      as.numeric(format(date, "%H")) + as.numeric(format(date, "%M"))/60,
      fill = log(total_illuminance)
    )
  ) +
  scale_fill_viridis_c(
    option = "B"
  ) +
  labs(
    title = "Diurnal cycles of total illuminance",
    subtitle = "(including the effect of moon contributions)",
    y = "Hour",
    x = ""
  ) +
  theme(
    legend.position = "bottom"
  )
```

Condensing this figure to only the mean half hourly value for a full moon captures the contributions of both sun and moon to the diurnal profile.

```{r echo = FALSE, warning=FALSE, message=FALSE, eval = process}
library(dplyr)
library(tidyr)

df <- input |>
  filter(
    moon_fraction == 100
  ) |>
  mutate(
    hour = as.numeric(format(date, "%H")),
    hour = ifelse(hour >= 12, hour - 24, hour)
  ) |>
  group_by(hour) |>
  summarize(
    sun_illuminance = mean(sun_illuminance),
    moon_illuminance = mean(moon_illuminance),
    total_illuminance = mean(total_illuminance),
    .groups = "drop"
  ) |>
  tidyr::pivot_longer(
    cols = ends_with("illuminance"),
    names_to = "illuminance"
  )

p <- ggplot(df) +
  geom_path(
    aes(
      hour,
      log(value),
      group = illuminance,
      colour = illuminance
    )
  ) +
  coord_cartesian(
    ylim = c(-4, 2),
    xlim = c(-5, 5)
  ) +
  labs(
    title = "Components of total illuminance",
    subtitle = "(30min increments during full moon)",
    x = "Hour (centered on midnight)"
  )

print(p)
```

### Piped data workflow

`skylight` supports piped data frames with appropriatedly named columns as input
to the function. This allows for fast processing of large data frames, with the
added advantage that input parameters are returned with the calculated data.

Note that you need a data frame with the three most basic parameters: **longitude**, 
**latitude**, **date**, named as such (all lower case). The function will complain 
if it doesn't find the required column names. Also note that due to the priority
of the piped construction over the other parameters **all parameters should be named**
when calling the function in a conventional way.

```{r eval = process}
# recreating the data frame with parameters
# as before
input <- data.frame(
  longitude = 0,
  latitude = 50,
  date =  as.POSIXct("2020-06-18 00:00:00", tz = "GMT") + seq(0, 1*24*3600, 1800),
  sky_condition = 1
)

# but now using the piped approach to calculate
# all values

df <- input |> skylight()

print(head(df))
```