From f31ec10c57a789147a1fee4915a5daf016d0396d Mon Sep 17 00:00:00 2001 From: Ben Haller Date: Sat, 22 Aug 2026 20:08:39 +0100 Subject: [PATCH 1/4] put binary derived state info in the mutation table metadata --- VERSIONS | 3 + core/slim_globals.cpp | 34 +++++++- core/species.cpp | 189 ++++++++++++++++++++---------------------- core/species.h | 4 +- 4 files changed, 127 insertions(+), 103 deletions(-) diff --git a/VERSIONS b/VERSIONS index 615f6c22..4e333750 100644 --- a/VERSIONS +++ b/VERSIONS @@ -208,6 +208,9 @@ multitrait branch: internally, split baselineOffset into a "from substitutions" component and a "from user" component, and expose those components in .trees metadata, to fix #661 the top-level JSON schema now has two keys, baselineOffsetFromUser and baselineOffsetFromSubstitutions, that provide these two components SLiM writes both of them, for use on the Python side; when reading a .trees, SLiM uses baselineOffsetFromUser and ignores baselineOffsetFromSubstitutions (calculating that itself) + shift to putting binary derived state info in the mutation table metadata, in addition to the ASCII derived state info in its own column, when on disk (i.e., in .trees files); see #664 + in-memory, the binary derived state info remains in the derived state column and the mutation table metadata is kept empty, as before + the schema for the mutation table metadata now encodes an array (of length determined by the metadata itself) of int64_t SLiM mutation IDs version 5.2 (Eidos version 4.2): diff --git a/core/slim_globals.cpp b/core/slim_globals.cpp index ed99ef68..a072ed81 100644 --- a/core/slim_globals.cpp +++ b/core/slim_globals.cpp @@ -2177,7 +2177,39 @@ const std::string gSLiM_tsk_metadata_binary_schema_FORMAT_SOURCE = R"V0G0N({ const std::string gSLiM_tsk_edge_metadata_schema_SOURCE = ""; const std::string gSLiM_tsk_site_metadata_schema_SOURCE = ""; -const std::string gSLiM_tsk_mutation_metadata_schema_SOURCE = ""; // this is now managed by gSLiM_tsk_metadata_binary_schema in top-level metadata + +// actual mutation metadata is now kept in the top-level metadata; see gSLiM_tsk_metadata_binary_schema +// the mutation metadata column now (on disk only) contains binary derived state info, in addition to the +// ASCII info kept in the derived state column; see DerivedStatesFromMetadata() and DerivedStatesToMetadata() +#pragma mark gSLiM_tsk_mutation_metadata_schema_SOURCE + +const std::string gSLiM_tsk_mutation_metadata_schema_SOURCE = R"V0G0N({ + "$schema": "http://json-schema.org/schema#", + "additionalProperties": false, + "codec": "struct", + "type": "object", + "description": "SLiM schema for representing binary derived state data in mutation metadata (the actual mutation metadata is stored in a table in top-level metadata now).", + "examples": [ + { + "derived_states": [0, 1, 17] + } + ], + "properties": { + "derived_states": { + "index": 1, + "type": "array", + "noLengthEncodingExhaustBuffer": true, + "description": "An array of SLiM mutation IDs (int64t), representing the (stacked) mutations contained by the derived state for the mutation. This information is also saved as a comma-separated ASCII string in the derived state column.", + "items": { + "binaryFormat": "q", + "type": "number" + } + } + }, + "required": [ + "derived_states" + ] +})V0G0N"; // BCH 12/10/2024: Removed the type field, and changed the treatment of is_vacant. We have a tricky problem // here, which is that is_vacant is now variable-length and there is no count. The number of byte (uint8_t) diff --git a/core/species.cpp b/core/species.cpp index 81923bef..090e111f 100644 --- a/core/species.cpp +++ b/core/species.cpp @@ -8771,115 +8771,101 @@ void Species::CheckAutoSimplification(void) } } -void Species::DerivedStatesFromAscii(tsk_table_collection_t *p_tables) +void Species::DerivedStatesFromMetadata(tsk_table_collection_t *p_tables) { - // This modifies p_tables in place, replacing the derived_state column of p_tables with a binary version. - tsk_mutation_table_t mutations_copy; - int ret = tsk_mutation_table_copy(&p_tables->mutations, &mutations_copy, 0); - if (ret < 0) handle_error("derived_from_ascii", ret); + // This is called when reading a .trees file. On disk, derived state information is kept in the mutation + // table's metadata column in binary, AND in the derived state column in ASCII (see #664 for discussion). + // Here we convert back to our in-memory format by swapping the derived state column's ASCII data with the + // metadata column's binary data, and then purging the ASCII data from the metadata column. - { - const char *derived_state = p_tables->mutations.derived_state; - tsk_size_t *derived_state_offset = p_tables->mutations.derived_state_offset; - std::vector binary_derived_state; - std::vector binary_derived_state_offset; - size_t derived_state_total_part_count = 0; - - binary_derived_state_offset.emplace_back(0); - - try { - for (size_t j = 0; j < p_tables->mutations.num_rows; j++) - { - std::string string_derived_state(derived_state + derived_state_offset[j], derived_state_offset[j+1] - derived_state_offset[j]); - - if (string_derived_state.size() == 0) - { - // nothing to do for an empty derived state - } - else if (string_derived_state.find(',') == std::string::npos) - { - // a single mutation can be handled more efficiently, and this is the common case so it's worth optimizing - binary_derived_state.emplace_back((slim_mutationid_t)std::stoll(string_derived_state)); - derived_state_total_part_count++; - } - else - { - // stacked mutations require that the derived state be separated to parse it - std::vector derived_state_parts = Eidos_string_split(string_derived_state, ","); - - for (std::string &derived_state_part : derived_state_parts) - binary_derived_state.emplace_back((slim_mutationid_t)std::stoll(derived_state_part)); - - derived_state_total_part_count += derived_state_parts.size(); - } - - binary_derived_state_offset.emplace_back((tsk_size_t)(derived_state_total_part_count * sizeof(slim_mutationid_t))); - } - } catch (...) { - EIDOS_TERMINATION << "ERROR (Species::DerivedStatesFromAscii): a mutation derived state was not convertible into an int64_t mutation id. The tree-sequence data may not be annotated for SLiM, or may be corrupted. If mutations were added in msprime, do you want to use the msprime.SLiMMutationModel?" << EidosTerminate(); - } - - if (binary_derived_state.size() == 0) - binary_derived_state.resize(1); - - ret = tsk_mutation_table_set_columns(&p_tables->mutations, - mutations_copy.num_rows, - mutations_copy.site, - mutations_copy.node, - mutations_copy.parent, - mutations_copy.time, - (char *)binary_derived_state.data(), - binary_derived_state_offset.data(), - mutations_copy.metadata, - mutations_copy.metadata_offset); - if (ret < 0) handle_error("derived_from_ascii", ret); - } + // We want to do this efficiently, without making copies of buffers, etc., so we munge around in tskit's + // structs. Maybe there is a way to do this more cleanly and safely? :-O + + assert(p_tables != nullptr); + assert(p_tables->mutations.derived_state != nullptr); + assert(p_tables->mutations.metadata != nullptr); + assert(p_tables->mutations.metadata_offset != nullptr); + + // Swap the derived state column and the metadata column; after this, binary derived states will be in the + // derived state column and ASCII derived states will be in the metadata column. + std::swap(p_tables->mutations.derived_state_length, p_tables->mutations.metadata_length); + std::swap(p_tables->mutations.max_derived_state_length, p_tables->mutations.max_metadata_length); + std::swap(p_tables->mutations.max_derived_state_length_increment, p_tables->mutations.max_metadata_length_increment); + std::swap(p_tables->mutations.derived_state, p_tables->mutations.metadata); + std::swap(p_tables->mutations.derived_state_offset, p_tables->mutations.metadata_offset); + + // Then empty out the metadata column, which we don't use while simulating; we want to free up the buffer + // for the metadata, rather than having the (perhaps large capacity) buffer hanging around forever. The + // only way to do that is by freeing the old buffer and mallocing a new minimal buffer. + free(p_tables->mutations.metadata); + p_tables->mutations.metadata = (char *)malloc(1); // avoids platform-dependencies on zero-length malloc + if (!p_tables->mutations.metadata) + EIDOS_TERMINATION << "ERROR (Species::DerivedStatesFromMetadata): allocation failed; you may need to raise the memory limit for SLiM." << EidosTerminate(nullptr); - tsk_mutation_table_free(&mutations_copy); + p_tables->mutations.metadata_length = 0; + p_tables->mutations.max_metadata_length = 1; + + // Zero out all the metadata offsets; there is an entry for each row, plus one + EIDOS_BZERO(p_tables->mutations.metadata_offset, (p_tables->mutations.num_rows + 1) * sizeof(tsk_size_t)); } -void Species::DerivedStatesToAscii(tsk_table_collection_t *p_tables) +void Species::DerivedStatesToMetadata(tsk_table_collection_t *p_tables) { - // This modifies p_tables in place, replacing the derived_state column of p_tables with an ASCII version. - tsk_mutation_table_t mutations_copy; - int ret = tsk_mutation_table_copy(&p_tables->mutations, &mutations_copy, 0); - if (ret < 0) handle_error("derived_to_ascii", ret); + // This is called when saving out a .trees file. On disk, derived state information is kept in the mutation + // table's metadata column in binary, AND in the derived state column in ASCII (see #664 for discussion). + // Here we convert from our in-memory format by swapping the derived state column's binary data with the + // empty metadata column, and then adding ASCII derived state information. + + // We want to do this efficiently, without making copies of buffers, etc., so we munge around in tskit's + // structs. Maybe there is a way to do this more cleanly and safely? :-O + + assert(p_tables != nullptr); + assert(p_tables->mutations.derived_state != nullptr); + assert(p_tables->mutations.metadata != nullptr); + assert(p_tables->mutations.metadata_offset != nullptr); + + // Generate the ASCII derived state column, for use on disk, and put it into the metadata column. + const char *derived_state = p_tables->mutations.derived_state; + tsk_size_t *derived_state_offset = p_tables->mutations.derived_state_offset; + std::string text_derived_state; + tsk_size_t *text_derived_state_offset = p_tables->mutations.metadata_offset; + + text_derived_state_offset[0] = 0; + for (size_t j = 0; j < p_tables->mutations.num_rows; j++) { - const char *derived_state = p_tables->mutations.derived_state; - tsk_size_t *derived_state_offset = p_tables->mutations.derived_state_offset; - std::string text_derived_state; - std::vector text_derived_state_offset; + slim_mutationid_t *int_derived_state = (slim_mutationid_t *)(derived_state + derived_state_offset[j]); + size_t cur_derived_state_length = (derived_state_offset[j+1] - derived_state_offset[j])/sizeof(slim_mutationid_t); - text_derived_state_offset.emplace_back(0); - - for (size_t j = 0; j < p_tables->mutations.num_rows; j++) + for (size_t i = 0; i < cur_derived_state_length; i++) { - slim_mutationid_t *int_derived_state = (slim_mutationid_t *)(derived_state + derived_state_offset[j]); - size_t cur_derived_state_length = (derived_state_offset[j+1] - derived_state_offset[j])/sizeof(slim_mutationid_t); - - for (size_t i = 0; i < cur_derived_state_length; i++) - { - if (i != 0) text_derived_state.append(","); - text_derived_state.append(std::to_string(int_derived_state[i])); - } - text_derived_state_offset.emplace_back((tsk_size_t)text_derived_state.size()); + if (i != 0) text_derived_state.append(","); + text_derived_state.append(std::to_string(int_derived_state[i])); } - - ret = tsk_mutation_table_set_columns(&p_tables->mutations, - mutations_copy.num_rows, - mutations_copy.site, - mutations_copy.node, - mutations_copy.parent, - mutations_copy.time, - text_derived_state.c_str(), - text_derived_state_offset.data(), - mutations_copy.metadata, - mutations_copy.metadata_offset); - if (ret < 0) handle_error("derived_to_ascii", ret); + text_derived_state_offset[j + 1] = (tsk_size_t)text_derived_state.size(); } - tsk_mutation_table_free(&mutations_copy); + // Copy the ASCII data into a new malloced block and replace any existing metadata with the ASCII. + // FIXME it'd be nice to do the work ourselves into a malloced buffer we own, to avoid the copy. + tsk_size_t metadata_size = text_derived_state.size() * sizeof(char); + char *new_metadata_buffer = (char *)malloc(metadata_size); + if (!new_metadata_buffer) + EIDOS_TERMINATION << "ERROR (Species::DerivedStatesToMetadata): allocation failed; you may need to raise the memory limit for SLiM." << EidosTerminate(nullptr); + + memcpy(new_metadata_buffer, text_derived_state.c_str(), metadata_size); + + free(p_tables->mutations.metadata); + p_tables->mutations.metadata = new_metadata_buffer; + p_tables->mutations.metadata_length = metadata_size; + p_tables->mutations.max_metadata_length = metadata_size; + + // Swap the derived state column and the metadata column; after this, binary derived states will + // be in the metadata column, and ASCII derived states will be in the derived state column. + std::swap(p_tables->mutations.derived_state_length, p_tables->mutations.metadata_length); + std::swap(p_tables->mutations.max_derived_state_length, p_tables->mutations.max_metadata_length); + std::swap(p_tables->mutations.max_derived_state_length_increment, p_tables->mutations.max_metadata_length_increment); + std::swap(p_tables->mutations.derived_state, p_tables->mutations.metadata); + std::swap(p_tables->mutations.derived_state_offset, p_tables->mutations.metadata_offset); } void Species::AddIndividualsToTable(Individual * const *p_individual, size_t p_num_individuals, tsk_table_collection_t *p_tables, INDIVIDUALS_HASH *p_individuals_hash, tsk_flags_t p_flags) @@ -10749,7 +10735,7 @@ void Species::WriteTreeSequence(std::string &p_recording_tree_path, bool p_simpl // Write out the copied tables { // derived state data must be in ASCII (or unicode) on disk, according to tskit policy - DerivedStatesToAscii(&output_tables); + DerivedStatesToMetadata(&output_tables); // In nucleotide-based models, put an ASCII representation of the reference sequence into the tables if (nucleotide_based_) @@ -13862,8 +13848,8 @@ slim_tick_t Species::_InitializePopulationFromTskitBinaryFile(const char *p_file ReadTreeSequenceMetadata(treeSeqInfo, &metadata_tick, &metadata_cycle, &file_model_type, &file_version, mut_metadata_table); - // convert ASCII derived-state data, which is the required format on disk, back to our in-memory binary format - DerivedStatesFromAscii(&treeSeqInfo.tables_); + // shift derived state information from the mutation metadata column to the derived state column + DerivedStatesFromMetadata(&treeSeqInfo.tables_); // in nucleotide-based models, read the ancestral sequence; we do this ourselves, directly from kastore, to avoid having // tskit make a full ASCII copy of the reference sequences from kastore into tables_; see tsk_table_collection_load() above @@ -14022,7 +14008,10 @@ slim_tick_t Species::_InitializePopulationFromTskitDirectory(std::string p_direc EIDOS_TERMINATION << "ERROR (Species::_InitializePopulationFromTskitDirectory): the .trees files for chromosomes have different file versions (" << file_version << " versus " << this_file_version << "). This must be consistent across all files." << EidosTerminate(); } - DerivedStatesFromAscii(&treeSeqInfo.tables_); + // shift derived state information from the mutation metadata column to the derived state column + DerivedStatesFromMetadata(&treeSeqInfo.tables_); + + // in nucleotide-based models, read the ancestral sequence _ReadAncestralSequence(expected_path.c_str(), *chromosome); // The first chromosome uses _InstantiateSLiMObjectsFromTables() and creates the subpopulations, etc., diff --git a/core/species.h b/core/species.h index 42a55d90..3b9c84f9 100644 --- a/core/species.h +++ b/core/species.h @@ -728,8 +728,8 @@ class Species : public EidosDictionaryUnretained void DisconnectCopiedSharedTables(tsk_table_collection_t &p_tables); // zeroes out the shared table copies in p_tables static void handle_error(const std::string &msg, int error) __attribute__((__noreturn__)) __attribute__((cold)) __attribute__((analyzer_noreturn)); - static void DerivedStatesFromAscii(tsk_table_collection_t *p_tables); - static void DerivedStatesToAscii(tsk_table_collection_t *p_tables); + static void DerivedStatesFromMetadata(tsk_table_collection_t *p_tables); + static void DerivedStatesToMetadata(tsk_table_collection_t *p_tables); bool _SubpopulationIDInUse(slim_objectid_t p_subpop_id) const; void RecordTablePosition(void); From fb2c56b2a2204906153652002e5fde4d94d75489 Mon Sep 17 00:00:00 2001 From: Ben Haller Date: Tue, 25 Aug 2026 12:23:06 +0100 Subject: [PATCH 2/4] suggestions from petrelharp code review Co-authored-by: Peter Ralph --- VERSIONS | 1 + core/slim_globals.cpp | 4 ++-- core/species.cpp | 4 +++- 3 files changed, 6 insertions(+), 3 deletions(-) diff --git a/VERSIONS b/VERSIONS index 4e333750..e50c2605 100644 --- a/VERSIONS +++ b/VERSIONS @@ -209,6 +209,7 @@ multitrait branch: the top-level JSON schema now has two keys, baselineOffsetFromUser and baselineOffsetFromSubstitutions, that provide these two components SLiM writes both of them, for use on the Python side; when reading a .trees, SLiM uses baselineOffsetFromUser and ignores baselineOffsetFromSubstitutions (calculating that itself) shift to putting binary derived state info in the mutation table metadata, in addition to the ASCII derived state info in its own column, when on disk (i.e., in .trees files); see #664 + however, we do not rely on the derived state being anything when reading in a .trees file in-memory, the binary derived state info remains in the derived state column and the mutation table metadata is kept empty, as before the schema for the mutation table metadata now encodes an array (of length determined by the metadata itself) of int64_t SLiM mutation IDs diff --git a/core/slim_globals.cpp b/core/slim_globals.cpp index a072ed81..9154908d 100644 --- a/core/slim_globals.cpp +++ b/core/slim_globals.cpp @@ -2188,7 +2188,7 @@ const std::string gSLiM_tsk_mutation_metadata_schema_SOURCE = R"V0G0N({ "additionalProperties": false, "codec": "struct", "type": "object", - "description": "SLiM schema for representing binary derived state data in mutation metadata (the actual mutation metadata is stored in a table in top-level metadata now).", + "description": "SLiM schema for representing binary derived state data in mutation metadata (the metadata for each unique SLiM mutation is stored in top-level metadata).", "examples": [ { "derived_states": [0, 1, 17] @@ -2199,7 +2199,7 @@ const std::string gSLiM_tsk_mutation_metadata_schema_SOURCE = R"V0G0N({ "index": 1, "type": "array", "noLengthEncodingExhaustBuffer": true, - "description": "An array of SLiM mutation IDs (int64t), representing the (stacked) mutations contained by the derived state for the mutation. This information is also saved as a comma-separated ASCII string in the derived state column.", + "description": "An array of SLiM mutation IDs (int64t), representing the (stacked) mutations contained by the derived state for the mutation.", "items": { "binaryFormat": "q", "type": "number" diff --git a/core/species.cpp b/core/species.cpp index 090e111f..e66aca21 100644 --- a/core/species.cpp +++ b/core/species.cpp @@ -8774,7 +8774,9 @@ void Species::CheckAutoSimplification(void) void Species::DerivedStatesFromMetadata(tsk_table_collection_t *p_tables) { // This is called when reading a .trees file. On disk, derived state information is kept in the mutation - // table's metadata column in binary, AND in the derived state column in ASCII (see #664 for discussion). + // table's metadata column in binary. When writing out the .trees file we also + // put this in the derived state column in ASCII (see #664 for discussion), but we + // specifically don't want to assume anything about the derived state column here. // Here we convert back to our in-memory format by swapping the derived state column's ASCII data with the // metadata column's binary data, and then purging the ASCII data from the metadata column. From ff92bb8236a3bdeb5ff21136459d28f381d00f17 Mon Sep 17 00:00:00 2001 From: Ben Haller Date: Tue, 25 Aug 2026 14:21:16 +0100 Subject: [PATCH 3/4] switch to using tsk_mutation_table_takeset_columns() --- core/species.cpp | 169 +++++++++++++++++++++++++++++++---------------- 1 file changed, 113 insertions(+), 56 deletions(-) diff --git a/core/species.cpp b/core/species.cpp index e66aca21..9aec81ae 100644 --- a/core/species.cpp +++ b/core/species.cpp @@ -8774,70 +8774,102 @@ void Species::CheckAutoSimplification(void) void Species::DerivedStatesFromMetadata(tsk_table_collection_t *p_tables) { // This is called when reading a .trees file. On disk, derived state information is kept in the mutation - // table's metadata column in binary. When writing out the .trees file we also - // put this in the derived state column in ASCII (see #664 for discussion), but we - // specifically don't want to assume anything about the derived state column here. - // Here we convert back to our in-memory format by swapping the derived state column's ASCII data with the - // metadata column's binary data, and then purging the ASCII data from the metadata column. - - // We want to do this efficiently, without making copies of buffers, etc., so we munge around in tskit's - // structs. Maybe there is a way to do this more cleanly and safely? :-O + // table's metadata column in binary. When writing out a .trees file we also put this in the derived state + // column in ASCII (see #664 for discussion), but we specifically don't want to assume anything about the + // derived state column here. Here we convert back to our in-memory format by swapping the derived state + // column's ASCII data with the metadata column's binary data, and then purging the ASCII data from the + // metadata column. To do this efficiently, without making copies of buffers, etc., we use a secondary + // mutation table with tsk_mutation_table_takeset_columns(), which is a bit tricky. assert(p_tables != nullptr); - assert(p_tables->mutations.derived_state != nullptr); - assert(p_tables->mutations.metadata != nullptr); - assert(p_tables->mutations.metadata_offset != nullptr); + + tsk_mutation_table_t &mutation_table = p_tables->mutations; + + assert(mutation_table.derived_state != nullptr); + assert(mutation_table.derived_state_offset != nullptr); + assert(mutation_table.metadata != nullptr); + assert(mutation_table.metadata_offset != nullptr); // Swap the derived state column and the metadata column; after this, binary derived states will be in the // derived state column and ASCII derived states will be in the metadata column. - std::swap(p_tables->mutations.derived_state_length, p_tables->mutations.metadata_length); - std::swap(p_tables->mutations.max_derived_state_length, p_tables->mutations.max_metadata_length); - std::swap(p_tables->mutations.max_derived_state_length_increment, p_tables->mutations.max_metadata_length_increment); - std::swap(p_tables->mutations.derived_state, p_tables->mutations.metadata); - std::swap(p_tables->mutations.derived_state_offset, p_tables->mutations.metadata_offset); + std::swap(mutation_table.derived_state_length, mutation_table.metadata_length); + std::swap(mutation_table.max_derived_state_length, mutation_table.max_metadata_length); + std::swap(mutation_table.max_derived_state_length_increment, mutation_table.max_metadata_length_increment); + std::swap(mutation_table.derived_state, mutation_table.metadata); + std::swap(mutation_table.derived_state_offset, mutation_table.metadata_offset); // Then empty out the metadata column, which we don't use while simulating; we want to free up the buffer // for the metadata, rather than having the (perhaps large capacity) buffer hanging around forever. The - // only way to do that is by freeing the old buffer and mallocing a new minimal buffer. - free(p_tables->mutations.metadata); - p_tables->mutations.metadata = (char *)malloc(1); // avoids platform-dependencies on zero-length malloc - if (!p_tables->mutations.metadata) - EIDOS_TERMINATION << "ERROR (Species::DerivedStatesFromMetadata): allocation failed; you may need to raise the memory limit for SLiM." << EidosTerminate(nullptr); - - p_tables->mutations.metadata_length = 0; - p_tables->mutations.max_metadata_length = 1; - - // Zero out all the metadata offsets; there is an entry for each row, plus one - EIDOS_BZERO(p_tables->mutations.metadata_offset, (p_tables->mutations.num_rows + 1) * sizeof(tsk_size_t)); + // best way to do this without using (much) private information is to (a) init a new mutation table, + // (b) use tsk_mutation_table_takeset_columns() to set up that table with our number of rows but with + // empty metadata (setting up whatever internal representation our version of tskit prefers for that), + // (c) swap that empty metadata column in to our mutation table, (d) free the metadata column info from the + // temporary mutation table (which now contains the old ASCII derived states), and then (e) let that mutation + // table go WITHOUT tsk_mutation_table_free(). (It will have pointers to all the same buffers as our own + // mutation table, so we don't want to free it!) + tsk_mutation_table_t temp_mutation_table; + + tsk_mutation_table_init(&temp_mutation_table, 0); + tsk_mutation_table_takeset_columns(&temp_mutation_table, + mutation_table.num_rows, + mutation_table.site, + mutation_table.node, + mutation_table.parent, + mutation_table.time, + mutation_table.derived_state, + mutation_table.derived_state_offset, + /* metadata */ NULL, + /* metadata_offset */ NULL); + + std::swap(temp_mutation_table.metadata_length, mutation_table.metadata_length); + std::swap(temp_mutation_table.max_metadata_length, mutation_table.max_metadata_length); + std::swap(temp_mutation_table.max_metadata_length_increment, mutation_table.max_metadata_length_increment); + std::swap(temp_mutation_table.metadata, mutation_table.metadata); + std::swap(temp_mutation_table.metadata_offset, mutation_table.metadata_offset); + + // WE DO NOT DO: + // tsk_mutation_table_free(&temp_mutation_table); + // INSTEAD WE DO: + tsk_safe_free(temp_mutation_table.metadata); // free the old ASCII derived state data + tsk_safe_free(temp_mutation_table.metadata_offset); // free the offsets for that old data + tsk_safe_free(temp_mutation_table.metadata_schema); // free the schema from tsk_mutation_table_init() + // AND THEN WE JUST LET IT GO OUT OF SCOPE AND DISAPPEAR. + + // Of course the above code will need to be updated if tskit's mutation table implementation changes! } void Species::DerivedStatesToMetadata(tsk_table_collection_t *p_tables) { // This is called when saving out a .trees file. On disk, derived state information is kept in the mutation // table's metadata column in binary, AND in the derived state column in ASCII (see #664 for discussion). - // Here we convert from our in-memory format by swapping the derived state column's binary data with the - // empty metadata column, and then adding ASCII derived state information. - - // We want to do this efficiently, without making copies of buffers, etc., so we munge around in tskit's - // structs. Maybe there is a way to do this more cleanly and safely? :-O + // Here we convert from our in-memory format by generating the ASCII column data, and then using a temporary + // mutation table with tsk_mutation_table_takeset_columns() to get the column in the right format for us. assert(p_tables != nullptr); - assert(p_tables->mutations.derived_state != nullptr); - assert(p_tables->mutations.metadata != nullptr); - assert(p_tables->mutations.metadata_offset != nullptr); - // Generate the ASCII derived state column, for use on disk, and put it into the metadata column. - const char *derived_state = p_tables->mutations.derived_state; - tsk_size_t *derived_state_offset = p_tables->mutations.derived_state_offset; + tsk_mutation_table_t &mutation_table = p_tables->mutations; + + assert(mutation_table.derived_state != nullptr); + assert(mutation_table.derived_state_offset != nullptr); + assert(mutation_table.metadata != nullptr); + assert(mutation_table.metadata_offset != nullptr); + assert(mutation_table.metadata_length == 0); // should be no existing metadata + + // Generate the ASCII derived state column, for use on disk, and put it into a malloced buffer. + // This could be made much faster, and with a lower memory high-water mark, by appending directly into + // a preallocated char buffer sufficiently large to hold all the values. Have a look at this repo: + // https://github.com/ramanawithu/fast_int_to_string/tree/master. But for now that is overkill. + const char *binary_derived_state = mutation_table.derived_state; + tsk_size_t *binary_derived_state_offset = mutation_table.derived_state_offset; std::string text_derived_state; - tsk_size_t *text_derived_state_offset = p_tables->mutations.metadata_offset; + tsk_size_t *text_derived_state_offset = mutation_table.metadata_offset; // note we put the ASCII offset data directly in here text_derived_state_offset[0] = 0; - for (size_t j = 0; j < p_tables->mutations.num_rows; j++) + for (size_t j = 0; j < mutation_table.num_rows; j++) { - slim_mutationid_t *int_derived_state = (slim_mutationid_t *)(derived_state + derived_state_offset[j]); - size_t cur_derived_state_length = (derived_state_offset[j+1] - derived_state_offset[j])/sizeof(slim_mutationid_t); + slim_mutationid_t *int_derived_state = (slim_mutationid_t *)(binary_derived_state + binary_derived_state_offset[j]); + size_t cur_derived_state_length = (binary_derived_state_offset[j+1] - binary_derived_state_offset[j])/sizeof(slim_mutationid_t); for (size_t i = 0; i < cur_derived_state_length; i++) { @@ -8847,8 +8879,6 @@ void Species::DerivedStatesToMetadata(tsk_table_collection_t *p_tables) text_derived_state_offset[j + 1] = (tsk_size_t)text_derived_state.size(); } - // Copy the ASCII data into a new malloced block and replace any existing metadata with the ASCII. - // FIXME it'd be nice to do the work ourselves into a malloced buffer we own, to avoid the copy. tsk_size_t metadata_size = text_derived_state.size() * sizeof(char); char *new_metadata_buffer = (char *)malloc(metadata_size); if (!new_metadata_buffer) @@ -8856,18 +8886,45 @@ void Species::DerivedStatesToMetadata(tsk_table_collection_t *p_tables) memcpy(new_metadata_buffer, text_derived_state.c_str(), metadata_size); - free(p_tables->mutations.metadata); - p_tables->mutations.metadata = new_metadata_buffer; - p_tables->mutations.metadata_length = metadata_size; - p_tables->mutations.max_metadata_length = metadata_size; - - // Swap the derived state column and the metadata column; after this, binary derived states will - // be in the metadata column, and ASCII derived states will be in the derived state column. - std::swap(p_tables->mutations.derived_state_length, p_tables->mutations.metadata_length); - std::swap(p_tables->mutations.max_derived_state_length, p_tables->mutations.max_metadata_length); - std::swap(p_tables->mutations.max_derived_state_length_increment, p_tables->mutations.max_metadata_length_increment); - std::swap(p_tables->mutations.derived_state, p_tables->mutations.metadata); - std::swap(p_tables->mutations.derived_state_offset, p_tables->mutations.metadata_offset); + // Now we want to toss any existing metadata column and put this new metadata column data in its place. The + // best way to do this without using (much) private information is to (a) init a new mutation table, (b) use + // tsk_mutation_table_takeset_columns() to set up that table with all our existing columns but with the new + // metadata column, (c) swap its representation of the new metadata column back into our own table, and + // (d) let the temporary mutation table go WITHOUT tsk_mutation_table_free(). (It will have pointers to + // all the same buffers as our own mutation table, so we don't want to free it!) + tsk_mutation_table_t temp_mutation_table; + + tsk_mutation_table_init(&temp_mutation_table, 0); + tsk_mutation_table_takeset_columns(&temp_mutation_table, + mutation_table.num_rows, + mutation_table.site, + mutation_table.node, + mutation_table.parent, + mutation_table.time, + mutation_table.derived_state, + mutation_table.derived_state_offset, + /* metadata */ new_metadata_buffer, + /* metadata_offset */ text_derived_state_offset); + + std::swap(temp_mutation_table.metadata_length, mutation_table.metadata_length); + std::swap(temp_mutation_table.max_metadata_length, mutation_table.max_metadata_length); + std::swap(temp_mutation_table.max_metadata_length_increment, mutation_table.max_metadata_length_increment); + std::swap(temp_mutation_table.metadata, mutation_table.metadata); + std::swap(temp_mutation_table.metadata_offset, mutation_table.metadata_offset); // identical already, actually + + // WE DO NOT DO: + // tsk_mutation_table_free(&temp_mutation_table); + // INSTEAD WE DO: + tsk_safe_free(temp_mutation_table.metadata_schema); // free the schema from tsk_mutation_table_init() + // AND THEN WE JUST LET IT GO OUT OF SCOPE AND DISAPPEAR. + + // Finally, swap the derived state column and the metadata column; after this, binary derived states + // will be in the metadata column, and ASCII derived states will be in the derived state column. + std::swap(mutation_table.derived_state_length, mutation_table.metadata_length); + std::swap(mutation_table.max_derived_state_length, mutation_table.max_metadata_length); + std::swap(mutation_table.max_derived_state_length_increment, mutation_table.max_metadata_length_increment); + std::swap(mutation_table.derived_state, mutation_table.metadata); + std::swap(mutation_table.derived_state_offset, mutation_table.metadata_offset); } void Species::AddIndividualsToTable(Individual * const *p_individual, size_t p_num_individuals, tsk_table_collection_t *p_tables, INDIVIDUALS_HASH *p_individuals_hash, tsk_flags_t p_flags) From c98f248ba1cd2a8ebe3fba7babedd298fb340693 Mon Sep 17 00:00:00 2001 From: Ben Haller Date: Tue, 25 Aug 2026 21:54:18 +0100 Subject: [PATCH 4/4] add some comments after discussion with Peter --- core/species.cpp | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/core/species.cpp b/core/species.cpp index 9aec81ae..fbc4076f 100644 --- a/core/species.cpp +++ b/core/species.cpp @@ -8781,6 +8781,15 @@ void Species::DerivedStatesFromMetadata(tsk_table_collection_t *p_tables) // metadata column. To do this efficiently, without making copies of buffers, etc., we use a secondary // mutation table with tsk_mutation_table_takeset_columns(), which is a bit tricky. + // This code might seem overly complex, but there are reasons. :-> Munging the tskit data structures more + // directly would not be safe, since tskit's internal implementation is not documented and is subject to + // change, which would risk breakage that might be silent. By using a temporary tsk_mutation_table_t with + // tsk_mutation_table_init() and tsk_mutation_table_takeset_columns(), we can get tskit to set things up for + // us. We still have to munge the internal data structures by swapping columns around with std::swap(), but + // that is relatively safe; it doesn't depend on how tskit is managing the internal state, it is just moving + // that internal state around from place to place. Still unsafe, but _less_ unsafe. To do this cleanly we + // would need new APIs added to tskit. + assert(p_tables != nullptr); tsk_mutation_table_t &mutation_table = p_tables->mutations; @@ -8845,6 +8854,15 @@ void Species::DerivedStatesToMetadata(tsk_table_collection_t *p_tables) // Here we convert from our in-memory format by generating the ASCII column data, and then using a temporary // mutation table with tsk_mutation_table_takeset_columns() to get the column in the right format for us. + // This code might seem overly complex, but there are reasons. :-> Munging the tskit data structures more + // directly would not be safe, since tskit's internal implementation is not documented and is subject to + // change, which would risk breakage that might be silent. By using a temporary tsk_mutation_table_t with + // tsk_mutation_table_init() and tsk_mutation_table_takeset_columns(), we can get tskit to set things up for + // us. We still have to munge the internal data structures by swapping columns around with std::swap(), but + // that is relatively safe; it doesn't depend on how tskit is managing the internal state, it is just moving + // that internal state around from place to place. Still unsafe, but _less_ unsafe. To do this cleanly we + // would need new APIs added to tskit. + assert(p_tables != nullptr); tsk_mutation_table_t &mutation_table = p_tables->mutations;