Commit d718d49b authored by Cecilia Nievas's avatar Cecilia Nievas
Browse files

Further adjustments of style to SERA_mapping_admin_units_to_cells.py

parent cdcfa89d
# Byte-compiled / optimized / DLL files
__pycache__/
.cache
*.py[cod]
# C extensions
......
......@@ -60,7 +60,8 @@ def run_this_file(config_dict):
DB_schema_name_grid = config_dict['Tiles Database']['db_tiles_schema_name']
DB_table_name_grid = config_dict['Tiles Database']['db_tiles_table_name']
# Countries to process:
countries = config_dict['SERA_mapping_admin_units_to_cells']['countries'].split(', ')
countries = \
config_dict['SERA_mapping_admin_units_to_cells']['countries'].split(', ')
####################################################
# SET UP LOGGER
####################################################
......@@ -74,16 +75,16 @@ def run_this_file(config_dict):
log_fn = 'log_file_mapping_adm_units_%s.log' % (now_str)
fh = logging.FileHandler(os.path.join(out_path, log_fn))
logger = logging.getLogger(__name__)
logging.basicConfig(level=logging.INFO)
####################################################
# START
####################################################
num_countries = len(countries)
for k,country_name in enumerate(countries):
for k, country_name in enumerate(countries):
for result in process_country(country_name, k, num_countries,
sera_models_path, sera_shp_path, logger):
# Unpack results:
cell_id, adm_level, admin_id, area_intersect, geom_intersect, cases_in_lev_grs, shp_not_found = result
cell_id, adm_level, admin_id, area_intersect, \
geom_intersect, cases_in_lev_grs, shp_nt_fnd = result
# Write entry to database:
err_msg = gdet_psql.insert_data_in_cell_admin_unit_mapping(
......@@ -99,10 +100,11 @@ def run_this_file(config_dict):
geom_intersect,
cases_in_lev_grs)
if len(err_msg)>0:
if len(err_msg) > 0:
logger.info(country_name+': '+err_msg)
logger.debug('\r %s' % (gdet_gral.add_to_log_shapefiles_found(country_name,
shp_not_found)))
logger.debug('\r %s'
% (gdet_gral.add_to_log_shapefiles_found(country_name,
shp_nt_fnd)))
def process_country(country_str, country_pos, number_countries,
......@@ -127,23 +129,24 @@ def process_country(country_str, country_pos, number_countries,
# definition (given by level_groups_levels):
for j in range(len(level_groups)):
try:
adm_boundary = gpd.GeoDataFrame.from_file(os.path.join(sera_shape_path,
'Adm%i_%s.shp'
% (level_groups_levels[j], country_str)))
adm_boundary = gpd.GeoDataFrame.from_file(
os.path.join(sera_shape_path,
'Adm%i_%s.shp'
% (level_groups_levels[j], country_str)))
except:
logger.info(' Shapefile not found for Admin Level %i'
% (level_groups_levels[j]))
# Keep track of admin IDs not found in the shapefile:
shp_not_found= 0
shp_not_found = 0
# Keep track of admin units that checked, because not all units might
# exist in each case (Res, Com, Ind):
admin_units_investigated= []
admin_units_investigated = []
# Go by case (Res, Com or Ind):
for case in level_groups[j]:
#unique_admins= np.unique(dfs_dict[case]['id_%i' % (level_groups_levels[j])].values)
unique_admins= np.unique(dfs_dict[case]['id_%i' % (level_groups_levels[j])].values[np.where(dfs_dict[case]['id_2']==4)[0]])
# unique_admins= np.unique(dfs_dict[case]['id_%i' % (level_groups_levels[j])].values)
unique_admins= np.unique(dfs_dict[case]['id_%i' % (level_groups_levels[j])].values[np.where(dfs_dict[case]['id_2'] == 4)[0]])
# Go one by one the admin units named in the exposure CSV file:
for m, admin_id in enumerate(unique_admins):
......@@ -151,29 +154,36 @@ def process_country(country_str, country_pos, number_countries,
continue
admin_units_investigated.append(admin_id)
# If the admin id cannot be found in the corresponding shapefile:
# If the admin id cannot be found in the corresponding shapef:
if admin_id not in adm_boundary['ID_%i' % (level_groups_levels[j])].values:
shp_not_found= shp_not_found + 1
shp_not_found = shp_not_found + 1
continue
# Identify row of this administrative unit in the adm_boundary DF:
# Identify row of this administrative unit in adm_boundary DF:
print(' %s, level group %i of %i, %s (%i associated admin IDs): Working on admin ID %i (position %i)'
% (country_str, j+1, len(level_groups), case, len(unique_admins), admin_id, m+1))
#logger.info(' %s, level group %i of %i, %s (%i associated admin IDs): Working on admin ID %i (position %i)'
# % (country_str, j+1, len(level_groups), case, len(unique_admins), admin_id, m+1))
which_row= np.where(adm_boundary['ID_%i' % (level_groups_levels[j])].values==admin_id)[0]
% (country_str, j+1, len(level_groups), case, len(unique_admins), admin_id, m+1))
# logger.info(' %s, level group %i of %i, %s (%i associated admin IDs): Working on admin ID %i (position %i)'
# % (country_str, j+1, len(level_groups), case, len(unique_admins), admin_id, m+1))
which_row = np.where(
adm_boundary['ID_%i' % (level_groups_levels[j])].values == admin_id)[0]
if len(which_row) != 1:
logger.info(gdet_gral.add_to_log_too_many_boundaries_found(country_str, len(which_row), admin_id))
logger.info(gdet_gral.add_to_log_too_many_boundaries_found(
country_str, len(which_row), admin_id))
continue
# Pass on adm_boundary and the row
# to process all cells associated with this admin unit:
for results_adm in process_admin_unit(adm_boundary, which_row, logger):
for results_adm in process_admin_unit(adm_boundary,
which_row, logger):
cell_id, area_intersect, geom_intersect = results_adm
yield cell_id, level_groups_levels[j], admin_id, area_intersect, geom_intersect, level_groups[j], shp_not_found
yield cell_id, level_groups_levels[j], admin_id, \
area_intersect, geom_intersect, level_groups[j], \
shp_not_found
break
print('\n')
def process_admin_unit(adm_boundary_gdf, which_row_of_gdf, logger):
def process_admin_unit(adm_boundary_gdf, which_row_of_gdf, logger=None):
"""
adm_boundary_gdf = GeoPandas DataFrame with admin boundaries
which_row_of_gdf = row of adm_boundary_gdf to be processed
......@@ -185,31 +195,38 @@ def process_admin_unit(adm_boundary_gdf, which_row_of_gdf, logger):
bbox[2], bbox[3])
# gdf_intersect has as many rows as cells do intersect this admin unit:
gdf_intersect = gpd.sjoin(gdf_cells,
adm_boundary_gdf.iloc[[which_row_of_gdf[0]],:],
adm_boundary_gdf.iloc[[which_row_of_gdf[0]], :],
how='right')
# Loop through all cells that intersect the admin boundary:
for i in range(gdf_intersect.shape[0]):
print('\r Intersection of cell '+str(gdf_intersect['cell_id'].values[i]), end='')
#logger.info('\r Intersection of cell %s'
# % (gdf_intersect['cell_id'].values[i]))
print('\r Intersection of cell '+str(gdf_intersect['cell_id'].values[i]),
end='')
# logger.info('\r Intersection of cell %s'
# % (gdf_intersect['cell_id'].values[i]))
geom_cell = gdf_cells['geometry'].values[gdf_intersect['index_left'].values[i]]
geom_adm_unit = adm_boundary_gdf['geometry'].values[which_row_of_gdf[0]]
geom_intersection = geom_adm_unit.intersection(geom_cell)
area_intersection = gdet_grid.get_area_polygon_on_Earth(geom_intersection,
units='km2')
yield gdf_intersect['cell_id'].values[i], area_intersection, geom_intersection
geom_intersect = geom_adm_unit.intersection(geom_cell)
area_intersect = gdet_grid.get_area_polygon_on_Earth(geom_intersect,
units='km2')
yield gdf_intersect['cell_id'].values[i], \
area_intersect, geom_intersect
if __name__ == '__main__':
logging.basicConfig(level=logging.INFO)
if __name__=='__main__':
# This code needs to be run from the command line as python3 namefile.py configfile.ini
# This code needs to be run from the command line as
# python3 namefile.py configfile.ini
# sys.argv retrieves all the commands entered in the command line;
# position [0] is this code, position [1] will be the config file name
config_filename = sys.argv[1]
section_names_to_validate = ['File Paths',
'Tiles Database',
os.path.basename(__file__).split('.')[0]]
config_dict = gdet_conf.read_config_parameters(os.path.join(os.getcwd(),
config_filename),
section_names_to_validate)
config_dict = gdet_conf.read_config_parameters(
os.path.join(os.getcwd(),
config_filename),
section_names_to_validate)
run_this_file(config_dict)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment