From 76af9086074929692e7f452f961b1013d5000042 Mon Sep 17 00:00:00 2001 From: jorgesicachanr Date: Thu, 8 Feb 2024 09:34:24 +0100 Subject: [PATCH 01/11] These are all the latest changes for refactoring co2 mass maps --- src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py | 63 ++++-- src/xtgeoapp_grd3dmaps/aggregate/_config.py | 3 +- .../aggregate/_grid_aggregation.py | 64 +++++- src/xtgeoapp_grd3dmaps/aggregate/_parser.py | 14 ++ .../aggregate/grid3d_aggregate_map.py | 24 +++ .../aggregate/grid3d_co2_mass.py | 201 ++++++++++++++++-- 6 files changed, 335 insertions(+), 34 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py index 4b3a0d4..a6305cb 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py @@ -65,10 +65,28 @@ def translate_co2data_to_property( grid_file, co2_mass_settings.unrst_source, properties_to_extract ) + print(dimensions) # Setting up the grid folder to store the gridproperties if not os.path.exists(grid_out_dir): os.makedirs(grid_out_dir) + #Patch + """ + zones = co2_mass_settings.zones + if zones is None: + zones = [] + elif isinstance(zones, str): + zones = [zones] + + zones = [zone_name.lower() for zone_name in zones] + """ + + #keep_zone_idx = [i for index,zone in enumerate(co2_data.zone) if zone in zones] + + total_mass_list = [] + dissolved_mass_list = [] + free_mass_list = [] + maps = co2_mass_settings.maps if maps is None: maps = [] @@ -82,25 +100,45 @@ def translate_co2data_to_property( store_all = "all" in maps or len(maps) == 0 for co2_at_date in co2_data.data_list: - date = str(co2_at_date.date) + date = str(co2_at_date.date)#Interesting here mass_as_grids = _convert_to_grid(co2_at_date, dimensions, triplets) + print("Printing mass_as_grids...") + print(mass_as_grids["MASS-TOTAL"].date) + if store_all or "total_co2" in maps: - mass_as_grids["mass-total"].to_file( - grid_out_dir + "/MASS_TOTAL_" + date + ".roff", fformat="roff" + mass_as_grids["MASS-TOTAL"].to_file( + grid_out_dir + "/CO2-MASS-TOTAL--" + date + ".roff", fformat="roff", ) - total_mass_list.append(mass_as_grids["mass-total"]) + total_mass_list.append(mass_as_grids["MASS-TOTAL"]) if store_all or "dissolved_co2" in maps: - mass_as_grids["mass-aqu-phase"].to_file( - grid_out_dir + "/MASS_AQU_PHASE_" + date + ".roff", + mass_as_grids["MASS-AQU-PHASE"].to_file( + grid_out_dir + "/CO2-MASS-AQU-PHASE--" + date + ".roff", fformat="roff", ) - dissolved_mass_list.append(mass_as_grids["mass-aqu-phase"]) + dissolved_mass_list.append(mass_as_grids["MASS-AQU-PHASE"]) if store_all or "free_co2" in maps: - mass_as_grids["mass-gas-phase"].to_file( - grid_out_dir + "/MASS_GAS_PHASE_" + date + ".roff", + mass_as_grids["MASS-GAS-PHASE"].to_file( + grid_out_dir + "/CO2-MASS-GAS-PHASE--" + date + ".roff", fformat="roff", ) - free_mass_list.append(mass_as_grids["mass-gas-phase"]) + free_mass_list.append(mass_as_grids["MASS-GAS-PHASE"]) + + + print("Printing property list...") + print(dissolved_mass_list[0].date) + #print(type(free_mass_list[0])) + print("Trying to export gridproperties") + grd_props = xtgeo.GridProperties(props=dissolved_mass_list) + print(grd_props) + print(dir(grd_props)) + print(grd_props.dates) + print(grd_props.props) + print(type(grd_props.props)) + + print("Done trying to export gridproperties") + + #print(type(total_mass_list[0])) + return [ free_mass_list, @@ -162,12 +200,13 @@ def _convert_to_grid( date = str(co2_at_date.date) for mass, name in zip( [co2_at_date.total_mass(), co2_at_date.aqu_phase, co2_at_date.gas_phase], - ["mass-total", "mass-aqu-phase", "mass-gas-phase"], + #["mass-total", "mass-aqu-phase", "mass-gas-phase"], + ["MASS-TOTAL", "MASS-AQU-PHASE", "MASS-GAS-PHASE"], ): mass_array = np.zeros(dimensions) for i, triplet in enumerate(triplets): mass_array[triplet] = mass[i] - mass_name = "co2-" + name + "--" + date + mass_name = "CO2-" + name #+ "--" + date grids[name] = xtgeo.grid3d.GridProperty( ncol=dimensions[0], nrow=dimensions[1], diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_config.py b/src/xtgeoapp_grd3dmaps/aggregate/_config.py index 6c0b7a0..f3368bd 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_config.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_config.py @@ -6,7 +6,7 @@ """ from dataclasses import dataclass, field from enum import Enum -from typing import Dict, List, Optional, Tuple +from typing import Dict, List, Optional, Tuple, Union class AggregationMethod(Enum): @@ -87,6 +87,7 @@ class CO2MassSettings: unrst_source: str init_source: str maps: Optional[List[str]] = None + zones: Optional[List[str]] = None def __post_init__(self): pass diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py b/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py index da057cc..1df40ae 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py @@ -43,19 +43,34 @@ def aggregate_maps( and the second index to `grid_props`. """ # pylint: disable=too-many-arguments + + print("----------------------------------------------") + print("Here we could have some problems ... this is how properties enter in aggregate_maps") + print(grid_props) + + print("And this is how inclusion_filters look like") + print(inclusion_filters) + props, active_cells, inclusion_filters = _read_properties_and_find_active_cells( grid, grid_props, inclusion_filters ) + + weights = grid.get_dz().values1d[active_cells] if weight_by_dz else None # Map nodes (pixel locations) and connections + print("Double checking active_cells") + print(active_cells) conn_data = _find_connections( grid, active_cells, map_template, ) + print("Done with find_connections") + # Iterate filters + print(props) results = _properties_to_maps( - inclusion_filters, + inclusion_filters, props, weights, method, @@ -71,11 +86,27 @@ def _read_properties_and_find_active_cells( ): active = grid.actnum_array.flatten().astype(bool) props = [p.values1d[active] for p in grid_props] + print("-----------------------------") + print("We are in _read_properties_and_find_active_cells and this how props look like") + print(props) + all_masked = np.all([p.mask for p in props], axis=0) + print("The logical condition not any(i is None for i in inclusion_filters) evaluates to: ") + print(not any(i is None for i in inclusion_filters)) + if not any(i is None for i in inclusion_filters): all_masked |= ~np.any(inclusion_filters, axis=0) active[active] = ~all_masked + + print("-----------------------------") + print("This is how all_masked looks like. This will impact props") + print(all_masked) + print("Props looks like this before applying all_masked") + print(props) + props = [p[~all_masked] for p in props] + print("Props looks like this after applying all_masked") + print(props) inclusion_filters = [ None if inc is None else inc[~all_masked] for inc in inclusion_filters ] @@ -108,6 +139,11 @@ def _find_connections( 0, map_template.nrow, dtype=float ) else: + print("Solving issues with mass maps") + print(footprints_x) + print(footprints_y) + print(map_template) + x_nodes, y_nodes = _derive_map_nodes( footprints_x, footprints_y, pixel_to_cell_size_ratio=map_template ) @@ -186,7 +222,10 @@ def _find_overlapped_nodes(nodes, cell_lower, cell_upper): def _cell_footprints(grid: xtgeo.Grid, active_cells): corners = grid.get_xyz_corners() + print("Printing the corners") + print(corners) xyz = [c.values1d[active_cells] for c in corners] + print(xyz) avg_xyz = [(xyz[i] + xyz[i + 12]) / 2 for i in range(12)] x_corners = avg_xyz[::3] y_corners = avg_xyz[1::3] @@ -271,6 +310,10 @@ def _properties_to_maps( map_ix = map_ix[~to_remove] grd_ix = grd_ix[~to_remove] results.append([]) + print("--------------------------------------------") + print("****We are in _properties_to_maps and this is how the properties object looks like****") + print(properties) + for prop in properties: results[-1].append( _property_to_map( @@ -295,11 +338,28 @@ def _property_to_map( ): rows = conn_data.node_indices cols = conn_data.grid_indices + + print("----------------------------") + + print("*** Here we have the problem ***") + + print("This is cols") + print(cols) + print("With size") + print(len(cols)) + + print("... and this is prop") + print(prop) + print("With size") + print(prop.shape) + assert rows.shape == cols.shape assert weights is None or weights.shape == prop.shape if weights is not None: assert method in [AggregationMethod.MEAN, AggregationMethod.SUM] - data = prop[0][cols] if len(prop) == 1 else prop[cols] + #data = prop[0][cols] if len(prop) == 1 else prop[cols] + data = prop[cols] + # Small hack due to a small difference between calculating mass and other properties weights = np.ones_like(data) if weights is None else weights[cols] if data.mask.any(): diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_parser.py b/src/xtgeoapp_grd3dmaps/aggregate/_parser.py index e36a090..0793793 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_parser.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_parser.py @@ -159,25 +159,38 @@ def extract_properties( Extract 3D grid properties based on the provided property specification """ properties: List[Property] = [] + print("Property_spec is: ") + print(property_spec) + if property_spec is None: return properties for spec in property_spec: + print(spec) try: names = "all" if spec.name is None else [spec.name] + print("names= " + str(names)) + print("dates or all means:" + str(dates or "all")) props = xtgeo.gridproperties_from_file( spec.source, names=names, grid=grid, dates=dates or "all", ).props + print(sum(props[0].values1d)) + print(dir(props[0])) except (RuntimeError, ValueError): props = [xtgeo.gridproperty_from_file(spec.source, name=spec.name)] + print("Line 183 and props is:") + print(props) + if spec.lower_threshold is not None: for prop in props: prop.values.mask[prop.values < spec.lower_threshold] = True # Check if any of the properties missing a date had date as part of the file # stem, separated by a "--" for prop in props: + print(prop.date) + print(spec.source) if prop.date is None and "--" in spec.source: date = pathlib.Path(spec.source.split("--")[-1]).stem try: @@ -186,6 +199,7 @@ def extract_properties( except ValueError: continue prop.date = date + print("The assigned prop.date is: "+ str(prop.date)) prop.name += f"--{date}" if len(dates) > 0: props = [p for p in props if p.date in dates] diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py index 7383349..38aa0b1 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py @@ -97,13 +97,29 @@ def generate_maps( """ _XTG.say("Reading grid, properties and zone(s)") grid = xtgeo.grid_from_file(input_.grid) + print(input_.properties) + print(grid) + print(input_.dates) + print("-------------------------") + properties = extract_properties(input_.properties, grid, input_.dates) + print(properties) + _filters = [] if computesettings.all: _filters.append(("all", None)) if computesettings.zone: _filters += extract_zonations(zonation, grid) + print(_filters) _XTG.say("Generating Property Maps") + print("These are the inputs of aggregate maps") + print("grid:") + print(grid) + print("properties") + print(properties) + print("[f[1] for f in _filters]") + print([f[1] for f in _filters]) + xn, yn, p_maps = _grid_aggregation.aggregate_maps( create_map_template(map_settings), grid, @@ -112,10 +128,12 @@ def generate_maps( computesettings.aggregation, computesettings.weight_by_dz, ) + print("a") prop_tags = [ _property_tag(p.name, computesettings.aggregation, output.aggregation_tag) for p in properties ] + print("b") surfs = _ndarray_to_regsurfs( [f[0] for f in _filters], prop_tags, @@ -124,6 +142,8 @@ def generate_maps( p_maps, output.lowercase, ) + print("c") + _write_surfaces(surfs, output.mapfolder, output.plotfolder, output.use_plotly) @@ -184,6 +204,10 @@ def generate_from_config(config: _config.RootConfig): """ Wrapper around `generate_maps` based on a configuration object (RootConfig) """ + print("Generating from config_") + print(config.zonation) + print(config.input) + generate_maps( config.input, config.zonation, diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 9d21713..3773681 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -1,16 +1,19 @@ #!/usr/bin/env python + +import time + import os import sys import tempfile import xtgeo -from typing import List +from typing import List, Optional from xtgeoapp_grd3dmaps.aggregate import ( _co2_mass, _config, _parser, grid3d_aggregate_map, ) -from xtgeoapp_grd3dmaps.aggregate._config import CO2MassSettings +from xtgeoapp_grd3dmaps.aggregate._config import (CO2MassSettings,Zonation) from ccs_scripts.co2_containment.co2_calculation import calculate_co2 PROPERTIES_TO_EXTRACT = [ @@ -41,12 +44,11 @@ # """ -def calculate_mass_property( - grid_file: str, - co2_mass_settings: CO2MassSettings, - out_folder: _config.Output, -) -> List[List[xtgeo.GridProperty]]: +def generate_co2_mass_maps(config_) : + """ + #-> List[List[xtgeo.GridProperty]]: + Calculates and exports 3D CO2 mass properties from the provided grid and config files Args: @@ -60,21 +62,106 @@ def calculate_mass_property( List[List[xtgeo.GridProperty]. """ - co2_data = calculate_co2(grid_file,co2_mass_settings.unrst_source,"mass",co2_mass_settings.init_source,None) + + _input = config_.input + co2_mass_settings = config_.co2_mass_settings + _output = config_.output + zonation = config_.zonation + + + print(zonation) + + zones = co2_mass_settings.zones + if zones is not None and isinstance(zones, str): + co2_mass_settings.zones = [zones] + + print(co2_mass_settings) + + grid_file = _input.grid + + dates = _input.dates + + co2_data = calculate_co2(grid_file,co2_mass_settings.unrst_source,"mass",co2_mass_settings.init_source,None) + + if len(dates)>0: + co2_data.data_list = [x for x in co2_data.data_list if x.date in dates] + + print("Done calculating co2_data") out_property_list = _co2_mass.translate_co2data_to_property( co2_data, grid_file, co2_mass_settings, PROPERTIES_TO_EXTRACT, - out_folder.mapfolder + "/grid", + _output.mapfolder + "/grid", ) + + + print(zonation.zranges) + print(co2_mass_settings) + + ## If all=Yes, get all the ranges from the zonation, good when zranges, missing the other part + if config_.computesettings.all: + if zonation.zranges is not None: + zranges_limits = [list(d.values())[0] for d in zonation.zranges] + elif zonation.Property is not None: + if zproperty.source.split(".")[-1] in ["yml", "yaml"]: + with open(zproperty.source, "r", encoding="utf8") as stream: + try: + zfile = yaml.safe_load(stream) + except yaml.YAMLError as exc: + print(exc) + sys.exit() + if "zranges" not in zfile: + error_text = "The yaml zone file must be in the format:\nzranges:\ + \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" + raise Exception(error_text) + zranges_limits = [list(d.values())[0] for d in zfile.zranges] + max_zvalue = max(sublist[-1] for sublist in zranges_limits) + min_zvalue = min(sublist[0] for sublist in zranges_limits) + all_zrange = [min_zvalue, max_zvalue] + + + config_.zonation.zranges = define_zones_to_plot(zonation,co2_mass_settings) + + print(config_.zonation.zranges)#.append({'all',all_zrange}) + + + if config_.computesettings.all: + tmp_zranges= config_.zonation.zranges + tmp_zranges.append({'all':all_zrange}) + config_.zonation.zranges = tmp_zranges + config_.computesettings.all = False + if not config_.computesettings.zone: + config_.computesettings.zone = True + #Remove whatever is not 'all' + config_.zonation.zranges = [zrange for zrange in config_.zonation.zranges if 'all' in zrange] + + + + print(config_.zonation.zranges)#.append({'all',all_zrange}) + + + print(config_) + + print("Getting into the last stage") + co2_mass_property_to_map(config_,out_property_list) + + + """ + config_.input.properties = [] + config_.computesettings.aggregation = _config.AggregationMethod.SUM + config_.output.aggregation_tag = False + _, temp_path = tempfile.mkstemp() + config_.input.properties.append(_config.Property(temp_path, None, None)) + print("Done here?") return out_property_list + """ def co2_mass_property_to_map( config_: _config.RootConfig, - t_prop: xtgeo.GridProperty, + property_list: List[xtgeo.GridProperty], ): """ Aggregates with SUM and writes a CO2 mass property to file using `grid3d_aggregate_map`. @@ -89,13 +176,63 @@ def co2_mass_property_to_map( config_.input.properties = [] config_.computesettings.aggregation = _config.AggregationMethod.SUM config_.output.aggregation_tag = False - _, temp_path = tempfile.mkstemp() - config_.input.properties.append(_config.Property(temp_path, t_prop.name, None)) - t_prop.to_file(temp_path) + + + for props in property_list: + if len(props)>0 : + for prop in props: + print(prop) + print(prop.date) + #co2_mass_property_to_map(config_, prop) + print("This is very relevant") + config_.input.properties.append(_config.Property(config_.output.mapfolder+"/grid/"+prop.name+"--"+prop.date+".roff", None, None)) + #print(config_) + + print("**** VERY IMPORTANT: printing config_ ****") + print(config_) + print("Trying a hack without affecting the existing code in grid3d_aggregate_map") + #config_.computesettings.all = False + #config_.computesettings.zone = True + print("This is how config_.zonation looks like") + print(config_.zonation) + #config_.zonation = ##Add all to the zonation with all the range in the required format:), remove everything else + grid3d_aggregate_map.generate_from_config(config_) - os.unlink(temp_path) +def define_zones_to_plot( + zonation: _config.Zonation, + co2_mass_settings: _config.CO2MassSettings, +): + """ + Based on the zones defined in CO2MassSettings determine for which zones maps are produced + """ + + if zonation.zranges is not None : + zone_names = [list(item.keys())[0] for item in zonation.zranges] + elif zonation.Property is not None : + if zproperty.source.split(".")[-1] in ["yml", "yaml"]: + with open(zproperty.source, "r", encoding="utf8") as stream: + try: + zfile = yaml.safe_load(stream) + except yaml.YAMLError as exc: + print(exc) + sys.exit() + if "zranges" not in zfile: + error_text = "The yaml zone file must be in the format:\nzranges:\ + \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" + raise Exception(error_text) + zone_names = [list(item.keys())[0] for item in zfile.zranges] + + + if co2_mass_settings.zones is not None: + zones_to_plot = [zone for zone in co2_mass_settings.zones if zone in zone_names] + if len(zones_to_plot)==0: + print("The zones specified in CO2 mass settings are not part of the zonation provided \n maps will be exported for all the existing zones") + else: + return( [item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot]) + + def main(arguments=None): """ Takes input arguments and calculates co2 mass as a property and aggregates it to a 2D map @@ -104,19 +241,45 @@ def main(arguments=None): if arguments is None: arguments = sys.argv[1:] config_ = _parser.process_arguments(arguments) + if config_.input.properties: raise ValueError("CO2 mass computation does not take a property as input") if config_.co2_mass_settings is None: raise ValueError("CO2 mass computation needs co2_mass_settings as input") - out_property_list = calculate_mass_property( - config_.input.grid, + #out_property_list = calculate_mass_property(config_) + + start_time = time.time() + generate_co2_mass_maps(config_) + end_time = time.time() + + execution_time = end_time - start_time + print("Execution time: ", execution_time, "seconds") + + """ + config_.input, config_.co2_mass_settings, - config_.output, + config_.output ) + define_zones_to_plot(config_) + print("Checking out some consistency here") + print("len(out_property_list)="+str(len(out_property_list))) + print("len(out_property_list[0])="+str(len(out_property_list[0]))) + print("len(out_property_list[1])="+str(len(out_property_list[1]))) + print("len(out_property_list[2])="+str(len(out_property_list[2]))) + print(out_property_list[1][0]) + print((out_property_list[1][0].date)) + + ##Until here we seem to be safe. However, there seems to be a more serious issue as when we export the properti(es). The date is not read afterwards + ##Idea: Improve the creation of the properties in calculate_mass_property and see if the feature "date" is recovered back when reading the property + for props in out_property_list: - for prop in props: - co2_mass_property_to_map(config_, prop) + if len(props)>0 : + for prop in props: + print(prop) + print(prop.date) + #co2_mass_property_to_map(config_, prop) + """ if __name__ == "__main__": From f480293c5d0c04f3524e384bbdc462189641c3e2 Mon Sep 17 00:00:00 2001 From: jorgesicachanr Date: Thu, 8 Feb 2024 14:33:36 +0100 Subject: [PATCH 02/11] Refactoring to improve various issues --- src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py | 43 +---- .../aggregate/_grid_aggregation.py | 61 ------- src/xtgeoapp_grd3dmaps/aggregate/_parser.py | 14 -- .../aggregate/grid3d_aggregate_map.py | 24 --- .../aggregate/grid3d_co2_mass.py | 162 +++++------------- 5 files changed, 48 insertions(+), 256 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py index a6305cb..921065d 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_co2_mass.py @@ -65,28 +65,10 @@ def translate_co2data_to_property( grid_file, co2_mass_settings.unrst_source, properties_to_extract ) - print(dimensions) # Setting up the grid folder to store the gridproperties if not os.path.exists(grid_out_dir): os.makedirs(grid_out_dir) - #Patch - """ - zones = co2_mass_settings.zones - if zones is None: - zones = [] - elif isinstance(zones, str): - zones = [zones] - - zones = [zone_name.lower() for zone_name in zones] - """ - - #keep_zone_idx = [i for index,zone in enumerate(co2_data.zone) if zone in zones] - - total_mass_list = [] - dissolved_mass_list = [] - free_mass_list = [] - maps = co2_mass_settings.maps if maps is None: maps = [] @@ -100,11 +82,8 @@ def translate_co2data_to_property( store_all = "all" in maps or len(maps) == 0 for co2_at_date in co2_data.data_list: - date = str(co2_at_date.date)#Interesting here + date = str(co2_at_date.date) mass_as_grids = _convert_to_grid(co2_at_date, dimensions, triplets) - print("Printing mass_as_grids...") - print(mass_as_grids["MASS-TOTAL"].date) - if store_all or "total_co2" in maps: mass_as_grids["MASS-TOTAL"].to_file( grid_out_dir + "/CO2-MASS-TOTAL--" + date + ".roff", fformat="roff", @@ -123,23 +102,6 @@ def translate_co2data_to_property( ) free_mass_list.append(mass_as_grids["MASS-GAS-PHASE"]) - - print("Printing property list...") - print(dissolved_mass_list[0].date) - #print(type(free_mass_list[0])) - print("Trying to export gridproperties") - grd_props = xtgeo.GridProperties(props=dissolved_mass_list) - print(grd_props) - print(dir(grd_props)) - print(grd_props.dates) - print(grd_props.props) - print(type(grd_props.props)) - - print("Done trying to export gridproperties") - - #print(type(total_mass_list[0])) - - return [ free_mass_list, dissolved_mass_list, @@ -200,13 +162,12 @@ def _convert_to_grid( date = str(co2_at_date.date) for mass, name in zip( [co2_at_date.total_mass(), co2_at_date.aqu_phase, co2_at_date.gas_phase], - #["mass-total", "mass-aqu-phase", "mass-gas-phase"], ["MASS-TOTAL", "MASS-AQU-PHASE", "MASS-GAS-PHASE"], ): mass_array = np.zeros(dimensions) for i, triplet in enumerate(triplets): mass_array[triplet] = mass[i] - mass_name = "CO2-" + name #+ "--" + date + mass_name = "CO2-" + name grids[name] = xtgeo.grid3d.GridProperty( ncol=dimensions[0], nrow=dimensions[1], diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py b/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py index 1df40ae..3e1525f 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_grid_aggregation.py @@ -43,32 +43,17 @@ def aggregate_maps( and the second index to `grid_props`. """ # pylint: disable=too-many-arguments - - print("----------------------------------------------") - print("Here we could have some problems ... this is how properties enter in aggregate_maps") - print(grid_props) - - print("And this is how inclusion_filters look like") - print(inclusion_filters) - props, active_cells, inclusion_filters = _read_properties_and_find_active_cells( grid, grid_props, inclusion_filters ) - - weights = grid.get_dz().values1d[active_cells] if weight_by_dz else None # Map nodes (pixel locations) and connections - print("Double checking active_cells") - print(active_cells) conn_data = _find_connections( grid, active_cells, map_template, ) - print("Done with find_connections") - # Iterate filters - print(props) results = _properties_to_maps( inclusion_filters, props, @@ -86,27 +71,11 @@ def _read_properties_and_find_active_cells( ): active = grid.actnum_array.flatten().astype(bool) props = [p.values1d[active] for p in grid_props] - print("-----------------------------") - print("We are in _read_properties_and_find_active_cells and this how props look like") - print(props) - all_masked = np.all([p.mask for p in props], axis=0) - print("The logical condition not any(i is None for i in inclusion_filters) evaluates to: ") - print(not any(i is None for i in inclusion_filters)) - if not any(i is None for i in inclusion_filters): all_masked |= ~np.any(inclusion_filters, axis=0) active[active] = ~all_masked - - print("-----------------------------") - print("This is how all_masked looks like. This will impact props") - print(all_masked) - print("Props looks like this before applying all_masked") - print(props) - props = [p[~all_masked] for p in props] - print("Props looks like this after applying all_masked") - print(props) inclusion_filters = [ None if inc is None else inc[~all_masked] for inc in inclusion_filters ] @@ -139,11 +108,6 @@ def _find_connections( 0, map_template.nrow, dtype=float ) else: - print("Solving issues with mass maps") - print(footprints_x) - print(footprints_y) - print(map_template) - x_nodes, y_nodes = _derive_map_nodes( footprints_x, footprints_y, pixel_to_cell_size_ratio=map_template ) @@ -222,10 +186,7 @@ def _find_overlapped_nodes(nodes, cell_lower, cell_upper): def _cell_footprints(grid: xtgeo.Grid, active_cells): corners = grid.get_xyz_corners() - print("Printing the corners") - print(corners) xyz = [c.values1d[active_cells] for c in corners] - print(xyz) avg_xyz = [(xyz[i] + xyz[i + 12]) / 2 for i in range(12)] x_corners = avg_xyz[::3] y_corners = avg_xyz[1::3] @@ -310,10 +271,6 @@ def _properties_to_maps( map_ix = map_ix[~to_remove] grd_ix = grd_ix[~to_remove] results.append([]) - print("--------------------------------------------") - print("****We are in _properties_to_maps and this is how the properties object looks like****") - print(properties) - for prop in properties: results[-1].append( _property_to_map( @@ -338,29 +295,11 @@ def _property_to_map( ): rows = conn_data.node_indices cols = conn_data.grid_indices - - print("----------------------------") - - print("*** Here we have the problem ***") - - print("This is cols") - print(cols) - print("With size") - print(len(cols)) - - print("... and this is prop") - print(prop) - print("With size") - print(prop.shape) - assert rows.shape == cols.shape assert weights is None or weights.shape == prop.shape if weights is not None: assert method in [AggregationMethod.MEAN, AggregationMethod.SUM] - #data = prop[0][cols] if len(prop) == 1 else prop[cols] data = prop[cols] - - # Small hack due to a small difference between calculating mass and other properties weights = np.ones_like(data) if weights is None else weights[cols] if data.mask.any(): invalid = data.mask diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_parser.py b/src/xtgeoapp_grd3dmaps/aggregate/_parser.py index 0793793..e36a090 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_parser.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_parser.py @@ -159,38 +159,25 @@ def extract_properties( Extract 3D grid properties based on the provided property specification """ properties: List[Property] = [] - print("Property_spec is: ") - print(property_spec) - if property_spec is None: return properties for spec in property_spec: - print(spec) try: names = "all" if spec.name is None else [spec.name] - print("names= " + str(names)) - print("dates or all means:" + str(dates or "all")) props = xtgeo.gridproperties_from_file( spec.source, names=names, grid=grid, dates=dates or "all", ).props - print(sum(props[0].values1d)) - print(dir(props[0])) except (RuntimeError, ValueError): props = [xtgeo.gridproperty_from_file(spec.source, name=spec.name)] - print("Line 183 and props is:") - print(props) - if spec.lower_threshold is not None: for prop in props: prop.values.mask[prop.values < spec.lower_threshold] = True # Check if any of the properties missing a date had date as part of the file # stem, separated by a "--" for prop in props: - print(prop.date) - print(spec.source) if prop.date is None and "--" in spec.source: date = pathlib.Path(spec.source.split("--")[-1]).stem try: @@ -199,7 +186,6 @@ def extract_properties( except ValueError: continue prop.date = date - print("The assigned prop.date is: "+ str(prop.date)) prop.name += f"--{date}" if len(dates) > 0: props = [p for p in props if p.date in dates] diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py index 38aa0b1..7383349 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_aggregate_map.py @@ -97,29 +97,13 @@ def generate_maps( """ _XTG.say("Reading grid, properties and zone(s)") grid = xtgeo.grid_from_file(input_.grid) - print(input_.properties) - print(grid) - print(input_.dates) - print("-------------------------") - properties = extract_properties(input_.properties, grid, input_.dates) - print(properties) - _filters = [] if computesettings.all: _filters.append(("all", None)) if computesettings.zone: _filters += extract_zonations(zonation, grid) - print(_filters) _XTG.say("Generating Property Maps") - print("These are the inputs of aggregate maps") - print("grid:") - print(grid) - print("properties") - print(properties) - print("[f[1] for f in _filters]") - print([f[1] for f in _filters]) - xn, yn, p_maps = _grid_aggregation.aggregate_maps( create_map_template(map_settings), grid, @@ -128,12 +112,10 @@ def generate_maps( computesettings.aggregation, computesettings.weight_by_dz, ) - print("a") prop_tags = [ _property_tag(p.name, computesettings.aggregation, output.aggregation_tag) for p in properties ] - print("b") surfs = _ndarray_to_regsurfs( [f[0] for f in _filters], prop_tags, @@ -142,8 +124,6 @@ def generate_maps( p_maps, output.lowercase, ) - print("c") - _write_surfaces(surfs, output.mapfolder, output.plotfolder, output.use_plotly) @@ -204,10 +184,6 @@ def generate_from_config(config: _config.RootConfig): """ Wrapper around `generate_maps` based on a configuration object (RootConfig) """ - print("Generating from config_") - print(config.zonation) - print(config.input) - generate_maps( config.input, config.zonation, diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 3773681..47dfde9 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -63,102 +63,45 @@ def generate_co2_mass_maps(config_) : """ - _input = config_.input co2_mass_settings = config_.co2_mass_settings - _output = config_.output zonation = config_.zonation - - - print(zonation) - zones = co2_mass_settings.zones if zones is not None and isinstance(zones, str): co2_mass_settings.zones = [zones] - print(co2_mass_settings) - - grid_file = _input.grid - - dates = _input.dates - + grid_file = config_.input.grid co2_data = calculate_co2(grid_file,co2_mass_settings.unrst_source,"mass",co2_mass_settings.init_source,None) + dates = config_.input.dates if len(dates)>0: co2_data.data_list = [x for x in co2_data.data_list if x.date in dates] - print("Done calculating co2_data") - out_property_list = _co2_mass.translate_co2data_to_property( co2_data, grid_file, co2_mass_settings, PROPERTIES_TO_EXTRACT, - _output.mapfolder + "/grid", + config_.output.mapfolder + "/grid", ) + if config_.computesettings.all: + if len(zonation.zranges)>0: + all_zrange = find_all_zrange(zonation=zonation) + else: + all_zrange = find_all_zrange(grid_file=grid_file) - print(zonation.zranges) - print(co2_mass_settings) - - ## If all=Yes, get all the ranges from the zonation, good when zranges, missing the other part - if config_.computesettings.all: - if zonation.zranges is not None: - zranges_limits = [list(d.values())[0] for d in zonation.zranges] - elif zonation.Property is not None: - if zproperty.source.split(".")[-1] in ["yml", "yaml"]: - with open(zproperty.source, "r", encoding="utf8") as stream: - try: - zfile = yaml.safe_load(stream) - except yaml.YAMLError as exc: - print(exc) - sys.exit() - if "zranges" not in zfile: - error_text = "The yaml zone file must be in the format:\nzranges:\ - \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" - raise Exception(error_text) - zranges_limits = [list(d.values())[0] for d in zfile.zranges] - max_zvalue = max(sublist[-1] for sublist in zranges_limits) - min_zvalue = min(sublist[0] for sublist in zranges_limits) - all_zrange = [min_zvalue, max_zvalue] - - - config_.zonation.zranges = define_zones_to_plot(zonation,co2_mass_settings) - - print(config_.zonation.zranges)#.append({'all',all_zrange}) - + if len(zonation.zranges)>0: + config_.zonation.zranges = define_zones_to_plot(zonation,co2_mass_settings) if config_.computesettings.all: - tmp_zranges= config_.zonation.zranges - tmp_zranges.append({'all':all_zrange}) - config_.zonation.zranges = tmp_zranges + config_.zonation.zranges.append({'all':all_zrange}) config_.computesettings.all = False if not config_.computesettings.zone: config_.computesettings.zone = True - #Remove whatever is not 'all' config_.zonation.zranges = [zrange for zrange in config_.zonation.zranges if 'all' in zrange] - - - - print(config_.zonation.zranges)#.append({'all',all_zrange}) - - - print(config_) - print("Getting into the last stage") co2_mass_property_to_map(config_,out_property_list) - - """ - config_.input.properties = [] - config_.computesettings.aggregation = _config.AggregationMethod.SUM - config_.output.aggregation_tag = False - _, temp_path = tempfile.mkstemp() - config_.input.properties.append(_config.Property(temp_path, None, None)) - print("Done here?") - return out_property_list - """ - - def co2_mass_property_to_map( config_: _config.RootConfig, property_list: List[xtgeo.GridProperty], @@ -177,26 +120,11 @@ def co2_mass_property_to_map( config_.computesettings.aggregation = _config.AggregationMethod.SUM config_.output.aggregation_tag = False - for props in property_list: if len(props)>0 : for prop in props: - print(prop) - print(prop.date) - #co2_mass_property_to_map(config_, prop) - print("This is very relevant") config_.input.properties.append(_config.Property(config_.output.mapfolder+"/grid/"+prop.name+"--"+prop.date+".roff", None, None)) - #print(config_) - print("**** VERY IMPORTANT: printing config_ ****") - print(config_) - print("Trying a hack without affecting the existing code in grid3d_aggregate_map") - #config_.computesettings.all = False - #config_.computesettings.zone = True - print("This is how config_.zonation looks like") - print(config_.zonation) - #config_.zonation = ##Add all to the zonation with all the range in the required format:), remove everything else - grid3d_aggregate_map.generate_from_config(config_) @@ -230,8 +158,42 @@ def define_zones_to_plot( if len(zones_to_plot)==0: print("The zones specified in CO2 mass settings are not part of the zonation provided \n maps will be exported for all the existing zones") else: - return( [item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot]) + return([item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot]) + else: + return(zonation.zranges) +def find_all_zrange( + zonation: Optional[_config.Zonation] = None, + grid_file: Optional[str] = None, +): + if zonation is not None: + if zonation.zranges is not None: + zranges_limits = [list(d.values())[0] for d in zonation.zranges] + elif zonation.Property is not None: + if zproperty.source.split(".")[-1] in ["yml", "yaml"]: + with open(zproperty.source, "r", encoding="utf8") as stream: + try: + zfile = yaml.safe_load(stream) + except yaml.YAMLError as exc: + print(exc) + sys.exit() + if "zranges" not in zfile: + error_text = "The yaml zone file must be in the format:\nzranges:\ + \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" + raise Exception(error_text) + zranges_limits = [list(d.values())[0] for d in zonation.zranges] + elif grid_file is not None: + grid_pf = xtgeo.grid_from_file(grid_file) + dimensions = (grid_pf.ncol, grid_pf.nrow, grid_pf.nlay) + zranges_limits = [[1,grid_pf.nlay]] + else: + error_text = "Either zonation or grid_file need to be provided" + raise Exception(error_text) + + max_zvalue = max(sublist[-1] for sublist in zranges_limits) + min_zvalue = min(sublist[0] for sublist in zranges_limits) + all_zrange = [min_zvalue, max_zvalue] + return(all_zrange) def main(arguments=None): """ @@ -246,40 +208,8 @@ def main(arguments=None): raise ValueError("CO2 mass computation does not take a property as input") if config_.co2_mass_settings is None: raise ValueError("CO2 mass computation needs co2_mass_settings as input") - #out_property_list = calculate_mass_property(config_) - start_time = time.time() generate_co2_mass_maps(config_) - end_time = time.time() - - execution_time = end_time - start_time - print("Execution time: ", execution_time, "seconds") - - """ - config_.input, - config_.co2_mass_settings, - config_.output - ) - define_zones_to_plot(config_) - print("Checking out some consistency here") - print("len(out_property_list)="+str(len(out_property_list))) - print("len(out_property_list[0])="+str(len(out_property_list[0]))) - print("len(out_property_list[1])="+str(len(out_property_list[1]))) - print("len(out_property_list[2])="+str(len(out_property_list[2]))) - print(out_property_list[1][0]) - print((out_property_list[1][0].date)) - - ##Until here we seem to be safe. However, there seems to be a more serious issue as when we export the properti(es). The date is not read afterwards - - ##Idea: Improve the creation of the properties in calculate_mass_property and see if the feature "date" is recovered back when reading the property - - for props in out_property_list: - if len(props)>0 : - for prop in props: - print(prop) - print(prop.date) - #co2_mass_property_to_map(config_, prop) - """ if __name__ == "__main__": From 769793730ce966ca7a0a6fe78e6d6bc2883e8e8b Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Thu, 8 Feb 2024 14:45:54 +0100 Subject: [PATCH 03/11] Update _config.py --- src/xtgeoapp_grd3dmaps/aggregate/_config.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/_config.py b/src/xtgeoapp_grd3dmaps/aggregate/_config.py index f3368bd..9de668f 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/_config.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/_config.py @@ -6,7 +6,7 @@ """ from dataclasses import dataclass, field from enum import Enum -from typing import Dict, List, Optional, Tuple, Union +from typing import Dict, List, Optional, Tuple class AggregationMethod(Enum): From d6d817c3b31ac37befecde3881c5d7fd91a3c9d7 Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Thu, 8 Feb 2024 14:46:21 +0100 Subject: [PATCH 04/11] Update grid3d_co2_mass.py --- src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 47dfde9..e23492f 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -1,7 +1,4 @@ #!/usr/bin/env python - -import time - import os import sys import tempfile From ef47ee89b21e24edde32cb8abefe87240eff926f Mon Sep 17 00:00:00 2001 From: Jorge Sicacha Date: Thu, 8 Feb 2024 14:53:05 +0100 Subject: [PATCH 05/11] Add sample yml with specific zones/dates --- tests/yaml/config_co2_mass_zones_dates.yml | 24 ++++++++++++++++++++++ 1 file changed, 24 insertions(+) create mode 100644 tests/yaml/config_co2_mass_zones_dates.yml diff --git a/tests/yaml/config_co2_mass_zones_dates.yml b/tests/yaml/config_co2_mass_zones_dates.yml new file mode 100644 index 0000000..d27182f --- /dev/null +++ b/tests/yaml/config_co2_mass_zones_dates.yml @@ -0,0 +1,24 @@ +input: + eclroot: tests/data/reek/REEK + grid: $eclroot.EGRID + dates: + - 21901113 + +co2_mass_settings: + unrst_source: $eclroot.UNRST + init_source: $eclroot.INIT + maps: "dissolved_co2" + zones: ["UPPER","ZERO"] + +zonation: + zranges: + - UPPER: [1, 6] + - LOWER: [8, 14] + - ZERO: [15, 15] + +output: + mapfolder: tmp + +computesettings: + zone: Yes + all: No From ac606aaf69f04c859e8702876e99cd75003386a6 Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Fri, 9 Feb 2024 13:09:53 +0100 Subject: [PATCH 06/11] Minor changes in grid3d_co2_mass.py --- .../aggregate/grid3d_co2_mass.py | 27 ++++++++++--------- 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index e23492f..ae54546 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -3,6 +3,7 @@ import sys import tempfile import xtgeo +import yaml from typing import List, Optional from xtgeoapp_grd3dmaps.aggregate import ( _co2_mass, @@ -87,8 +88,9 @@ def generate_co2_mass_maps(config_) : else: all_zrange = find_all_zrange(grid_file=grid_file) - if len(zonation.zranges)>0: + if len(zonation.zranges)>0 or zonation.zproperty is not None: config_.zonation.zranges = define_zones_to_plot(zonation,co2_mass_settings) + config_.zonation.zproperty = None if config_.computesettings.all: config_.zonation.zranges.append({'all':all_zrange}) @@ -133,11 +135,11 @@ def define_zones_to_plot( Based on the zones defined in CO2MassSettings determine for which zones maps are produced """ - if zonation.zranges is not None : + if len(zonation.zranges) > 0 : zone_names = [list(item.keys())[0] for item in zonation.zranges] - elif zonation.Property is not None : - if zproperty.source.split(".")[-1] in ["yml", "yaml"]: - with open(zproperty.source, "r", encoding="utf8") as stream: + elif zonation.zproperty is not None : + if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: + with open(zonation.zproperty.source, "r", encoding="utf8") as stream: try: zfile = yaml.safe_load(stream) except yaml.YAMLError as exc: @@ -146,8 +148,9 @@ def define_zones_to_plot( if "zranges" not in zfile: error_text = "The yaml zone file must be in the format:\nzranges:\ \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" - raise Exception(error_text) - zone_names = [list(item.keys())[0] for item in zfile.zranges] + raise Exception(error_text) + zonation.zranges = zfile['ranges'] + zone_names = [list(item.keys())[0] for item in zonation.zranges] if co2_mass_settings.zones is not None: @@ -164,11 +167,11 @@ def find_all_zrange( grid_file: Optional[str] = None, ): if zonation is not None: - if zonation.zranges is not None: + if len(zonation.zranges) > 0: zranges_limits = [list(d.values())[0] for d in zonation.zranges] - elif zonation.Property is not None: - if zproperty.source.split(".")[-1] in ["yml", "yaml"]: - with open(zproperty.source, "r", encoding="utf8") as stream: + elif zonation.zproperty is not None: + if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: + with open(zonation.zproperty.source, "r", encoding="utf8") as stream: try: zfile = yaml.safe_load(stream) except yaml.YAMLError as exc: @@ -178,7 +181,7 @@ def find_all_zrange( error_text = "The yaml zone file must be in the format:\nzranges:\ \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" raise Exception(error_text) - zranges_limits = [list(d.values())[0] for d in zonation.zranges] + zranges_limits = [list(d.values())[0] for d in zfile.zranges] elif grid_file is not None: grid_pf = xtgeo.grid_from_file(grid_file) dimensions = (grid_pf.ncol, grid_pf.nrow, grid_pf.nlay) From 169d7ccdcb735d2a2abb91f29a88537f3392fb45 Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Fri, 9 Feb 2024 13:20:15 +0100 Subject: [PATCH 07/11] Solving typo --- src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index ae54546..7e78691 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -149,7 +149,7 @@ def define_zones_to_plot( error_text = "The yaml zone file must be in the format:\nzranges:\ \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" raise Exception(error_text) - zonation.zranges = zfile['ranges'] + zonation.zranges = zfile['zranges'] zone_names = [list(item.keys())[0] for item in zonation.zranges] From 7d5d9b039e323e62c639511437cb40c7d48117e8 Mon Sep 17 00:00:00 2001 From: Jorge Sicacha Date: Mon, 12 Feb 2024 15:31:34 +0100 Subject: [PATCH 08/11] Remove duplicate code in zonation processing --- .../aggregate/grid3d_co2_mass.py | 69 ++++++++++++++++--- 1 file changed, 60 insertions(+), 9 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 7e78691..6cd3109 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -82,16 +82,21 @@ def generate_co2_mass_maps(config_) : config_.output.mapfolder + "/grid", ) - if config_.computesettings.all: + config_.zonation.zranges, all_zrange = process_zonation(co2_mass_settings,zonation,grid_file) + + """ + if config_.computesettings.all : if len(zonation.zranges)>0: all_zrange = find_all_zrange(zonation=zonation) else: all_zrange = find_all_zrange(grid_file=grid_file) + """ - if len(zonation.zranges)>0 or zonation.zproperty is not None: - config_.zonation.zranges = define_zones_to_plot(zonation,co2_mass_settings) + if len(config_.zonation.zranges)>0: config_.zonation.zproperty = None + + if config_.computesettings.all: config_.zonation.zranges.append({'all':all_zrange}) config_.computesettings.all = False @@ -126,14 +131,11 @@ def co2_mass_property_to_map( grid3d_aggregate_map.generate_from_config(config_) - +""" def define_zones_to_plot( zonation: _config.Zonation, co2_mass_settings: _config.CO2MassSettings, ): - """ - Based on the zones defined in CO2MassSettings determine for which zones maps are produced - """ if len(zonation.zranges) > 0 : zone_names = [list(item.keys())[0] for item in zonation.zranges] @@ -161,7 +163,8 @@ def define_zones_to_plot( return([item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot]) else: return(zonation.zranges) - +""" +""" def find_all_zrange( zonation: Optional[_config.Zonation] = None, grid_file: Optional[str] = None, @@ -194,7 +197,55 @@ def find_all_zrange( min_zvalue = min(sublist[0] for sublist in zranges_limits) all_zrange = [min_zvalue, max_zvalue] return(all_zrange) - +""" +def process_zonation(co2_mass_settings: _config.CO2MassSettings, + zonation: Optional[_config.Zonation]=None, + grid_file: Optional[str] = None): + if zonation is not None: + if zonation.zproperty is not None: + if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: + zfile = read_yml_file(zonation.zproperty.source) + zonation.zranges = zfile['ranges'] + if len(zonation.zranges) > 0: + zone_names = [list(item.keys())[0] for item in zonation.zranges] + zranges_limits = [list(d.values())[0] for d in zonation.zranges] + elif grid_file is not None: + grid_pf = xtgeo.grid_from_file(grid_file) + zranges_limits = [[1,grid_pf.nlay]] + else: + error_text = "Either zonation or grid_file need to be provided" + raise Exception(error_text) + max_zvalue = max(sublist[-1] for sublist in zranges_limits) + min_zvalue = min(sublist[0] for sublist in zranges_limits) + all_zrange = [min_zvalue, max_zvalue] + + if zone_names is not None: + if co2_mass_settings.zones is not None: + zones_to_plot = [zone for zone in co2_mass_settings.zones if zone in zone_names] + if len(zones_to_plot) == 0: + print( + "The zones specified in CO2 mass settings are not part of the zonation provided \n maps will be exported for all the existing zones") + return zonation.zranges,all_zrange + else: + return [item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot],all_zrange + else: + return zonation.zranges,all_zrange + else: + return [], all_zrange + +def read_yml_file(file_path: str): + with open(file_path, "r", encoding="utf8") as stream: + try: + zfile = yaml.safe_load(stream) + except yaml.YAMLError as exc: + print(exc) + sys.exit() + if "zranges" not in zfile: + error_text = "The yaml zone file must be in the format:\nzranges:\ + \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" + raise Exception(error_text) + return zfile + def main(arguments=None): """ Takes input arguments and calculates co2 mass as a property and aggregates it to a 2D map From 5b7bd877a375636219fc8dd293c18fdd956698a3 Mon Sep 17 00:00:00 2001 From: Jorge Sicacha Date: Mon, 12 Feb 2024 16:50:19 +0100 Subject: [PATCH 09/11] Few more changes to grid3d_co2_mass.py --- .../aggregate/grid3d_co2_mass.py | 148 ++++-------------- 1 file changed, 33 insertions(+), 115 deletions(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 6cd3109..b0711d1 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -4,7 +4,7 @@ import tempfile import xtgeo import yaml -from typing import List, Optional +from typing import List, Optional, Dict, Tuple from xtgeoapp_grd3dmaps.aggregate import ( _co2_mass, _config, @@ -45,20 +45,10 @@ def generate_co2_mass_maps(config_) : """ - #-> List[List[xtgeo.GridProperty]]: - - Calculates and exports 3D CO2 mass properties from the provided grid and config files + Calculates and exports 2D and 3D CO2 mass properties from the provided config file Args: - grid_file (str): Path to EGRID-file - co2_mass_settings (CO2MassSettings): Settings from config file for calculation - of CO2 mass maps. - out_folder (str): Path to store the produced 3D GridProperties. - - - Returns: - List[List[xtgeo.GridProperty]. - + config_: Arguments in the config file """ co2_mass_settings = config_.co2_mass_settings @@ -66,14 +56,11 @@ def generate_co2_mass_maps(config_) : zones = co2_mass_settings.zones if zones is not None and isinstance(zones, str): co2_mass_settings.zones = [zones] - grid_file = config_.input.grid - co2_data = calculate_co2(grid_file,co2_mass_settings.unrst_source,"mass",co2_mass_settings.init_source,None) - + co2_data = calculate_co2(grid_file,co2_mass_settings.unrst_source,"mass",co2_mass_settings.init_source,None) dates = config_.input.dates if len(dates)>0: co2_data.data_list = [x for x in co2_data.data_list if x.date in dates] - out_property_list = _co2_mass.translate_co2data_to_property( co2_data, grid_file, @@ -81,29 +68,15 @@ def generate_co2_mass_maps(config_) : PROPERTIES_TO_EXTRACT, config_.output.mapfolder + "/grid", ) - config_.zonation.zranges, all_zrange = process_zonation(co2_mass_settings,zonation,grid_file) - - """ - if config_.computesettings.all : - if len(zonation.zranges)>0: - all_zrange = find_all_zrange(zonation=zonation) - else: - all_zrange = find_all_zrange(grid_file=grid_file) - """ - if len(config_.zonation.zranges)>0: config_.zonation.zproperty = None - - - if config_.computesettings.all: config_.zonation.zranges.append({'all':all_zrange}) config_.computesettings.all = False if not config_.computesettings.zone: config_.computesettings.zone = True config_.zonation.zranges = [zrange for zrange in config_.zonation.zranges if 'all' in zrange] - co2_mass_property_to_map(config_,out_property_list) def co2_mass_property_to_map( @@ -111,114 +84,56 @@ def co2_mass_property_to_map( property_list: List[xtgeo.GridProperty], ): """ - Aggregates with SUM and writes a CO2 mass property to file using `grid3d_aggregate_map`. - The property is written to a temporary file while performing the - aggregation. + Aggregates with SUM and writes a list of CO2 mass property to files + using `grid3d_aggregate_map`. Args: config_: Arguments in the config file - t_prop: Grid property to be aggregated + property_list: List of Grid property objects to be aggregated """ config_.input.properties = [] config_.computesettings.aggregation = _config.AggregationMethod.SUM config_.output.aggregation_tag = False - for props in property_list: if len(props)>0 : for prop in props: - config_.input.properties.append(_config.Property(config_.output.mapfolder+"/grid/"+prop.name+"--"+prop.date+".roff", None, None)) - + config_.input.properties.append(_config.Property(config_.output.mapfolder+ + "/grid/"+prop.name+"--"+ + prop.date+".roff", None, None)) grid3d_aggregate_map.generate_from_config(config_) -""" -def define_zones_to_plot( - zonation: _config.Zonation, - co2_mass_settings: _config.CO2MassSettings, -): - - if len(zonation.zranges) > 0 : - zone_names = [list(item.keys())[0] for item in zonation.zranges] - elif zonation.zproperty is not None : - if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: - with open(zonation.zproperty.source, "r", encoding="utf8") as stream: - try: - zfile = yaml.safe_load(stream) - except yaml.YAMLError as exc: - print(exc) - sys.exit() - if "zranges" not in zfile: - error_text = "The yaml zone file must be in the format:\nzranges:\ - \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" - raise Exception(error_text) - zonation.zranges = zfile['zranges'] - zone_names = [list(item.keys())[0] for item in zonation.zranges] - +def process_zonation(co2_mass_settings: _config.CO2MassSettings, + grid_file: str, + zonation: Optional[_config.Zonation]=None + ) -> Tuple[List,List]: + """ + Processes a zonation file, if existing, and extracts both zranges per zone + and the complete range in the zaxis. Otherwise, uses the grid_file. - if co2_mass_settings.zones is not None: - zones_to_plot = [zone for zone in co2_mass_settings.zones if zone in zone_names] - if len(zones_to_plot)==0: - print("The zones specified in CO2 mass settings are not part of the zonation provided \n maps will be exported for all the existing zones") - else: - return([item for item in zonation.zranges if list(item.keys())[0] in zones_to_plot]) - else: - return(zonation.zranges) -""" -""" -def find_all_zrange( - zonation: Optional[_config.Zonation] = None, - grid_file: Optional[str] = None, -): - if zonation is not None: - if len(zonation.zranges) > 0: - zranges_limits = [list(d.values())[0] for d in zonation.zranges] - elif zonation.zproperty is not None: - if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: - with open(zonation.zproperty.source, "r", encoding="utf8") as stream: - try: - zfile = yaml.safe_load(stream) - except yaml.YAMLError as exc: - print(exc) - sys.exit() - if "zranges" not in zfile: - error_text = "The yaml zone file must be in the format:\nzranges:\ - \n - Zone1: [1, 5]\n - Zone2: [6, 10]\n - Zone3: [11, 14])" - raise Exception(error_text) - zranges_limits = [list(d.values())[0] for d in zfile.zranges] - elif grid_file is not None: - grid_pf = xtgeo.grid_from_file(grid_file) - dimensions = (grid_pf.ncol, grid_pf.nrow, grid_pf.nlay) - zranges_limits = [[1,grid_pf.nlay]] - else: - error_text = "Either zonation or grid_file need to be provided" - raise Exception(error_text) + Args: + co2_mass_settings: Arguments in CO2 mass settings + grid_file: Path to grid file + zonation: Arguments in zonation - max_zvalue = max(sublist[-1] for sublist in zranges_limits) - min_zvalue = min(sublist[0] for sublist in zranges_limits) - all_zrange = [min_zvalue, max_zvalue] - return(all_zrange) -""" -def process_zonation(co2_mass_settings: _config.CO2MassSettings, - zonation: Optional[_config.Zonation]=None, - grid_file: Optional[str] = None): - if zonation is not None: + Returns: + Tuple[List,List] + """ + if zonation.zproperty is not None or len(zonation.zranges)>0: if zonation.zproperty is not None: if zonation.zproperty.source.split(".")[-1] in ["yml", "yaml"]: zfile = read_yml_file(zonation.zproperty.source) - zonation.zranges = zfile['ranges'] + zonation.zranges = zfile['zranges'] if len(zonation.zranges) > 0: zone_names = [list(item.keys())[0] for item in zonation.zranges] zranges_limits = [list(d.values())[0] for d in zonation.zranges] - elif grid_file is not None: + else: grid_pf = xtgeo.grid_from_file(grid_file) zranges_limits = [[1,grid_pf.nlay]] - else: - error_text = "Either zonation or grid_file need to be provided" - raise Exception(error_text) + zone_names = None max_zvalue = max(sublist[-1] for sublist in zranges_limits) min_zvalue = min(sublist[0] for sublist in zranges_limits) all_zrange = [min_zvalue, max_zvalue] - if zone_names is not None: if co2_mass_settings.zones is not None: zones_to_plot = [zone for zone in co2_mass_settings.zones if zone in zone_names] @@ -233,7 +148,10 @@ def process_zonation(co2_mass_settings: _config.CO2MassSettings, else: return [], all_zrange -def read_yml_file(file_path: str): +def read_yml_file(file_path: str) -> Dict[List]: + """ + Reads a yml from a given path in file_path argument + """ with open(file_path, "r", encoding="utf8") as stream: try: zfile = yaml.safe_load(stream) @@ -249,7 +167,7 @@ def read_yml_file(file_path: str): def main(arguments=None): """ Takes input arguments and calculates co2 mass as a property and aggregates it to a 2D map - at each time step, divided into different phases and locations(TODO). + at each time step, divided into different phases and locations. """ if arguments is None: arguments = sys.argv[1:] From 6e0e8cf0c0a142574f3c5d9961ff9bcc0edb0133 Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Tue, 13 Feb 2024 13:05:11 +0100 Subject: [PATCH 10/11] Solve issue in read_yml_file --- src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index b0711d1..3eab02b 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -148,7 +148,7 @@ def process_zonation(co2_mass_settings: _config.CO2MassSettings, else: return [], all_zrange -def read_yml_file(file_path: str) -> Dict[List]: +def read_yml_file(file_path: str) -> Dict[str,List]: """ Reads a yml from a given path in file_path argument """ From d5b9beaf3d6cdacd8cae756dacbe2ed35f403e28 Mon Sep 17 00:00:00 2001 From: jorgesicachanr <114475076+jorgesicachanr@users.noreply.github.com> Date: Tue, 13 Feb 2024 13:15:49 +0100 Subject: [PATCH 11/11] Update grid3d_co2_mass.py --- src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py index 3eab02b..4b9d5b7 100644 --- a/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py +++ b/src/xtgeoapp_grd3dmaps/aggregate/grid3d_co2_mass.py @@ -68,7 +68,7 @@ def generate_co2_mass_maps(config_) : PROPERTIES_TO_EXTRACT, config_.output.mapfolder + "/grid", ) - config_.zonation.zranges, all_zrange = process_zonation(co2_mass_settings,zonation,grid_file) + config_.zonation.zranges, all_zrange = process_zonation(co2_mass_settings,grid_file,zonation) if len(config_.zonation.zranges)>0: config_.zonation.zproperty = None if config_.computesettings.all: