From c6930cab692850151fe9d95bb1963bb8aaf8f20d Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sat, 14 May 2022 20:42:06 +0100 Subject: [PATCH 1/8] First attempt to implement a tool that patches trees for halo finder instabilities, a little like a simple version of Behroozi's consistent-trees --- tangos/tools/__init__.py | 2 +- tangos/tools/merger_tree_patcher.py | 155 ++++++++++++++++++++++++++++ 2 files changed, 156 insertions(+), 1 deletion(-) create mode 100644 tangos/tools/merger_tree_patcher.py diff --git a/tangos/tools/__init__.py b/tangos/tools/__init__.py index 81c4ffd7..b669b54c 100644 --- a/tangos/tools/__init__.py +++ b/tangos/tools/__init__.py @@ -65,4 +65,4 @@ def add_tools(cls, subparse): from . import (add_simulation, ahf_merger_tree_importer, consistent_trees_importer, crosslink, property_deleter, - property_importer, property_writer) + property_importer, property_writer, merger_tree_patcher) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py new file mode 100644 index 00000000..6c4d36e5 --- /dev/null +++ b/tangos/tools/merger_tree_patcher.py @@ -0,0 +1,155 @@ +import os +import re + +import tangos as db + +from .. import config +from .. import core +from ..core import get_or_create_dictionary_item +from ..core.halo_data import HaloLink, HaloProperty +from ..input_handlers import ahf_trees as at +from ..log import logger +from . import GenericTangosTool +import numpy as np + +class MergerTreePatcher(GenericTangosTool): + tool_name = 'patch-trees' + tool_description = "Attempt to patch up merger trees" + parallel = False + + @classmethod + def add_parser_arguments(self, parser): + parser.add_argument('--sims', '--for', action='store', nargs='*', + metavar='simulation_name', + help='Specify a simulation (or multiple simulations) to run on') + parser.add_argument('--include-only', action='store', type=str, nargs='*', default=['NDM()>100'], + help="Specify a filter that describes which objects the calculation should be executed for. " + "Multiple filters may be specified, in which case they must all evaluate to true for " + "the object to be included.") + parser.add_argument('--relative-match', action='store', type=str, default=['NStar()','NDM()'], + help='Variables to minimise when finding a halo match. If multiple expressions' + 'are provided, their relative differences are added in quadrature') + parser.add_argument('--max-for-match', action='store', type=float, default=0.02, + help='The maximum fractional offset in values of variables specified in relative-match' + 'that qualifies as a match') + parser.add_argument('--position-variable', action='store', type=str, default='shrink_center/a()', + help='The name of the variable which represents the 3d position') + parser.add_argument('--max-position-offset', action='store', type=float, default=10.0, + help='The maximum value of the position offset that is allowed') + parser.add_argument('--max-timesteps', action='store', type=int, default=3, + help='The maximum number of steps back in time to scan for a match') + + + def process_options(self, options): + self.options = options + + + @staticmethod + def filename_to_snapnum(filename): + match = re.match(".*snapdir_([0-9]{3})/?", filename) + if match: + return int(match.group(1)) + else: + raise ValueError("Unable to convert %s to snapshot number"%filename) + + def create_timestep_halo_dictionary(self, ts): + halos = ts.halos.all() + out = {} + for h in halos: + out[h.finder_id] = h + return out + + def create_links(self, ts, ts_next, link_dictionary): + session = db.get_default_session() + d_id = get_or_create_dictionary_item(session, "ahf_tree_link") + objs_this = self.create_timestep_halo_dictionary(ts) + objs_next = self.create_timestep_halo_dictionary(ts_next) + links = [] + for this_id, (next_id, merger_ratio) in link_dictionary: + this_obj = objs_this.get(this_id, None) + next_obj = objs_next.get(next_id, None) + if this_obj is not None and next_obj is not None: + links.append(HaloLink(this_obj, next_obj, d_id, merger_ratio)) + links.append(HaloLink(next_obj, this_obj, d_id, merger_ratio)) + session.add_all(links) + session.commit() + logger.info("%d links created between %s and %s",len(links), ts, ts_next) + + def create_link(self, halo, matched): + session = db.get_default_session() + core.creator.get_creator(session) # just to prevent any commits happening later when not expected + + d_id = get_or_create_dictionary_item(session, "patch-trees") + logger.info(f"Creating a link from {halo} to {matched}") + + next = halo + + ts = halo.timestep.previous + + while ts!=matched.timestep: + phantom = core.halo.PhantomHalo(ts, halo.halo_number, 0) + session.add(phantom) + session.add(core.HaloLink(next, phantom, d_id, 1.0)) + session.add(core.HaloLink(phantom, next, d_id, 1.0)) + + next = phantom + ts = ts.previous + + session.add(core.HaloLink(next, matched, d_id, 1.0)) + session.add(core.HaloLink(matched, next, d_id, 1.0)) + session.commit() + + + def fixup(self, halo): + success = False + ts_candidate = halo.timestep.previous + source_values = [halo.calculate(m) for m in self.options.relative_match] + source_position = halo.calculate(self.options.position_variable) + + for i in range(self.options.max_timesteps): + + dbid, candidate_positions, *candidate_values = ts_candidate.calculate_all("dbid()", + self.options.position_variable, + *self.options.relative_match) + + offsets = [cv-sv for sv, cv in zip(source_values, candidate_values)] + rel_offsets = [(o/sv)**2 for o, sv in zip(offsets, source_values)] + rel_offsets = np.sqrt(sum(rel_offsets)) + + pos_offset = np.linalg.norm(source_position - candidate_positions, axis=1) + mask = pos_offset < self.options.max_position_offset + + dbid = dbid[mask] + rel_offsets = rel_offsets[mask] + + if len(rel_offsets)>0: + if np.min(rel_offsets) < self.options.max_for_match: + match_dbid = dbid[np.argmin(rel_offsets)] + match = db.get_halo(match_dbid) + self.create_link(halo,match) + success = True + break + + ts_candidate = ts_candidate.previous + + if not success: + logger.info(f"No luck in finding a match for {halo}") + + + def run_calculation_loop(self): + simulations = db.sim_query_from_name_list(self.options.sims) + logger.info("This tool is experimental. Use at your own risk!") + for simulation in simulations: + logger.info("Processing %s",simulation) + logger.info(self.options.include_only) + include_criterion = "(!has_link(earlier(1))) & " + (" & ".join(self.options.include_only)) + logger.info("Query for missing links is %s",include_criterion) + for ts in simulation.timesteps[::-1]: + dbids, flag = ts.calculate_all("dbid()", include_criterion) + dbids = dbids[flag] + logger.info(f"Timestep {ts} has {len(dbids)} broken links to consider") + for dbid in dbids: + obj = db.get_halo(dbid) + self.fixup(obj) + + From 0448586c998174668c85799bb62c841669f53763 Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sat, 14 May 2022 20:49:43 +0100 Subject: [PATCH 2/8] Add TODO messages --- tangos/live_calculation/builtin_functions/__init__.py | 3 +-- tangos/tools/merger_tree_patcher.py | 7 +++++++ 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/tangos/live_calculation/builtin_functions/__init__.py b/tangos/live_calculation/builtin_functions/__init__.py index cf1b841a..da31b4e9 100644 --- a/tangos/live_calculation/builtin_functions/__init__.py +++ b/tangos/live_calculation/builtin_functions/__init__.py @@ -11,13 +11,12 @@ @BuiltinFunction.register def match(source_halos, target): - timestep = consistent_collection.ConsistentCollection(source_halos).timestep if target is None: results = [None]*len(source_halos) else: from ... import relation_finding if not isinstance(target, core.Base): - target = tangos.get_item(target, core.Session.object_session(timestep)) + target = tangos.get_item(target, core.Session.object_session(source_halos[0])) results = relation_finding.MultiSourceMultiHopStrategy(source_halos, target).all() # if following assert fails, it might be duplicate links in the database which the # current MultiSourceMultiHop implementation cannot de-duplicate: diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index 6c4d36e5..e8c23ad2 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -88,6 +88,9 @@ def create_link(self, halo, matched): while ts!=matched.timestep: phantom = core.halo.PhantomHalo(ts, halo.halo_number, 0) + # TODO: shouldn't assume this halo number is "free" for use, probably phantom halo IDs should be + # allocated sequentially + session.add(phantom) session.add(core.HaloLink(next, phantom, d_id, 1.0)) session.add(core.HaloLink(phantom, next, d_id, 1.0)) @@ -143,6 +146,10 @@ def run_calculation_loop(self): logger.info("Processing %s",simulation) logger.info(self.options.include_only) include_criterion = "(!has_link(earlier(1))) & " + (" & ".join(self.options.include_only)) + # TODO: this is all very well when there actually isn't a link, but if a "bad" link is stored, + # it's not much use (e.g. if a tiny progenitor is included but not the major progenitor, the + # tool currently won't fix the major progenitor branch) + logger.info("Query for missing links is %s",include_criterion) for ts in simulation.timesteps[::-1]: dbids, flag = ts.calculate_all("dbid()", include_criterion) From bab1ca15997267618beff252f0baa3b0e73b5e21 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Sat, 14 May 2022 21:18:02 +0000 Subject: [PATCH 3/8] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tangos/tools/__init__.py | 4 ++-- tangos/tools/merger_tree_patcher.py | 9 ++++----- 2 files changed, 6 insertions(+), 7 deletions(-) diff --git a/tangos/tools/__init__.py b/tangos/tools/__init__.py index b669b54c..a8d456c4 100644 --- a/tangos/tools/__init__.py +++ b/tangos/tools/__init__.py @@ -64,5 +64,5 @@ def add_tools(cls, subparse): c.add_tools(subparse) from . import (add_simulation, ahf_merger_tree_importer, - consistent_trees_importer, crosslink, property_deleter, - property_importer, property_writer, merger_tree_patcher) + consistent_trees_importer, crosslink, merger_tree_patcher, + property_deleter, property_importer, property_writer) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index e8c23ad2..3822f397 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -1,16 +1,17 @@ import os import re +import numpy as np + import tangos as db -from .. import config -from .. import core +from .. import config, core from ..core import get_or_create_dictionary_item from ..core.halo_data import HaloLink, HaloProperty from ..input_handlers import ahf_trees as at from ..log import logger from . import GenericTangosTool -import numpy as np + class MergerTreePatcher(GenericTangosTool): tool_name = 'patch-trees' @@ -158,5 +159,3 @@ def run_calculation_loop(self): for dbid in dbids: obj = db.get_halo(dbid) self.fixup(obj) - - From 407c12cb8823b2832b795f6354294d90713d0e4f Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sun, 15 May 2022 14:13:34 +0100 Subject: [PATCH 4/8] Add tool to 'prune' trees, i.e. remove questionable links which otherwise may get in the way of patching trees --- tangos/tools/merger_tree_patcher.py | 123 ++++++++++++++++++---------- 1 file changed, 80 insertions(+), 43 deletions(-) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index 3822f397..34a26533 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -1,21 +1,89 @@ -import os -import re +import argparse import numpy as np +from sqlalchemy import and_, or_, orm -import tangos as db - -from .. import config, core +from .. import core, query from ..core import get_or_create_dictionary_item -from ..core.halo_data import HaloLink, HaloProperty -from ..input_handlers import ahf_trees as at from ..log import logger from . import GenericTangosTool +class MergerTreePruner(GenericTangosTool): + tool_name = 'prune-trees' + tool_description = '[Experimental] Remove questionable links from merger trees' + parallel = False + + @classmethod + def add_parser_arguments(self, parser): + parser.add_argument('--sims', '--for', action='store', nargs='*', + metavar='simulation_name', + help='Specify a simulation (or multiple simulations) to run on') + + parser.add_argument('--min-weight-progenitor', action='store', type=float, default=0, + help='The minimum weight that should be associated with a progenitor for retention. ' + 'To retain minor progenitors, this should normally be set to zero.') + + parser.add_argument('--min-weight-descendant', action='store', type=float, default=0.05, + help='The minimum weight that should be associated with a descendant for retention. ' + 'Weights close to zero are normally suspect, since they indicate the halo' + 'transfers only a small fraction of its mass into the descendant.') + # TO IMPLEMENT: + parser.add_argument('--max-NDM-progenitor', action='store', type=float, default=1.0, + help='The maximum relative increase in number of dark matter particles in the progenitor' + 'for retention. For example, if this is set to 1.0 (default), the link is dropped' + 'if there are more than twice as many DM particles in the progenitor as in the' + 'descendant, since this suggests a mis-identification.') + + + def process_options(self, options: argparse.ArgumentParser): + self.options = options + + def process_timestep(self, ts: core.TimeStep): + session = core.get_default_session() + + next_ts = ts.next + if next_ts is None: + return + + + progenitor_from_alias = orm.aliased(core.halo.SimulationObjectBase) + progenitor_to_alias = orm.aliased(core.halo.SimulationObjectBase) + progenitor_halolink = orm.aliased(core.HaloLink) + descendant_halolink = orm.aliased(core.HaloLink) + + progenitor_links = session.query(progenitor_halolink.id, descendant_halolink.id)\ + .join(progenitor_from_alias, and_(progenitor_halolink.halo_from_id == progenitor_from_alias.id, progenitor_from_alias.timestep_id == ts.id))\ + .join(progenitor_to_alias, and_(progenitor_halolink.halo_to_id == progenitor_to_alias.id, progenitor_to_alias.timestep_id == next_ts.id))\ + .join(descendant_halolink, and_(progenitor_halolink.halo_to_id == descendant_halolink.halo_from_id, + progenitor_halolink.halo_from_id == descendant_halolink.halo_to_id))\ + .filter(or_(progenitor_halolink.weight < self.options.min_weight_progenitor, + descendant_halolink.weight < self.options.min_weight_descendant, + progenitor_to_alias.NDM > (1.0 + self.options.max_NDM_progenitor) * progenitor_from_alias.NDM)) + + delete_links = [] + progenitors_and_descendants = progenitor_links.all() + for p,d in progenitors_and_descendants: + delete_links+=[p,d] + + row_count = session.query(core.HaloLink).filter(core.HaloLink.id.in_(delete_links)).delete() + session.commit() + if row_count>0: + logger.info(f"Deleted {row_count} links between timesteps {ts} and {next_ts}") + + def run_calculation_loop(self): + simulations = core.sim_query_from_name_list(self.options.sims) + logger.info("This tool is experimental. Use at your own risk!") + for simulation in simulations: + logger.info("Processing %s", simulation) + + for ts in simulation.timesteps[::-1]: + self.process_timestep(ts) + + class MergerTreePatcher(GenericTangosTool): tool_name = 'patch-trees' - tool_description = "Attempt to patch up merger trees" + tool_description = "[Experimental] Attempt to patch up merger trees" parallel = False @classmethod @@ -45,39 +113,8 @@ def process_options(self, options): self.options = options - @staticmethod - def filename_to_snapnum(filename): - match = re.match(".*snapdir_([0-9]{3})/?", filename) - if match: - return int(match.group(1)) - else: - raise ValueError("Unable to convert %s to snapshot number"%filename) - - def create_timestep_halo_dictionary(self, ts): - halos = ts.halos.all() - out = {} - for h in halos: - out[h.finder_id] = h - return out - - def create_links(self, ts, ts_next, link_dictionary): - session = db.get_default_session() - d_id = get_or_create_dictionary_item(session, "ahf_tree_link") - objs_this = self.create_timestep_halo_dictionary(ts) - objs_next = self.create_timestep_halo_dictionary(ts_next) - links = [] - for this_id, (next_id, merger_ratio) in link_dictionary: - this_obj = objs_this.get(this_id, None) - next_obj = objs_next.get(next_id, None) - if this_obj is not None and next_obj is not None: - links.append(HaloLink(this_obj, next_obj, d_id, merger_ratio)) - links.append(HaloLink(next_obj, this_obj, d_id, merger_ratio)) - session.add_all(links) - session.commit() - logger.info("%d links created between %s and %s",len(links), ts, ts_next) - def create_link(self, halo, matched): - session = db.get_default_session() + session = core.get_default_session() core.creator.get_creator(session) # just to prevent any commits happening later when not expected d_id = get_or_create_dictionary_item(session, "patch-trees") @@ -129,7 +166,7 @@ def fixup(self, halo): if len(rel_offsets)>0: if np.min(rel_offsets) < self.options.max_for_match: match_dbid = dbid[np.argmin(rel_offsets)] - match = db.get_halo(match_dbid) + match = query.get_halo(match_dbid) self.create_link(halo,match) success = True break @@ -141,7 +178,7 @@ def fixup(self, halo): def run_calculation_loop(self): - simulations = db.sim_query_from_name_list(self.options.sims) + simulations = core.sim_query_from_name_list(self.options.sims) logger.info("This tool is experimental. Use at your own risk!") for simulation in simulations: logger.info("Processing %s",simulation) @@ -157,5 +194,5 @@ def run_calculation_loop(self): dbids = dbids[flag] logger.info(f"Timestep {ts} has {len(dbids)} broken links to consider") for dbid in dbids: - obj = db.get_halo(dbid) + obj = query.get_halo(dbid) self.fixup(obj) From c4d672c4b6f873ed4aab43bde5082e9dc21fe857 Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sun, 15 May 2022 15:38:09 +0100 Subject: [PATCH 5/8] Improvements in merger tree tools: * prune-trees: Add minimum (as well as maximum) bound on NDM of a progenitor * prune-trees: Improve clarity of query * patch-trees: Make position offset sphere a comoving radius, to improve identifications at high redshift * patch-trees: Search over more timesteps by default * patch-trees: Cache gathered information about halos in timestep to improve performance --- tangos/tools/merger_tree_patcher.py | 69 +++++++++++++++++++---------- 1 file changed, 45 insertions(+), 24 deletions(-) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index 34a26533..d2b4032e 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -28,48 +28,58 @@ def add_parser_arguments(self, parser): help='The minimum weight that should be associated with a descendant for retention. ' 'Weights close to zero are normally suspect, since they indicate the halo' 'transfers only a small fraction of its mass into the descendant.') - # TO IMPLEMENT: - parser.add_argument('--max-NDM-progenitor', action='store', type=float, default=1.0, + + parser.add_argument('--max-NDM-progenitor', action='store', type=float, default=2.0, help='The maximum relative increase in number of dark matter particles in the progenitor' - 'for retention. For example, if this is set to 1.0 (default), the link is dropped' + 'for retention. For example, if this is set to 2.0 (default), the link is dropped' 'if there are more than twice as many DM particles in the progenitor as in the' 'descendant, since this suggests a mis-identification.') + parser.add_argument('--min-NDM-progenitor', action='store', type=float, default=0.33, + help='The maximum relative number of dark matter particles in the progenitor' + 'for retention. For example, if this is set to 0.33 (default), the link is dropped' + 'if there are less than a third as many DM particles in the progenitor as in the' + 'descendant, since this suggests a mis-identification.') + def process_options(self, options: argparse.ArgumentParser): self.options = options - def process_timestep(self, ts: core.TimeStep): + def process_timestep(self, early_timestep: core.TimeStep): session = core.get_default_session() - next_ts = ts.next - if next_ts is None: + late_timestep = early_timestep.next + if late_timestep is None: return - progenitor_from_alias = orm.aliased(core.halo.SimulationObjectBase) - progenitor_to_alias = orm.aliased(core.halo.SimulationObjectBase) - progenitor_halolink = orm.aliased(core.HaloLink) - descendant_halolink = orm.aliased(core.HaloLink) + descendant_alias = orm.aliased(core.halo.SimulationObjectBase, name="descendant") + progenitor_alias = orm.aliased(core.halo.SimulationObjectBase, name="progenitor") + progenitor_halolink = orm.aliased(core.HaloLink, name="descendant_to_progenitor_link") + descendant_halolink = orm.aliased(core.HaloLink, name="progenitor_to_descendant_link") progenitor_links = session.query(progenitor_halolink.id, descendant_halolink.id)\ - .join(progenitor_from_alias, and_(progenitor_halolink.halo_from_id == progenitor_from_alias.id, progenitor_from_alias.timestep_id == ts.id))\ - .join(progenitor_to_alias, and_(progenitor_halolink.halo_to_id == progenitor_to_alias.id, progenitor_to_alias.timestep_id == next_ts.id))\ - .join(descendant_halolink, and_(progenitor_halolink.halo_to_id == descendant_halolink.halo_from_id, + .join(descendant_alias, and_(progenitor_halolink.halo_from_id == descendant_alias.id, descendant_alias.timestep_id == late_timestep.id))\ + .join(progenitor_alias, and_(progenitor_halolink.halo_to_id == progenitor_alias.id, progenitor_alias.timestep_id == early_timestep.id))\ + .outerjoin(descendant_halolink, and_(progenitor_halolink.halo_to_id == descendant_halolink.halo_from_id, progenitor_halolink.halo_from_id == descendant_halolink.halo_to_id))\ .filter(or_(progenitor_halolink.weight < self.options.min_weight_progenitor, descendant_halolink.weight < self.options.min_weight_descendant, - progenitor_to_alias.NDM > (1.0 + self.options.max_NDM_progenitor) * progenitor_from_alias.NDM)) + progenitor_alias.NDM > self.options.max_NDM_progenitor * descendant_alias.NDM, + progenitor_alias.NDM < self.options.min_NDM_progenitor * descendant_alias.NDM)) + delete_links = [] progenitors_and_descendants = progenitor_links.all() for p,d in progenitors_and_descendants: - delete_links+=[p,d] + delete_links.append(p) + if d is not None: + delete_links.append(d) row_count = session.query(core.HaloLink).filter(core.HaloLink.id.in_(delete_links)).delete() session.commit() if row_count>0: - logger.info(f"Deleted {row_count} links between timesteps {ts} and {next_ts}") + logger.info(f"Deleted {row_count} links between timesteps {early_timestep} and {late_timestep}") def run_calculation_loop(self): simulations = core.sim_query_from_name_list(self.options.sims) @@ -98,14 +108,15 @@ def add_parser_arguments(self, parser): parser.add_argument('--relative-match', action='store', type=str, default=['NStar()','NDM()'], help='Variables to minimise when finding a halo match. If multiple expressions' 'are provided, their relative differences are added in quadrature') - parser.add_argument('--max-for-match', action='store', type=float, default=0.02, + parser.add_argument('--max-for-match', action='store', type=float, default=0.1, help='The maximum fractional offset in values of variables specified in relative-match' 'that qualifies as a match') parser.add_argument('--position-variable', action='store', type=str, default='shrink_center/a()', help='The name of the variable which represents the 3d position') - parser.add_argument('--max-position-offset', action='store', type=float, default=10.0, - help='The maximum value of the position offset that is allowed') - parser.add_argument('--max-timesteps', action='store', type=int, default=3, + parser.add_argument('--max-position-offset', action='store', type=float, default=100.0, + help='The maximum value of the position offset that is allowed in comoving units (i.e.' + 'will be divided by scalefactor)') + parser.add_argument('--max-timesteps', action='store', type=int, default=6, help='The maximum number of steps back in time to scan for a match') @@ -140,25 +151,34 @@ def create_link(self, halo, matched): session.add(core.HaloLink(matched, next, d_id, 1.0)) session.commit() + _candidates_cache = {} + + def _get_candidate_id_position_values(self, ts_candidate): + if ts_candidate.id not in self._candidates_cache: + self._candidates_cache[ts_candidate.id] = ts_candidate.calculate_all("dbid()", + self.options.position_variable, + *self.options.relative_match) + return self._candidates_cache[ts_candidate.id] def fixup(self, halo): success = False ts_candidate = halo.timestep.previous + if ts_candidate is None: + return + source_values = [halo.calculate(m) for m in self.options.relative_match] source_position = halo.calculate(self.options.position_variable) for i in range(self.options.max_timesteps): - dbid, candidate_positions, *candidate_values = ts_candidate.calculate_all("dbid()", - self.options.position_variable, - *self.options.relative_match) + dbid, candidate_positions, *candidate_values = self._get_candidate_id_position_values(ts_candidate) offsets = [cv-sv for sv, cv in zip(source_values, candidate_values)] rel_offsets = [(o/sv)**2 for o, sv in zip(offsets, source_values)] rel_offsets = np.sqrt(sum(rel_offsets)) pos_offset = np.linalg.norm(source_position - candidate_positions, axis=1) - mask = pos_offset < self.options.max_position_offset + mask = pos_offset < self.options.max_position_offset*(1.+ts_candidate.redshift) dbid = dbid[mask] rel_offsets = rel_offsets[mask] @@ -180,6 +200,7 @@ def fixup(self, halo): def run_calculation_loop(self): simulations = core.sim_query_from_name_list(self.options.sims) logger.info("This tool is experimental. Use at your own risk!") + for simulation in simulations: logger.info("Processing %s",simulation) logger.info(self.options.include_only) From 69d133bf2a796580552288f920fe4c2eccd2af46 Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sun, 15 May 2022 15:41:41 +0100 Subject: [PATCH 6/8] Fix bug where patch-trees would crash at early timesteps --- tangos/tools/merger_tree_patcher.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index d2b4032e..4024e843 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -192,6 +192,8 @@ def fixup(self, halo): break ts_candidate = ts_candidate.previous + if ts_candidate is None: + break if not success: logger.info(f"No luck in finding a match for {halo}") From de76160c5faa8d62d06b4915c284bf45cc5b11c0 Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Sun, 15 May 2022 15:45:26 +0100 Subject: [PATCH 7/8] Remove outdated comment --- tangos/tools/merger_tree_patcher.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/tangos/tools/merger_tree_patcher.py b/tangos/tools/merger_tree_patcher.py index 4024e843..8040ee61 100644 --- a/tangos/tools/merger_tree_patcher.py +++ b/tangos/tools/merger_tree_patcher.py @@ -207,9 +207,6 @@ def run_calculation_loop(self): logger.info("Processing %s",simulation) logger.info(self.options.include_only) include_criterion = "(!has_link(earlier(1))) & " + (" & ".join(self.options.include_only)) - # TODO: this is all very well when there actually isn't a link, but if a "bad" link is stored, - # it's not much use (e.g. if a tiny progenitor is included but not the major progenitor, the - # tool currently won't fix the major progenitor branch) logger.info("Query for missing links is %s",include_criterion) for ts in simulation.timesteps[::-1]: From 31e0c06e4d16d7590e63c5f1fc7eae9ae2ff8e99 Mon Sep 17 00:00:00 2001 From: Andrew Pontzen Date: Tue, 26 Jul 2022 09:48:35 +0100 Subject: [PATCH 8/8] Fix failing pre-commit hooks --- tangos/tools/__init__.py | 1 - 1 file changed, 1 deletion(-) diff --git a/tangos/tools/__init__.py b/tangos/tools/__init__.py index 9102b029..35a25043 100644 --- a/tangos/tools/__init__.py +++ b/tangos/tools/__init__.py @@ -67,4 +67,3 @@ def add_tools(cls, subparse): consistent_trees_importer, crosslink, merger_tree_patcher, property_deleter, property_importer, property_writer, subfind_merger_tree_importer) -