reshaping-meteo.Rmd
library(meteoland)
#> Package 'meteoland' [ver. 2.2.1]
library(stars)
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUE
library(dplyr)
#>
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#>
#> filter, lag
#> The following objects are masked from 'package:base':
#>
#> intersect, setdiff, setequal, unionFor the interpolation of meteorological variables on our target
locations, we need meteorological data for a reference set of locations.
In meteoland this is a sf object with the
spatial coordinates of our reference locations (usually meteorological
stations) and daily values of the meteorological variables needed to
perform the interpolation, i.e.:
meteoland_meteo_example
#> Simple feature collection with 5652 features and 18 fields
#> Geometry type: POINT
#> Dimension: XY
#> Bounding box: xmin: 0.30565 ymin: 40.55786 xmax: 3.18165 ymax: 42.77011
#> Geodetic CRS: WGS 84
#> # A tibble: 5,652 × 19
#> dates service stationID station_name station_province elevation
#> * <dttm> <chr> <chr> <chr> <chr> <dbl>
#> 1 2022-04-01 00:00:00 meteoc… C6 Castellnou … Lleida 264
#> 2 2022-04-01 00:00:00 meteoc… C7 Tàrrega Lleida 427
#> 3 2022-04-01 00:00:00 meteoc… C8 Cervera Lleida 554
#> 4 2022-04-01 00:00:00 meteoc… C9 Mas de Barb… Tarragona 240
#> 5 2022-04-01 00:00:00 meteoc… CC Orís Barcelona 626
#> 6 2022-04-01 00:00:00 meteoc… CD la Seu d'Ur… Lleida 849
#> 7 2022-04-01 00:00:00 meteoc… CE els Hostale… Barcelona 316
#> 8 2022-04-01 00:00:00 meteoc… CG Molló - Fab… Girona 1405
#> 9 2022-04-01 00:00:00 meteoc… CI Sant Pau de… Girona 852
#> 10 2022-04-01 00:00:00 meteoc… CJ Organyà Lleida 566.
#> # ℹ 5,642 more rows
#> # ℹ 13 more variables: MeanTemperature <dbl>, MinTemperature <dbl>,
#> # MaxTemperature <dbl>, MeanRelativeHumidity <dbl>,
#> # MinRelativeHumidity <dbl>, MaxRelativeHumidity <dbl>, Precipitation <dbl>,
#> # WindDirection <dbl>, WindSpeed <dbl>, Radiation <dbl>, geom <POINT [°]>,
#> # aspect <dbl>, slope <dbl>meteoland expects names to be as in the example:
names(meteoland_meteo_example)
#> [1] "dates" "service" "stationID"
#> [4] "station_name" "station_province" "elevation"
#> [7] "MeanTemperature" "MinTemperature" "MaxTemperature"
#> [10] "MeanRelativeHumidity" "MinRelativeHumidity" "MaxRelativeHumidity"
#> [13] "Precipitation" "WindDirection" "WindSpeed"
#> [16] "Radiation" "geom" "aspect"
#> [19] "slope"The only mandatory variables are MinTemperature and
MaxTemperature. Other variables
(Precipitation, WindSpeed…), when present,
allow for a more complete interpolation.
meteoland format
Meteorological data can come in many formats and with different names
for the same variables. As we saw above, we need to convert it to a
meteoland compliant format (sf object with
correct names).
If we have our meteorological data in a data.frame (i.e. we obtained it from our own weather stations, or download it from other source in this format) we can just simply transform it to the desired format:
unformatted_meteo
#> # A tibble: 15 × 7
#> date station latitude longitude min_temp max_temp rh
#> <date> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 2022-12-01 a 41.4 -0.33 1.86 16.8 20
#> 2 2022-12-02 a 41.4 -0.33 7.68 12.5 68.1
#> 3 2022-12-03 a 41.4 -0.33 7.31 19.4 72.1
#> 4 2022-12-04 a 41.4 -0.33 14.6 24.5 54.4
#> 5 2022-12-05 a 41.4 -0.33 5.71 16.7 47.4
#> 6 2022-12-01 b 40.1 0.12 3.62 10.1 61.1
#> 7 2022-12-02 b 40.1 0.12 8.77 17.5 46.8
#> 8 2022-12-03 b 40.1 0.12 16.6 21.4 20
#> 9 2022-12-04 b 40.1 0.12 13.8 30.3 20
#> 10 2022-12-05 b 40.1 0.12 -0.782 7.45 57.9
#> 11 2022-12-01 c 42 1.12 23.1 32.0 30.0
#> 12 2022-12-02 c 42 1.12 14.5 28.0 30.3
#> 13 2022-12-03 c 42 1.12 7.66 16.9 60.1
#> 14 2022-12-04 c 42 1.12 8.71 18.2 20
#> 15 2022-12-05 c 42 1.12 12.7 24.4 22.6
ready_meteo <- unformatted_meteo |>
# convert names to correct ones
dplyr::mutate(
MinTemperature = min_temp,
MaxTemperature = max_temp,
MeanRelativeHumidity = rh
) |>
# transform to sf (WGS84)
sf::st_as_sf(
coords = c("longitude", "latitude"),
crs = sf::st_crs(4326)
)
ready_meteo
#> Simple feature collection with 15 features and 8 fields
#> Geometry type: POINT
#> Dimension: XY
#> Bounding box: xmin: -0.33 ymin: 40.11 xmax: 1.12 ymax: 42
#> Geodetic CRS: WGS 84
#> # A tibble: 15 × 9
#> date station min_temp max_temp rh MinTemperature MaxTemperature
#> * <date> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 2022-12-01 a 1.86 16.8 20 1.86 16.8
#> 2 2022-12-02 a 7.68 12.5 68.1 7.68 12.5
#> 3 2022-12-03 a 7.31 19.4 72.1 7.31 19.4
#> 4 2022-12-04 a 14.6 24.5 54.4 14.6 24.5
#> 5 2022-12-05 a 5.71 16.7 47.4 5.71 16.7
#> 6 2022-12-01 b 3.62 10.1 61.1 3.62 10.1
#> 7 2022-12-02 b 8.77 17.5 46.8 8.77 17.5
#> 8 2022-12-03 b 16.6 21.4 20 16.6 21.4
#> 9 2022-12-04 b 13.8 30.3 20 13.8 30.3
#> 10 2022-12-05 b -0.782 7.45 57.9 -0.782 7.45
#> 11 2022-12-01 c 23.1 32.0 30.0 23.1 32.0
#> 12 2022-12-02 c 14.5 28.0 30.3 14.5 28.0
#> 13 2022-12-03 c 7.66 16.9 60.1 7.66 16.9
#> 14 2022-12-04 c 8.71 18.2 20 8.71 18.2
#> 15 2022-12-05 c 12.7 24.4 22.6 12.7 24.4
#> # ℹ 2 more variables: MeanRelativeHumidity <dbl>, geometry <POINT [°]>And voilà, we have our meteo data in the correct format
meteoland offers transformation functions for
meteorological data downloaded from the meteospain
and worldmet
R packages.
meteospain data
For data coming from meteospain package we have the
meteospain2meteoland function that transforming the data
for us:
library(meteospain)
get_meteo_from(
"meteogalicia",
meteogalicia_options('daily', as.Date("2022-12-01"), as.Date("2022-12-05"))
) |>
meteospain2meteoland()worldmet data
For data coming from worldmet package we have the
worldmet2meteoland function that do the reshaping for
us:
library(worldmet)
worldmet::importNOAA("081120-99999", year = 2022) |>
worldmet2meteoland()As we have seen, we need the meteorological reference data in a
sf object (points). To be able to create an interpolator
from a raster, we need to transform the cell values to points
(i.e. using the cell center coordinates) for each variable and
day and use it to create the interpolator.
So if we have a multi-layered raster with several dates of meteorological data:
> Remember that the raster, besides the meteo data, needs to contain also the topographic
(elevation, aspect and slope) data for the interpolator to work.
raster_meteo_reference
#> stars object with 3 dimensions and 14 attributes
#> attribute(s):
#> Min. 1st Qu. Median Mean
#> MeanTemperature 3.47453413 11.4535839 13.940206 13.471788
#> MinTemperature -3.27100714 4.3824425 6.969667 5.860040
#> MaxTemperature 7.16265820 14.7798101 18.901768 18.420680
#> Precipitation 0.00000000 0.0000000 0.000000 1.217155
#> MeanRelativeHumidity 35.08268335 56.7720829 64.269391 65.520741
#> MinRelativeHumidity 26.85783275 38.7257795 45.058448 48.591400
#> MaxRelativeHumidity 53.84205927 100.0000000 100.000000 96.878881
#> Radiation 7.54372475 15.9837124 20.690395 19.595816
#> WindSpeed 0.02164994 0.9119314 1.252240 1.385929
#> WindDirection 0.18621266 71.4293342 198.955399 181.615512
#> PET 1.09027039 2.4105288 3.176816 3.087621
#> elevation 240.00000000 370.0000000 447.000000 460.322314
#> slope 1.43209624 5.7204332 11.348120 13.073426
#> aspect 5.19442749 74.7448807 174.369324 181.679232
#> 3rd Qu. Max. NA's
#> MeanTemperature 16.1442989 20.781408 0
#> MinTemperature 8.2924049 11.071295 0
#> MaxTemperature 22.0316730 28.801911 0
#> Precipitation 0.2922964 21.000409 0
#> MeanRelativeHumidity 75.5124486 100.000000 0
#> MinRelativeHumidity 55.7307897 90.243329 0
#> MaxRelativeHumidity 100.0000000 100.000000 0
#> Radiation 23.7299248 27.982274 0
#> WindSpeed 1.7444530 5.811866 30
#> WindDirection 264.0404014 359.908196 1350
#> PET 3.7736818 5.625855 0
#> elevation 525.0000000 786.000000 0
#> slope 19.7585106 31.071896 0
#> aspect 291.0375061 360.000000 0
#> dimension(s):
#> from to offset delta refsys point x/y
#> x 1 11 1.671 0.01058 WGS 84 FALSE [x]
#> y 1 11 41.76 -0.01058 WGS 84 FALSE [y]
#> date 1 30 2022-04-01 UTC 1 days POSIXct FALSEwe need to convert it to points:
points_meteo_reference <- names(raster_meteo_reference) |>
# for each variable
purrr::map(
# take the variable raster
~ raster_meteo_reference[.x] |>
# convert to sf
sf::st_as_sf(as_points = TRUE, na.rm = FALSE) |>
# pivot the data for dates to be in one column
tidyr::pivot_longer(cols = -geometry, names_to = "dates", values_to = .x) |>
# convert to tibble to fasten the process
dplyr::as_tibble() |>
# convert to date and create stationID
dplyr::mutate(
dates = as.Date(dates),
stationID = as.character(geometry)
)
) |>
# join all variables
purrr::reduce(dplyr::left_join) |>
# create the points sf object
sf::st_as_sf()
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
#> Joining with `by = join_by(geometry, dates, stationID)`
points_meteo_reference
#> Simple feature collection with 3630 features and 16 fields
#> Geometry type: POINT
#> Dimension: XY
#> Bounding box: xmin: 1.6762 ymin: 41.65133 xmax: 1.781988 ymax: 41.75712
#> Geodetic CRS: WGS 84
#> # A tibble: 3,630 × 17
#> geometry dates MeanTemperature stationID MinTemperature
#> <POINT [°]> <date> <dbl> <chr> <dbl>
#> 1 (1.6762 41.75712) 2022-04-01 6.24 c(1.676199866843… 0.541
#> 2 (1.6762 41.75712) 2022-04-02 6.64 c(1.676199866843… -1.50
#> 3 (1.6762 41.75712) 2022-04-03 4.82 c(1.676199866843… -3.13
#> 4 (1.6762 41.75712) 2022-04-04 6.53 c(1.676199866843… -1.59
#> 5 (1.6762 41.75712) 2022-04-05 9.05 c(1.676199866843… -1.51
#> 6 (1.6762 41.75712) 2022-04-06 11.8 c(1.676199866843… 3.70
#> 7 (1.6762 41.75712) 2022-04-07 14.1 c(1.676199866843… 3.67
#> 8 (1.6762 41.75712) 2022-04-08 15.1 c(1.676199866843… 7.71
#> 9 (1.6762 41.75712) 2022-04-09 13.6 c(1.676199866843… 6.79
#> 10 (1.6762 41.75712) 2022-04-10 11.6 c(1.676199866843… 4.94
#> # ℹ 3,620 more rows
#> # ℹ 12 more variables: MaxTemperature <dbl>, Precipitation <dbl>,
#> # MeanRelativeHumidity <dbl>, MinRelativeHumidity <dbl>,
#> # MaxRelativeHumidity <dbl>, Radiation <dbl>, WindSpeed <dbl>,
#> # WindDirection <dbl>, PET <dbl>, elevation <dbl>, slope <dbl>, aspect <dbl>And now we can use it to build an interpolator object:
with_meteo(points_meteo_reference) |>
create_meteo_interpolator()
#> ℹ Checking meteorology object...
#> ✔ meteorology object ok
#> ℹ Creating interpolator...
#> Warning: No interpolation parameters provided, using defaults
#> ℹ Set the `params` argument to modify parameter default values
#> • Calculating smoothed variables...
#> • Updating intial_Rp parameter with the actual stations mean distance...
#> ✔ Interpolator created.
#> stars object with 2 dimensions and 13 attributes
#> attribute(s):
#> Min. 1st Qu. Median Mean
#> Temperature 3.47453413 11.4535839 13.940206 13.471788
#> MinTemperature -3.27100714 4.3824425 6.969667 5.860040
#> MaxTemperature 7.16265820 14.7798101 18.901768 18.420680
#> RelativeHumidity 35.08268335 56.7720829 64.269391 65.520741
#> Precipitation 0.00000000 0.0000000 0.000000 1.217155
#> Radiation 7.54372475 15.9837124 20.690395 19.595816
#> WindDirection 0.18621266 71.4293342 198.955399 181.615512
#> WindSpeed 0.02164994 0.9119314 1.252240 1.385929
#> elevation 240.00000000 370.0000000 447.000000 460.322314
#> aspect 5.19442749 74.7448807 174.369324 181.679232
#> slope 1.43209624 5.7204332 11.348120 13.073426
#> SmoothedPrecipitation 0.25336109 1.5179424 3.874980 3.862412
#> SmoothedTemperatureRange 9.41552220 11.8383211 12.587418 12.458785
#> 3rd Qu. Max. NA's
#> Temperature 16.1442989 20.781408 0
#> MinTemperature 8.2924049 11.071295 0
#> MaxTemperature 22.0316730 28.801911 0
#> RelativeHumidity 75.5124486 100.000000 0
#> Precipitation 0.2922964 21.000409 0
#> Radiation 23.7299248 27.982274 0
#> WindDirection 264.0404014 359.908196 1350
#> WindSpeed 1.7444530 5.811866 30
#> elevation 525.0000000 786.000000 0
#> aspect 291.0375061 360.000000 0
#> slope 19.7585106 31.071896 0
#> SmoothedPrecipitation 6.2789283 11.753856 505
#> SmoothedTemperatureRange 13.2623789 15.014616 0
#> dimension(s):
#> from to offset delta refsys point
#> date 1 30 2022-04-01 1 days Date FALSE
#> station 1 121 NA NA WGS 84 TRUE
#> values
#> date NULL
#> station POINT (1.6762 41.65133),...,POINT (1.781988 41.75712)