diff --git a/hermesv3_gr/modules/emision_inventories/emission_inventory.py b/hermesv3_gr/modules/emision_inventories/emission_inventory.py
index 0c3df52470998c04e3c611ed282cff29992e0426..3101a2fccdef07220c15f2b448d56b2971c302d3 100644
--- a/hermesv3_gr/modules/emision_inventories/emission_inventory.py
+++ b/hermesv3_gr/modules/emision_inventories/emission_inventory.py
@@ -106,7 +106,6 @@ class EmissionInventory(object):
world_mask_file=os.path.join(os.path.dirname(options.auxiliar_files_path),
'{0}_WorldMask.nc'.format(inventory_name)))
- self.input_pollutants = pollutants
self.pollutant_dicts = self.create_pollutants_dicts(pollutants)
self.masking.check_regrid_mask(self.pollutant_dicts[0]['path'])
@@ -208,12 +207,8 @@ class EmissionInventory(object):
if self.source_type == 'area':
extension = 'nc'
-
elif self.source_type == 'point':
- if self.inventory_name == 'GFASv12':
- extension = 'nc'
- else:
- extension = 'csv'
+ extension = 'csv'
else:
settings.write_log('ERROR: Check the .err file to get more info.')
if settings.rank == 0:
@@ -232,9 +227,8 @@ class EmissionInventory(object):
elif self.input_frequency == 'monthly':
file_name = "{0}_{1}{2}.{3}".format(pollutant, self.reference_year, self.date.strftime("%m"), extension)
elif self.input_frequency == 'daily':
- print pollutant, self.reference_year, self.date.strftime("%m"), self.date.strftime("%d"), extension
file_name = "{0}_{1}{2}{3}.{4}".format(pollutant, self.reference_year, self.date.strftime("%m"),
- self.date.strftime("%d"), extension)
+ self.date.strftime("%D"), extension)
else:
settings.write_log('ERROR: Check the .err file to get more info.')
if settings.rank == 0:
@@ -287,7 +281,6 @@ class EmissionInventory(object):
"""
import pandas as pd
import re
- from point_gfas_emission_inventory import PointGfasEmissionInventory
from gfas_emission_inventory import GfasEmissionInventory
from point_source_emission_inventory import PointSourceEmissionInventory
@@ -352,30 +345,19 @@ class EmissionInventory(object):
p_hour=emission_inventory.p_hour,
p_speciation=emission_inventory.p_speciation))
elif emission_inventory.source_type == 'point':
- if emission_inventory.ei == 'GFASv12':
- emission_inventory_list.append(
- PointGfasEmissionInventory(
- options, grid, date, emission_inventory.ei, emission_inventory.source_type,
- emission_inventory.sector, pollutants, emission_inventory_path,
- emission_inventory.frequency, vertical_output_profile,
- reference_year=emission_inventory.ref_year, factors=emission_inventory.factor_mask,
- regrid_mask=emission_inventory.regrid_mask, p_vertical=emission_inventory.p_vertical,
- p_month=p_month, p_day=emission_inventory.p_day, p_hour=emission_inventory.p_hour,
- p_speciation=emission_inventory.p_speciation))
- else:
- emission_inventory_list.append(
- PointSourceEmissionInventory(options, grid, date, emission_inventory.ei,
- emission_inventory.source_type, emission_inventory.sector,
- pollutants, emission_inventory_path,
- emission_inventory.frequency, vertical_output_profile,
- reference_year=emission_inventory.ref_year,
- factors=emission_inventory.factor_mask,
- regrid_mask=emission_inventory.regrid_mask,
- p_vertical=emission_inventory.p_vertical,
- p_month=p_month,
- p_day=emission_inventory.p_day,
- p_hour=emission_inventory.p_hour,
- p_speciation=emission_inventory.p_speciation))
+ emission_inventory_list.append(
+ PointSourceEmissionInventory(options, grid, date, emission_inventory.ei,
+ emission_inventory.source_type, emission_inventory.sector, pollutants,
+ emission_inventory_path,
+ emission_inventory.frequency, vertical_output_profile,
+ reference_year=emission_inventory.ref_year,
+ factors=emission_inventory.factor_mask,
+ regrid_mask=emission_inventory.regrid_mask,
+ p_vertical=emission_inventory.p_vertical,
+ p_month=p_month,
+ p_day=emission_inventory.p_day,
+ p_hour=emission_inventory.p_hour,
+ p_speciation=emission_inventory.p_speciation))
else:
settings.write_log('ERROR: Check the .err file to get more info.')
if settings.rank == 0:
diff --git a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py b/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py
deleted file mode 100755
index 35663df73328033b9d423ba3195146e731750c49..0000000000000000000000000000000000000000
--- a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py
+++ /dev/null
@@ -1,428 +0,0 @@
-#!/usr/bin/env python
-
-# Copyright 2018 Earth Sciences Department, BSC-CNS
-#
-# This file is part of HERMESv3_GR.
-#
-# HERMESv3_GR is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# HERMESv3_GR is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with HERMESv3_GR. If not, see .
-
-import sys
-import os
-import timeit
-import warnings
-
-import hermesv3_gr.config.settings as settings
-from emission_inventory import EmissionInventory
-
-
-class PointGfasEmissionInventory(EmissionInventory):
- """
- Class that defines the content and the methodology for the GFAS emission inventories
-
- :param current_date: Date of the day to simulate.
- :type current_date: datetime.datetime
-
- :param inventory_name: Name of the inventory to use.
- :type inventory_name: str
-
- :param sector: Name of the sector of the inventory to use.
- :type sector: str
-
- :param pollutants: List of the pollutant name to take into account.
- :type pollutants: list of str
-
- :param frequency: Frequency of the inputs. [yearly, monthly, daily]
- :type frequency: str
-
- :param reference_year: year of reference of the information of the dataset.
- :type reference_year: int
-
- :param factors: NOT IMPLEMENTED YET
- :type factors: NOT IMPLEMENTED YET
-
- :param regrid_mask: NOT IMPLEMENTED YET
- :type regrid_mask: NOT IMPLEMENTED YET
-
- :param p_vertical: ID of the vertical profile to use.
- :type p_vertical: str
-
- :param p_month: ID of the temporal monthly profile to use.
- :type p_month: str
-
- :param p_day: ID of the temporal daily profile to use.
- :type p_day: str
-
- :param p_hour: ID of the temporal hourly profile to use.
- :type p_hour: str
-
- :param p_speciation: ID of the speciation profile to use.
- :type p_speciation: str
- """
-
- def __init__(self, options, grid, current_date, inventory_name, source_type, sector, pollutants, inputs_path,
- frequency, vertical_output_profile,
- reference_year=2010, factors=None, regrid_mask=None, p_vertical=None, p_month=None, p_day=None,
- p_hour=None, p_speciation=None):
- from hermesv3_gr.modules.vertical.vertical_gfas import GfasVerticalDistribution
-
- st_time = timeit.default_timer()
- settings.write_log('\t\tCreating GFAS emission inventory as point source.', level=3)
- super(PointGfasEmissionInventory, self).__init__(
- options, grid, current_date, inventory_name, source_type, sector, pollutants, inputs_path, frequency,
- vertical_output_profile,
- reference_year=reference_year, factors=factors, regrid_mask=regrid_mask, p_vertical=None,
- p_month=p_month, p_day=p_day, p_hour=p_hour, p_speciation=p_speciation)
-
- # self.approach = self.get_approach(p_vertical)
- self.method = self.get_method(p_vertical)
-
- # self.altitude = self.get_altitude()
-
- self.vertical = GfasVerticalDistribution(vertical_output_profile, self.get_approach(p_vertical),
- self.get_altitude())
-
- settings.write_time('PointGfasEmissionInventory', 'Init', timeit.default_timer() - st_time, level=3)
-
- def __str__(self):
- string = "PointGfasEmissionInventory:\n"
- string += "\t self.method = {0}\n".format(self.method)
- string += "\t self.vertical = {0}\n".format(self.vertical)
-
- return string
-
- def get_input_path(self, pollutant=None, extension='nc'):
- """
- Completes the path of the NetCDF that contains the needed information of the given pollutant.
-
- :param pollutant: Name of the pollutant of the NetCDF.
- :type pollutant: str
-
- :param extension: Extension of the input file.
- :type: str
-
- :return: Full path of the needed NetCDF.
- :rtype: str
- """
- st_time = timeit.default_timer()
-
- file_path = os.path.join(self.inputs_path, 'multivar', 'ga_{0}.{1}'.format(
- self.date.strftime('%Y%m%d'), extension))
-
- # Checking input file
- if not os.path.exists(file_path):
- settings.write_log('ERROR: Check the .err file to get more info.')
- if settings.rank == 0:
- raise IOError('ERROR: File {0} not found.'.format(file_path))
- sys.exit(1)
-
- settings.write_time('PointGfasEmissionInventory', 'get_input_path', timeit.default_timer() - st_time, level=3)
-
- return file_path
-
- def get_altitude(self):
- """
- Extract the altitude values depending on the choosen method.
-
- :return: Array with the alittude of each fire.
- :rtype: numpy.array
- """
- from hermesv3_gr.tools.netcdf_tools import extract_vars
-
- st_time = timeit.default_timer()
-
- if self.method == 'sovief':
- alt_var = 'apt'
- elif self.method == 'prm':
- alt_var = 'mami'
- else:
- alt_var = None
-
- print "ERROR: Only 'sovief' and 'prm' methods are accepted."
-
- [alt] = extract_vars(self.get_input_path(), [alt_var])
-
- alt = alt['data']
-
- settings.write_time('PointGfasEmissionInventory', 'get_altitude', timeit.default_timer() - st_time, level=3)
- return alt
-
- @ staticmethod
- def get_approach(p_vertical):
- """
- Extract the given approach value.
-
- :return: Approach value
- :rtype: str
- """
- import re
-
- st_time = timeit.default_timer()
-
- return_value = None
- aux_list = re.split(', |,| , | ,', p_vertical)
- for element in aux_list:
- aux_value = re.split('=| =|= | = ', element)
- if aux_value[0] == 'approach':
- return_value = aux_value[1]
-
- settings.write_time('PointGfasEmissionInventory', 'get_approach', timeit.default_timer() - st_time, level=3)
-
- return return_value
-
- @ staticmethod
- def get_method(p_vertical):
- """
- Extract the given method value.
-
- :return: Method value
- :rtype: str
- """
- import re
-
- st_time = timeit.default_timer()
-
- return_value = None
- aux_list = re.split(', |,| , | ,', p_vertical)
- for element in aux_list:
- aux_value = re.split('=| =|= | = ', element)
- if aux_value[0] == 'method':
- return_value = aux_value[1]
-
- settings.write_time('PointGfasEmissionInventory', 'get_method', timeit.default_timer() - st_time, level=3)
-
- return return_value
-
- def do_vertical_allocation(self, values):
- """
- Allocates the fire emissions on their top level.
-
- :param values: 2D array with the fire emissions
- :type values: numpy.array
-
- :return: Emissions already allocated on the top altitude of each fire.
- :rtype: numpy.array
- """
- print 'do_vertical_allocation'
- sys.exit()
- st_time = timeit.default_timer()
-
- return_value = self.vertical.do_vertical_interpolation_allocation(values, self.altitude)
-
- settings.write_time('PointGfasEmissionInventory', 'do_vertical_allocation', timeit.default_timer() - st_time,
- level=3)
-
- return return_value
-
- def do_regrid(self):
- import pandas as pd
- import geopandas as gpd
- from shapely.geometry import Point
- import numpy as np
- from netCDF4 import Dataset
-
- st_time = timeit.default_timer()
- settings.write_log("\tAllocating GFAS as point sources on grid:", level=2)
-
- netcdf = Dataset(self.pollutant_dicts[0]['path'], mode='r')
- gdf = gpd.GeoDataFrame(crs={'init': 'epsg:4326'})
-
- # gdf['src_index'] = np.where(self.vertical.altitude.flatten() > 0)[0]
- # print self.input_pollutants[0]
- first_var = netcdf.variables[self.input_pollutants[0]][:]
- # first_var = netcdf.variables['bc'][:]
-
- # print first_var
- gdf['src_index'] = np.where(first_var.flatten() > 0)[0]
-
- # Masking
- if self.masking.regrid_mask is not None:
- gdf = gdf.loc[gdf['src_index'].isin(np.where(self.masking.regrid_mask.flatten() > 0)[0]), ]
-
- gdf['altitude'] = self.vertical.altitude.flatten()[gdf['src_index']]
-
- lat_1d = netcdf.variables['lat'][:]
- lon_1d = netcdf.variables['lon'][:]
- lon_1d[lon_1d > 180] -= 360
- # 1D to 2D
- lats = np.array([lat_1d] * len(lon_1d)).T
- lons = np.array([lon_1d] * len(lat_1d))
-
- gdf['geometry'] = [Point(xy) for xy in zip(lons.flatten()[gdf['src_index']],
- lats.flatten()[gdf['src_index']])]
- gdf['src_area'] = netcdf.variables['cell_area'][:].flatten()[gdf['src_index']]
- grid_shp = self.grid.to_shapefile(full_grid=False)
-
- temp_coords = Dataset(os.path.join(self.grid.temporal_path, 'temporal_coords.nc'), mode='r')
- cell_area = temp_coords.variables['cell_area'][:].flatten()
- grid_shp['dst_area'] = cell_area[grid_shp['FID']]
- # grid_shp['dst_area'] = grid_shp.to_crs({'init': 'epsg:3857'}).area
- # print grid_shp.crs['units']
- # sys.exit()
-
- gdf = gpd.sjoin(gdf.to_crs(grid_shp.crs), grid_shp, how='inner')
-
- for num, pollutant in enumerate(self.pollutant_dicts):
- settings.write_log('\t\tPollutant {0} ({1}/{2})'.format(
- pollutant['name'], num + 1, len(self.pollutant_dicts)), level=3)
-
- aux = netcdf.variables[pollutant['name']][:].flatten()[gdf['src_index']]
-
- if self.masking.scale_mask is not None:
- aux = aux * self.masking.scale_mask.flatten()[gdf['src_index']]
-
- gdf[pollutant['name']] = aux * gdf['src_area']
- # print 'masa {0}: {1} '.format(pollutant['name'], gdf[pollutant['name']].sum())
- gdf[pollutant['name']] = (aux / gdf['dst_area']) * netcdf.variables['cell_area'][:].flatten()[
- gdf['src_index']]
- # print netcdf.variables['bc'][:].sum()
-
- netcdf.close()
-
- settings.write_time('PointGfasEmissionInventory', 'do_regrid', timeit.default_timer() - st_time, level=2)
- del gdf['src_index'], gdf['index_right']
- self.emissions = gdf
-
- return True
-
- def do_regrid_balanced_NOT_USED(self):
- import pandas as pd
- import geopandas as gpd
- from shapely.geometry import Point
- import numpy as np
- from netCDF4 import Dataset
-
- st_time = timeit.default_timer()
- settings.write_log("\tAllocating GFAS as point sources on grid:", level=2)
-
- netcdf = Dataset(self.pollutant_dicts[0]['path'], mode='r')
- if settings.rank == 0:
- gdf = gpd.GeoDataFrame(crs={'init': 'epsg:4326'})
-
- first_var = netcdf.variables[self.input_pollutants[0]][:]
-
- gdf['src_index'] = np.where(first_var.flatten() > 0)[0]
-
- gdf['altitude'] = self.vertical.altitude.flatten()[gdf['src_index']]
-
- lat_1d = netcdf.variables['lat'][:]
- lon_1d = netcdf.variables['lon'][:]
- lon_1d[lon_1d > 180] -= 360
- # 1D to 2D
- lats = np.array([lat_1d] * len(lon_1d)).T
- lons = np.array([lon_1d] * len(lat_1d))
-
- gdf['geometry'] = [Point(xy) for xy in zip(lons.flatten()[gdf['src_index']],
- lats.flatten()[gdf['src_index']])]
- gdf['src_area'] = netcdf.variables['cell_area'][:].flatten()[gdf['src_index']]
- grid_shp = self.grid.to_shapefile()
-
- temp_coords = Dataset(os.path.join(self.grid.temporal_path, 'temporal_coords.nc'), mode='r')
- cell_area = temp_coords.variables['cell_area'][:].flatten()
- grid_shp['dst_area'] = cell_area[grid_shp['FID']]
- # grid_shp['dst_area'] = grid_shp.to_crs({'init': 'epsg:3857'}).area
- # print grid_shp.crs['units']
- # sys.exit()
-
- gdf = gpd.sjoin(gdf.to_crs(grid_shp.crs), grid_shp, how='inner')
- print gdf
- gdf = np.array_split(gdf, settings.size)
- else:
- gdf = None
-
- gdf = settings.comm.scatter(gdf, root=0)
-
- for num, pollutant in enumerate(self.pollutant_dicts):
- settings.write_log('\t\tPollutant {0} ({1}/{2})'.format(
- pollutant['name'], num + 1, len(self.pollutant_dicts)), level=3)
- print ('\t\tPollutant {0} ({1}/{2})'.format(pollutant['name'], num + 1, len(self.pollutant_dicts)))
- aux = netcdf.variables[pollutant['name']][:].flatten()[gdf['src_index']]
- gdf[pollutant['name']] = aux * gdf['src_area']
- # print 'masa {0}: {1} '.format(pollutant['name'], gdf[pollutant['name']].sum())
- gdf[pollutant['name']] = (aux / gdf['dst_area']) * netcdf.variables['cell_area'][:].flatten()[
- gdf['src_index']]
- # print netcdf.variables['bc'][:].sum()
-
- netcdf.close()
-
- settings.write_time('PointGfasEmissionInventory', 'do_regrid', timeit.default_timer() - st_time, level=2)
- print 'regrid done'
- del gdf['src_index'], gdf['index_right']
- self.emissions = gdf
-
- return True
-
- def calculate_altitudes(self, vertical_description_path):
- """
- Calculate the number layer to allocate the point source.
-
- :param vertical_description_path: Path to the file that contains the vertical description
- :type vertical_description_path: str
-
- :return: True
- :rtype: bool
- """
- import pandas as pd
-
- st_time = timeit.default_timer()
- settings.write_log("\t\tCalculating vertical allocation.", level=3)
- df = pd.read_csv(vertical_description_path, sep=';')
- if self.vertical.approach == 'surface':
- self.emissions['altitude'] = 0
-
- self.emissions['layer'] = None
- for i, line in df.iterrows():
- layer = line['Ilayer'] - 1
- self.emissions.loc[self.emissions['altitude'] <= line['height_magl'], 'layer'] = layer
- self.emissions.loc[self.emissions['altitude'] <= line['height_magl'], 'altitude'] = None
-
- # Fires with higher altitudes than the max layer limit goes to the last layer
- if len(self.emissions[~self.emissions['altitude'].isna()]) > 0:
- self.emissions.loc[~self.emissions['altitude'].isna(), 'layer'] = layer
- warnings.warn('WARNING: One or more fires have an altitude of fire emission injection higher than the top' +
- ' layer of the model defined in the {0} file'.format(vertical_description_path))
-
- print self.emissions.loc[~self.emissions['altitude'].isna()]
- del self.emissions['altitude']
-
- self.emissions = self.emissions.groupby(['FID', 'layer']).sum()
- self.emissions.reset_index(inplace=True)
- # print '---> Rank {0}, {1}, {2}'.format(settings.rank, self.emissions, len(self.emissions))
- # sys.exit()
- self.emissions = self.vertical.distribute_vertically(self.emissions, self.input_pollutants)
-
- settings.write_time('PointGfasEmissionInventory', 'calculate_altitudes', timeit.default_timer() - st_time,
- level=2)
- return True
-
- def point_source_by_cell(self):
- """
- Sums the different emissions that are allocated in the same cell and layer.
-
- :return: None
- """
- emissions = self.emissions
- self.location = emissions.loc[:, ['FID', 'layer']]
-
- self.emissions = []
- for input_pollutant in self.input_pollutants:
- dict_aux = {
- 'name': input_pollutant,
- 'units': '...',
- 'data': emissions.loc[:, input_pollutant].values
- }
- self.emissions.append(dict_aux)
-
-
-if __name__ == "__main__":
- pass
diff --git a/hermesv3_gr/modules/emision_inventories/point_source_emission_inventory.py b/hermesv3_gr/modules/emision_inventories/point_source_emission_inventory.py
index a66f65b638c0665d6d174262d298f27a46eb6cf4..8c22e114aa539794cd8a57109f9cb83820ab6e0d 100755
--- a/hermesv3_gr/modules/emision_inventories/point_source_emission_inventory.py
+++ b/hermesv3_gr/modules/emision_inventories/point_source_emission_inventory.py
@@ -121,9 +121,9 @@ class PointSourceEmissionInventory(EmissionInventory):
df = gpd.GeoDataFrame(df.loc[:, ['Emis', 'Alt_Injection']], crs=self.crs, geometry=geometry)
df = df.to_crs(grid_shape.crs)
-
df = gpd.sjoin(df, grid_shape, how="inner", op='intersects')
+ # Drops duplicates when the point source is on the boundary of the cell
df = df[~df.index.duplicated(keep='first')]
if self.location is None:
diff --git a/hermesv3_gr/modules/grids/grid.py b/hermesv3_gr/modules/grids/grid.py
index fe15d965030915faaf3fe781bb0df1b9964a2d3f..0c424d1ecd90aa27f068164d514c93216535a649 100644
--- a/hermesv3_gr/modules/grids/grid.py
+++ b/hermesv3_gr/modules/grids/grid.py
@@ -250,7 +250,6 @@ class Grid(object):
regular_latlon=True)
# Calculate the cell area of the auxiliary NetCDF file
-
self.cell_area = self.get_cell_area()
# Re-writes the NetCDF adding the cell area
diff --git a/hermesv3_gr/modules/vertical/vertical_gfas.py b/hermesv3_gr/modules/vertical/vertical_gfas.py
index b845208ad25749d2392741840b10322b7c9f2796..600ab7f6fef3559a064f0f37ea209847ada65260 100644
--- a/hermesv3_gr/modules/vertical/vertical_gfas.py
+++ b/hermesv3_gr/modules/vertical/vertical_gfas.py
@@ -17,7 +17,7 @@
# You should have received a copy of the GNU General Public License
# along with HERMESv3_GR. If not, see .
-import sys
+
import timeit
import hermesv3_gr.config.settings as settings
from vertical import VerticalDistribution
@@ -42,14 +42,6 @@ class GfasVerticalDistribution(VerticalDistribution):
settings.write_time('GfasVerticalDistribution', 'Init', timeit.default_timer() - st_time, level=3)
- def __str__(self):
- string = "GFAS Vertical distribution:\n"
- string += "\t self.approach = {0}\n".format(self.approach)
- string += "\t self.output_heights = {0}\n".format(self.output_heights)
- string += "\t self.altitude = {0}\n".format(self.altitude)
-
- return string
-
@staticmethod
def calculate_widths(heights_list):
"""
@@ -180,50 +172,6 @@ class GfasVerticalDistribution(VerticalDistribution):
level=3)
return fire_list
- def calculate_weight_layer_dict(self, layer):
- weight_layer_dict = {x: None for x in xrange(layer + 1)}
- if self.approach == '50_top':
- weight_layer_dict[layer] = 0.5
- to_distribute = 0.5
- layer = layer - 1
- elif self.approach == 'uniform':
- to_distribute = 1.
-
- total_width = self.output_heights[layer]
-
- previous_height = 0
- for i, height in enumerate(self.output_heights[0:layer + 1]):
- partial_width = height - previous_height
- weight_layer_dict[i] = to_distribute * partial_width / total_width
-
- previous_height = height
-
- return weight_layer_dict
-
- def distribute_vertically(self, emissions_df, input_pollutant_list):
- import pandas as pd
-
- vert_emissions = []
- for layer, emis in emissions_df.groupby('layer'):
- if layer == 0:
- vert_emissions.append(emis)
- else:
- weight_layer_dict = self.calculate_weight_layer_dict(layer)
- for layer, weight in weight_layer_dict.iteritems():
- aux_emis = emis.copy()
- aux_emis.loc[:, 'layer'] = layer
- aux_emis.loc[:, input_pollutant_list] = aux_emis[input_pollutant_list].multiply(weight)
-
- vert_emissions.append(aux_emis)
- if len(vert_emissions) == 0:
- vert_emissions = emissions_df
- else:
- vert_emissions = pd.concat(vert_emissions)
-
- vert_emissions = vert_emissions.groupby(['FID', 'layer']).sum()
- vert_emissions.reset_index(inplace=True)
- return vert_emissions
-
if __name__ == '__main__':
pass
diff --git a/hermesv3_gr/modules/writing/writer_monarch.py b/hermesv3_gr/modules/writing/writer_monarch.py
index e5019262657d9e68302d2e049212a64ed015dafc..3321b06c0f0945da74668f8e8d50d6c6bb503663 100644
--- a/hermesv3_gr/modules/writing/writer_monarch.py
+++ b/hermesv3_gr/modules/writing/writer_monarch.py
@@ -369,8 +369,8 @@ class WriterMonarch(Writer):
# Rotated pole
mapping = netcdf.createVariable('rotated_pole', 'c')
mapping.grid_mapping_name = 'rotated_latitude_longitude'
- mapping.grid_north_pole_latitude = 90 - self.grid.new_pole_latitude_degrees
- mapping.grid_north_pole_longitude = self.grid.new_pole_longitude_degrees
+ mapping.grid_north_pole_latitude = self.grid.new_pole_latitude_degrees
+ mapping.grid_north_pole_longitude = 90 - self.grid.new_pole_longitude_degrees
elif LambertConformalConic:
# CRS
mapping = netcdf.createVariable('Lambert_conformal', 'i')
diff --git a/preproc/GFAS_hourly_Parameters.csv b/preproc/GFAS_hourly_Parameters.csv
deleted file mode 100644
index 7b8f414d2ed167ab10262cdfea3a90e8826b63b4..0000000000000000000000000000000000000000
--- a/preproc/GFAS_hourly_Parameters.csv
+++ /dev/null
@@ -1,27 +0,0 @@
-Name;Short_Name;Units;id
-Carbon Monoxide;co;kg m-2 s-1;param81.210.192
-Non-Methane Hydro-Carbons;nmhc;kg m-2 s-1;param83.210.192
-Nitrogen Oxides expressed as monoxide nitrogen;nox_no;kg m-2 s-1;param85.210.192
-Particulate Matter PM2.5;pm25;kg m-2 s-1;param87.210.192
-Organic Carbon;oc;kg m-2 s-1;param90.210.192
-Black Carbon;bc;kg m-2 s-1;param91.210.192
-Sulfur Dioxide;so2;kg m-2 s-1;param102.210.192
-Methanol;ch3oh;kg m-2 s-1;param103.210.192
-Higher Alkanes;hialkanes;kg m-2 s-1;param112.210.192
-Formaldehyde;ch2o;kg m-2 s-1;param113.210.192
-Acetaldehyde;c2h4o;kg m-2 s-1;param114.210.192
-Acetone;c3h6o;kg m-2 s-1;param115.210.192
-Ammonia;nh3;kg m-2 s-1;param116.210.192
-Dimethyl Sulfide;c2h6s;kg m-2 s-1;param117.210.192
-Ethane;c2h6;kg m-2 s-1;param118.210.192
-Toluene;c7h8;kg m-2 s-1;param231.210.192
-Benzene;c6h6;kg m-2 s-1;param232.210.192
-Xylene;c8h10;kg m-2 s-1;param233.210.192
-Butenes;c4h8;kg m-2 s-1;param234.210.192
-Pentenes;c5h10;kg m-2 s-1;param235.210.192
-Hexene;c6h12;kg m-2 s-1;param236.210.192
-Octene;c8h16;kg m-2 s-1;param237.210.192
-Butanes;c4h10;kg m-2 s-1;param238.210.192
-Pentanes;c5h12;kg m-2 s-1;param239.210.192
-Hexanes;c6h14;kg m-2 s-1;param240.210.192
-Heptane;c7h16;kg m-2 s-1;param241.210.192
diff --git a/preproc/gfas12_h_preproc.py b/preproc/gfas12_h_preproc.py
deleted file mode 100755
index 9cc93fef153e21e8666355ab350c815145f1cf56..0000000000000000000000000000000000000000
--- a/preproc/gfas12_h_preproc.py
+++ /dev/null
@@ -1,282 +0,0 @@
-#!/usr/bin/env python
-
-# Copyright 2018 Earth Sciences Department, BSC-CNS
-#
-# This file is part of HERMESv3_GR.
-#
-# HERMESv3_GR is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# HERMESv3_GR is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with HERMESv3_GR. If not, see .
-
-
-import os
-from netCDF4 import Dataset
-import cf_units
-import pandas as pd
-import datetime
-from datetime import datetime, timedelta
-
-# ============== README ======================
-"""
-downloading website: http://apps.ecmwf.int/datasets/data/cams-gfas/
-reference: https://www.biogeosciences.net/9/527/2012/
-Besides citing HERMESv3_GR, users must also acknowledge the use of the corresponding emission inventories in their works
-"""
-
-# ============== CONFIGURATION PARAMETERS ======================
-INPUT_PATH = '/esarchive/recon/ecmwf/gfas/original_files/ga_mc_sfc_gfas_ecmf/'
-INPUT_NAME = 'gfas_hourly_.grb'
-OUTPUT_PATH = '/esarchive/recon/ecmwf/gfas'
-
-STARTING_DATE = datetime(year=2018, month=11, day=01)
-ENDIND_DATE = datetime(year=2018, month=11, day=01)
-
-PARAMETERS_FILE = '/esarchive/recon/ecmwf/gfas/original_files/ga_mc_sfc_gfas_ecmf/GFAS_hourly_Parameters.csv'
-# ==============================================================
-
-
-def create_bounds(coords, number_vertices=2):
- """
- Calculate the vertices coordinates.
-
- :param coords: Coordinates in degrees (latitude or longitude)
- :type coords: numpy.array
-
- :param number_vertices: Non mandatory parameter that informs the number of vertices that must have the boundaries.
- (by default 2)
- :type number_vertices: int
-
- :return: Array with as many elements as vertices for each value of coords.
- :rtype: numpy.array
- """
- import numpy as np
-
- interval = coords[1] - coords[0]
-
- coords_left = coords - interval / 2
- coords_right = coords + interval / 2
- if number_vertices == 2:
- bound_coords = np.dstack((coords_left, coords_right))
- elif number_vertices == 4:
- bound_coords = np.dstack((coords_left, coords_right, coords_right, coords_left))
- else:
- raise ValueError('The number of vertices of the boudaries must be 2 or 4')
-
- return bound_coords
-
-
-def get_grid_area(filename):
- """
- Calculate the area for each cell of the grid using CDO
-
- :param filename: Path to the file to calculate the cell area
- :type filename: str
-
- :return: Area of each cell of the grid.
- :rtype: numpy.array
- """
- from cdo import Cdo
- from netCDF4 import Dataset
-
- cdo = Cdo()
- s = cdo.gridarea(input=filename)
- nc_aux = Dataset(s, mode='r')
- grid_area = nc_aux.variables['cell_area'][:]
- nc_aux.close()
-
- return grid_area
-
-
-def write_netcdf(output_name_path, data_list, center_lats, center_lons, grid_cell_area, date):
- """
- Write a NetCDF with the given information.
-
- :param output_name_path: Complete path to the output NetCDF to be stored.
- :type output_name_path: str
-
- :param data_list
-
- :param center_lats: Latitudes of the center of each cell.
- :type center_lats: numpy.array
-
- :param center_lons: Longitudes of the center of each cell.
- :type center_lons: numpy.array
-
- :param grid_cell_area: Area of each cell of the grid.
- :type grid_cell_area: numpy.array
-
- :param date: Date of the current netCDF.
- :type date: datetime.datetime
-
- """
- print output_name_path
- # Creating NetCDF & Dimensions
- nc_output = Dataset(output_name_path, mode='w', format="NETCDF4")
- nc_output.createDimension('nv', 2)
- nc_output.createDimension('lon', center_lons.shape[0])
- nc_output.createDimension('lat', center_lats.shape[0])
- nc_output.createDimension('time', None)
-
- # TIME
- time = nc_output.createVariable('time', 'd', ('time',), zlib=True)
- # time.units = "{0} since {1}".format(tstep_units, global_atts['Start_DateTime'].strftime('%Y-%m-%d %H:%M:%S'))
- time.units = "hours since {0}".format(date.strftime('%Y-%m-%d %H:%M:%S'))
- time.standard_name = "time"
- time.calendar = "gregorian"
- time.long_name = "time"
- time[:] = [0]
-
- # LATITUDE
- lat = nc_output.createVariable('lat', 'f', ('lat',), zlib=True)
- lat.bounds = "lat_bnds"
- lat.units = "degrees_north"
- lat.axis = "Y"
- lat.long_name = "latitude"
- lat.standard_name = "latitude"
- lat[:] = center_lats
-
- lat_bnds = nc_output.createVariable('lat_bnds', 'f', ('lat', 'nv',), zlib=True)
- lat_bnds[:] = create_bounds(center_lats)
-
- # LONGITUDE
- lon = nc_output.createVariable('lon', 'f', ('lon',), zlib=True)
- lon.bounds = "lon_bnds"
- lon.units = "degrees_east"
- lon.axis = "X"
- lon.long_name = "longitude"
- lon.standard_name = "longitude"
- lon[:] = center_lons
-
- lon_bnds = nc_output.createVariable('lon_bnds', 'f', ('lon', 'nv',), zlib=True)
- lon_bnds[:] = create_bounds(center_lons)
-
- for var in data_list:
- # VARIABLE
- nc_var = nc_output.createVariable(var['name'], 'f', ('time', 'lat', 'lon',), zlib=True)
- nc_var.units = var['units'].symbol
- nc_var.long_name = var['long_name']
- nc_var.coordinates = 'lat lon'
- nc_var.grid_mapping = 'crs'
- nc_var.cell_measures = 'area: cell_area'
- nc_var[:] = var['data']
-
- # CELL AREA
- cell_area = nc_output.createVariable('cell_area', 'f', ('lat', 'lon',))
- cell_area.long_name = "area of the grid cell"
- cell_area.standard_name = "area"
- cell_area.units = "m2"
- cell_area[:] = grid_cell_area
-
- # CRS
- crs = nc_output.createVariable('crs', 'i')
- crs.grid_mapping_name = "latitude_longitude"
- crs.semi_major_axis = 6371000.0
- crs.inverse_flattening = 0
-
- nc_output.setncattr('title', 'GFASv1.2 inventory')
- nc_output.setncattr('Conventions', 'CF-1.6', )
- nc_output.setncattr('institution', 'ECMWF', )
- nc_output.setncattr('source', 'GFAS', )
- nc_output.setncattr('history', 'Re-writing of the GFAS input to follow the CF-1.6 conventions;\n' +
- '2019-01-02: Added boundaries;\n' +
- '2019-01-02: Added global attributes;\n' +
- '2019-01-02: Re-naming pollutant;\n' +
- '2019-01-02: Added cell_area variable;\n' +
- '2019-01-02: Added new varaibles;\n')
- nc_output.setncattr('references', 'downloading website: http://apps.ecmwf.int/datasets/data/cams-gfas/\n' +
- 'reference: https://www.biogeosciences.net/9/527/2012/', )
- nc_output.setncattr('comment', 'Re-writing done by Carles Tena (carles.tena@bsc.es) from the BSC-CNS ' +
- '(Barcelona Supercomputing Center)', )
-
- nc_output.close()
- return True
-
-
-def do_transformation(input_file, date, output_dir, variables_list):
- """
- Transform the original file into a NEtCDF file that follows the conventions.
-
- :param input_file:
- :type input_file: str
-
- :param date: Date of the file to do the transformation.
- :type date: datetime.datetime
-
- :param output_dir: Path where have to be stored the output file.
- :type output_dir: str
-
- :param variables_list: LIst of dictionaries with the information of each variable of the output files.
- :type variables_list: list
- """
- from cdo import Cdo
- cdo = Cdo()
-
- nc_temp = cdo.copy(input=input_file, options='-R -r -f nc4c -z zip_4')
-
- nc_in = Dataset(nc_temp, mode='r')
-
- cell_area = get_grid_area(nc_temp)
-
- lats = nc_in.variables['lat'][:]
- lons = nc_in.variables['lon'][:]
-
- for variable in variables_list:
- variable['data'] = nc_in.variables[variable['original_name']][:]
-
- nc_in.close()
-
- out_path_aux = os.path.join(output_dir, '1hourly', 'multivar')
- if not os.path.exists(out_path_aux):
- os.makedirs(out_path_aux)
- out_path_aux = os.path.join(out_path_aux, 'ga_{0}.nc'.format(date.strftime('%Y%m%d')))
- write_netcdf(out_path_aux, variables_list, lats, lons, cell_area, date)
-
- return True
-
-
-def do_var_list(variables_file):
- """
- Create the List of dictionaries
-
- :param variables_file: CSV file with the information of each variable
- :type variables_file: str
-
- :return: Dictionaries list with the information of each variable.
- :rtype: list
- """
- df = pd.read_csv(variables_file, sep=';')
- list_aux = []
- for i, element in df.iterrows():
- # print element
- dict_aux = {
- 'original_name': str(element.id),
- 'name': element['Short_Name'],
- 'long_name': element['Name'],
- 'units': cf_units.Unit(element['Units']),
- }
- list_aux.append(dict_aux)
- return list_aux
-
-
-if __name__ == '__main__':
-
- var_list = do_var_list(PARAMETERS_FILE)
-
- date_aux = STARTING_DATE
- while date_aux <= ENDIND_DATE:
- f = os.path.join(INPUT_PATH, INPUT_NAME.replace('', date_aux.strftime('%Y%m%d')))
- if os.path.isfile(f):
- do_transformation(f, date_aux, OUTPUT_PATH, var_list)
- else:
- print 'ERROR: file {0} not found'.format(f)
-
- date_aux = date_aux + timedelta(days=1)
diff --git a/run_test.py b/run_test.py
index 711c425bbcb719bc46c277b18d1f4f754326274b..6eb2926793f63a9c62eee79c06273d12f95fab82 100644
--- a/run_test.py
+++ b/run_test.py
@@ -1,22 +1,4 @@
-#!/usr/bin/env python
# coding=utf-8
-
-# Copyright 2018 Earth Sciences Department, BSC-CNS
-#
-# This file is part of HERMESv3_GR.
-#
-# HERMESv3_GR is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# HERMESv3_GR is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with HERMESv3_GR. If not, see .
"""Script to run the tests for EarthDiagnostics and generate the code coverage report"""
import os