diff --git a/dggrid4py/dggrid_runner.py b/dggrid4py/dggrid_runner.py index b45b068..654bc63 100644 --- a/dggrid4py/dggrid_runner.py +++ b/dggrid4py/dggrid_runner.py @@ -30,12 +30,10 @@ import geopandas as gpd from .interrupt import crosses_interruption, interrupt_cell, get_geom_coords -if TYPE_CHECKING: - # place here to avoid imposing 'shapely' if not installed, but supported interface - from geopandas.base import GeometryArray - from shapely.geometry.base import BaseGeometry # noqa +from geopandas.base import GeometryArray +from shapely.geometry.base import BaseGeometry # noqa - AnyGeometry = BaseGeometry | GeometryArray +AnyGeometry = BaseGeometry | GeometryArray fiona_drivers = fiona.supported_drivers @@ -818,42 +816,46 @@ def dgapi_grid_gen( metafile.append(cmd) # clip_subset_types - if subset_conf['clip_subset_type'] == 'WHOLE_EARTH': - metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) - elif subset_conf['clip_subset_type'] in [ 'SHAPEFILE' , 'AIGEN', 'GDAL'] and not subset_conf['clip_region_files'] is None: - metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) - metafile.append("clip_region_files " + subset_conf['clip_region_files']) - elif subset_conf['clip_subset_type'] in [ 'SEQNUMS', 'COARSE_CELLS', 'INPUT_ADDRESS_TYPE'] and not subset_conf['clip_region_files'] is None: - if not dggs.dggs_type in ['PLANETRISK']: - # metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) - # metafile.append("clip_region_files " + subset_conf['clip_region_files']) - # COARSE_CELLS needs clip_cell_res - for elem in filter(lambda x: x.startswith('clip_') , subset_conf.keys()): - metafile.append(f"{elem} " + str(subset_conf[elem])) - - else: - # if dggs_aperture_type would be SEQUENCE - # have to reset to WHOLE_EARTH and clip based on - # output_first_seqnum and output_last_seqnum - # this should now almost never happen - subset_conf['clip_subset_type'] = 'WHOLE_EARTH' + if subset_conf.get("clip_subset_type"): + if subset_conf['clip_subset_type'] == 'WHOLE_EARTH': metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) - # loading seqnums - files = subset_conf['clip_region_files'].split(' ') - seqnums = [] - for file in files: - seqs = pd.read_csv( file, header=None)[0].values - seqnums.append(seqs.min()) - seqnums.append(seqs.max()) - - first_seqnum = np.array(seqnums).min() - last_seqnum = np.array(seqnums).max() - subset_conf['output_first_seqnum'] = first_seqnum - metafile.append("output_first_seqnum " + str(subset_conf['output_first_seqnum'])) - subset_conf['output_last_seqnum'] = last_seqnum - metafile.append("output_last_seqnum " + str(subset_conf['output_last_seqnum'])) - else: - raise ValueError('something is not correct in subset_conf') + elif subset_conf['clip_subset_type'] in [ 'SHAPEFILE' , 'AIGEN', 'GDAL'] and not subset_conf['clip_region_files'] is None: + metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) + metafile.append("clip_region_files " + subset_conf['clip_region_files']) + elif subset_conf['clip_subset_type'] in [ 'SEQNUMS', 'COARSE_CELLS', 'INPUT_ADDRESS_TYPE'] and not subset_conf['clip_region_files'] is None: + if not dggs.dggs_type in ['PLANETRISK']: + # metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) + # metafile.append("clip_region_files " + subset_conf['clip_region_files']) + # COARSE_CELLS needs clip_cell_res + for elem in filter(lambda x: x.startswith('clip_') , subset_conf.keys()): + metafile.append(f"{elem} " + str(subset_conf[elem])) + + else: + # if dggs_aperture_type would be SEQUENCE + # have to reset to WHOLE_EARTH and clip based on + # output_first_seqnum and output_last_seqnum + # this should now almost never happen + subset_conf['clip_subset_type'] = 'WHOLE_EARTH' + metafile.append("clip_subset_type " + subset_conf['clip_subset_type']) + # loading seqnums + files = subset_conf['clip_region_files'].split(' ') + seqnums = [] + for file in files: + seqs = pd.read_csv( file, header=None)[0].values + seqnums.append(seqs.min()) + seqnums.append(seqs.max()) + + first_seqnum = np.array(seqnums).min() + last_seqnum = np.array(seqnums).max() + subset_conf['output_first_seqnum'] = first_seqnum + metafile.append("output_first_seqnum " + str(subset_conf['output_first_seqnum'])) + subset_conf['output_last_seqnum'] = last_seqnum + metafile.append("output_last_seqnum " + str(subset_conf['output_last_seqnum'])) + else: + raise ValueError('something is not correct in subset_conf') + # clipper_scale_factor + if subset_conf.get("clipper_scale_factor"): + metafile.append("clipper_scale_factor " + subset_conf['clipper_scale_factor']) if 'input_address_type' in subset_conf.keys() and subset_conf.get('input_address_type', 'NOPE') in self.input_address_types: metafile.append("input_address_type " + subset_conf['input_address_type']) @@ -1307,7 +1309,7 @@ def grid_cell_polygons_for_extent( dggs_type: DggsTypeT, resolution: int, mixed_aperture_level: int = None, - clip_geom: "AnyGeometry" = None, + clip_geom: AnyGeometry = None, split_dateline: bool = False, output_address_type: DggsOutputAddressTypeT | None = None, **conf_extra: DggridMetaConfigParameterT, @@ -1322,24 +1324,11 @@ def grid_cell_polygons_for_extent( dggs = dgselect(dggs_type = dggs_type, res= resolution, mixed_aperture_level=mixed_aperture_level) dggs.update(**conf_extra, strict=True) - subset_conf: DggridMetaConfigT = { 'update_frequency': 100000, 'clip_subset_type': 'WHOLE_EARTH' } - - if not clip_geom is None and clip_geom.area > 0: - - clip_gdf = gpd.GeoDataFrame(pd.DataFrame({'id' : [1], 'geometry': [clip_geom]}), geometry='geometry', crs=4326) - clip_path = Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}" - clip_gdf.to_file(str(clip_path), driver=self.tmp_geo_out['driver']) - - subset_conf.update({ - 'clip_subset_type': 'GDAL', - 'clip_region_files': str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}").resolve()), - }) - - if not self.has_gdal: - subset_conf.update({ - 'clip_subset_type': cast(DggsClipSubsetTypeT, self.tmp_geo_out['legacy_driver'].upper()) - }) - + subset_conf: DggridMetaConfigT = {} + clipper_metafile_settings, _ = specify_clip_settings(tmp_dir=tmp_dir, tmp_id=tmp_id, tmp_geo_out=self.tmp_geo_out, + clip_geom=clip_geom, has_gdal=self.has_gdal, resolution=resolution, + **conf_extra) + subset_conf.update(clipper_metafile_settings) subset_conf.update(specify_resolution(**conf_extra)) subset_conf.update(specify_orient_type_args(**conf_extra)) subset_conf.update(specify_topo_aperture(**conf_extra)) @@ -1368,7 +1357,6 @@ def grid_cell_polygons_for_extent( else: if self.debug: print(f"ignoring unknown output_address_type: {output_address_type}") - dggs_ops = self.dgapi_grid_gen(dggs, subset_conf, output_conf) if self.debug: print(dggs_ops) @@ -1406,7 +1394,7 @@ def grid_cell_centroids_for_extent( dggs_type: DggsTypeT, resolution: int, mixed_aperture_level: int = None, - clip_geom: "AnyGeometry" = None, + clip_geom: AnyGeometry = None, output_address_type: DggsOutputAddressTypeT | None = None, **conf_extra: DggridMetaConfigParameterT, ) -> gpd.GeoDataFrame: @@ -1420,24 +1408,11 @@ def grid_cell_centroids_for_extent( dggs = dgselect(dggs_type = dggs_type, res= resolution, mixed_aperture_level=mixed_aperture_level) dggs.update(**conf_extra, strict=True) - subset_conf: DggridMetaConfigT = { 'update_frequency': 100000, 'clip_subset_type': 'WHOLE_EARTH' } - - if not clip_geom is None and clip_geom.area > 0: - - clip_gdf = gpd.GeoDataFrame(pd.DataFrame({'id' : [1], 'geometry': [clip_geom]}), geometry='geometry', crs=4326) - clip_path = Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}" - clip_gdf.to_file(str(clip_path), driver=self.tmp_geo_out['driver'] ) - - subset_conf.update({ - 'clip_subset_type': 'GDAL', - 'clip_region_files': str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}").resolve()), - }) - - if not self.has_gdal: - subset_conf.update({ - 'clip_subset_type': cast(DggsClipSubsetTypeT, self.tmp_geo_out['legacy_driver'].upper()) - }) - + subset_conf: DggridMetaConfigT = {} + clipper_metafile_settings, _ = specify_clip_settings(tmp_dir=tmp_dir, tmp_id=tmp_id, tmp_geo_out=self.tmp_geo_out, + clip_geom=clip_geom, has_gdal=self.has_gdal, resolution=resolution, + **conf_extra) + subset_conf.update(clipper_metafile_settings) subset_conf.update(specify_resolution(**conf_extra)) subset_conf.update(specify_orient_type_args(**conf_extra)) subset_conf.update(specify_topo_aperture(**conf_extra)) @@ -1519,48 +1494,12 @@ def grid_cell_polygons_from_cellids( dggs = dgselect(dggs_type = dggs_type, res= resolution, mixed_aperture_level=mixed_aperture_level) dggs.update(**conf_extra, strict=True) - subset_conf: DggridMetaConfigT = { 'update_frequency': 100000, 'clip_subset_type': clip_subset_type } + subset_conf: DggridMetaConfigT = {} seq_df = None - - if not cell_id_list is None and len(cell_id_list) > 0: - - seq_df = pd.DataFrame({ input_address_type: cell_id_list}) - seq_df.to_csv( str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve()) , header=False, index=False, columns=[input_address_type], sep=' ') - - subset_conf.update({ - 'clip_subset_type': 'SEQNUMS', - 'clip_region_files': str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve()), - }) - - # TODO, for Z3, Z7, ZORDER can potentially also be COARSE_CELLS / aka parent cells? - # clip_subset_type should INPUT_ADDRESS_TYPE for the equivalent of SEQNUM (tp use input_address_type Z3 ...), or COARSE_CELLS as an actual parent cell type clip (also for Z3 ..) - if ( - input_address_type in self.input_address_types - and output_address_type in self.output_address_types - ): - subset_conf.update( - { - 'clip_subset_type': 'INPUT_ADDRESS_TYPE', - 'input_address_type': input_address_type - } - ) - - if not clip_subset_type is None and clip_subset_type in ['COARSE_CELLS']: - subset_conf.update( - { - 'clip_subset_type': clip_subset_type, - 'input_address_type': input_address_type, - 'clip_cell_res' : clip_cell_res, - 'clip_cell_addresses' : " ".join([str(address) for address in cell_id_list]) - } - ) - if "clip_cell_densification" in conf_extra: - subset_conf.update( - { - 'clip_cell_densification' : conf_extra['clip_cell_densification'] - } - ) - + clip_metafile_settings, seq_df = specify_clip_settings(clip_subset_type, tmp_dir, tmp_id, input_address_type=input_address_type, + clip_cell_res=clip_cell_res, cell_id_list=cell_id_list, + resolution=resolution, **conf_extra) + subset_conf.update(clip_metafile_settings) subset_conf.update(specify_resolution(**conf_extra)) subset_conf.update(specify_orient_type_args(**conf_extra)) subset_conf.update(specify_topo_aperture(**conf_extra)) @@ -1666,43 +1605,11 @@ def grid_cell_centroids_from_cellids( dggs = dgselect(dggs_type = dggs_type, res= resolution, mixed_aperture_level=mixed_aperture_level) dggs.update(**conf_extra, strict=True) - subset_conf: DggridMetaConfigT = { 'update_frequency': 100000, 'clip_subset_type': clip_subset_type } + subset_conf: DggridMetaConfigT = {} seq_df = None - - if not cell_id_list is None and len(cell_id_list) > 0: - - seq_df = pd.DataFrame({ input_address_type: cell_id_list}) - seq_df.to_csv( str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve()) , header=False, index=False, columns=[input_address_type], sep=' ') - - # TODO, for Z3, Z7, ZORDER can potentially also be COARSE_CELLS / aka parent cells? - subset_conf.update({ - 'clip_subset_type': 'SEQNUMS', - 'clip_region_files': str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve()), - }) - - # TODO, for Z3, Z7, ZORDER can potentially also be COARSE_CELLS / aka parent cells? - # clip_subset_type should INPUT_ADDRESS_TYPE for the equivalent of SEQNUM (tp use input_address_type Z3 ...), or COARSE_CELLS as an actual parent cell type clip (also for Z3 ..) - if ( - input_address_type in self.input_address_types - and output_address_type in self.output_address_types - ): - subset_conf.update( - { - 'clip_subset_type': 'INPUT_ADDRESS_TYPE', - 'input_address_type': input_address_type - } - ) - - if not clip_subset_type is None and clip_subset_type in ['COARSE_CELLS']: - subset_conf.update( - { - 'clip_subset_type': clip_subset_type, - 'input_address_type': input_address_type, - 'clip_cell_res' : clip_cell_res, - 'clip_cell_addresses' : " ".join([str(address) for address in cell_id_list]) - } - ) - + clip_metafile_settings, seq_df = specify_clip_settings(clip_subset_type, tmp_dir, tmp_id, input_address_type=input_address_type, + clip_cell_res=clip_cell_res, cell_id_list=cell_id_list, **conf_extra) + subset_conf.update(clip_metafile_settings) subset_conf.update(specify_resolution(**conf_extra)) subset_conf.update(specify_orient_type_args(**conf_extra)) subset_conf.update(specify_topo_aperture(**conf_extra)) @@ -1778,7 +1685,7 @@ def grid_cellids_for_extent( dggs_type: DggsTypeT, resolution: int, mixed_aperture_level: int = None, - clip_geom: "AnyGeometry" = None, + clip_geom: AnyGeometry = None, output_address_type: DggsOutputAddressTypeT | None = None, **conf_extra: DggridMetaConfigParameterT, ) -> pd.DataFrame: @@ -1793,24 +1700,11 @@ def grid_cellids_for_extent( dggs = dgselect(dggs_type = dggs_type, res= resolution, mixed_aperture_level=mixed_aperture_level) dggs.update(**conf_extra, strict=True) - subset_conf = { 'update_frequency': 100000, 'clip_subset_type': 'WHOLE_EARTH' } - - if not clip_geom is None and clip_geom.area > 0: - - clip_gdf = gpd.GeoDataFrame(pd.DataFrame({'id' : [1], 'geometry': [clip_geom]}), geometry='geometry', crs=4326) - clip_path = Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}" - clip_gdf.to_file(str(clip_path), driver=self.tmp_geo_out['driver'] ) - - subset_conf.update({ - 'clip_subset_type': 'GDAL', - 'clip_region_files': str( (Path(tmp_dir) / f"temp_clip_{tmp_id}.{self.tmp_geo_out['ext']}").resolve()), - }) - - if not self.has_gdal: - subset_conf.update({ - 'clip_subset_type': cast(DggsClipSubsetTypeT, self.tmp_geo_out['legacy_driver'].upper()) - }) - + subset_conf: DggridMetaConfigT = {} + clipper_metafile_settings, _ = specify_clip_settings(tmp_dir=tmp_dir, tmp_id=tmp_id, tmp_geo_out=self.tmp_geo_out, + clip_geom=clip_geom, has_gdal=self.has_gdal, resolution=resolution, + **conf_extra) + subset_conf.update(clipper_metafile_settings) subset_conf.update(specify_resolution(**conf_extra)) subset_conf.update(specify_orient_type_args(**conf_extra)) subset_conf.update(specify_topo_aperture(**conf_extra)) @@ -1937,12 +1831,13 @@ def cells_for_geo_points( else: # grid_gen from seqnums dggs_conf = dggs.to_dict() + new_input_address_type = output_address_type if (output_address_type is not None) else "SEQNUM" gdf = self.grid_cell_polygons_from_cellids( cell_id_list=cell_id_list, # dggs_type=dggs_type, passed via dggs_conf resolution=resolution, mixed_aperture_level=mixed_aperture_level, - input_address_type=output_address_type, + input_address_type=new_input_address_type, output_address_type=output_address_type, **dggs_conf, # ensure any extra parameters are passed on ) @@ -2229,6 +2124,72 @@ def specify_resolution( return {} +def specify_clip_settings( + clip_subset_type: DggsClipSubsetTypeT = 'WHOLE_EARTH', + tmp_dir: str = None, + tmp_id: uuid = None, + tmp_geo_out: dict = None, + clip_geom: AnyGeometry = None, + input_address_type: str = "SEQNUM", + has_gdal: bool = True, + clip_cell_res: int = None, + cell_id_list: list = None, + resolution: int = 9, + **conf_extra: DggridMetaConfigParameterT, +) -> DggridMetaConfigT: + + clip_metafile_settings = {"update_frequency": 100000, "clip_subset_type": clip_subset_type} + # Adjust the clipper_scale_factor for refinement level >= 16 + # https://github.com/sahrk/DGGRID/issues/97#issuecomment-3642871117 + rf_level_scalefactor = {17: "10000000", 18: "100000000", + 19: "1000000000", 20: "10000000000"} + # handling grid_cell_polygons/centroids_from_cellids where clip_subset_type != COARSE_CELLS + seq_df = None + if (cell_id_list is not None and len(cell_id_list) > 0): + if (tmp_id is None or tmp_dir is None): + raise ValueError("Either tmp dir or tmp id is not set") + seq_df = pd.DataFrame({input_address_type: cell_id_list}) + seq_df.to_csv(str((Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve()), + header=False, index=False, columns=[input_address_type], sep=' ') + clip_metafile_settings.update({'clip_region_files': str((Path(tmp_dir) / f"temp_clip_{tmp_id}.txt").resolve())}) + if (input_address_type == "SEQNUM"): + clip_metafile_settings.update({'clip_subset_type': 'SEQNUMS'}) + if (input_address_type != "SEQNUM"): + clip_metafile_settings.update({'clip_subset_type': 'INPUT_ADDRESS_TYPE'}) + clip_metafile_settings.update({'input_address_type': input_address_type}) + if (clip_subset_type == "COARSE_CELLS"): + clip_metafile_settings.update({'clip_subset_type': 'COARSE_CELLS'}) + clip_metafile_settings.update({'clip_cell_res': clip_cell_res}) + clip_metafile_settings.update({'clip_cell_addresses': " ".join([str(address) for address in cell_id_list])}) + if "clip_cell_densification" in conf_extra: + clip_metafile_settings.update({'clip_cell_densification': conf_extra['clip_cell_densification']}) + if (resolution >= 17 and resolution <= 20): + clip_metafile_settings.update({'clipper_scale_factor': rf_level_scalefactor[resolution]}) + return clip_metafile_settings, seq_df + + # handling grid_cell_polygons/centroids_for_extent, grid_cellids_for_extent + if clip_geom is not None and clip_geom.area > 0: + if (tmp_id is None or tmp_dir is None or tmp_geo_out is None): + raise ValueError("Either tmp dir or tmp id is not set") + if has_gdal: + clip_metafile_settings.update({'clip_subset_type': 'GDAL'}) + else: + clip_metafile_settings.update({'clip_subset_type': cast(DggsClipSubsetTypeT, tmp_geo_out['legacy_driver'].upper())}) + clip_gdf = gpd.GeoDataFrame(pd.DataFrame({'id': [1], 'geometry': [clip_geom]}), geometry='geometry', crs=4326) + clip_path = Path(tmp_dir) / f"temp_clip_{tmp_id}.{tmp_geo_out['ext']}" + clip_gdf.to_file(str(clip_path), driver=tmp_geo_out['driver']) + clip_metafile_settings.update({'clip_region_files': + str((Path(tmp_dir) / f"temp_clip_{tmp_id}.{tmp_geo_out['ext']}").resolve())}) + if "clipper_scale_factor" in conf_extra: + clip_metafile_settings.update({'clipper_scale_factor': str(conf_extra["clipper_scale_factor"])}) + elif (resolution >= 17 and resolution <= 20): + clip_metafile_settings.update({'clipper_scale_factor': rf_level_scalefactor[resolution]}) + return clip_metafile_settings, seq_df + + return clip_metafile_settings, seq_df + + + def dgconstruct(dggs_type: DggsTypeT = "CUSTOM", # dggs_type projection: DggsProjectionT = 'ISEA', # dggs_projection aperture: DggsApertureT = 3, # dggs_aperture_type / dggs_aperture diff --git a/dggrid4py/igeo7.py b/dggrid4py/igeo7.py index 251a2a9..dc92934 100644 --- a/dggrid4py/igeo7.py +++ b/dggrid4py/igeo7.py @@ -16,14 +16,14 @@ def decode_z7hex_index(z7_hex_str): """ Decode a Z7 hexadecimal index (provided as hex string). - + Format: - First 4 bits: base cell number (0-11) - Remaining 60 bits: 20 groups of 3 bits each for resolution digits (0-6, 7 for beyond resolution) - + Args: z7_hex_str (str): Hexadecimal string representing the Z7 cell index - + Returns: tuple: (base_cell, resolution_digits) - base_cell: integer 0-11 @@ -31,17 +31,17 @@ def decode_z7hex_index(z7_hex_str): """ # Convert hex to binary string, ensuring we have full 64 bits binary = bin(int(z7_hex_str, 16))[2:].zfill(64) - + # Extract base cell (first 4 bits) base_cell = int(binary[:4], 2) - + # Extract resolution digits (20 groups of 3 bits each) resolution_digits = [] for i in range(20): # 60 remaining bits = 20 groups of 3 bits start = 4 + (i * 3) # Start after the first 4 bits value = int(binary[start:start + 3], 2) resolution_digits.append(value) - + return base_cell, resolution_digits @@ -62,7 +62,7 @@ def z7hex_to_z7string(z7_hex_str): return "".join(str_rep) -def z7hex_to_z7int( z7_hex_str): +def z7hex_to_z7int(z7_hex_str): # From hex string to integer value = int(z7_hex_str, 16) return value @@ -74,6 +74,25 @@ def z7int_to_z7hex(z7_int): return hex_back +def z7string_to_z7int(z7string): + """ + Convert z7string to z7int. Inverse of decode_z7hex_index + + Args: + z7string (str): The z7 digits representation of zone ID + + Returns: + int: integer representation of the z7 string. + """ + base, digit = z7string[:2], z7string[2:] + digit = digit.ljust(20, '7') + binary_repr = [np.binary_repr(int(base), width=4)] + binary_digit = [np.binary_repr(int(d), width=3) for d in digit] + binary_repr += binary_digit + binary_repr = ''.join(binary_repr) + return int(binary_repr, 2) + + def get_z7hex_resolution(z7_hex_str): """ Get the resolution of a Z7 cell. @@ -81,7 +100,7 @@ def get_z7hex_resolution(z7_hex_str): Args: z7_hex_str (str): Z7 hexadecimal string represention of the cell index """ - base_cell, resolution_digits = decode_z7hex_index(z7_hex_str) + base_cell, resolution_digits = decode_z7hex_index(z7_hex_str) return len( list( filter(lambda x: x >= 0 and x < 7, resolution_digits) ) ) @@ -168,5 +187,5 @@ def apply_convert_z7hex_to_z7string(z7_hex_str): z7_res = get_z7hex_resolution(z7_hex_str) parent, local_pos, is_center = get_z7hex_local_pos(z7_hex_str) return pd.Series([z7_string, z7_res, parent, local_pos, is_center]) - + diff --git a/tests/test_clipper_scale_factor.py b/tests/test_clipper_scale_factor.py new file mode 100644 index 0000000..fdcdfed --- /dev/null +++ b/tests/test_clipper_scale_factor.py @@ -0,0 +1,151 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +import decimal +import inspect +import os +import tempfile + +import pytest +import shapely +import geopandas as gpd +from geopandas.testing import assert_geodataframe_equal + +from dggrid4py import DGGRIDv8, Dggs + +dggrid_path = os.getenv("DGGRID_PATH") +dggrid = DGGRIDv8(executable=dggrid_path, debug=True) + + +def mock_dggrid_run(__metafile): + return 0 + + +def test_grid_cell_polygons_for_extent(monkeypatch): + metafile = [] + + def mock_dggrid_grid_gen_run(__metafile): + metafile[:] = __metafile + return -1 # cause grid_gen to early-exit in error + + monkeypatch.setattr(dggrid, "run", mock_dggrid_grid_gen_run) + + clip_bound = shapely.geometry.box(27.2, 57.5, 29.3, 59.2) + with pytest.raises(ValueError): # catch and ignore (early-abort "run error") + dggrid.grid_cell_polygons_for_extent( + dggs_type="IGEO7", + resolution=17, + densification=5, + geodetic_densify=0.01, + clip_geom=clip_bound, + # use string to preserve precision and training zeros explicitly + dggs_vert0_azimuth=0.0, + dggs_vert0_lat="58.282525588538994675786", # default: 58.28252559 + dggs_vert0_lon="11.20", # default: 11.25 + ) + + # pre-check temp file paths to ignore in check of specific values + meta_args = dict([line.split(" ") for line in metafile]) + assert meta_args["clip_region_files"].startswith("/tmp/dggrid") + assert meta_args["cell_output_file_name"].startswith("/tmp/dggrid") + meta_args.pop("clip_region_files") + meta_args.pop("cell_output_file_name") + metafile_patched = [f"{key} {val}" for key, val in meta_args.items()] + assert set(metafile_patched) == { + "dggrid_operation GENERATE_GRID", + "dggs_type IGEO7", + "dggs_proj ISEA", + "dggs_aperture 7", + "dggs_topology HEXAGON", + # NOTE: 'dggs_res_spec' to be inferred from Dggs() attribute, not passing it explicitly to 'specify_resolution' + "dggs_res_spec 17", + "precision 7", + "densification 5", + "geodetic_densify 0.01", + "clip_subset_type GDAL", + "clipper_scale_factor 10000000", + # "clip_region_files /tmp/dggrid/...", + # "cell_output_file_name /tmp/dggrid/...", + "cell_output_type GDAL", + "cell_output_gdal_format FlatGeobuf", + # following set explicitly by input parameters + "dggs_orient_specify_type SPECIFIED", + "dggs_vert0_azimuth 0.0", + "dggs_vert0_lat 58.282525588538994675786", + "dggs_vert0_lon 11.20", + # WARNING: following technically not set by Dggs(), though it probably should for 'IGEO7' ? + # "output_cell_label_type OUTPUT_ADDRESS_TYPE", + # "output_address_type HIERNDX", + # "output_hier_ndx_system Z7", + # "output_hier_ndx_form DIGIT_STRING", + "point_output_type NONE" + } + +def test_grid_cell_polygons_from_cellids(monkeypatch): + metafile = [] + + def mock_dggrid_grid_gen_run(__metafile): + metafile[:] = __metafile + return -1 # cause grid_gen to early-exit in error + + monkeypatch.setattr(dggrid, "run", mock_dggrid_grid_gen_run) + + with pytest.raises(ValueError): # catch and ignore (early-abort "run error") + dggrid.grid_cell_polygons_from_cellids( + dggs_type="IGEO7", + resolution=18, + cell_id_list=["023255620345"], + # clip_cell_densification=5, + clip_subset_type="COARSE_CELLS", # required for densification to take effect + clip_cell_res=10, + input_address_type="HIERNDX", + input_hier_ndx_forms="DIGIT_STRING", + input_hier_ndx_systems="Z7", + output_cell_label_type="OUTPUT_ADDRESS_TYPE", + output_address_type="HIERNDX", + output_hier_ndx_forms="DIGIT_STRING", + output_hier_ndx_systems="Z7", + # use string to preserve precision and training zeros explicitly + dggs_vert0_azimuth=0.0, + dggs_vert0_lat="58.282525588538994675786", # default: 58.28252559 + dggs_vert0_lon="11.20", # default: 11.25 + ) + + # pre-check temp file paths to ignore in check of specific values + meta_args = dict([line.split(" ") for line in metafile]) + assert meta_args["clip_region_files"].startswith("/tmp/dggrid") + assert meta_args["cell_output_file_name"].startswith("/tmp/dggrid") + meta_args.pop("clip_region_files") + meta_args.pop("cell_output_file_name") + metafile_patched = [f"{key} {val}" for key, val in meta_args.items()] + + assert set(metafile_patched) == { + "dggrid_operation GENERATE_GRID", + "dggs_type IGEO7", + "dggs_proj ISEA", + "dggs_aperture 7", + "dggs_topology HEXAGON", + # NOTE: 'dggs_res_spec' to be inferred from Dggs() attribute, not passing it explicitly to 'specify_resolution' + "dggs_res_spec 18", + "precision 7", + # "clip_region_files /tmp/dggrid/...", + # "cell_output_file_name /tmp/dggrid/...", + "cell_output_type GDAL", + "cell_output_gdal_format FlatGeobuf", + "clip_cell_addresses 023255620345", + "clip_subset_type COARSE_CELLS", + "clip_cell_res 10", + "clipper_scale_factor 100000000", + # following set explicitly by input parameters + "dggs_orient_specify_type SPECIFIED", + "dggs_vert0_azimuth 0.0", + "dggs_vert0_lat 58.282525588538994675786", + "dggs_vert0_lon 11.20", + "input_address_type HIERNDX", + "output_cell_label_type OUTPUT_ADDRESS_TYPE", + "output_address_type HIERNDX", + # WARNING: following technically not set by Dggs(), though it probably should for 'IGEO7' ? + # "output_hier_ndx_system Z7", + # "output_hier_ndx_form DIGIT_STRING", + "point_output_type NONE" + } +