In lhmet-ped/fuse.prep: Preparation Of Input Data For Fuse

options(rmarkdown.html_vignette.check_title = FALSE)
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

Introdução

O pacote {fuse.prep} prepara os dados de elevação do terreno e de forçantes hidrometeorológicas necessários para criação dos arquivos NetCDF de entrada do FUSE para uma bacia hidrográfica (BH).

Nesta vinheta descreve-se o pré-processamento dos dados necessários para geração do arquivo NetCDF de bandas de elevação do terreno. Este processo é ilustrado na Figura \@ref(fig:fluxograma).

## Colar texto abaixo em
# https://mermaid-js.github.io/mermaid-live-editor
## e Download svg

## referência
# https://mermaid-js.github.io/mermaid/diagrams-and-syntax-and-examples/flowchart.html

#library(DiagrammeR)
#univar_flowchart = "
graph TD

 B1[(Dados Meteorológicos BR)] -->|brick-prec-25km-19800101-20151231.nc| B2[import_nc] 
 A1[(CON Hydrosheds)] -->|sa_con_3s_hydrosheds.gr*| A2

  subgraph fuse.prep

    B2 --> PB
    A2 --> |con_dem| Z[elev_bands_nc]

     subgraph HEgis
      C1[/Parte do nome do posto/] -->|name_regex| C2[info_station]
      C2 --> IP[/info_posto/]
      IP --> |posto| C3[extract_poly]
      %%C2 --> |station| C3[extract_poly]
      PS[/poly_posto/]     
      A2[extract_condem] 
     end

    %% poly_station unico
    C3 --> PS
    PS --> |poly_station|A2
    PS --> |poly_station|F[annual_climatology]
    PS --> |poly_station|C4[centroids]

    C4 --> |ccords| Z

    PB[/prec_brick/]

    PB --> |meteo_brick| F
    F --> |meteo_raster| Z

  end

style C4 stroke:#333,stroke-width:4px
style F stroke:#333,stroke-width:4px
style A2 stroke:#333,stroke-width:4px
style fuse.prep fill:#FFFFFF,stroke:#333,stroke-width:2px

#"
#mermaid(diagram = univar_flowchart, width = 800, height = 800)

``r. Os demais são feitos com o {fuse.prep}**.", out.width='95%'} knitr::include_graphics("mermaid-diagram-20201019171247.svg")

# Pré-requisitos

## Pacotes 

```r
library(raster)
library(HEgis)
library(lhmetools)
library(fuse.prep)
library(dplyr)

Base de dados

Para demonstrar o pré-processamento dos dados são usadas duas (grandes) base de dados:

Modelo de elevação do terreno hidrologicamente consistente (CON) da base de dados Hydrosheds para América do Sul com resolução espacial de aproximadamente 90 m.
Dados diários gradeados de variáveis meteorológicas para o Brasil (Xavier et al 2017, Xavier et al. 2016)

A primeira base de dados é distribuída em arquivos divididos em tiles que cobrem quase todo o globo. Estes arquivos foram combinados em único arquivo (15 GB) para Amérca do Sul no formato nativo do pacote raster (sa_con_3s_hydrosheds.grd) e pode ser baixado aqui [^1]. Este arquivo pode ser baixado manualmente no diretório onde está sendo realizada a preparação de dados para o FUSE. Recomenda-se armazenar os dados em um sub-diretório nomeado input.

A segunda base de dados é distribuída em arquivos NetCDF divididos por décadas. Os arquivos dos dados diários de precipitação e evapotranspiração de referência foram combinados em único arquivo NetCDF (1.5 GB para cada variável) e podem ser baixados aqui e aqui, respectivamente. O trecho de código abaixo utiliza a função import_nc que faz o download dos arquivos NetCDF de cada variável e retorna o caminho de acesso aos arquivos baixados.

[^1]: Certifique-se da disponibilidade de espaço para armazenamento dos dados em disco.

(nc_files <- meteo_nc(c("prec", "et0"), dest_dir = "input"))

prec 
"input/brick-prec-25km-19800101-20151231.nc" 
et0 
"input/brick-ETo-25km-19800101-20170731.nc"

Exemplo de dados pré-processados

O pré-processamento dos dados requeridos pela função que cria o arquivo NetCDF de bandas de elevação (elev_bands_nc()) consiste nos processos representados pelos retângulos com linha mais grossa na Figura \@ref(fig:fluxograma). Os dados de saída daquelas funções são o resultado do pré-processamento. Neste exemplo, os dados referem-se ao posto 74 do ONS, e incluem:

con_dem: raster do modelo de elevação do terreno recortado e com máscara para BH. Um exemplo de dados é o condem74 no pacote {HEgis} e resulta da função extract_condem().

condem74

ccoords: tibble ou data frame com as coordenadas do centróide do polígono da BH. Estas coordenadas resultam da função centroids() que requer o polígono da BH (p.ex: poly74 no pacote HEgis).

poly74

raster da climatologia anual de uma variável meteorológica (e.g. precipitação, temperatura, evapotranspiração potencial). Um exemplo de dados é o precclim74 no pacote {HEgis} e resulta da função annual_climatology().

precclim74

Pré-processamento

Agora veremos como obter cada um dos dados listados na seção anterior.

Centróide da BH (`ccoords`)

As coordenadas do centróide da BH usado de entrada da função elev_bands_nc(..., ccoords) são obtidas no código abaixo seguindo o fluxograma desta etapa (coluna central da Figura \@ref(fig:fluxograma)).

# parte do nome do posto da UHE
parte_nome <- "MUNHOZ" 
# polígono do posto
poly_posto <- parte_nome %>%
  # informações do posto
  info_station(name_regex = .) %>%
  # acessa valor da coluna posto
  pull(posto) %>%
  # extrai plígono do posto
  extract_poly(station = .)
# centróides do polígono do posto
(ccoords <- centroids(poly_station = poly_posto))

Elevação do terreno da BH (`con_dem`)

O raster de entrada da função elev_bands_nc(con_dem, ...) é o obtido no código abaixo seguindo o fluxograma desta etapa (primeira coluna da Figura \@ref(fig:fluxograma)).

# local onde foi salvo o 
path_hs <- "input/sa_con_3s_hydrosheds.grd"

path_hs <- "~/Dropbox/datasets/GIS/hydrosheds/sa_con_3s_hydrosheds.grd"

con_hs <- raster(path_hs)
con_dem <- extract_condem(condem = con_hs, poly_station = poly_posto)
con_dem

class      : RasterLayer 
dimensions : 1957, 3329, 6514853  (nrow, ncol, ncell)
resolution : 0.0008333333, 0.0008333333  (x, y)
extent     : -51.7225, -48.94833, -26.855, -25.22417  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : memory
names      : layer 
values     : 588, 1506  (min, max)

Climatologia da precipitação anual (`meteo_raster`)

O raster de entrada da função elev_bands_nc(..., meteo_raster) é o obtido no código abaixo seguindo o fluxograma desta etapa (última coluna da Figura \@ref(fig:fluxograma)).

# lê arquivo de prec diária ou faz download se ele não estiver no dest_dir
prec_brick <- import_nc(
  varnc = "prec",
  dest_dir = "input"
)
# climatologia anual na BH
meteo_raster <- annual_climatology(
  meteo_brick = prec_brick,
  poly_station = poly_posto, 
  cutoff = 0
)
meteo_raster

class      : RasterLayer 
dimensions : 6, 11, 66  (nrow, ncol, ncell)
resolution : 0.25, 0.25  (x, y)
extent     : -51.75, -49, -26.75, -25.25  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : memory
names      : layer 
values     : 1449.61, 2070.507  (min, max)

Criação do NetCDF de bandas de elevação

Com os dados pré-processados temos os principais dados de entrada da função elev_bands_nc(). Então podemos gerar o arquivo NetCDF de bandas de elevação para o posto 74, usando 14 bandas de elevação usando o código:

# altere aqui se desejar
suffix <- paste0("posto", poly_posto$codONS)
dest_dir <- "output"
elev_bands_ncname <- file.path(
  dest_dir, 
  paste(suffix, "elev_bands.nc", sep = "_")
)
elev_bands_ncname

suffix <- paste0("posto", poly_posto$codONS)
dest_dir <- tempdir()
elev_bands_ncname <- file.path(
  dest_dir, 
  paste(suffix, "elev_bands.nc", sep = "_")
)
print("output/posto74_elev_bands.nc")

elev_bands_file <- elev_bands_nc(
  con_dem,
  meteo_raster,
  dz = NULL,
  nbands = 14,
  ccoords,
  file_nc = elev_bands_ncname,
  na = -9999,
  force_v4 = TRUE,
  quiet = TRUE
)