Monthly Values, Evapotranspiration and Derived Variables

This vignette goes from hourly data to quantities used in agriculture, hydrology and climatology: derived hourly variables, reference evapotranspiration and monthly summaries. It uses the one-station example dataset shipped with the package, so every chunk runs offline. With your own data, replace rmet_example by the result of inmet_get().

library(rmet)
hourly <- rmet_example
hourly[, c("station_code", "datetime", "temp_dry_c", "humid_rel_pct", "precip_mm")]
#> <rmet hourly> 8,663 rows x 5 columns, 1 station
#> Period: 2023-01-01 to 2023-12-31 23:00:00 (UTC)
#>    station_code            datetime temp_dry_c humid_rel_pct precip_mm
#> 1          A801 2023-01-01 00:00:00       26.5            55         0
#> 2          A801 2023-01-01 01:00:00       26.8            49         0
#> 3          A801 2023-01-01 02:00:00       26.3            50         0
#> 4          A801 2023-01-01 03:00:00       25.2            59         0
#> 5          A801 2023-01-01 04:00:00       24.1            65         0
#> 6          A801 2023-01-01 05:00:00       23.1            73         0
#> 7          A801 2023-01-01 06:00:00       22.0            81         0
#> 8          A801 2023-01-01 07:00:00       21.5            89         0
#> 9          A801 2023-01-01 08:00:00       21.0            91         0
#> 10         A801 2023-01-01 09:00:00       21.1            91         0
#> ... with 8,653 more rows

Hourly data: class and summary

Results of inmet_read() and inmet_get() are rmet_hourly data frames. They are ordinary data frames, with a compact print() and a summary() that shows, per station, the period and how complete the measurements are:

summary(hourly)
#> <summary of rmet hourly data>
#>   station_code n_rows      first                last complete_pct
#> 1         A801   8663 2023-01-01 2023-12-31 23:00:00         97.3
#> complete_pct: mean percentage of non-missing values in the measurement columns

1. Derived variables

inmet_derive() adds vapour pressure, vapour pressure deficit (VPD), wind components and the heat index. VPD, in kPa, is the usual measure of the drying power of the air:

hourly <- inmet_derive(hourly)
summary(as.data.frame(hourly[, c("vapor_pressure_kpa", "vpd_kpa", "heat_index_c")]))
#>  vapor_pressure_kpa    vpd_kpa         heat_index_c  
#>  Min.   :0.4933     Min.   :0.02905   Min.   : 3.02  
#>  1st Qu.:1.5274     1st Qu.:0.20492   1st Qu.:16.53  
#>  Median :1.8534     Median :0.44492   Median :20.67  
#>  Mean   :1.8670     Mean   :0.67634   Mean   :20.86  
#>  3rd Qu.:2.2669     3rd Qu.:0.93302   3rd Qu.:24.45  
#>  Max.   :3.4621     Max.   :5.10102   Max.   :45.60  
#>  NA's   :3          NA's   :3         NA's   :3

The wind components follow the meteorological convention: the direction is where the wind comes from, so a wind from the west (270 degrees) has a positive wind_u_ms (it blows towards the east).

Values flagged by inmet_qc() are ignored, so run the quality control first:

hourly <- inmet_qc(hourly, quiet = TRUE)
hourly <- inmet_derive(hourly)

2. Daily values and reference evapotranspiration

inmet_daily() builds daily summaries in each station’s local time. inmet_et0() adds the FAO-56 Penman-Monteith grass reference evapotranspiration, et0_mm, in mm per day:

daily <- inmet_daily(hourly)
daily <- inmet_et0(daily)
head(daily[, c("date", "temp_max_c", "temp_min_c", "radiation_mj_m2", "et0_mm")])
#> <rmet daily> 6 rows x 5 columns
#> Period: 2022-12-31 to 2023-01-05
#>         date temp_max_c temp_min_c radiation_mj_m2   et0_mm
#> 1 2022-12-31         NA         NA              NA       NA
#> 2 2023-01-01       37.0       20.8         28.9801 6.653436
#> 3 2023-01-02       33.6       23.6         19.3385 4.735852
#> 4 2023-01-03       30.0       20.8         28.6707 5.870368
#> 5 2023-01-04       32.9       20.0         29.5361 6.263259
#> 6 2023-01-05       33.0       19.3         30.7514 6.600444

Two points deserve care. First, et0_mm is NA on days when any input is missing, so the min_hours rule of inmet_daily() also applies to it. Second, daily radiation is the sum of the hourly values, so it underestimates the true value on days with missing daytime hours. Compare n_hours and the inputs before using ET0 for decisions such as irrigation scheduling.

3. Monthly values

inmet_monthly() aggregates the daily values by calendar month. A monthly value is NA unless at least min_days days (20 by default) have a valid value:

monthly <- inmet_monthly(daily)
monthly[, c("month", "n_days", "temp_mean_c", "precip_mm", "rain_days", "et0_mm")]
#> <rmet monthly> 13 rows x 6 columns
#> Period: 2022-12-01 to 2023-12-01
#>         month n_days temp_mean_c precip_mm rain_days    et0_mm
#> 1  2022-12-01      1          NA        NA        NA        NA
#> 2  2023-01-01     31    26.10725      53.2         5 172.28004
#> 3  2023-02-01     28    24.73563     116.4         8 128.35298
#> 4  2023-03-01     31    25.46325      89.8         9 119.69306
#> 5  2023-04-01     30    20.59958      35.2         3  78.73280
#> 6  2023-05-01     31    18.35235      98.8         5  46.28680
#> 7  2023-06-01     30    16.29611     247.6         9  39.75671
#> 8  2023-07-01     31    15.46237     210.8         9  44.28155
#> 9  2023-08-01     31    16.69462      77.0         7  68.55122
#> 10 2023-09-01     30    18.31847     394.8        17  62.83033
#> ... with 3 more rows

Months are local calendar months, so the first month of a year of UTC data can contain a single day (here, the evening of 31 December). Such partial months are NA because of min_days.

A simple water balance indicator, rainfall minus reference evapotranspiration, follows directly:

monthly$balance_mm <- monthly$precip_mm - monthly$et0_mm
round(monthly$balance_mm, 1)
#>  [1]     NA -119.1  -12.0  -29.9  -43.5   52.5  207.8  166.5    8.4  332.0
#> [11]   15.9  235.9  -58.7

4. Potential evapotranspiration and water balance

When only temperature is reliable, inmet_et0(method = "hargreaves") needs no radiation, wind or humidity. For monthly data, inmet_pet_thornthwaite() gives the potential evapotranspiration, and inmet_water_balance() follows the soil water month after month (Thornthwaite and Mather):

monthly <- inmet_pet_thornthwaite(monthly)
balance <- inmet_water_balance(monthly, awc = 100)
round(balance[, c("precip_mm", "pet_mm", "storage_mm", "deficit_mm", "excess_mm")][2:13, ], 1)
#> <rmet monthly> 12 rows x 5 columns
#>    precip_mm pet_mm storage_mm deficit_mm excess_mm
#> 2       53.2  154.5       36.3       37.6       0.0
#> 3      116.4  117.2       36.0        0.5       0.0
#> 4       89.8  129.8       24.1       28.1       0.0
#> 5       35.2   71.6       16.8       29.0       0.0
#> 6       98.8   53.0       62.6        0.0       0.0
#> 7      247.6   37.7      100.0        0.0     172.5
#> 8      210.8   35.1      100.0        0.0     175.7
#> 9       77.0   44.6      100.0        0.0      32.4
#> 10     394.8   57.7      100.0        0.0     337.1
#> 11     110.0   69.7      100.0        0.0      40.3
#> ... with 2 more rows

awc is the available water capacity of the soil, in mm. The first months depend on the initial storage (initial), so use a series that starts before the period of interest; a single year, as here, only illustrates the output.

5. Climate indices

inmet_indices() computes ETCCDI-style indices per station and year (or month) from the daily values: extremes, counts of hot days and frost days, growing degree days, total rainfall on wet days and the longest dry and wet spells. An index is NA when too many days are missing (min_valid):

indices <- inmet_indices(daily)
t(indices[, c("txx", "tnn", "summer_days", "gdd", "prcptot", "rx1day", "cdd", "cwd")])
#>             [,1]   [,2]
#> txx           NA   38.8
#> tnn           NA    3.9
#> summer_days   NA  210.0
#> gdd           NA 4091.6
#> prcptot       NA 1851.2
#> rx1day        NA   89.4
#> cdd           NA   19.0
#> cwd           NA    6.0

The first column is the evening of 31 December 2022 in local time, the start of the UTC year of the example data, which is too short to give any index.

With 20 or more years of monthly data, inmet_spi() gives the Standardized Precipitation Index, and inmet_normals() and inmet_anomaly() the climatological normals and the anomalies relative to them.

6. Short gaps

inmet_fill() fills short gaps by linear interpolation and, if there are several years, with the typical value of the day of the year. Every filled value is marked, so the result can be traced back to measurements. Rainfall is never filled:

holes <- hourly
holes$temp_dry_c[100:102] <- NA
filled <- inmet_fill(holes, "temp_dry_c", method = "linear")
table(filled$temp_dry_c_filled, useNA = "ifany")
#> 
#> linear   <NA> 
#>      3   8660

7. Stations and locations

inmet_nearest() accepts many locations at once, inmet_stations_within() selects the stations inside a box or polygon, and inmet_interpolate() estimates a variable at other places from the stations by inverse-distance weighting:

inmet_nearest(c(-30.03, -23.55), c(-51.23, -46.63), n = 2, source = "bundled")[
  , c("location", "rank", "station_code", "station_name", "distance_km")]
#>   location rank station_code                   station_name distance_km
#> 1        1    1         A801 PORTO ALEGRE - JARDIM BOTANICO    5.933358
#> 2        1    2         B808                   CACHOEIRINHA   13.427377
#> 3        2    1         A701            SAO PAULO - MIRANTE    6.057680
#> 4        2    2         A771         SAO PAULO - INTERLAGOS   19.991750

8. Finding stations with data

Before downloading and reading, inmet_availability() lists the stations and periods in archives that are already in the cache, using only the file names:

inmet_availability(dest_dir = file.path(tempdir(), "rmet-vignette-cache"))
#> [1] station_code state        year         start        end         
#> [6] n_files     
#> <0 linhas> (ou row.names de comprimento 0)

With a populated cache it has one row per station and year. To keep the cache somewhere else for a whole session, set options(rmet.cache_dir = "path") or the environment variable RMET_CACHE_DIR; inmet_cache_dir() shows the directory in use.