From 028bc7f5fa8887b91a7294085dad486774a9cc61 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 13:32:13 +0200 Subject: [PATCH 01/11] Added custom COMM Added first time Added erase auxiliary files --- conf/hermes.conf | 2 ++ hermesv3_gr/config/config.py | 33 ++++++++++++++++++------ hermesv3_gr/config/settings.py | 6 +++-- hermesv3_gr/hermes.py | 4 +-- hermesv3_gr/modules/temporal/temporal.py | 2 +- hermesv3_gr/tools/netcdf_tools.py | 5 ---- 6 files changed, 34 insertions(+), 18 deletions(-) diff --git a/conf/hermes.conf b/conf/hermes.conf index a9d1393..8797122 100755 --- a/conf/hermes.conf +++ b/conf/hermes.conf @@ -11,6 +11,8 @@ start_date = 2017/01/01 00:00:00 output_timestep_type = hourly output_timestep_num = 24 output_timestep_freq = 1 +first_time = 1 +erase_auxiliary_files = 1 [DOMAIN] diff --git a/hermesv3_gr/config/config.py b/hermesv3_gr/config/config.py index f8d1d4a..9199455 100755 --- a/hermesv3_gr/config/config.py +++ b/hermesv3_gr/config/config.py @@ -20,13 +20,20 @@ from configargparse import ArgParser import sys +import os +from shutil import rmtree +from mpi4py import MPI class Config(ArgParser): """ Initialization of the arguments that the parser can handle. """ - def __init__(self): + def __init__(self, comm=None): + if comm is None: + comm = MPI.COMM_WORLD + self.comm = comm + super(Config, self).__init__() self.options = self.read_options() @@ -61,6 +68,10 @@ class Config(ArgParser): type=str, choices=['hourly', 'daily', 'monthly', 'yearly']) p.add_argument('--output_timestep_num', required=True, help='Number of timesteps to simulate.', type=int) p.add_argument('--output_timestep_freq', required=True, help='Frequency between timesteps.', type=int) + p.add_argument('--first_time', required=False, default='False', type=str, + help='Indicates if you want to run it for first time (only create auxiliary files).') + p.add_argument('--erase_auxiliary_files', required=False, default='False', type=str, + help='Indicates if you want to start from scratch removing the auxiliary files already created.') p.add_argument('--output_model', required=True, help='Name of the output model.', choices=['MONARCH', 'CMAQ', 'WRF_CHEM']) @@ -177,8 +188,16 @@ class Config(ArgParser): arguments.start_date = self._parse_start_date(arguments.start_date) arguments.end_date = self._parse_end_date(arguments.end_date, arguments.start_date) + arguments.first_time = self._parse_bool(arguments.first_time) + arguments.erase_auxiliary_files = self._parse_bool(arguments.erase_auxiliary_files) + self.create_dir(arguments.output_dir) - self.create_dir(arguments.auxiliar_files_path) + if arguments.erase_auxiliary_files: + if os.path.exists(arguments.auxiliary_files_path): + if self.comm.Get_rank() == 0: + rmtree(arguments.auxiliary_files_path) + self.comm.Barrier() + self.create_dir(arguments.auxiliary_files_path) return arguments @@ -207,8 +226,7 @@ class Config(ArgParser): full_path = os.path.join(self.options.output_dir, file_name) return full_path - @staticmethod - def create_dir(path): + def create_dir(self, path): """ Create the given folder if it is not created yet. @@ -217,8 +235,7 @@ class Config(ArgParser): """ import os from mpi4py import MPI - icomm = MPI.COMM_WORLD - comm = icomm.Split(color=0, key=0) + comm = self.comm.Split(color=0, key=0) rank = comm.Get_rank() if rank == 0: @@ -304,12 +321,12 @@ class Config(ArgParser): else: return self._parse_start_date(end_date) - def set_log_level(self): + def set_log_level(self, comm): """ Defines the log_level using the common script settings. """ from . import settings - settings.define_global_vars(self.options.log_level) + settings.define_global_vars(self.options.log_level, comm) if __name__ == '__main__': diff --git a/hermesv3_gr/config/settings.py b/hermesv3_gr/config/settings.py index 18b6a6a..46b10a6 100755 --- a/hermesv3_gr/config/settings.py +++ b/hermesv3_gr/config/settings.py @@ -46,7 +46,7 @@ global log_file global df_times -def define_global_vars(in_log_level): +def define_global_vars(in_log_level, icomm): # TODO Documentation from mpi4py import MPI @@ -55,7 +55,9 @@ def define_global_vars(in_log_level): global rank global size - icomm = MPI.COMM_WORLD + if icomm is None: + icomm = MPI.COMM_WORLD + comm = icomm.Split(color=0, key=0) rank = comm.Get_rank() size = comm.Get_size() diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 52d46cb..00f6977 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -34,7 +34,7 @@ class Hermes(object): """ Interface class for HERMESv3. """ - def __init__(self, config, new_date=None): + def __init__(self, config, new_date=None, comm=None): from hermesv3_gr.modules.grids.grid import Grid from hermesv3_gr.modules.temporal.temporal import TemporalDistribution from hermesv3_gr.modules.writing.writer import Writer @@ -48,7 +48,7 @@ class Hermes(object): if new_date is not None: self.options.start_date = new_date - config.set_log_level() + config.set_log_level(comm) settings.define_log_file(self.options.output_dir, self.options.start_date) settings.define_times_file() diff --git a/hermesv3_gr/modules/temporal/temporal.py b/hermesv3_gr/modules/temporal/temporal.py index 457d0c8..5606513 100755 --- a/hermesv3_gr/modules/temporal/temporal.py +++ b/hermesv3_gr/modules/temporal/temporal.py @@ -688,7 +688,7 @@ class TemporalDistribution(object): nc_in = Dataset(self.netcdf_timezones) timezones = nc_in.variables['timezone_id'][:, self.grid.x_lower_bound:self.grid.x_upper_bound, - self.grid.y_lower_bound:self.grid.y_upper_bound].astype(int) + self.grid.y_lower_bound:self.grid.y_upper_bound].astype(int) nc_in.close() tz_list = np.chararray(timezones.shape, itemsize=32) diff --git a/hermesv3_gr/tools/netcdf_tools.py b/hermesv3_gr/tools/netcdf_tools.py index 48869a1..1709435 100755 --- a/hermesv3_gr/tools/netcdf_tools.py +++ b/hermesv3_gr/tools/netcdf_tools.py @@ -23,11 +23,6 @@ from netCDF4 import Dataset from mpi4py import MPI from functools import reduce -ICOMM = MPI.COMM_WORLD -COMM = ICOMM.Split(color=0, key=0) -RANK = COMM.Get_rank() -SIZE = COMM.Get_size() - def open_netcdf(netcdf_path): """ -- GitLab From 8485dd7b5ee097ab6892f5488fe205068a871aac Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 13:37:31 +0200 Subject: [PATCH 02/11] WIP --- hermesv3_gr/hermes.py | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 00f6977..c2a92c3 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -41,6 +41,9 @@ class Hermes(object): global full_time st_time = full_time = timeit.default_timer() + if comm is None: + comm = MPI.COMM_WORLD + self.comm = comm self.config = config self.options = config.options @@ -48,7 +51,7 @@ class Hermes(object): if new_date is not None: self.options.start_date = new_date - config.set_log_level(comm) + config.set_log_level(self.comm) settings.define_log_file(self.options.output_dir, self.options.start_date) settings.define_times_file() @@ -95,6 +98,9 @@ class Hermes(object): """ from datetime import timedelta + if self.options.fisrt_time: + self.comm.Abort(0) + st_time = timeit.default_timer() settings.write_log('') settings.write_log('***** Starting HERMESv3 *****') @@ -139,10 +145,10 @@ class Hermes(object): return None -def run(): +def run(comm=None): date = Hermes(Config()).main() while date is not None: - date = Hermes(Config(), new_date=date).main() + date = Hermes(Config(comm), new_date=date, comm=comm).main() sys.exit(0) -- GitLab From 57abf198c060d962146714aff5c37a594a504ad6 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 13:43:17 +0200 Subject: [PATCH 03/11] auxiliar to auxiliary --- conf/hermes.conf | 4 ++-- hermesv3_gr/config/config.py | 2 +- hermesv3_gr/config/settings.py | 6 ++++-- hermesv3_gr/hermes.py | 2 +- .../modules/emision_inventories/emission_inventory.py | 6 +++--- 5 files changed, 11 insertions(+), 9 deletions(-) diff --git a/conf/hermes.conf b/conf/hermes.conf index 8797122..58dae1a 100755 --- a/conf/hermes.conf +++ b/conf/hermes.conf @@ -12,7 +12,7 @@ output_timestep_type = hourly output_timestep_num = 24 output_timestep_freq = 1 first_time = 1 -erase_auxiliary_files = 1 +erase_auxiliary_files = 0 [DOMAIN] @@ -29,7 +29,7 @@ domain_type = global # domain_type = rotated # domain_type = mercator vertical_description = /data/profiles/vertical/Benchmark_15layers_vertical_description.csv -auxiliar_files_path = /data/auxiliar_files/_ +auxiliary_files_path = /data/auxiliar_files/_ # if domain_type == global: inc_lat = 1. diff --git a/hermesv3_gr/config/config.py b/hermesv3_gr/config/config.py index 9199455..7775770 100755 --- a/hermesv3_gr/config/config.py +++ b/hermesv3_gr/config/config.py @@ -80,7 +80,7 @@ class Config(ArgParser): p.add_argument('--domain_type', required=True, help='Type of domain to simulate.', choices=['global', 'lcc', 'rotated', 'mercator', 'regular']) - p.add_argument('--auxiliar_files_path', required=True, + p.add_argument('--auxiliary_files_path', required=True, help='Path to the directory where the necessary auxiliary files will be created if them are ' + 'not created yet.') diff --git a/hermesv3_gr/config/settings.py b/hermesv3_gr/config/settings.py index 46b10a6..038f004 100755 --- a/hermesv3_gr/config/settings.py +++ b/hermesv3_gr/config/settings.py @@ -46,7 +46,7 @@ global log_file global df_times -def define_global_vars(in_log_level, icomm): +def define_global_vars(in_log_level, mycomm): # TODO Documentation from mpi4py import MPI @@ -55,8 +55,10 @@ def define_global_vars(in_log_level, icomm): global rank global size - if icomm is None: + if mycomm is None: icomm = MPI.COMM_WORLD + else: + icomm = mycomm comm = icomm.Split(color=0, key=0) rank = comm.Get_rank() diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index c2a92c3..8d6ddda 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -68,7 +68,7 @@ class Hermes(object): self.grid = Grid.select_grid( self.options.domain_type, self.options.vertical_description, self.options.output_timestep_num, - self.options.auxiliar_files_path, self.options.inc_lat, self.options.inc_lon, self.options.lat_orig, + self.options.auxiliary_files_path, self.options.inc_lat, self.options.inc_lon, self.options.lat_orig, self.options.lon_orig, self.options.n_lat, self.options.n_lon, self.options.centre_lat, self.options.centre_lon, self.options.west_boundary, self.options.south_boundary, self.options.inc_rlat, self.options.inc_rlon, self.options.lat_1, self.options.lat_2, self.options.lon_0, self.options.lat_0, diff --git a/hermesv3_gr/modules/emision_inventories/emission_inventory.py b/hermesv3_gr/modules/emision_inventories/emission_inventory.py index 43e0db7..e03a3bd 100755 --- a/hermesv3_gr/modules/emision_inventories/emission_inventory.py +++ b/hermesv3_gr/modules/emision_inventories/emission_inventory.py @@ -103,7 +103,7 @@ class EmissionInventory(object): # It will also create the WoldMasks necessaries self.masking = Masking( options.world_info, factors, regrid_mask, grid, - world_mask_file=os.path.join(os.path.dirname(options.auxiliar_files_path), + world_mask_file=os.path.join(os.path.dirname(options.auxiliary_files_path), '{0}_WorldMask.nc'.format(inventory_name))) self.input_pollutants = pollutants @@ -116,7 +116,7 @@ class EmissionInventory(object): if self.source_type == 'area': self.regrid = ConservativeRegrid( self.pollutant_dicts, - os.path.join(options.auxiliar_files_path, + os.path.join(options.auxiliary_files_path, "Weight_Matrix_{0}_{1}.nc".format(self.inventory_name, settings.size)), grid, masking=self.masking) @@ -136,7 +136,7 @@ class EmissionInventory(object): self.temporal = TemporalDistribution( current_date, options.output_timestep_type, options.output_timestep_num, options.output_timestep_freq, options.p_month, p_month, options.p_week, p_week, options.p_day, p_day, options.p_hour, p_hour, - options.world_info, options.auxiliar_files_path, grid) + options.world_info, options.auxiliary_files_path, grid) else: self.temporal = None settings.write_log('\t\tNone temporal profile set.', level=2) -- GitLab From dfa880540200a0420ab381521166579766b8bd8f Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 13:53:22 +0200 Subject: [PATCH 04/11] WIP --- .codacy.yml | 2 +- .gitlab-ci.yml | 4 ++-- CHANGELOG | 5 +++++ environment.yml | 2 +- hermesv3_gr/hermes.py | 4 +++- run_test.py | 1 - 6 files changed, 12 insertions(+), 6 deletions(-) diff --git a/.codacy.yml b/.codacy.yml index bf12fc2..0bc0d39 100755 --- a/.codacy.yml +++ b/.codacy.yml @@ -13,7 +13,7 @@ engines: enabled: true pylint: enabled: true - python_version: 2 + python_version: 3 exclude_paths: [ 'doc/**', diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 08c1355..d4676c9 100755 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -20,11 +20,11 @@ prepare: test_python2: stage: test script: - - conda env update -f environment.yml -n hermesv3_gr python=2.7 + - conda env update -f environment.yml -n hermesv3_gr python=3.7 - source activate hermesv3_gr - python run_test.py - pip install codacy-coverage --upgrade - - python-codacy-coverage -r tests/report/python2/coverage.xml + - python-codacy-coverage -r tests/report/python3/coverage.xml #test_python3: # stage: test diff --git a/CHANGELOG b/CHANGELOG index 7d4bbaf..947c8d8 100755 --- a/CHANGELOG +++ b/CHANGELOG @@ -1,3 +1,8 @@ +2.0.0 + - Python 3 + - first_time option + - erase_auxiliary_files option + 1.0.4 2019/07/19 - Specified version on timezonefinder to be less than 4.0.0 (not further support to python 2.7.X) diff --git a/environment.yml b/environment.yml index bebcce3..dc46ab7 100755 --- a/environment.yml +++ b/environment.yml @@ -6,7 +6,7 @@ channels: - conda-forge dependencies: - - python = 2 + - python = 3 - numpy - netcdf4 >= 1.3.1 - python-cdo >= 1.3.6 diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 8d6ddda..bfc8aec 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -98,7 +98,9 @@ class Hermes(object): """ from datetime import timedelta - if self.options.fisrt_time: + if self.options.first_time: + settings.write_log('FINISHED FIRST TIME RUN') + self.comm.Barrier() self.comm.Abort(0) st_time = timeit.default_timer() diff --git a/run_test.py b/run_test.py index 6eb2926..655465f 100755 --- a/run_test.py +++ b/run_test.py @@ -1,5 +1,4 @@ # coding=utf-8 -"""Script to run the tests for EarthDiagnostics and generate the code coverage report""" import os import sys -- GitLab From f9dfbfe62bb4f4de19fa1a74666806af435c3160 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 14:02:41 +0200 Subject: [PATCH 05/11] PEP8 corrections --- hermesv3_gr/modules/speciation/speciation.py | 2 +- hermesv3_gr/tools/netcdf_tools.py | 5 ++--- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/hermesv3_gr/modules/speciation/speciation.py b/hermesv3_gr/modules/speciation/speciation.py index 9d0c1b8..efb3061 100755 --- a/hermesv3_gr/modules/speciation/speciation.py +++ b/hermesv3_gr/modules/speciation/speciation.py @@ -146,7 +146,7 @@ class Speciation(object): raise KeyError('ERROR: {0} pollutant is not in the molecular weights file {1} .'.format( emission['name'], self.molecular_weights_path)) sys.exit(1) - exec ("{0} = np.array(emission['data'], dtype=settings.precision)".format(emission['name'])) + exec("{0} = np.array(emission['data'], dtype=settings.precision)".format(emission['name'])) emission['units'] = '' input_pollutants.append(emission['name']) diff --git a/hermesv3_gr/tools/netcdf_tools.py b/hermesv3_gr/tools/netcdf_tools.py index 1709435..8678e9f 100755 --- a/hermesv3_gr/tools/netcdf_tools.py +++ b/hermesv3_gr/tools/netcdf_tools.py @@ -601,9 +601,8 @@ def create_netcdf(netcdf_path, center_latitudes, center_longitudes, data_list, # if variable['data'] is not 0: # print var[:].shape, variable['data'].shape, variable['data'].max() shape = tuple() - exec ("shape = (len(hours), {0}.size, {1}.size, {2}.size)".format(var_dim[0], var_dim[1], var_dim[2])) - # exit() - print(shape) + exec("shape = (len(hours), {0}.size, {1}.size, {2}.size)".format(var_dim[0], var_dim[1], var_dim[2])) + var[:] = np.zeros(shape) # Grid mapping -- GitLab From 3b558b8534242b03c74266b4dbaa5e9c035794d9 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Tue, 22 Oct 2019 17:48:01 +0200 Subject: [PATCH 06/11] WIP --- conf/hermes.conf | 4 +- hermesv3_gr/hermes.py | 6 +-- .../emision_inventories/emission_inventory.py | 36 ++++++------- .../gfas_emission_inventory.py | 10 ++-- .../point_gfas_emission_inventory.py | 11 ++-- hermesv3_gr/modules/grids/grid.py | 45 +++++++++++++---- hermesv3_gr/modules/grids/grid_global.py | 12 ++--- hermesv3_gr/modules/grids/grid_latlon.py | 12 ++--- hermesv3_gr/modules/grids/grid_lcc.py | 12 ++--- hermesv3_gr/modules/grids/grid_mercator.py | 8 +-- hermesv3_gr/modules/grids/grid_rotated.py | 16 ++---- .../modules/regrid/regrid_conservative.py | 2 +- hermesv3_gr/modules/speciation/speciation.py | 2 +- hermesv3_gr/modules/temporal/temporal.py | 50 ++++++++----------- 14 files changed, 103 insertions(+), 123 deletions(-) diff --git a/conf/hermes.conf b/conf/hermes.conf index 58dae1a..7a12597 100755 --- a/conf/hermes.conf +++ b/conf/hermes.conf @@ -11,8 +11,8 @@ start_date = 2017/01/01 00:00:00 output_timestep_type = hourly output_timestep_num = 24 output_timestep_freq = 1 -first_time = 1 -erase_auxiliary_files = 0 +first_time = 0 +erase_auxiliary_files = 1 [DOMAIN] diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index bfc8aec..8f67bef 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -67,6 +67,7 @@ class Hermes(object): self.levels = VerticalDistribution.get_vertical_output_profile(self.options.vertical_description) self.grid = Grid.select_grid( + self.comm, self.options.domain_type, self.options.vertical_description, self.options.output_timestep_num, self.options.auxiliary_files_path, self.options.inc_lat, self.options.inc_lon, self.options.lat_orig, self.options.lon_orig, self.options.n_lat, self.options.n_lon, self.options.centre_lat, @@ -98,11 +99,6 @@ class Hermes(object): """ from datetime import timedelta - if self.options.first_time: - settings.write_log('FINISHED FIRST TIME RUN') - self.comm.Barrier() - self.comm.Abort(0) - st_time = timeit.default_timer() settings.write_log('') settings.write_log('***** Starting HERMESv3 *****') diff --git a/hermesv3_gr/modules/emision_inventories/emission_inventory.py b/hermesv3_gr/modules/emision_inventories/emission_inventory.py index e03a3bd..1670424 100755 --- a/hermesv3_gr/modules/emision_inventories/emission_inventory.py +++ b/hermesv3_gr/modules/emision_inventories/emission_inventory.py @@ -119,16 +119,16 @@ class EmissionInventory(object): os.path.join(options.auxiliary_files_path, "Weight_Matrix_{0}_{1}.nc".format(self.inventory_name, settings.size)), grid, masking=self.masking) - - # Creating Vertical Object - if p_vertical is not None: - self.vertical = VerticalDistribution( - self.get_profile(p_vertical), vertical_profile_path=options.p_vertical, - vertical_output_profile=vertical_output_profile) - else: - self.vertical = None - settings.write_log('\t\tNone vertical profile set.', level=2) - self.vertical_factors = None + if not options.first_time: + # Creating Vertical Object + if p_vertical is not None: + self.vertical = VerticalDistribution( + self.get_profile(p_vertical), vertical_profile_path=options.p_vertical, + vertical_output_profile=vertical_output_profile) + else: + self.vertical = None + settings.write_log('\t\tNone vertical profile set.', level=2) + self.vertical_factors = None # Creating Temporal Object # It will also create the necessaries timezone files @@ -136,18 +136,18 @@ class EmissionInventory(object): self.temporal = TemporalDistribution( current_date, options.output_timestep_type, options.output_timestep_num, options.output_timestep_freq, options.p_month, p_month, options.p_week, p_week, options.p_day, p_day, options.p_hour, p_hour, - options.world_info, options.auxiliary_files_path, grid) + options.world_info, options.auxiliary_files_path, grid, options.first_time) else: self.temporal = None settings.write_log('\t\tNone temporal profile set.', level=2) self.temporal_factors = None - - # Creating Speciation Object - if p_speciation is not None: - self.speciation = Speciation(p_speciation, options.p_speciation, options.molecular_weights) - else: - self.speciation = None - settings.write_log('\t\tNone speciation profile set.', level=2) + if not options.first_time: + # Creating Speciation Object + if p_speciation is not None: + self.speciation = Speciation(p_speciation, options.p_speciation, options.molecular_weights) + else: + self.speciation = None + settings.write_log('\t\tNone speciation profile set.', level=2) self.vertical_weights = None diff --git a/hermesv3_gr/modules/emision_inventories/gfas_emission_inventory.py b/hermesv3_gr/modules/emision_inventories/gfas_emission_inventory.py index 33a5328..c062615 100755 --- a/hermesv3_gr/modules/emision_inventories/gfas_emission_inventory.py +++ b/hermesv3_gr/modules/emision_inventories/gfas_emission_inventory.py @@ -82,13 +82,13 @@ class GfasEmissionInventory(EmissionInventory): vertical_output_profile, reference_year=reference_year, factors=factors, regrid_mask=regrid_mask, p_vertical=None, p_month=p_month, p_week=p_week, p_day=p_day, p_hour=p_hour, p_speciation=p_speciation) + if not options.first_time: + self.approach = self.get_approach(p_vertical) + self.method = self.get_method(p_vertical) - self.approach = self.get_approach(p_vertical) - self.method = self.get_method(p_vertical) + self.altitude = self.get_altitude() - self.altitude = self.get_altitude() - - self.vertical = GfasVerticalDistribution(vertical_output_profile, self.approach, self.get_altitude()) + self.vertical = GfasVerticalDistribution(vertical_output_profile, self.approach, self.get_altitude()) settings.write_time('GFAS_EmissionInventory', 'Init', timeit.default_timer() - st_time, level=3) diff --git a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py b/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py index 78167bb..10ca3a1 100755 --- a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py +++ b/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py @@ -84,13 +84,14 @@ class PointGfasEmissionInventory(EmissionInventory): reference_year=reference_year, factors=factors, regrid_mask=regrid_mask, p_vertical=None, p_month=p_month, p_week=p_week, 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) + if not options.first_time: + # self.approach = self.get_approach(p_vertical) + self.method = self.get_method(p_vertical) - # self.altitude = self.get_altitude() + # self.altitude = self.get_altitude() - self.vertical = GfasVerticalDistribution(vertical_output_profile, self.get_approach(p_vertical), - 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) diff --git a/hermesv3_gr/modules/grids/grid.py b/hermesv3_gr/modules/grids/grid.py index dca2013..b3acdb9 100755 --- a/hermesv3_gr/modules/grids/grid.py +++ b/hermesv3_gr/modules/grids/grid.py @@ -40,10 +40,10 @@ class Grid(object): :type temporal_path: str """ - def __init__(self, grid_type, vertical_description_path, temporal_path): + def __init__(self, grid_type, vertical_description_path, temporal_path, comm): st_time = timeit.default_timer() # settings.write_log('Creating Grid...', level=1) - + self.comm = comm # Defining class atributes self.procs_array = None self.nrows = 0 @@ -66,7 +66,7 @@ class Grid(object): self.temporal_path = temporal_path self.shapefile_path = None - self.esmf_grid = None + # self.esmf_grid = None self.x_lower_bound = None self.x_upper_bound = None self.y_lower_bound = None @@ -78,6 +78,27 @@ class Grid(object): settings.write_time('Grid', 'Init', timeit.default_timer() - st_time, level=1) + def calculate_bounds(self): + x_min = 0 + y_min = 0 + if len(self.center_latitudes.shape) == 1: + x_max = self.center_latitudes.shape[0] + y_max = self.center_longitudes.shape[0] + else: + x_max = self.center_latitudes.shape[0] + y_max = self.center_longitudes.shape[1] + y_len = y_max // self.comm.Get_size() + exceeded_ranks = y_max % self.comm.Get_size() + self.x_lower_bound = x_min + self.x_upper_bound = x_max + if self.comm.Get_rank() < exceeded_ranks: + y_len = y_len + 1 + self.y_lower_bound = self.comm.Get_rank() * y_len + self.y_upper_bound = (self.comm.Get_rank() + 1) * y_len + else: + self.y_lower_bound = (self.comm.Get_rank() * y_len) + exceeded_ranks + self.y_upper_bound = ((self.comm.Get_rank() + 1) * y_len) + exceeded_ranks + @staticmethod def create_esmf_grid_from_file(file_name, sphere=True): import ESMF @@ -94,7 +115,7 @@ class Grid(object): return grid @staticmethod - def select_grid(grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, lat_orig, + def select_grid(comm, grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, lat_orig, lon_orig, n_lat, n_lon, centre_lat, centre_lon, west_boundary, south_boundary, inc_rlat, inc_rlon, lat_1, lat_2, lon_0, lat_0, nx, ny, inc_x, inc_y, x_0, y_0, lat_ts): # TODO describe better the rotated parameters @@ -180,26 +201,27 @@ class Grid(object): # Creating a different object depending on the grid type if grid_type == 'global': from hermesv3_gr.modules.grids.grid_global import GlobalGrid - grid = GlobalGrid(grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon) + grid = GlobalGrid(grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, + comm=comm) elif grid_type == 'regular': from hermesv3_gr.modules.grids.grid_latlon import LatLonGrid grid = LatLonGrid(grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, - lat_orig, lon_orig, n_lat, n_lon) + lat_orig, lon_orig, n_lat, n_lon, comm=comm) elif grid_type == 'rotated': from hermesv3_gr.modules.grids.grid_rotated import RotatedGrid grid = RotatedGrid(grid_type, vertical_description_path, timestep_num, temporal_path, - centre_lat, centre_lon, west_boundary, south_boundary, inc_rlat, inc_rlon) + centre_lat, centre_lon, west_boundary, south_boundary, inc_rlat, inc_rlon, comm=comm) elif grid_type == 'lcc': from hermesv3_gr.modules.grids.grid_lcc import LccGrid grid = LccGrid(grid_type, vertical_description_path, timestep_num, temporal_path, lat_1, lat_2, lon_0, - lat_0, nx, ny, inc_x, inc_y, x_0, y_0) + lat_0, nx, ny, inc_x, inc_y, x_0, y_0, comm=comm) elif grid_type == 'mercator': from hermesv3_gr.modules.grids.grid_mercator import MercatorGrid grid = MercatorGrid(grid_type, vertical_description_path, timestep_num, temporal_path, lat_ts, lon_0, - nx, ny, inc_x, inc_y, x_0, y_0) + nx, ny, inc_x, inc_y, x_0, y_0, comm=comm) else: settings.write_log('ERROR: Check the .err file to get more info.') if settings.rank == 0: @@ -389,8 +411,9 @@ class Grid(object): st_time = timeit.default_timer() settings.write_log('\t\tGetting 2D coordinates from ESMPy Grid', level=3) - lat = self.esmf_grid.get_coords(1, ESMF.StaggerLoc.CENTER).T - lon = self.esmf_grid.get_coords(0, ESMF.StaggerLoc.CENTER).T + esmf_grid = self.create_esmf_grid_from_file(self.coords_netcdf_file) + lat = esmf_grid.get_coords(1, ESMF.StaggerLoc.CENTER).T + lon = esmf_grid.get_coords(0, ESMF.StaggerLoc.CENTER).T settings.write_time('Grid', 'get_coordinates_2d', timeit.default_timer() - st_time, level=3) diff --git a/hermesv3_gr/modules/grids/grid_global.py b/hermesv3_gr/modules/grids/grid_global.py index a3ca236..2f63c3c 100755 --- a/hermesv3_gr/modules/grids/grid_global.py +++ b/hermesv3_gr/modules/grids/grid_global.py @@ -54,14 +54,13 @@ class GlobalGrid(Grid): """ def __init__(self, grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, - center_longitude=float(0)): - import ESMF + center_longitude=float(0), comm=None): st_time = timeit.default_timer() settings.write_log('\tCreating Global grid.', level=2) # Initialize the class using parent - super(GlobalGrid, self).__init__(grid_type, vertical_description_path, temporal_path) + super(GlobalGrid, self).__init__(grid_type, vertical_description_path, temporal_path, comm) self.center_lat = float(0) self.center_lon = center_longitude @@ -76,12 +75,7 @@ class GlobalGrid(Grid): super(GlobalGrid, self).write_coords_netcdf() settings.comm.Barrier() - self.esmf_grid = super(GlobalGrid, self).create_esmf_grid_from_file(self.coords_netcdf_file) - - self.x_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][1] - self.x_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][1] - self.y_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][0] - self.y_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][0] + self.calculate_bounds() self.shape = (timestep_num, len(self.vertical_description), self.x_upper_bound-self.x_lower_bound, self.y_upper_bound-self.y_lower_bound) diff --git a/hermesv3_gr/modules/grids/grid_latlon.py b/hermesv3_gr/modules/grids/grid_latlon.py index 0409d7b..df0a676 100755 --- a/hermesv3_gr/modules/grids/grid_latlon.py +++ b/hermesv3_gr/modules/grids/grid_latlon.py @@ -63,14 +63,13 @@ class LatLonGrid(Grid): """ def __init__(self, grid_type, vertical_description_path, timestep_num, temporal_path, inc_lat, inc_lon, lat_orig, - lon_orig, n_lat, n_lon): - import ESMF + lon_orig, n_lat, n_lon, comm=None): st_time = timeit.default_timer() settings.write_log('\tCreating Regional regular lat-lon grid.', level=2) # Initialize the class using parent - super(LatLonGrid, self).__init__(grid_type, vertical_description_path, temporal_path) + super(LatLonGrid, self).__init__(grid_type, vertical_description_path, temporal_path, comm) self.inc_lat = inc_lat self.inc_lon = inc_lon @@ -87,12 +86,7 @@ class LatLonGrid(Grid): self.write_coords_netcdf() settings.comm.Barrier() - self.esmf_grid = self.create_esmf_grid_from_file(self.coords_netcdf_file, sphere=False) - - self.x_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][1] - self.x_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][1] - self.y_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][0] - self.y_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][0] + self.calculate_bounds() self.shape = (timestep_num, len(self.vertical_description), self.x_upper_bound-self.x_lower_bound, self.y_upper_bound-self.y_lower_bound) diff --git a/hermesv3_gr/modules/grids/grid_lcc.py b/hermesv3_gr/modules/grids/grid_lcc.py index d360d74..42aaebc 100755 --- a/hermesv3_gr/modules/grids/grid_lcc.py +++ b/hermesv3_gr/modules/grids/grid_lcc.py @@ -77,13 +77,12 @@ class LccGrid(Grid): """ def __init__(self, grid_type, vertical_description_path, timestep_num, temporal_path, lat_1, lat_2, lon_0, lat_0, - nx, ny, inc_x, inc_y, x_0, y_0, earth_radius=6370000.000): - import ESMF + nx, ny, inc_x, inc_y, x_0, y_0, earth_radius=6370000.000, comm=None): st_time = timeit.default_timer() settings.write_log('\tCreating Lambert Conformal Conic (LCC) grid.', level=2) # Initialises with parent class - super(LccGrid, self).__init__(grid_type, vertical_description_path, temporal_path) + super(LccGrid, self).__init__(grid_type, vertical_description_path, temporal_path, comm) # Setting parameters self.lat_1 = lat_1 @@ -113,12 +112,7 @@ class LccGrid(Grid): self.write_coords_netcdf() settings.comm.Barrier() - self.esmf_grid = super(LccGrid, self).create_esmf_grid_from_file(self.coords_netcdf_file, sphere=False) - # - self.x_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][1] - self.x_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][1] - self.y_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][0] - self.y_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][0] + self.calculate_bounds() self.shape = (timestep_num, len(self.vertical_description), self.x_upper_bound-self.x_lower_bound, self.y_upper_bound-self.y_lower_bound) diff --git a/hermesv3_gr/modules/grids/grid_mercator.py b/hermesv3_gr/modules/grids/grid_mercator.py index 38bb722..d18c32d 100755 --- a/hermesv3_gr/modules/grids/grid_mercator.py +++ b/hermesv3_gr/modules/grids/grid_mercator.py @@ -69,7 +69,6 @@ class MercatorGrid(Grid): def __init__(self, grid_type, vertical_description_path, timestep_num, temporal_path, lat_ts, lon_0, nx, ny, inc_x, inc_y, x_0, y_0, earth_radius=6370000.000): - import ESMF st_time = timeit.default_timer() settings.write_log('\tCreating Mercator grid.', level=2) @@ -101,12 +100,7 @@ class MercatorGrid(Grid): self.write_coords_netcdf() settings.comm.Barrier() - self.esmf_grid = super(MercatorGrid, self).create_esmf_grid_from_file(self.coords_netcdf_file, sphere=False) - # - self.x_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][1] - self.x_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][1] - self.y_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][0] - self.y_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][0] + self.calculate_bounds() self.shape = (timestep_num, len(self.vertical_description), self.x_upper_bound-self.x_lower_bound, self.y_upper_bound-self.y_lower_bound) diff --git a/hermesv3_gr/modules/grids/grid_rotated.py b/hermesv3_gr/modules/grids/grid_rotated.py index 24c0f9a..8c99036 100755 --- a/hermesv3_gr/modules/grids/grid_rotated.py +++ b/hermesv3_gr/modules/grids/grid_rotated.py @@ -40,14 +40,12 @@ class RotatedGrid(Grid): """ def __init__(self, grid_type, vertical_description_path, timestep_num, temporal_path, centre_lat, centre_lon, - west_boundary, south_boundary, inc_rlat, inc_rlon): - import ESMF - + west_boundary, south_boundary, inc_rlat, inc_rlon, comm=None): st_time = timeit.default_timer() settings.write_log('\tCreating Rotated grid.', level=2) # Initialises with parent class - super(RotatedGrid, self).__init__(grid_type, vertical_description_path, temporal_path) + super(RotatedGrid, self).__init__(grid_type, vertical_description_path, temporal_path, comm) # Setting parameters self.new_pole_longitude_degrees = -180 + centre_lon @@ -71,18 +69,10 @@ class RotatedGrid(Grid): if not os.path.exists(self.coords_netcdf_file): if settings.rank == 0: - # super(RotatedGrid, self).write_coords_netcdf() self.write_coords_netcdf() settings.comm.Barrier() - # self.write_coords_netcdf() - - self.esmf_grid = super(RotatedGrid, self).create_esmf_grid_from_file(self.coords_netcdf_file, sphere=False) - - self.x_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][1] - self.x_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][1] - self.y_lower_bound = self.esmf_grid.lower_bounds[ESMF.StaggerLoc.CENTER][0] - self.y_upper_bound = self.esmf_grid.upper_bounds[ESMF.StaggerLoc.CENTER][0] + self.calculate_bounds() self.shape = (timestep_num, len(self.vertical_description), self.x_upper_bound-self.x_lower_bound, self.y_upper_bound-self.y_lower_bound) diff --git a/hermesv3_gr/modules/regrid/regrid_conservative.py b/hermesv3_gr/modules/regrid/regrid_conservative.py index 35c718f..74c219e 100755 --- a/hermesv3_gr/modules/regrid/regrid_conservative.py +++ b/hermesv3_gr/modules/regrid/regrid_conservative.py @@ -48,7 +48,7 @@ class ConservativeRegrid(Regrid): src_field.read(filename=self.pollutant_dicts[0]['path'], variable=self.pollutant_dicts[0]['name'], timeslice=0) - dst_grid = self.grid.esmf_grid + dst_grid = self.grid.create_esmf_grid_from_file(self.grid.coords_netcdf_file) dst_field = ESMF.Field(dst_grid, name='my outut field') ESMF.Regrid(src_field, dst_field, filename=self.weight_matrix_file, regrid_method=ESMF.RegridMethod.CONSERVE,) diff --git a/hermesv3_gr/modules/speciation/speciation.py b/hermesv3_gr/modules/speciation/speciation.py index efb3061..713a484 100755 --- a/hermesv3_gr/modules/speciation/speciation.py +++ b/hermesv3_gr/modules/speciation/speciation.py @@ -182,7 +182,7 @@ class Speciation(object): raise AttributeError( "Error in speciation profile {0}: ".format(self.id) + "The output specie {0} cannot be calculated ".format(pollutant['name']) + - "with the expression {0} because{1}".format(formula, e.message)) + "with the expression {0} because{1}".format(formula, str(e))) else: sys.exit(1) speciated_emissions.append(dict_aux) diff --git a/hermesv3_gr/modules/temporal/temporal.py b/hermesv3_gr/modules/temporal/temporal.py index 5606513..6f2c457 100755 --- a/hermesv3_gr/modules/temporal/temporal.py +++ b/hermesv3_gr/modules/temporal/temporal.py @@ -70,7 +70,7 @@ class TemporalDistribution(object): """ def __init__(self, starting_date, timestep_type, timestep_num, timestep_freq, monthly_profile_path, monthly_profile_id, weekly_profile_path, weekly_profile_id, daily_profile_path, daily_profile_id, - hourly_profile_path, hourly_profile_id, world_info_path, auxiliar_files_dir, grid): + hourly_profile_path, hourly_profile_id, world_info_path, auxiliar_files_dir, grid, first_time=False): from timezonefinder import TimezoneFinder import pandas as pd @@ -80,40 +80,34 @@ class TemporalDistribution(object): self.grid = grid - # self.starting_date = starting_date + if not first_time: + self.date_array = self.caclulate_date_array(starting_date, timestep_type, timestep_num, timestep_freq) - self.date_array = self.caclulate_date_array(starting_date, timestep_type, timestep_num, timestep_freq) + if timestep_type in ['daily', 'hourly']: + self.daily_profile = self.get_daily_profile(daily_profile_id, daily_profile_path) + else: + self.daily_profile = None - # self.timestep_type = timestep_type - # self.timestep_num = timestep_num - # self.timestep_freq = timestep_freq - # - # self.ending_date = self.calculate_ending_date() - if timestep_type in ['daily', 'hourly']: - self.daily_profile = self.get_daily_profile(daily_profile_id, daily_profile_path) - else: - self.daily_profile = None + if self.daily_profile is None: + if timestep_type in ['monthly', 'daily', 'hourly']: + self.monthly_profile = self.get_monthly_profile(monthly_profile_id, monthly_profile_path) + else: + self.monthly_profile = None - if self.daily_profile is None: - if timestep_type in ['monthly', 'daily', 'hourly']: - self.monthly_profile = self.get_monthly_profile(monthly_profile_id, monthly_profile_path) + if timestep_type in ['daily', 'hourly']: + self.weekly_profile = self.get_weekly_profile(weekly_profile_id, weekly_profile_path) + else: + self.weekly_profile = None else: self.monthly_profile = None - - if timestep_type in ['daily', 'hourly']: - self.weekly_profile = self.get_weekly_profile(weekly_profile_id, weekly_profile_path) - else: self.weekly_profile = None - else: - self.monthly_profile = None - self.weekly_profile = None - if monthly_profile_id is not None or weekly_profile_id is not None: - warn('WARNING: Daily profile is set. Monthly or Weekly profiles will be ignored.') + if monthly_profile_id is not None or weekly_profile_id is not None: + warn('WARNING: Daily profile is set. Monthly or Weekly profiles will be ignored.') - if timestep_type in ['hourly']: - self.hourly_profile = self.get_hourly_profile(hourly_profile_id, hourly_profile_path) - else: - self.hourly_profile = None + if timestep_type in ['hourly']: + self.hourly_profile = self.get_hourly_profile(hourly_profile_id, hourly_profile_path) + else: + self.hourly_profile = None self.world_info = world_info_path self.netcdf_timezones = os.path.join(auxiliar_files_dir, 'timezones.nc') -- GitLab From a9e802fed48c4e64a3cd96d5ae1e48aa6938316e Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Wed, 23 Oct 2019 10:17:48 +0200 Subject: [PATCH 07/11] drop pandas warning --- hermesv3_gr/modules/grids/grid.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/hermesv3_gr/modules/grids/grid.py b/hermesv3_gr/modules/grids/grid.py index b3acdb9..59f42be 100755 --- a/hermesv3_gr/modules/grids/grid.py +++ b/hermesv3_gr/modules/grids/grid.py @@ -534,7 +534,7 @@ class Grid(object): # Make a list of list of tuples # [[(point_1.1), (point_1.2), (point_1.3), (point_1.4)], # [(point_2.1), (point_2.2), (point_2.3), (point_2.4)], ...] - list_points = df.as_matrix() + list_points = df.values del df['p1'], df['p2'], df['p3'], df['p4'] # List of polygons from the list of points -- GitLab From 17d596f304ae8bf75f5f42e21c8236c466c8ab83 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Wed, 23 Oct 2019 11:52:55 +0200 Subject: [PATCH 08/11] bug fix --- hermesv3_gr/hermes.py | 74 ++++++++++--------- .../point_gfas_emission_inventory.py | 11 ++- hermesv3_gr/modules/grids/grid.py | 42 +---------- 3 files changed, 45 insertions(+), 82 deletions(-) diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 8f67bef..8718644 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -99,41 +99,45 @@ class Hermes(object): """ from datetime import timedelta - st_time = timeit.default_timer() - settings.write_log('') - settings.write_log('***** Starting HERMESv3 *****') - num = 1 - for ei in self.emission_list: - settings.write_log('Processing emission inventory {0} for the sector {1} ({2}/{3}):'.format( - ei.inventory_name, ei.sector, num, len(self.emission_list))) - num += 1 - - ei.do_regrid() - - if ei.vertical is not None: - settings.write_log("\tCalculating vertical distribution.", level=2) - if ei.source_type == 'area': - ei.vertical_factors = ei.vertical.calculate_weights() - elif ei.source_type == 'point': - ei.calculate_altitudes(self.options.vertical_description) - ei.point_source_by_cell() - # To avoid use point source as area source when is going to apply vertical factors while writing. - ei.vertical = None - else: - settings.write_log('ERROR: Check the .err file to get more info.') - if settings.rank == 0: - raise AttributeError('Unrecognized emission source type {0}'.format(ei.source_type)) - sys.exit(1) - - if ei.temporal is not None: - ei.temporal_factors = ei.temporal.calculate_3d_temporal_factors() - if ei.speciation is not None: - ei.emissions = ei.speciation.do_speciation(ei.emissions) - - self.writer.write(self.emission_list) - - settings.write_log("***** HERMESv3 simulation finished successfully *****") - settings.write_time('HERMES', 'main', timeit.default_timer() - st_time) + if self.options.first_time: + # Stop run + settings.write_log('***** HERMESv3_BU First Time finished successfully *****') + else: + st_time = timeit.default_timer() + settings.write_log('') + settings.write_log('***** Starting HERMESv3_GR *****') + num = 1 + for ei in self.emission_list: + settings.write_log('Processing emission inventory {0} for the sector {1} ({2}/{3}):'.format( + ei.inventory_name, ei.sector, num, len(self.emission_list))) + num += 1 + + ei.do_regrid() + + if ei.vertical is not None: + settings.write_log("\tCalculating vertical distribution.", level=2) + if ei.source_type == 'area': + ei.vertical_factors = ei.vertical.calculate_weights() + elif ei.source_type == 'point': + ei.calculate_altitudes(self.options.vertical_description) + ei.point_source_by_cell() + # To avoid use point source as area source when is going to apply vertical factors while writing. + ei.vertical = None + else: + settings.write_log('ERROR: Check the .err file to get more info.') + if settings.rank == 0: + raise AttributeError('Unrecognized emission source type {0}'.format(ei.source_type)) + sys.exit(1) + + if ei.temporal is not None: + ei.temporal_factors = ei.temporal.calculate_3d_temporal_factors() + if ei.speciation is not None: + ei.emissions = ei.speciation.do_speciation(ei.emissions) + + self.writer.write(self.emission_list) + + settings.write_log("***** HERMESv3_GR simulation finished successfully *****") + settings.write_time('HERMES', 'main', timeit.default_timer() - st_time) settings.write_time('HERMES', 'TOTAL', timeit.default_timer() - full_time) settings.finish_logs(self.options.output_dir, self.options.start_date) diff --git a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py b/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py index 10ca3a1..78167bb 100755 --- a/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py +++ b/hermesv3_gr/modules/emision_inventories/point_gfas_emission_inventory.py @@ -84,14 +84,13 @@ class PointGfasEmissionInventory(EmissionInventory): reference_year=reference_year, factors=factors, regrid_mask=regrid_mask, p_vertical=None, p_month=p_month, p_week=p_week, p_day=p_day, p_hour=p_hour, p_speciation=p_speciation) - if not options.first_time: - # self.approach = self.get_approach(p_vertical) - self.method = self.get_method(p_vertical) + # self.approach = self.get_approach(p_vertical) + self.method = self.get_method(p_vertical) - # self.altitude = self.get_altitude() + # self.altitude = self.get_altitude() - self.vertical = GfasVerticalDistribution(vertical_output_profile, self.get_approach(p_vertical), - 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) diff --git a/hermesv3_gr/modules/grids/grid.py b/hermesv3_gr/modules/grids/grid.py index 59f42be..520b56b 100755 --- a/hermesv3_gr/modules/grids/grid.py +++ b/hermesv3_gr/modules/grids/grid.py @@ -450,21 +450,13 @@ class Grid(object): settings.write_log('\t\tGrid shapefile not done. Lets try to create it.', level=3) # Create Shapefile - # Use the meters coordiantes to create the shapefile - y = self.boundary_latitudes x = self.boundary_longitudes - # sys.exit() if self.grid_type in ['global', 'regular']: x = x.reshape((x.shape[1], x.shape[2])) y = y.reshape((y.shape[1], y.shape[2])) - # x_aux = np.empty((x.shape[0], y.shape[0], 4)) - # x_aux[:, :, 0] = x[:, np.newaxis, 0] - # x_aux[:, :, 1] = x[:, np.newaxis, 1] - # x_aux[:, :, 2] = x[:, np.newaxis, 1] - # x_aux[:, :, 3] = x[:, np.newaxis, 0] aux_shape = (y.shape[0], x.shape[0], 4) x_aux = np.empty(aux_shape) x_aux[:, :, 0] = x[np.newaxis, :, 0] @@ -475,12 +467,6 @@ class Grid(object): x = x_aux del x_aux - # y_aux = np.empty((x.shape[0], y.shape[0], 4)) - # y_aux[:, :, 0] = y[np.newaxis, :, 0] - # y_aux[:, :, 1] = y[np.newaxis, :, 0] - # y_aux[:, :, 2] = y[np.newaxis, :, 1] - # y_aux[:, :, 3] = y[np.newaxis, :, 1] - y_aux = np.empty(aux_shape) y_aux[:, :, 0] = y[:, np.newaxis, 0] y_aux[:, :, 1] = y[:, np.newaxis, 0] @@ -490,8 +476,6 @@ class Grid(object): y = y_aux del y_aux - # exit() - if not full_grid: y = y[self.x_lower_bound:self.x_upper_bound, self.y_lower_bound:self.y_upper_bound, :] x = x[self.x_lower_bound:self.x_upper_bound, self.y_lower_bound:self.y_upper_bound, :] @@ -499,29 +483,12 @@ class Grid(object): aux_b_lats = y.reshape((y.shape[0] * y.shape[1], y.shape[2])) aux_b_lons = x.reshape((x.shape[0] * x.shape[1], x.shape[2])) - # The regular lat-lon projection has only 2 (laterals) points for each cell instead of 4 (corners) - # if aux_b_lats.shape[1] == 2: - # aux_b = np.empty((aux_b_lats.shape[0], 4)) - # aux_b[:, 0] = aux_b_lats[:, 0] - # aux_b[:, 1] = aux_b_lats[:, 0] - # aux_b[:, 2] = aux_b_lats[:, 1] - # aux_b[:, 3] = aux_b_lats[:, 1] - # aux_b_lats = aux_b - # - # if aux_b_lons.shape[1] == 2: - # aux_b = np.empty((aux_b_lons.shape[0], 4)) - # aux_b[:, 0] = aux_b_lons[:, 0] - # aux_b[:, 1] = aux_b_lons[:, 1] - # aux_b[:, 2] = aux_b_lons[:, 1] - # aux_b[:, 3] = aux_b_lons[:, 0] - # aux_b_lons = aux_b - # Create one dataframe with 8 columns, 4 points with two coordinates each one df_lats = pd.DataFrame(aux_b_lats, columns=['b_lat_1', 'b_lat_2', 'b_lat_3', 'b_lat_4']) df_lons = pd.DataFrame(aux_b_lons, columns=['b_lon_1', 'b_lon_2', 'b_lon_3', 'b_lon_4']) df = pd.concat([df_lats, df_lons], axis=1) - # Substituate 8 columns by 4 with the two coordinates + # Substitute 8 columns by 4 with the two coordinates df['p1'] = list(zip(df.b_lon_1, df.b_lat_1)) del df['b_lat_1'], df['b_lon_1'] df['p2'] = list(zip(df.b_lon_2, df.b_lat_2)) @@ -539,13 +506,6 @@ class Grid(object): # List of polygons from the list of points geometry = [Polygon(list(points)) for points in list_points] - # geometry = [] - # for point in list_points: - # print point - # geometry.append(Polygon(list(point))) - # print geometry[0] - # sys.exit() - # print len(geometry), len(df), gdf = gpd.GeoDataFrame(df, crs={'init': 'epsg:4326'}, geometry=geometry) gdf = gdf.to_crs(self.crs) -- GitLab From 71cdd3d9b8200bc667070a1bfc32292c9ec40220 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Wed, 23 Oct 2019 17:59:21 +0200 Subject: [PATCH 09/11] wip --- hermesv3_gr/hermes.py | 8 +++++--- hermesv3_gr/modules/writing/writer.py | 25 ++++++++++++++++++++++--- 2 files changed, 27 insertions(+), 6 deletions(-) diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 8718644..0dcf0e0 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -93,7 +93,7 @@ class Hermes(object): settings.write_time('HERMES', 'Init', timeit.default_timer() - st_time, level=1) # @profile - def main(self): + def main(self, return_emis=True): """ Main functionality of the model. """ @@ -133,8 +133,10 @@ class Hermes(object): ei.temporal_factors = ei.temporal.calculate_3d_temporal_factors() if ei.speciation is not None: ei.emissions = ei.speciation.do_speciation(ei.emissions) - - self.writer.write(self.emission_list) + if return_emis: + return self.writer.get_emis(self.emission_list) + else: + self.writer.write(self.emission_list) settings.write_log("***** HERMESv3_GR simulation finished successfully *****") settings.write_time('HERMES', 'main', timeit.default_timer() - st_time) diff --git a/hermesv3_gr/modules/writing/writer.py b/hermesv3_gr/modules/writing/writer.py index 8cb46f5..e12df93 100755 --- a/hermesv3_gr/modules/writing/writer.py +++ b/hermesv3_gr/modules/writing/writer.py @@ -25,6 +25,7 @@ from mpi4py import MPI from netCDF4 import Dataset from hermesv3_gr.config import settings from functools import reduce +from pandas import DataFrame, MultiIndex class Writer(object): @@ -186,7 +187,7 @@ class Writer(object): return True - def calculate_data_by_var(self, variable, inventory_list, shape): + def calculate_data_by_var(self, variable, inventory_list, shape, change_units=True): """ Calculate the date of the given variable throw the inventory list. @@ -243,7 +244,8 @@ class Writer(object): settings.write_time('TemporalDistribution', 'calculate_data_by_var', timeit.default_timer() - temporal_time, level=2) # Unit changes - data = self.unit_change(variable, data) + if change_units: + data = self.unit_change(variable, data) if data is not None: data[data < 0] = 0 settings.write_time('Writer', 'calculate_data_by_var', timeit.default_timer() - st_time, level=3) @@ -253,7 +255,7 @@ class Writer(object): """ Implement on inner class """ - return np.array([0]) + return data @staticmethod def calculate_displacements(counts): @@ -603,3 +605,20 @@ class Writer(object): netcdf.setncatts(global_attributes) netcdf.close() + + def get_emis(self, emission_list): + self.set_variable_attributes(emission_list) + self.change_variable_attributes() + data = DataFrame(columns=self.variables_attributes.keys(), + index=['FID', 'layer', 'tstep']) + for var_name in self.variables_attributes.keys(): + print(var_name) + var_data = self.calculate_data_by_var(var_name, emission_list, self.grid.shape, change_units=False) + if var_data is None: + data[var_name] = 0 + else: + for tstep in range(data.shape[0]): + for layer in range(data.shape[1]): + + print(data.shape) + return None -- GitLab From 7ef41b3f8b266862140dbbba64ce82084ff8e4a8 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Fri, 25 Oct 2019 17:00:05 +0200 Subject: [PATCH 10/11] return_emis in memory done --- hermesv3_gr/hermes.py | 6 ++--- hermesv3_gr/modules/writing/writer.py | 38 ++++++++++++++++++++------- 2 files changed, 32 insertions(+), 12 deletions(-) diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index 0dcf0e0..c7f6053 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -93,7 +93,7 @@ class Hermes(object): settings.write_time('HERMES', 'Init', timeit.default_timer() - st_time, level=1) # @profile - def main(self, return_emis=True): + def main(self, return_emis=False): """ Main functionality of the model. """ @@ -149,10 +149,10 @@ class Hermes(object): return None -def run(comm=None): +def run(): date = Hermes(Config()).main() while date is not None: - date = Hermes(Config(comm), new_date=date, comm=comm).main() + date = Hermes(Config(), new_date=date).main() sys.exit(0) diff --git a/hermesv3_gr/modules/writing/writer.py b/hermesv3_gr/modules/writing/writer.py index e12df93..535d50e 100755 --- a/hermesv3_gr/modules/writing/writer.py +++ b/hermesv3_gr/modules/writing/writer.py @@ -25,6 +25,7 @@ from mpi4py import MPI from netCDF4 import Dataset from hermesv3_gr.config import settings from functools import reduce +import pandas as pd from pandas import DataFrame, MultiIndex @@ -609,16 +610,35 @@ class Writer(object): def get_emis(self, emission_list): self.set_variable_attributes(emission_list) self.change_variable_attributes() - data = DataFrame(columns=self.variables_attributes.keys(), - index=['FID', 'layer', 'tstep']) + + data_list = [] + + fid_array = np.array(range(self.grid.full_shape[-2] * self.grid.full_shape[-1])) + fid_array = fid_array.reshape((self.grid.full_shape[-2], self.grid.full_shape[-1])) + fid_array = fid_array[self.grid.x_lower_bound:self.grid.x_upper_bound, + self.grid.y_lower_bound:self.grid.y_upper_bound].flatten() + + zero_vars = [] for var_name in self.variables_attributes.keys(): - print(var_name) var_data = self.calculate_data_by_var(var_name, emission_list, self.grid.shape, change_units=False) if var_data is None: - data[var_name] = 0 + zero_vars.append(var_name) else: - for tstep in range(data.shape[0]): - for layer in range(data.shape[1]): - - print(data.shape) - return None + var_data_list = [] + for tstep in range(var_data.shape[0]): + for layer in range(var_data.shape[1]): + aux_data = DataFrame(data=var_data[tstep, layer, :, :].flatten(), columns=[var_name], + index=MultiIndex.from_product([fid_array, [layer], [tstep]], + names=['FID', 'layer', 'tstep'])) + aux_data = aux_data.replace(0, np.nan) + aux_data.dropna(how='all', inplace=True) + var_data_list.append(aux_data) + data_list.append(pd.concat(var_data_list, sort=False)) + + data = pd.concat(data_list, axis=1, sort=False) + data = data.replace(np.nan, 0) + + for var_name in zero_vars: + data[var_name] = 0 + + return data -- GitLab From 37a1d5f72dfedf73387bf08d72f2955e52095844 Mon Sep 17 00:00:00 2001 From: Carles Tena Date: Fri, 25 Oct 2019 17:12:53 +0200 Subject: [PATCH 11/11] PEP8 corrections --- hermesv3_gr/hermes.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/hermesv3_gr/hermes.py b/hermesv3_gr/hermes.py index c7f6053..bb047b9 100755 --- a/hermesv3_gr/hermes.py +++ b/hermesv3_gr/hermes.py @@ -121,7 +121,7 @@ class Hermes(object): elif ei.source_type == 'point': ei.calculate_altitudes(self.options.vertical_description) ei.point_source_by_cell() - # To avoid use point source as area source when is going to apply vertical factors while writing. + # To avoid use point source as area source when is going to apply vertical factors while writing ei.vertical = None else: settings.write_log('ERROR: Check the .err file to get more info.') -- GitLab