,
Sergio Ibarra-Espinosa [aut] | Title: | Tools for Create Emissions for Air Quality Models |
|---|---|
| Description: | Processing tools to create emissions for use in numerical air quality models. Emissions can be calculated both using emission factors and activity data (Schuch et al 2018) <doi:10.21105/joss.00662> or using pollutant inventories (Schuch et al., 2018) <doi:10.30564/jasr.v1i1.347>. Functions to process individual point emissions, line emissions and area emissions of pollutants are available as well as methods to incorporate alternative data for Spatial distribution of emissions such as satellite images (Gavidia-Calderon et. al, 2018) <doi:10.1016/j.atmosenv.2018.09.026> or openstreetmap data (Andrade et al, 2015) <doi:10.3389/fenvs.2015.00009>. |
| Authors: | Daniel Schuch [aut, cre] ,
Sergio Ibarra-Espinosa [aut] |
| Maintainer: | Daniel Schuch <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.665.8.0 |
| Built: | 2024-11-22 04:28:59 UTC |
| Source: | https://github.com/atmoschem/emissv |
Calculate the spatial distribution by a raster masked by shape/model grid information.
areaSource(s, r, grid = NA, name = "", as_frac = FALSE, verbose = TRUE)areaSource(s, r, grid = NA, name = "", as_frac = FALSE, verbose = TRUE)
s |
input shape object |
r |
input raster object |
grid |
grid with the output format |
name |
area name |
as_frac |
return a fraction instead of the raster value |
verbose |
display additional data |
About the DMSP and example data https://en.wikipedia.org/wiki/Defense_Meteorological_Satellite_Program
a raster object containing the spatial distribution of emissions
Data avaliable https://www.nesdis.noaa.gov/current-satellite-missions/currently-flying/defense-meteorological-satellite-program
shape <- raster::shapefile(paste(system.file("extdata", package = "EmissV"), "/BR.shp",sep=""),verbose = FALSE) shape <- shape[22,1] # subset for Sao Paulo - BR raster <- raster::raster(paste(system.file("extdata", package = "EmissV"), "/dmsp.tiff",sep="")) grid <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d02",sep="")) SP <- areaSource(shape,raster,grid,name = "SPMA") raster::plot(SP,ylab="Lat",xlab="Lon", main ="Spatial Distribution by Lights for Sao Paulo - Brazil")shape <- raster::shapefile(paste(system.file("extdata", package = "EmissV"), "/BR.shp",sep=""),verbose = FALSE) shape <- shape[22,1] # subset for Sao Paulo - BR raster <- raster::raster(paste(system.file("extdata", package = "EmissV"), "/dmsp.tiff",sep="")) grid <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d02",sep="")) SP <- areaSource(shape,raster,grid,name = "SPMA") raster::plot(SP,ylab="Lat",xlab="Lon", main ="Spatial Distribution by Lights for Sao Paulo - Brazil")
Combine area sources and total emissions to model output
emission( inventory = NULL, grid, mm = 1, aerosol = FALSE, check = TRUE, total, pol, area, plot = FALSE, verbose = TRUE )emission( inventory = NULL, grid, mm = 1, aerosol = FALSE, check = TRUE, total, pol, area, plot = FALSE, verbose = TRUE )
inventory |
a inventory raster from read |
grid |
grid information |
mm |
pollutant molar mass |
aerosol |
TRUE for aerosols and FALSE (defoult) for gazes |
check |
TRUE (defoult) to check negative and NA values and replace it for zero |
total |
list of total emission |
pol |
pollutant name |
area |
list of area sources or matrix with a spatial distribution |
plot |
TRUE for plot the final emissions |
verbose |
display additional information |
matrix of emission
a vector of emissions in MOL / mk2 h for gases and ug / m2 s for aerosols.
if Inventory is provided, the firsts tree arguments are not be used by the function.
Is a good practice use the set_units(fe,your_unity), where fe is your emission factory and your_unity is usually g/km on your emission factory
the list of area must be in the same order as defined in vehicles and total emission.
just WRF-Chem is suported by now
fleet <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(fleet,EmissionFactors,pol = c("CO"),verbose = TRUE) grid <- gridInfo(paste0(system.file("extdata", package = "EmissV"),"/wrfinput_d01")) shape <- raster::shapefile(paste0(system.file("extdata", package = "EmissV"),"/BR.shp")) raster <- raster::raster(paste0(system.file("extdata", package = "EmissV"),"/dmsp.tiff")) SP <- areaSource(shape[22,1],raster,grid,name = "SP") RJ <- areaSource(shape[17,1],raster,grid,name = "RJ") e_CO <- emission(total = TOTAL, pol = "CO", area = list(SP = SP, RJ = RJ), grid = grid, mm = 28)fleet <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(fleet,EmissionFactors,pol = c("CO"),verbose = TRUE) grid <- gridInfo(paste0(system.file("extdata", package = "EmissV"),"/wrfinput_d01")) shape <- raster::shapefile(paste0(system.file("extdata", package = "EmissV"),"/BR.shp")) raster <- raster::raster(paste0(system.file("extdata", package = "EmissV"),"/dmsp.tiff")) SP <- areaSource(shape[22,1],raster,grid,name = "SP") RJ <- areaSource(shape[17,1],raster,grid,name = "RJ") e_CO <- emission(total = TOTAL, pol = "CO", area = list(SP = SP, RJ = RJ), grid = grid, mm = 28)
Return a data frame for emission for multiple pollutants.
emissionFactor( ef, poluttant = names(ef), vnames = NA, unit = "g/km", example = FALSE, verbose = TRUE )emissionFactor( ef, poluttant = names(ef), vnames = NA, unit = "g/km", example = FALSE, verbose = TRUE )
ef |
list with emission factors |
poluttant |
poluttant names |
vnames |
name of each vehicle categoy (optional) |
unit |
tring with unit from unit package, for default is "g/km" |
example |
TRUE to diaplay a simple example |
verbose |
display additional information |
a emission factor data frame
a emission factor data.frame for totalEmission function
EF <- emissionFactor(example = TRUE) # or the code for the same result EF <- emissionFactor(ef = list(CO = c(1.75,10.04,0.39,0.45,0.77,1.48,1.61,0.75), PM = c(0.0013,0.0,0.0010,0.0612,0.1052,0.1693,0.0,0.0)), vnames = c("Light Duty Vehicles Gasohol","Light Duty Vehicles Ethanol", "Light Duty Vehicles Flex","Diesel Trucks","Diesel Urban Busses", "Diesel Intercity Busses","Gasohol Motorcycles", "Flex Motorcycles"))EF <- emissionFactor(example = TRUE) # or the code for the same result EF <- emissionFactor(ef = list(CO = c(1.75,10.04,0.39,0.45,0.77,1.48,1.61,0.75), PM = c(0.0013,0.0,0.0010,0.0612,0.1052,0.1693,0.0,0.0)), vnames = c("Light Duty Vehicles Gasohol","Light Duty Vehicles Ethanol", "Light Duty Vehicles Flex","Diesel Trucks","Diesel Urban Busses", "Diesel Intercity Busses","Gasohol Motorcycles", "Flex Motorcycles"))
Return a list containing information of a regular grid / domain
gridInfo( file = file.choose(), z = FALSE, missing_time = "1984-03-10", verbose = TRUE )gridInfo( file = file.choose(), z = FALSE, missing_time = "1984-03-10", verbose = TRUE )
file |
file name/path to a wrfinput, wrfchemi or geog_em file |
z |
TRUE for read wrfinput vertical coordinades |
missing_time |
time if the variable Times is missing |
verbose |
display additional information |
a list with grid information from air quality model
just WRF-Chem is suported by now
grid_d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d01",sep="")) grid_d2 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d02",sep="")) grid_d3 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d03",sep="")) names(grid_d1) # for plot the shapes shape <- raster::shapefile(paste0(system.file("extdata", package = "EmissV"), "/BR.shp")) raster::plot(shape,xlim = c(-55,-40),ylim = c(-30,-15), main="3 nested domains") axis(1); axis(2); box(); grid() lines(grid_d1$boundary, col = "red") text(grid_d1$xlim[2],grid_d1$Ylim[1],"d1",pos=4, offset = 0.5) lines(grid_d2$boundary, col = "red") text(grid_d2$xlim[2],grid_d2$Ylim[1],"d2",pos=4, offset = 0.5) lines(grid_d3$boundary, col = "red") text(grid_d3$xlim[1],grid_d3$Ylim[2],"d3",pos=2, offset = 0.0)grid_d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d01",sep="")) grid_d2 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d02",sep="")) grid_d3 <- gridInfo(paste(system.file("extdata", package = "EmissV"), "/wrfinput_d03",sep="")) names(grid_d1) # for plot the shapes shape <- raster::shapefile(paste0(system.file("extdata", package = "EmissV"), "/BR.shp")) raster::plot(shape,xlim = c(-55,-40),ylim = c(-30,-15), main="3 nested domains") axis(1); axis(2); box(); grid() lines(grid_d1$boundary, col = "red") text(grid_d1$xlim[2],grid_d1$Ylim[1],"d1",pos=4, offset = 0.5) lines(grid_d2$boundary, col = "red") text(grid_d2$xlim[2],grid_d2$Ylim[1],"d2",pos=4, offset = 0.5) lines(grid_d3$boundary, col = "red") text(grid_d3$xlim[1],grid_d3$Ylim[2],"d3",pos=2, offset = 0.0)
Create a emission distribution from 'sp' or 'sf' spatial lines data.frame or spatial lines.
There 3 modes available to create the emission grid: - using gridInfo function output (defoult) - using the patch to "wrfinput" (output from real.exe) file or "geo" for (output from geog.exe) - "sf" (and "sp") uses a grid in SpatialPolygons format
The variable is the column of the data.frame with contains the variable to be used as emissions, by defoult the idstribution taken into acount the lench distribution of lines into each grid cell and the output is normalized.
lineSource( s, grid, as_raster = FALSE, type = "info", gcol = 100, grow = 100, variable = "length", verbose = TRUE )lineSource( s, grid, as_raster = FALSE, type = "info", gcol = 100, grow = 100, variable = "length", verbose = TRUE )
s |
SpatialLinesDataFrame of SpatialLines object |
grid |
grid object with the grid information or filename |
as_raster |
output format, TRUE for raster, FALSE for matrix |
type |
"info" (default), "wrfinput", "geo", "sp" or "sf" for grid type |
gcol |
grid points for a "sp" or "sf" type |
grow |
grid points for a "sp" or "sf" type |
variable |
variable to use, default is line length |
verbose |
display additional information |
a raster object containing the spatial distribution of emissions
OpenstreetMap data avaliable https://www.openstreetmap.org/ and https://download.geofabrik.de/
gridInfo and rasterSource
# loading a shapefile with osm data for sao paulo metropolitan area roads <- sf::st_read(paste0(system.file("extdata",package="EmissV"),"/streets.shp")) d3 <- gridInfo(paste0(system.file("extdata", package = "EmissV"),"/wrfinput_d03")) # calculating only for 2 streets roadLength <- lineSource(roads[1:2,],d3,as_raster=TRUE) # to generate a quick plot raster::plot(roadLength,ylab="Lat", xlab="Lon",main="Length of roads") # lines(road_lines)# loading a shapefile with osm data for sao paulo metropolitan area roads <- sf::st_read(paste0(system.file("extdata",package="EmissV"),"/streets.shp")) d3 <- gridInfo(paste0(system.file("extdata", package = "EmissV"),"/wrfinput_d03")) # calculating only for 2 streets roadLength <- lineSource(roads[1:2,],d3,as_raster=TRUE) # to generate a quick plot raster::plot(roadLength,ylab="Lat", xlab="Lon",main="Length of roads") # lines(road_lines)
Set of hourly profiles that represents the mean activity for each hour (local time) of the week.
Light Duty vehicles
Heavy Duty vehicles
passenger cars counted in June 2012
passenger cars counted in June 2013
passenger cars counted in June 2014
light comercial vehicles counted in June 2012
light comercial vehicles counted in June 2013
light comercial vehicles counted in June 2014
motorcycles counted in June 2012
motorcycles counted in June 2013
motorcycles counted in June 2014
Heavy good vehicles counted in June 2012
Heavy good vehicles counted in June 2013
Heavy good vehicles counted in June 2014
passenger cars counted in january 2012
passenger cars counted in january 2013
passenger cars counted in january 2014
light comercial vehicles counted in january 2012
light comercial vehicles counted in january 2013
light comercial vehicles counted in january 2014
Motorcycles counted in january 2012
Motorcycles counted in january 2014
Heavy good vehicles counted in january 2012
Heavy good vehicles counted in january 2013
Heavy good vehicles counted in january 2014
Power generation emission profile
Industrial emission profile
Residencial emission profile
Transport emission profile
Agriculture emission profile
Emission profile for ships
Solvent use emission constant profile
Waste burning emisssion constant profile
passenger cars at Janio Quadros on November 2018
heavy good vehicles at Janio Quadros on November 2018
total vehicle at Janio Quadros on November 2018
passenger cars at Anhanguera-Castello Branco on October 2018
Motorcycles cars at Anhanguera-Castello Branco on October 2018
data(perfil)data(perfil)
A list of data frames with activity by hour and weekday.
- Profiles 1 to 2 are from traffic count at São Paulo city from Perez Martínez et al (2014).
- Profiles 3 to 25 comes from traffic counted of toll stations located in São Paulo city, for summer and winters of 2012, 2013 and 2014.
- Profiles 26 to 33 are for different sectors from Oliver et al (2003).
- Profiles 34 to 36 are for volumetric mechanized traffic count at Janio Quadros tunnel on November 2018.
- Profiles 37 to 38 are for volumetric mechanized traffic count at Anhanguera-Castello Branco on October 2018.
The profile is normalized by days (but is balanced for a complete week) it means diary_emission x profile = hourly_emission.
Pérez-Martínez, P. J., Miranda, R. M., Nogueira, T., Guardani, M. L., Fornaro, A., Ynoue, R., & Andrade, M. F. (2014). Emission factors of air pollutants from vehicles measured inside road tunnels in São Paulo: case study comparison. International Journal of Environmental Science and Technology, 11(8), 2155-2168.
Olivier, J., J. Peters, C. Granier, G. Pétron, J.F. Müller, and S. Wallens, Present and future surface emissions of atmospheric compounds, POET Report #2, EU project EVK2-1999-00011, 2003.
# load the data data(perfil) # function to simple view plot.perfil <- function(per = perfil$LDV, text="", color = "#0000FFBB"){ plot(per[,1],ty = "l", ylim = range(per),axe = FALSE, xlab = "hour",ylab = "Intensity",main = text,col=color) for(i in 2:7){ lines(per[,i],col = color) } for(i in 1:7){ points(per[,i],col = "black", pch = 20) } axis(1,at=0.5+c(0,6,12,18,24),labels = c("00:00","06:00","12:00","18:00","00:00")) axis(2) box() } # view all profiles in perfil data for(i in 1:length(names(perfil))){ cat(paste("profile",i,names(perfil)[i],"\n")) plot.perfil(perfil[[i]],names(perfil)[i]) }# load the data data(perfil) # function to simple view plot.perfil <- function(per = perfil$LDV, text="", color = "#0000FFBB"){ plot(per[,1],ty = "l", ylim = range(per),axe = FALSE, xlab = "hour",ylab = "Intensity",main = text,col=color) for(i in 2:7){ lines(per[,i],col = color) } for(i in 1:7){ points(per[,i],col = "black", pch = 20) } axis(1,at=0.5+c(0,6,12,18,24),labels = c("00:00","06:00","12:00","18:00","00:00")) axis(2) box() } # view all profiles in perfil data for(i in 1:length(names(perfil))){ cat(paste("profile",i,names(perfil)[i],"\n")) plot.perfil(perfil[[i]],names(perfil)[i]) }
Calculate the maximum height of rise based on Brigs (1975), the height is calculated using different formulations depending on stability and wind conditions.
plumeRise(df, imax = 10, ermax = 1/100, Hmax = TRUE, verbose = TRUE)plumeRise(df, imax = 10, ermax = 1/100, Hmax = TRUE, verbose = TRUE)
df |
data.frame with micrometeorological and emission data |
imax |
maximum number of iteractions |
ermax |
maximum error |
Hmax |
use weil limit for plume rise, see details |
verbose |
display additional information |
data.frame with the input, rise (m) and effective higt (m)
The input data.frame must contains the folloging colluns:
- z: height of the emission (m)
- r: source raius (m)
- Ve: emission velocity (m/s)
- Te: emission temperature (K)
- ws: wind speed (m/s)
- Temp: ambient temperature (K)
- h: height of the Atmospheric Boundary Layer-ABL (m)
- L: Monin-Obuhkov Lench (m)
- dtdz: lapse ration of potential temperature, used only for stable ABL (K/m)
- Ustar: atriction velocity, used only for neutral ABL (m/s)
- Wstar: scale of convectie velocity, used only for convective ABL (m/s)
Addcitionaly some combination of ws, Wstar and Ustar can produce inacurate results, Weil (1979) propose a geometric limit of 0.62 * (h - Hs) for the rise value.
a data.frame with effective height of emissions for pointSource function
The plume rise formulas are from Brigs (1975):"Brigs, G. A. Plume rise predictions, Lectures on Air Pollution and Environmental Impact Analyses. Amer. Meteor. Soc. p. 59-111, 1975." and Arya 1999: "Arya, S.P., 1999, Air Pollution Meteorology and Dispersion, Oxford University Press, New York, 310 p."
The limits are from Weil (1979): "WEIL, J.C. Assessmet of plume rise and dispersion models using LIDAR data, PPSP-MP-24. Prepared by Environmental Center, Martin Marietta Corporation, for Maryland Department of Natural Resources. 1979."
The example is data from a chimney of the Candiota thermoelectric powerplant from Arabage et al (2006) "Arabage, M. C.; Degrazia, G. A.; Moraes O. L. Simulação euleriana da dispersão local da pluma de poluente atmosférico de Candiota-RS. Revista Brasileira de Meteorologia, v.21, n.2, p. 153-160, 2006."
candiota <- matrix(c(150,1,20,420,3.11,273.15 + 3.16,200,-34.86,3.11,0.33, 150,1,20,420,3.81,273.15 + 4.69,300,-34.83,3.81,0.40, 150,1,20,420,3.23,273.15 + 5.53,400,-24.43,3.23,0.48, 150,1,20,420,3.47,273.15 + 6.41,500,-15.15,3.48,0.52, 150,1,20,420,3.37,273.15 + 6.35,600, -8.85,3.37,2.30, 150,1,20,420,3.69,273.15 + 5.93,800,-10.08,3.69,2.80, 150,1,20,420,3.59,273.15 + 6.08,800, -7.23,3.49,1.57, 150,1,20,420,4.14,273.15 + 6.53,900,-28.12,4.14,0.97), ncol = 10, byrow = TRUE) candiota <- data.frame(candiota) names(candiota) <- c("z","r","Ve","Te","ws","Temp","h","L","Ustar","Wstar") row.names(candiota) <- c("08:00","09:00",paste(10:15,":00",sep="")) candiota <- plumeRise(candiota,Hmax = TRUE) print(candiota)candiota <- matrix(c(150,1,20,420,3.11,273.15 + 3.16,200,-34.86,3.11,0.33, 150,1,20,420,3.81,273.15 + 4.69,300,-34.83,3.81,0.40, 150,1,20,420,3.23,273.15 + 5.53,400,-24.43,3.23,0.48, 150,1,20,420,3.47,273.15 + 6.41,500,-15.15,3.48,0.52, 150,1,20,420,3.37,273.15 + 6.35,600, -8.85,3.37,2.30, 150,1,20,420,3.69,273.15 + 5.93,800,-10.08,3.69,2.80, 150,1,20,420,3.59,273.15 + 6.08,800, -7.23,3.49,1.57, 150,1,20,420,4.14,273.15 + 6.53,900,-28.12,4.14,0.97), ncol = 10, byrow = TRUE) candiota <- data.frame(candiota) names(candiota) <- c("z","r","Ve","Te","ws","Temp","h","L","Ustar","Wstar") row.names(candiota) <- c("08:00","09:00",paste(10:15,":00",sep="")) candiota <- plumeRise(candiota,Hmax = TRUE) print(candiota)
Transform a set of points into a grinded output
pointSource(emissions, grid, verbose = TRUE)pointSource(emissions, grid, verbose = TRUE)
emissions |
list of points |
grid |
grid object with the grid information |
verbose |
display additional information |
a raster
gridInfo and rasterSource
d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) p = data.frame(lat = c(-22,-22,-23.5), lon = c(-46,-48,-47 ), z = c(0 , 0, 0 ), emission = c(666,444,111 ) ) p_emissions <- pointSource(emissions = p, grid = d1) raster::plot(p_emissions,ylab="Lat", xlab="Lon", main = "3 point sources for domain d1")d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) p = data.frame(lat = c(-22,-22,-23.5), lon = c(-46,-48,-47 ), z = c(0 , 0, 0 ), emission = c(666,444,111 ) ) p_emissions <- pointSource(emissions = p, grid = d1) raster::plot(p_emissions,ylab="Lat", xlab="Lon", main = "3 point sources for domain d1")
Calculate the spatial distribution by a raster
rasterSource(r, grid, nlevels = "all", conservative = TRUE, verbose = TRUE)rasterSource(r, grid, nlevels = "all", conservative = TRUE, verbose = TRUE)
r |
input raster object |
grid |
grid object with the grid information |
nlevels |
number of vertical levels off the emission array |
conservative |
TRUE (default) to conserve total mass, FALSE to conserve flux |
verbose |
display additional information |
About the DMSP and example data https://en.wikipedia.org/wiki/Defense_Meteorological_Satellite_Program
Returns a matrix
Exemple data is a low resolution cutting from image of persistent lights of the Defense Meteorological Satellite Program (DMSP) https://pt.wikipedia.org/wiki/Defense_Meteorological_Satellite_Program
Data avaliable https://www.nesdis.noaa.gov/current-satellite-missions/currently-flying/defense-meteorological-satellite-program
gridInfo and lineSource
grid <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) x <- raster::raster(paste(system.file("extdata", package = "EmissV"),"/dmsp.tiff",sep="")) test <- rasterSource(x, grid) image(test, axe = FALSE, main = "Spatial distribution by Persistent Nocturnal Lights from DMSP")grid <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) x <- raster::raster(paste(system.file("extdata", package = "EmissV"),"/dmsp.tiff",sep="")) test <- rasterSource(x, grid) image(test, axe = FALSE, main = "Spatial distribution by Persistent Nocturnal Lights from DMSP")
Read data from global inventories. Several files can be read to produce one emission output and/or can be splitted into several species
read( file = file.choose(), version = NA, coef = rep(1, length(file)), spec = NULL, year = 1, month = 1, hour = 1, categories, reproject = TRUE, as_raster = TRUE, skip_missing = FALSE, verbose = TRUE )read( file = file.choose(), version = NA, coef = rep(1, length(file)), spec = NULL, year = 1, month = 1, hour = 1, categories, reproject = TRUE, as_raster = TRUE, skip_missing = FALSE, verbose = TRUE )
file |
file name or names (variables are summed) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
version |
Character; One of of the following:
|
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
coef |
coefficients to merge different sources (file) into one emission |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
spec |
numeric speciation vector to split emission into different species |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
year |
scenario index (only for GAINS and VULCAN-y) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
month |
the desired month of the inventory (MACCITY and ODIAC) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
hour |
hour of the emission (only for ACES and VULCAN-h) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
categories |
considered categories (for MACCITY/GAINS variable names), empty for use all |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
reproject |
to project the output to "+proj=longlat" needed for emission function (only for VULCAN and ACES) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
as_raster |
return a raster (default) or matrix (with units) |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
skip_missing |
return a zero emission and a warning for missing files/variables |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
verbose |
display additional information |
Matrix or raster
for EDGAR (all versions), GAINS, RCP and MACCTITY, please use flux (kg m-2 s-1) NetCDF file.
Read abbout EDGAR at http://edgar.jrc.ec.europa.eu and MACCITY at http://accent.aero.jussieu.fr/MACC_metadata.php
More info for EDGARv8.1 https://edgar.jrc.ec.europa.eu/dataset_ap81 for short live species and https://edgar.jrc.ec.europa.eu/dataset_ghg80 for GHG
Janssens-Maenhout, G., Dentener, F., Van Aardenne, J., Monni, S., Pagliari, V., Orlandini, L., ... & Wankmüller, R. (2012). EDGAR-HTAP: a harmonized gridded air pollution emission dataset based on national inventories. European Commission Joint Research Centre Institute for Environment and Sustainability. JRC 68434 UR 25229 EUR 25229, ISBN 978-92-79-23123-0.
Lamarque, J.-F., Bond, T. C., Eyring, V., Granier, C., Heil, A., Klimont, Z., Lee, D., Liousse, C., Mieville, A., Owen, B., Schultz, M. G., Shindell, D., Smith, S. J., Stehfest, E., Van Aardenne, J., Cooper, O. R., Kainuma, M., Mahowald, N., McConnell, J. R., Naik, V., Riahi, K., and van Vuuren, D. P.: Historical (1850-2000) gridded anthropogenic and biomass burning emissions of reactive gases and aerosols: methodology and application, Atmos. Chem. Phys., 10, 7017-7039, doi:10.5194/acp-10-7017-2010, 2010.
Z Klimont, S. J. Smith and J Cofala The last decade of global anthropogenic sulfur dioxide: 2000–2011 emissions Environmental Research Letters 8, 014003, 2013
Gurney, Kevin R., Jianming Liang, Risa Patarasuk, Yang Song, Jianhua Huang, and Geoffrey Roest (2019) The Vulcan Version 3.0 High-Resolution Fossil Fuel CO2 Emissions for the United States. Nature Scientific Data.
rasterSource and gridInfo
dir.create(file.path(tempdir(), "EDGARv432")) folder <- setwd(file.path(tempdir(), "EDGARv432")) url <- "http://jeodpp.jrc.ec.europa.eu/ftp/jrc-opendata/EDGAR/datasets/v432_AP/NOx" file <- 'v432_NOx_2012.0.1x0.1.zip' download.file(paste0(url,'/TOTALS/',file), file) unzip('v432_NOx_2012.0.1x0.1.zip') nox <- read(file = dir(pattern = '.nc'), version = 'EDGAR', spec = c(E_NO = 0.9 , # 90% of NOx is NO E_NO2 = 0.1 )) # 10% of NOx is NO2 setwd(folder) # creating a color scale cor <- colorRampPalette(colors = c(c("#031057", "#0522FC", "#7E0AFA", "#EF0AFF", "#FFA530", "#FFF957"))) raster::plot(nox$E_NO,xlab="Lat", ylab="Lon", col = cor(12),zlim = c(-6.5e-7,1.4e-5), main="NO emissions from EDGAR (in g / m2 s)") d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) NO <- emission(grid = d1, inventory = nox$E_NO, pol = "NO", mm = 30.01, plot = TRUE)dir.create(file.path(tempdir(), "EDGARv432")) folder <- setwd(file.path(tempdir(), "EDGARv432")) url <- "http://jeodpp.jrc.ec.europa.eu/ftp/jrc-opendata/EDGAR/datasets/v432_AP/NOx" file <- 'v432_NOx_2012.0.1x0.1.zip' download.file(paste0(url,'/TOTALS/',file), file) unzip('v432_NOx_2012.0.1x0.1.zip') nox <- read(file = dir(pattern = '.nc'), version = 'EDGAR', spec = c(E_NO = 0.9 , # 90% of NOx is NO E_NO2 = 0.1 )) # 10% of NOx is NO2 setwd(folder) # creating a color scale cor <- colorRampPalette(colors = c(c("#031057", "#0522FC", "#7E0AFA", "#EF0AFF", "#FFA530", "#FFF957"))) raster::plot(nox$E_NO,xlab="Lat", ylab="Lon", col = cor(12),zlim = c(-6.5e-7,1.4e-5), main="NO emissions from EDGAR (in g / m2 s)") d1 <- gridInfo(paste(system.file("extdata", package = "EmissV"),"/wrfinput_d01",sep="")) NO <- emission(grid = d1, inventory = nox$E_NO, pol = "NO", mm = 30.01, plot = TRUE)
Distribute the total mass of estimated emissions into model species.
speciation(total, spec = NULL, verbose = TRUE)speciation(total, spec = NULL, verbose = TRUE)
total |
emissions from totalEmissions |
spec |
numeric speciation vector of species |
verbose |
display additional information |
Return a list with the daily total emission by interest area (cityes, states, countries, etc).
veic <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(veic,EmissionFactors,pol = "PM") pm_iag <- c(E_PM25I = 0.0509200, E_PM25J = 0.1527600, E_ECI = 0.1196620, E_ECJ = 0.0076380, E_ORGI = 0.0534660, E_ORGJ = 0.2279340, E_SO4I = 0.0063784, E_SO4J = 0.0405216, E_NO3J = 0.0024656, E_NO3I = 0.0082544, E_PM10 = 0.3300000) PM <- speciation(TOTAL,pm_iag)veic <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(veic,EmissionFactors,pol = "PM") pm_iag <- c(E_PM25I = 0.0509200, E_PM25J = 0.1527600, E_ECI = 0.1196620, E_ECJ = 0.0076380, E_ORGI = 0.0534660, E_ORGJ = 0.2279340, E_SO4I = 0.0063784, E_SO4J = 0.0405216, E_NO3J = 0.0024656, E_NO3I = 0.0082544, E_PM10 = 0.3300000) PM <- speciation(TOTAL,pm_iag)
Set of tables for speciation:
Volatile organic compounds for RADM2
Volatile organic compounds for CBMZ
Volatile organic compounds for MOZART
volatile organic compounds for SAPRC99
Vehicular volatile organic compounds for RADM2 (MOL)
Vehicular volatile organic compounds for CBMZ (MOL)
Vehicular volatile organic compounds for MOZART (MOL)
Vehicular volatile organic compounds for SAPRC99 (MOL)
Particulate matter for made/sorgan
Fine particulate matter for made/sorgan
Nox split Perez Martínez et al (2014)
Nox split usin Ntziachristos and Zamaras (2016)
Volatile organic compounds species from EDGAR 4.3.2 for RADM2 (MOL)
Volatile organic compounds species from EDGAR 4.3.2 for MOZART (MOL)
- Volatile organic compounds species map from 1 to 4 are from Li et al (2014) taken into account several sources of pollutants.
- Volatile organic compounds from vehicular activity species map 5 to 8 is a by fuel and emission process from USP-IAG tunel experiments (Rafee et al., 2017) emited by the process of exhaust (through the exhaust pipe), liquid (carter and evaporative) and vapor (fuel transfer operations).
- Particulate matter speciesmap for made/sorgan emissions 9 and 10.
- Nox split using Perez Martínez et al (2014) data (11).
- Nox split using mean of Ntziachristos and Zamaras (2016) data (12).
- Volatile organic compounds species map 13 and 14 are the corespondence from EDGAR 4.3.2 VOC specialization to RADM2 and MOZART.
data(species)data(species)
List of numeric vectors with the 'names()' of the species and the values of each species.
iag-voc: After estimating all the emissions of NMHC, it was used the speciation presented in (RAFEE et al., 2017). This speciation is based on tunnel measurements in São Paulo, depends on the type of fuel (E25, E100 and B5) and provides the mass of each chemical compound as mol/g. This speciation splits the NMHC from evaporative, liquid and exhaust emissions of E25, E100 and B5, into minimum compounds required for the Carbon Bond Mechanism (CBMZ) (ZAVERI; PETERS, 1999). Atmospheric simulations using the same pollutants in Brazil have resulted in good agreement with observations (ANDRADE et al., 2015).
iag-pm: data tunnel experiments at São Paulo in Perez Martínez et al (2014)
iag-nox: common NOx split for São Paulo Metropolitan area.
bcom-nox: mean of Ntziachristos and Zamaras (2016) data.
mic: from Li et al (2014).
edgar: Edgar 4.3.2 emissions Crippa et al. (2018).
The units are mass ratio (mass/mass) or MOL (MOL), this last case do not change the default 'mm' into 'emission()' function.
Li, M., Zhang, Q., Streets, D. G., He, K. B., Cheng, Y. F., Emmons, L. K., ... & Su, H. (2014). Mapping Asian anthropogenic emissions of non-methane volatile organic compounds to multiple chemical mechanisms. Atmos. Chem. Phys, 14(11), 5617-5638.
Huang, G., Brook, R., Crippa, M., Janssens-Maenhout, G., Schieberle, C., Dore, C., ... & Friedrich, R. (2017). Speciation of anthropogenic emissions of non-methane volatile organic compounds: a global gridded data set for 1970–2012. Atmospheric Chemistry and Physics, 17(12), 7683.
Abou Rafee, S. A., Martins, L. D., Kawashima, A. B., Almeida, D. S., Morais, M. V. B., Souza, R. V. A., Oliveira, M. B. L., Souza, R. A. F., Medeiros, A. S. S., Urbina, V., Freitas, E. D., Martin, S. T., and Martins, J. A.: Contributions of mobile, stationary and biogenic sources to air pollution in the Amazon rainforest: a numerical study with the WRF-Chem model, Atmos. Chem. Phys., 17, 7977-7995, https://doi.org/10.5194/acp-17-7977-2017, 2017.
Martins, L. D., Andrade, M. F. D., Freitas, E., Pretto, A., Gatti, L. V., Junior, O. M. A., et al. (2006). Emission factors for gas-powered vehicles traveling through road tunnels in Sao Paulo, Brazil. Environ. Sci. Technol. 40, 6722–6729. doi: 10.1021/es052441u
Pérez-Martínez, P. J., Miranda, R. M., Nogueira, T., Guardani, M. L., Fornaro, A., Ynoue, R., & Andrade, M. F. (2014). Emission factors of air pollutants from vehicles measured inside road tunnels in São Paulo: case study comparison. International Journal of Environmental Science and Technology, 11(8), 2155-2168.
ANDRADE, M. d. F. et al. Air quality forecasting system for southeastern brazil. Frontiers in Environmental Science, Frontiers, v. 3, p. 1–12, 2015.
Crippa, M., Guizzardi, D., Muntean, M., Schaaf, E., Dentener, F., Aardenne, J. A. V., ... & Janssens-Maenhout, G. (2018). Gridded emissions of air pollutants for the period 1970–2012 within EDGAR v4.3.2. Earth System Science Data, 10(4), 1987-2013.
speciation and read
# load the mapping tables data(species) # names of eath mapping tables for(i in 1:length(names(species))) cat(paste0("specie map ",i," ",names(species)[i],"\n")) # names of species contained in the (first) mapping table names(species[[1]]) # The first mapping table / species and values species[1]# load the mapping tables data(species) # names of eath mapping tables for(i in 1:length(names(species))) cat(paste0("specie map ",i," ",names(species)[i],"\n")) # names of species contained in the (first) mapping table names(species[[1]]) # The first mapping table / species and values species[1]
Caculate the total emission with:
Emission(pollutant) = sum( Vehicles(n) * Km_day_use(n) * Emission_Factor(n,pollutant) )
where n is the type of the veicle
totalEmission(v, ef, pol, verbose = TRUE)totalEmission(v, ef, pol, verbose = TRUE)
v |
dataframe with the vehicle data |
ef |
emission factor |
pol |
pollutant name in ef |
verbose |
display additional information |
Return a list with the daily total emission by interest area (cityes, states, countries, etc).
the units (set_units("value",unit) where the recomended unit is g/d) must be used to make the ef data.frame
rasterSource, lineSource and emission
veic <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(veic,EmissionFactors,pol = c("CO","PM"))veic <- vehicles(example = TRUE) EmissionFactors <- emissionFactor(example = TRUE) TOTAL <- totalEmission(veic,EmissionFactors,pol = c("CO","PM"))
Return a data frame with 4 columns (vehicle category, type, fuel and avarage kilometers driven) and an aditional column with the number of vehicles for each interest area (cityes, states, countries, etc).
Average daily kilometres driven are defined by vehicle type:
- LDV (Light duty Vehicles) 41 km / day
- TRUCKS (Trucks) 110 km / day
- BUS (Busses) 165 km / day
- MOTO (motorcycles and other vehicles) 140 km / day
The number of vehicles are defined by the distribution of vehicles by vehicle classs and the total number of vehicles by area.
vehicles( total_v, area_name = names(total_v), distribution, type, category = NA, fuel = NA, vnames = NA, example = FALSE, verbose = TRUE )vehicles( total_v, area_name = names(total_v), distribution, type, category = NA, fuel = NA, vnames = NA, example = FALSE, verbose = TRUE )
total_v |
total of vehicles by area (area length) |
area_name |
area names (area length) |
distribution |
distribution of vehicles by vehicle class |
type |
type of vehicle by vehicle class (distribution length) |
category |
category name (distribution length / NA) |
fuel |
fuel type by vehicle class (distribution length / NA) |
vnames |
name of each vehicle class (distribution length / NA) |
example |
a simple example |
verbose |
display additional information |
a fleet distribution data.frame for totalEmission function
total_v and area_name must have the same length.
distribution, type, category (if used), fuel (if used) and vnames (if used) must have the same length.
fleet <- vehicles(example = TRUE) # or the code bellow for the same result # DETRAN 2016 data for total number of vehicles for 5 Brazilian states (Sao Paulo, # Rio de Janeiro, Minas Gerais, Parana and Santa Catarina) # vahicle distribution of Sao Paulo fleet <- vehicles(total_v = c(27332101, 6377484, 10277988, 7140439, 4772160), area_name = c("SP", "RJ", "MG", "PR", "SC"), distribution = c( 0.4253, 0.0320, 0.3602, 0.0260, 0.0290, 0.0008, 0.1181, 0.0086), category = c("LDV_E25","LDV_E100","LDV_F","TRUCKS_B5", "CBUS_B5","MBUS_B5","MOTO_E25","MOTO_F"), type = c("LDV", "LDV", "LDV","TRUCKS", "BUS","BUS","MOTO", "MOTO"), fuel = c("E25", "E100", "FLEX","B5", "B5","B5","E25", "FLEX"), vnames = c("Light duty Vehicles Gasohol","Light Duty Vehicles Ethanol", "Light Duty Vehicles Flex","Diesel trucks","Diesel urban busses", "Diesel intercity busses","Gasohol motorcycles", "Flex motorcycles"))fleet <- vehicles(example = TRUE) # or the code bellow for the same result # DETRAN 2016 data for total number of vehicles for 5 Brazilian states (Sao Paulo, # Rio de Janeiro, Minas Gerais, Parana and Santa Catarina) # vahicle distribution of Sao Paulo fleet <- vehicles(total_v = c(27332101, 6377484, 10277988, 7140439, 4772160), area_name = c("SP", "RJ", "MG", "PR", "SC"), distribution = c( 0.4253, 0.0320, 0.3602, 0.0260, 0.0290, 0.0008, 0.1181, 0.0086), category = c("LDV_E25","LDV_E100","LDV_F","TRUCKS_B5", "CBUS_B5","MBUS_B5","MOTO_E25","MOTO_F"), type = c("LDV", "LDV", "LDV","TRUCKS", "BUS","BUS","MOTO", "MOTO"), fuel = c("E25", "E100", "FLEX","B5", "B5","B5","E25", "FLEX"), vnames = c("Light duty Vehicles Gasohol","Light Duty Vehicles Ethanol", "Light Duty Vehicles Flex","Diesel trucks","Diesel urban busses", "Diesel intercity busses","Gasohol motorcycles", "Flex motorcycles"))